多序列比对,进化树分析,motif预测,密码子偏好性分析(2)
- 前言:假如类鼻疽杆菌中发现有一个感兴趣的gene,想查看此gene在各类鼻疽菌株甚至原核生物和真核生物中是否也表达此gene。若表达,各菌株之间序列差异大不大,进而看其保守性,密码子的使用是否有所偏好。
以BopA蛋白为例
第一部分:下载所需物种的BopA蛋白序列(多物种)
1下载BopA蛋白序列
1.1 第一种方式
打开https://www.ncbi.nlm.nih.gov/protein,输入BopA,search
然后直接复制,或Run BLAST,或send to
1.2 下载BopA蛋白序列:第二种方式
- 1 打开https://www.ncbi.nlm.nih.gov/gene/
-
2 输入BopA,search,结果中选中第一个k96243标准株
-
3 点击下图中的GeneBank
-
4 点击YP_111530.1或直接复制序列
-
5 点击进入,Run BLAST。或复制序列。
2 BLAST(nr protein sequeces database)
2.1 条件设置如下
2.2 结果如下
一般情况下,按下面这个步骤1234导出即可。
但是现在这个情况,按上面这样是不可能的,因为完全看不出物种信息,做进化树没什么意义。
所以需要稍微麻烦一些,具体如下:
点击image9中的[Taxonomy reports]
上图中的Accesion即为蛋白的登录号,左边中括号内为物种,可见非常多重复,继续往下看,可见
- 可以看到物种及对应的accession ID,我们的目的是得到有物种名和蛋白序列的fasta文件。
对上图来说,只需要得到ACCESSION ID,再搜索一次protein database即可 - 需要注意的是,只需要得到每一个物种的第一个Accession ID
- 这些ID,应该可以用python爬取相应数据。但是假如你需要的物种少,可以直接复制这些ID。
- 假如需要几十甚至100多个物种,逐个复制粘贴会很麻烦,这里我采取excel来做,大概需要几分钟时间。 具体步骤为,excel导入数据-Organism Report部分-替换部分数据-分列-去重复。最终得到如图image14(部分)。
2.3 提取上述Accession ID的protein fasta序列
注意,NCBI一次提取不能超过100个accession ID,上述138,所以分两次提取
或者
最终得到下列文件
下面可以进行进化树制作了,但是之前需要进行一下改变,那就是fasta文件说明部分只保留物种名即可,也就是说,只保留[ ]里面的内容。内容少的话可以手动编辑,我这里100多,所以用一个脚本进行。
$ awk -F '[][]' '{if ($0~/>/) {print ">"$2} else {print $0}}' BopAsequence.fasta > BopAoutput.fasta
$ less BopAoutput.fasta
结果如下,这样就可以进行多序列比对及进化树分析了
第二部分 多序列比对
MEGAX软件下载地址,安装
找到上面文件,打开
默认参数,OK,多序列比对完成。下面site输入氨基酸位置可以选择查看自己想看的序列,保存文件。
第三部分进化树分析
默认参数,关于gap和extention罚分方法,可以自行google。
物种太多,我选择环形展示。明显看出thailandensis,mallei,peudomallei在其中的分布,各自成一cluster。也就说BopA具有进化意义。
备注:
- 1 通过blast可以看到,BopA在真核生物中并不存在