-1.仅作笔记使用
我所用到的所有脚本可以+Q 840842906索
简单介绍一下方向也可以，说不定以后能互相帮忙～
0.
功能注释这一块大体上就是用自己的序列去和各种数据库进行比对的一个过程
比对算法有Diamond(blast)和HMMe
为什么要用dbCAN2呢，因为他一个软件可以同时跑两种算法，四舍五入少做了一次比对，血赚
1.
dbCAN2的安装与使用参考如下
https://zhuanlan.zhihu.com/p/403142931
https://github.com/linnabrown/run_dbcan
这里主要是下游的处理
我的命令行是：
~~~
run_dbcan.py --db_dir /data/01/user164/workspace_for_whole_object/database_metagenome/dbCAN2 \
--hmm_cov 0.35 --hmm_eval 1e-15 --hmm_cpu 20 \
--dia_eval 1e-102 --dia_cpu 20 \
--out_dir dbcan2 --out_pre genemark prodigal.pro.GreatThan2.fa protein
~~~
2.
会输出 xxdiamond,out 和 xxhmmer.out，即两种算法的输出结果
下面以 diamond的结果为例
3. 先处理Level1
先看一眼输出文件长啥样
第二列 xxxxxx | GH31 xxxxx | GH3 管道符后面的这些字母数字就是Level2，把数字去掉就是Level1
偶尔还能看到管道符后面有纯数字的，这就是level3
~~~
awk '{print $1,$2}' genemarkdiamond.out | uniq | awk '{print $2}' | sed 's/|/ /g' | awk '{print $2}' | sed 's/[|_]/ /g' | awk '{print $1}' | sed 's/[0123456789]//g' > tmp
awk '{print $1}' genemarkdiamond.out > tmp2
paste tmp2 tmp > outfile
~~~
~~~
python CaZyAbun.py out | awk 'BEGIN{ FS=",";OFS="," }{a[$1]+=$2}{b[$1]+=$3}{c[$1]+=$4}{d[$1]+=$5}{e[$1]+=$6}{f[$1]+=$7}{g[$1]+=$8}{h[$1]+=$9}{m[$1]+=$10}{j[$1]+=$11}{k[$1]+=$12}{l[$1]+=$13}END{for(i in a){print i,a[i],b[i],c[i],d[i],e[i],f[i],g[i],h[i],m[i],j[i],k[i],l[i]}}' > CaZyLv1.Abun.csv
~~~
到这可以做堆叠柱状图，热图，下面就很简单了
比如在R里面
a = read.csv("CaZyLv1.Abun.csv",header = T,row.names = 1)
pheatmap(a)
4.Level2
2022.9更新：具体的CaZy比如GH1，可以在CaZy查询
关于这个的处理我只找到了《Seasonal cycling in the gut microbiome of the Hadza hunter-gatherers of Tanzania》
DOI: 10.1126/science.aan4834 (Science封面文章
下面是这篇文章中关于CaZy注释的处理
作者是想通过两大组(旱季、雨季时候人的粪便)样品间，动植物碳水相关的有显著差异，对应旱季、雨季的时候这些人饮食有显著差异
再和前文发现的肠道微生物的季节循环相对应起来
作者通过自己收集文献，把CaZy Level2的注释和具体的物质对应了起来，Table S9从上到下分别是：
植物(细胞壁)碳水、动物碳水、果聚糖/蔗糖、粘多糖(越多表示吃的植物越少,作者原文是这样写的)，就很清晰的分成了吃肉和吃素(或者果子)，然后和前文照应
还注意到图中这几个箱线图的纵坐标是 RPM，这是转录组里面常用的概念
基因A在样本1的RPM就是，样本1中比对到基因A的reads数/样本1比对上的reads数量
我这里用的非冗余基因集当然只能算相对丰度了~
关于如何做还有待更新。。。。
2022.3更新
像上文那样处理是存在问题的，所以写了两个脚本
最开始是处理eggnogmapper的输出结果，但是dbcan2是一模一样的，把前两列提出来就行
用到的脚本
