问题
在转录组或基因组注释文件中有基因的GO注释，为了做GO富集分析，需要从中提取每个有GO注释的基因对应的GO号。
Interpro_output.gff3文件内容如下：
$ head CIA_Interpro_output.gff3
##gff-version 3
##feature-ontology http://song.cvs.sourceforge.net/viewvc/song/ontology/sofa.obo?revision=1.269
##interproscan-version 5.35-74.0
##sequence-region CIA_ORF9411 1 466
CIA_ORF9411 . polypeptide 1 466 . + . ID=CIA_ORF9411;md5=2f764379a092e5b22282e2c0f506f89f
CIA_ORF9411 Pfam protein_match 41 133 1.1E-20 + . date=17-11-2019;Target=CIA_ORF9411 41 133;ID=match$1_41_133;signature_desc=Domain associated at C-terminal with AAA;Name=PF14363;status=T;Dbxref="InterPro:IPR025753","KEGG:00230+3.6.1.3"
CIA_ORF9411 Gene3D protein_match 211 401 2.7E-44 + . date=17-11-2019;Target=CIA_ORF9411 211 401;ID=match$2_211_401;Name=G3DSA:3.40.50.300;status=T
CIA_ORF9411 Pfam protein_match 262 400 7.7E-17 + . date=17-11-2019;Target=CIA_ORF9411 262 400;Ontology_term="GO:0005524";ID=match$3_262_400;signature_desc=ATPase family associated with various cellular activities (AAA);Name=PF00004;status=T;Dbxref="InterPro:IPR003959"
CIA_ORF9411 SUPERFAMILY protein_match 164 444 3.18E-49 + . date=17-11-2019;Target=CIA_ORF9411 164 444;ID=match$4_164_444;Name=SSF52540;status=T;Dbxref="InterPro:IPR027417"
CIA_ORF9411 ProSitePatterns protein_match 371 389 . + . date=17-11-2019;Target=CIA_ORF9411 371 389;Ontology_term="GO:0005524";ID=match$5_371_389;signature_desc=AAA-protein family signature.;Name=PS00674;status=T;Dbxref="InterPro:IPR003960"
想要获得这种格式的GO注释信息
$ head CIA_GO1.txt
CIA_ORF100395 GO:0016747
CIA_ORF100980 GO:0015267 GO:0016020 GO:0055085
CIA_ORF102278 GO:0005200 GO:0005525 GO:0005874 GO:0007017
CIA_ORF105663 GO:0003824 GO:0046872
CIA_ORF105663 GO:0004222 GO:0006508
CIA_ORF12198 GO:0005524 GO:0009058
CIA_ORF12200 GO:0005524 GO:0009058
CIA_ORF12201 GO:0005524 GO:0009058
CIA_ORF12202 GO:0005524 GO:0009058
CIA_ORF12203 GO:0005524 GO:0009058
或者这种格式的：
$ head CIA_GO2.txt
CIA_ORF100395 GO:0016747
CIA_ORF100980 GO:0015267
CIA_ORF100980 GO:0016020
CIA_ORF100980 GO:0055085
CIA_ORF105663 GO:0003824
CIA_ORF105663 GO:0004222
CIA_ORF105663 GO:0006508
CIA_ORF105663 GO:0046872
CIA_ORF12198 GO:0005524
CIA_ORF12198 GO:0009058
python代码实现
import re
import os
f = open('CIA_Interpro_output.gff3', 'r')
result1 = open('CIA_GO1.txt', 'w')
result2 = open('CIA_GO2.txt', 'w')
lines = f.readlines()
for line in lines:
result = re.findall('GO:\d+', line)
if result:
print(line.split('\t')[0], "\t".join(result), file=result1)
for GO in result:
print(line.split('\t')[0], GO, file=result2)
#排序去重
os.system('cat CIA_GO1.txt | sort | uniq > tmp ; mv tmp CIA_GO1.txt')
os.system('cat CIA_GO2.txt | sort | uniq > tmp ; mv tmp CIA_GO2.txt')
Interpro_output.gff3文件下载
