Rosalind工具库：使用EMBOSS处理生物数据

序列全局联配：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)。开放阅读框指的是以起始密码子开始，到终止密码子结束的一段序列

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