太久没有更新了，最近忙着组装好几个基因组，后面我也会吧基因组组装的一些东西分享出来哦~
那就步入正题吧，我们在前面准备好了数据之后，就可以开始鉴定所需要的基因家族成员了，这样才能开始后续的分析
先从，最最最最基础的linux的一些使用开始吧，要是不知道在自己电脑上怎么安装linux系统可以参考下面网址哟~
https://www.jianshu.com/p/821e1b409e95?v=1689247563616
1、最最最最基础的就是需要在linux中安装软件了，就像买了一个新手机一样；安装软件当然需要一个软件商店喽
就需要首先安装一个linux中的软件商店，叫做conda，有Anaconda和Miniconda两种，区别的话.......感兴趣的可以了解一下，自己电脑上推荐Miniconda，不占内存，哈哈哈哈
#首先用mkdir software && cd software 创建一个文件夹
#安装Anaconda
#下载
wget -c http://mirrors.aliyun.com/anaconda/archive/Anaconda3-5.3.1-Linux-x86_64.sh?spm=a2c6h.25603864.0.0.49f16323pglbMi -P ~/software
#运行安装包（一路yes）
bash ~/software/Anaconda3-5.3.1-Linux-x86_64.sh
#启动conda程序
source ~/.bashrc
#安装Miniaconda
##Download Miniaconda3
wget https://repo.anaconda.com/miniconda/Miniconda3-py310_23.1.0-1-Linux-x86_64.sh -P ~/software
bash ~/software/Anaconda3-5.3.1-Linux-x86_64.sh ##一路yes
source ~/.bashrc
#生成.condarc文件
conda config --set show_channel_urls yes
#添加conda镜像源
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/free/
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/main/
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/conda-forge
conda config --add channels https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/bioconda
conda config --set show_channel_urls yes
2、安装完conda就需要借助conda安装基因家族分析需要的软件（工具）
#基因家族成员鉴定的软件，也就几个，三下五除二就能搞定！！！
conda install seqtk -y
conda install seqkit -y
conda install blast -y
conda install hmmer -y
没啦，是不是几分钟！
3、那我们就开始正式鉴定家族成员，继续以CYP超家族为例吧
#你们是不是进ubuntn系统里面去找不到自己电脑上的文件呀？？很简单，一个命令!!
cd /mnt/ #就可以看到你的c盘d盘f盘啥的
############################################################
###废话少说，解释一下下面那些带着$的符号。
#对氨基酸序列进行结构域搜索
hmmsearch --cut_tc --domtblout CYP.domtblout -o CYP.hmmout CYP.hmm genome.pep.fa
#基于E-value值进行过滤
awk '$7<0.01 && $1 !~ /^#/ {print $0 }' CYP.domtblout > CYP.domtblout.filter
#提取预测的目标基因ID
awk '{print $1}' CYP.domtblout.filter | sort -u > CYP.domtblout.filter.id
#构建blast比对数据库
makeblastdb -in Ath.CYP.pep.fa -dbtype prot
#blastp比对
blastp -query $fasta -db Ath.CYP.pep.fa -evalue 1e-10 -outfmt '6 std qlen slen ' -out CYP.blastout
#blast结果基于identity30%过滤
awk '$3 > 30 && $4> 155 {print $1}' CYP.blastout | sort -u > CYP.blastout.filter.id
###具体情况具体分析，二选一即可
#合并hmm和blast结果,也就是取交集（避免缺漏）
cat CYP.domtblout.filter.id CYP.blastout.filter.id | sort -u > CYP.geneID
##当然也可以取交集（更加精确）
cat CYP.domtblout.filter.id CYP.blastout.filter.id | sort | uniq -C | awk '$1 == 2 {print $2}' > CYP.geneID
#提取目标基因的pep和cds序列
seqtk subseq genome.pep.fa CYP.geneID >CYP.pep.fasta
seqtk subseq genome.cds.fa CYP.geneID > CYP.cds.fasta
现在就得到家族基因的基因ID，蛋白序列和cds序列了！很简单吧！
还有其他复杂一点的家族基因，比如多结构域的、结构域嵌合的等，下次再说吧......
我也要去搬砖去了，呜呜呜~~~
网友评论