KOG 注释简明指南

午间，收到朋友求助，问及如何进行序列的 KOG 注释...Emmm，这个可就太久远了。事实上，无论是 COG 还是 KOG 数据库，基本都没怎么更新过。其中 KOG 某种程度属于 COG 的拓展，后者主要用于原核生物，前者用于真核生物。尽管COG在2014年（我刚接触生信的时候）更新了一版，但KOG就没更新过，停留在 20年前。这个数据库主要是基于 7 个主要真核生物构建，其中植物应该就拟南芥... 收到朋友的求助，我还是非常惊讶，怎么这个年代还有人做KOG注释。然而，无论如何，还是看看。
下述为操作记录，比较简单

wget https://ftp.ncbi.nih.gov/pub/COG/KOG/fun.txtwget https://ftp.ncbi.nih.gov/pub/COG/KOG/kogwget https://ftp.ncbi.nih.gov/pub/COG/KOG/kyva# wget https://ftp.ncbi.nih.gov/pub/COG/KOG/kyva=gb# wget https://ftp.ncbi.nih.gov/pub/COG/KOG/lse# wget https://ftp.ncbi.nih.gov/pub/COG/KOG/pa# wget https://ftp.ncbi.nih.gov/pub/COG/KOG/readmewget https://ftp.ncbi.nih.gov/pub/COG/KOG/twog

建库

# makeblastdb -in kyva -input_type fasta -dbtype prot

比对（逻辑上可以用 ghostz 或者 diamond 替代）

# blastx -db kyva -query transcripts.fa -out transcripts.2.kog.blast.tab -evalue 1e-5 -qcov_hsp_perc 33 -outfmt 6 -max_target_seqs 20 -num_threads 40

实在太慢，我用 diamond

diamond makedb --in kyva -d kyva.diamonddiamond blastx --db kyva.diamond --query transcripts.fa --out transcripts.2.kog.diamond.tab --threads 20

过滤

perl -F'\t' -lane 'next if $seen{$F[0]}++;next if $F[-2] 1e-5;print' transcripts.2.kog.diamond.tab transcripts.2.kog.diamond.tab.filtered

perl -lane 'next unless $_;if(/^\[([A-Z]+)\]\s*(.*?)$/){$anno=;@flag=split "",;}else{print join qq{\t},$F[1],,}' kog kog.parsed.tab

表格合并

perl -F'\t' -lane 'if(!$flag){$anno{$F[0]}=q{[}.$F[1].q{]}.$F[2]}else{print join qq{\t},$F[0],qq{$F[1]/$F[-2]/$anno{$F[1]}}}$flag=1 if eof ARGV' kog.parsed.tab transcripts.2.kog.diamond.tab.filtered kog.anno.tab

结果示例

对应的后续可以做一下结果可视化，比如大家都喜欢进行归类信息，做个柱形图云云，具体参考「fun.txt」的分类整理，然后绘图就可以了（PS: 注意，有一些Group同时归属于两个大字母，比如 T 或者 U）

啊，怀念小时候。2014年我刚接触生物信息的时候，这些东西看起来都是高大上，只要会整，那就.... 现在来看...Emmm