快速处理和分析基因家族的蛋白质结构域

yy951010

1. 使用InterProScan
将蛋白质序列提交给InterProScan进行分析，得到每个蛋白质的结构域位置。
2. InterProScan获得的输出格式如下：
gene_id start_position end_position
gene1   10             50
gene1   100            150
gene2   30             90
gene3   50             120
gene3   200            250
3.Python代码处理结构域注释文件
import csv
from Bio import SeqIO

def load_domain_positions(domain_file):
    """
    读取结构域注释文件
    :param domain_file: 结构域位置文件（TSV格式）
    :return: 结构域位置信息的字典 {序列ID: [(起始, 结束), ...]}
    """
    domain_positions = {}
   
    # 读取结构域位置信息文件
    with open(domain_file, "r") as file:
        reader = csv.reader(file, delimiter="\t")
        for row in reader:
            gene_id = row[0]
            start = int(row[1])
            end = int(row[2])
            if gene_id not in domain_positions:
                domain_positions[gene_id] = []
            domain_positions[gene_id].append((start, end))
   
    return domain_positions

def extract_domain_sequences(input_fasta, domain_positions, output_fasta):
    """
    提取蛋白结构域序列
    :param input_fasta: 输入的基因家族氨基酸序列FASTA文件
    :param domain_positions: 结构域的位置列表，格式为 {序列ID: [(起始, 结束), ...]}
    :param output_fasta: 输出的结构域序列FASTA文件
    """
    with open(output_fasta, "w") as output_handle:
        for record in SeqIO.parse(input_fasta, "fasta"):
            if record.id in domain_positions:
                positions = domain_positions[record.id]
                for i, (start, end) in enumerate(positions):
                    # 提取结构域序列
                    domain_seq = record.seq[start-1:end]  # 因为序列索引从0开始，FASTA的氨基酸编号从1开始
                    domain_id = f"{record.id}_domain_{i+1}"  # 给结构域序列命名
                    new_record = record[:0]  # 创建一个新的空记录
                    new_record.id = domain_id
                    new_record.description = f"Domain from {start} to {end}"
                    new_record.seq = domain_seq
                    
                    # 将结构域序列写入输出文件
                    SeqIO.write(new_record, output_handle, "fasta")
                    
if __name__ == "__main__":
    # 输入的基因家族FASTA文件
    input_fasta = "gene_family_sequences.fasta"
   
    # 结构域位置文件（由Pfam或InterProScan生成）
    domain_file = "domain_positions.tsv"
   
    # 输出的结构域序列FASTA文件
    output_fasta = "gene_family_domain_sequences.fasta"
   
    # 加载结构域位置信息
    domain_positions = load_domain_positions(domain_file)
   
    # 提取结构域序列
    extract_domain_sequences(input_fasta, domain_positions, output_fasta)
   
    print("结构域序列已提取并导出。")
代码说明

• load_domain_positions: 读取Pfam或InterProScan生成的结构域注释文件，并将其保存为字典格式，字典的键为基因ID，值为对应的结构域位置列表。
• extract_domain_sequences: 读取输入的FASTA文件，按照结构域的位置提取序列，并生成一个新的FASTA文件。
• domain_file: 是一个结构域注释文件（例如TSV格式），你可以从Pfam或InterProScan结果中获取这些位置信息并保存为文件。