根据单独文件中的条目从FASTA文件中提取序列

Question

我有两个文件.

文件1:带有基因序列的FASTA文件,如下例所示:

>PITG_00002 | Phytophthora infestans T30-4 conserved hypothetical protein (426 nt)
ATGCATCGCTCGGGTTCCGCACGGAAAGCCCAAGGTCTGGGATTACGGGGTGGTGGTCGG
TTACACTTGGAATAACCTCGCAAATTCAGAATCTCTACAGGCTACGTTCGCGGATGGAAC
>PITG_00003 | Phytophthora infestans T30-4 protein kinase (297 nt)
ATGACGGCTGGGGTCGGTACGCCCTACTGGATCGCACCGGAGATTCTTGAAGGCAAACGG
TACACTGAGCAAGCGGATATTTACTCGTTCGGAGTGGTTTTATCCGAGCTGGACACGTGC
AAGATGCCGTTCTCTGACGTCGTTACGGCAGAGGGAAAGAAACCCAAACCAGTTCAGATC
>PITG_00004 | Phytophthora infestans T30-4 protein kinase, putative (1969 nt)
ATGCGCGTGTCTGGTCTCCTTTCAATTCTTGCAGCCACTTTGACCACGGCCCAAGACTAC

文件2:一个简单的文本文件,其中包含基因的登录标识.像这样.

PITG_00003
PITG_00005
PITG_00023

文件2中的每个条目都在文件1中的某个位置,但不是文件1中的每个条目都在文件2中.我需要从文件1中删除不在文件2中的所有条目.我觉得在biopython中必须有一些东西可以帮助我的模块,我只是不知道是什么.例如,我原本以为我可以使用该SeqIO.parse函数从我的FASTA文件中提取加入,但这真的让我有两个加入号文件.我不知道如何有选择地提取其他文件中的加入.也许喜欢将文件2中的所有条目读入字典,然后将该条目与文件1中的匹配条目相关联并用于SeqIO.parse提取整个序列...但我真的不知道....任何人都可以给予任何帮助我非常感谢!

Answer 1

试试这个:

f2 = open('accessionids.txt','r')
f1 = open('fasta.txt','r')
f3 = open('fasta_parsed.txt','w')

AI_DICT = {}
for line in f2:
    AI_DICT[line[:-1]] = 1

skip = 0
for line in f1:
    if line[0] == '>':
        _splitline = line.split('|')
        accessorIDWithArrow = _splitline[0]
        accessorID = accessorIDWithArrow[1:-1]
        # print accessorID
        if accessorID in AI_DICT:
            f3.write(line)
            skip = 0
        else:
            skip = 1
    else:
        if not skip:
            f3.write(line)

f1.close()
f2.close()
f3.close()

简要解释这里发生了什么... accessionids.txt是你的文件2,而fasta.txt你的文件1.显然,您需要将这些文件名替换为代码中的实际文件名.

首先,我们创建一个字典(有时被称为散列或关联数组)和为每个登录ID在文件2我们创建了一个条目,其中关键是保藏ID和值被设置为1(不该值真正重要在这种情况下).

接下来我们查看文件1并再次查看该文件中的每一行.如果文件中的行开头,>那么我们知道它包含一个入藏ID.我们接受该行并将其拆分,|因为具有入藏ID的每一行都将|在字符串中包含a .接下来,按照指定的方式获取拆分的第一部分_splitline[0].我们使用accessorIDWithArrow[1:-1]砍掉字符串中的第一个和最后一个字符,这是>在前面的符号和后面的空格.

此时,accessorID现在包含我们期望从文件2格式的Accession ID .

接下来,我们检查我们之前创建并填充的字典是否将此Accession ID定义为键.如果是,我们立即将具有登录ID的行写入新文件fasta_parsed.txt,并将skip'flag'变量设置/重置为0.然后,else包含该if not skip段的语句将允许与我们发现打印到该fasta_parsed.txt文件的登录ID相关联的后续行.

对于保藏ID从文件1不在字典中(未在发现文件2),我们没有行写入fasta_parsed.txt和我们设置skip标志为0.因此,直到另一个登录ID是在发现文件1中存在的文件2,将跳过所有后续行.