序列全局联配:Pairwise Global Alignment
序列相似性意味着遗传的同源性是遗传和进化研究上的一个假说。
给定不完全相同的两个序列,如果允许错配和开口(gap)的话,会有无数多的联配(alignment)方式,为了找到最有可能的联配方法,我们需要制定一套评价标准,这就是打分机制,打分机制分为三个部分:
- 得分矩阵:用来计算正确配对(match)和错误配对(mismatch)的得分
- 开口惩罚(gap opening penalty): 如果联配中有开口被打开,那么要会对起做相应的惩罚
- 开口延伸惩罚(gap extension penaly): 开口打开后在原先的基础上进行延伸,也会被惩罚,并且惩罚力度比较小。这是基于如果碱基缺失是独立时间,如果有两个碱基缺失,它更有可能是2个碱基同时缺失,而不是两个独立时间,也就是先缺失一个,然后在缺失一个
GARFIELDTHEVERYFA-TCAT
:::|||||||| || ||||
WINFIELDTHE----FASTCAT
如上是一个联配方式,如果的:
表示错配,得分为-1, |
为正确联配得分为1, -
gap打开为-5,延伸为-1。那么得分就是14-3-(2*5+3) = -2. 我们可以自己写程序做这样的计算,但是目前来说就用别人写好的工具进行,在Rosalind兵工厂系列结束后再说。
问题:给定两个GenBank ID,请基于Needl网页工具https://www.ebi.ac.uk/Tools/psa/emboss_needle/nucleotide.html, 对序列进行联配。要求
- 使用DNAfull得分矩阵, 参考https://en.wikipedia.org/wiki/Nucleic_acid_notation
- 开口是10分
- 延伸是1分
如果用的是网页工具,则是先去Genbank查找两个序列,然后网页工具不要用默认选项,且EN GAP PENALTY一定设置为True。
如果用的EMBOSS命令行,保存两个序列,设置gap的得分
~/opt/biosoft/emboss-6.6.0/bin/needle a.fasta b.fasta \
-gapopen 10.0 -gapextend 1.0 \
-endweight true -endopen 10 -endextend 1.0
寻找开放阅读框: Finding Genes with ORFs
在DNA序列上寻找基因的第一步是找到序列上的开放阅读框(ORF)。开放阅读框指的是以起始密码子开始,到终止密码子结束的一段序列
对于一条DNA序列,有六种方式来寻找ORF,先对当前DNA链从第一个碱基、第二个碱基、第三个碱基开始翻译密码子,然后再从反向互补链里也从第一个、第二个、第三个碱基进行翻译。对于真核生物而言,这个方法可能太过简单了,因为真核生物还有内含子。对于原核生物,由于基因组比较简单,因此可以试一试。
问题:给你一条DNA序列,找到最长的ORF。
思路: EMBOSS的getorf
能以6种方式寻找给定序列的ORF,注意两个参数:
-
-table
: 翻译密码表,0表示标准翻译密码表 -
-find
: 结果输出哪些氨基酸,0表示从到终止密码子的所有氨基酸,1表示起始密码子到终止密码子中间的氨基酸,2表示到终止密码子的所有碱基,3示起始密码子到终止密码子中间的碱基,4 包含起始密码子前几个碱基,5包含第一个终止密码子后几个碱基,65包含最后一个终止密码子后几个碱基
cat rosalind_orfr.txt | ~/opt/biosoft/emboss-6.6.0/bin/getorf -filter -find 1 \
| bioawk -c fastx '{print length($seq) "\t" $seq}' \
| sort -nrk1,1