什么是Genbank文件？如何利用它？
当前生物信息学分析中对Genbank文件的使用好像并不是很多。（至少我的研究领域中很少有用到Genbank文件）但是这并不能说明Genbank文件没有用处。我第一次用到Genbank文件是在2017年6月我大三的时候参加南开大学生科院的夏令营，生科院生物信息学的高山老师让我用R语言写一个程序：从人类线粒体基因组中提取所有的gene序列和CDS序列，并且将每个序列单独用一条FASTA序列表示，最后写入一个文本文件中。那个时候我并没有用R语言做，当时我的同学钟权（现在在华中农业大学信息学院读研）和我一起想这个程序该怎么写，最后用biopython解决了。大致思路是将Genbank文件中的注释信息打断成每个class，每个class放入字典中，最后根据class信息从Genbank的序列部分提取序列。大概代码如下
不过高老师的要求是用R语言写，所以不久之后我想到了用R语言结合GFF3文件和fasta文件提取全部gene和CDS，因为GFF3文件每行是一个Gene或者CDS的信息，包括了在基因组上的起始位点，结束位点，在正链还是在负链，是否有功能注释。但是需要用awk命令先把GFF3要提取的序列位置信息提出。具体代码如下：
以上则是我分别利用Genbank，FASTA，GFF3文件提取基因序列及信息的代码。后来的后来我才知道R语言有一个叫做genbankr的bioconductor包，仅仅使用简单的几个函数就可以提取基因序列和信息。
不过话说回来，什么是Genbank文件？Genbank文件又叫GBFF文件，这个文件储存了基因的基因信息，特征注释信息和序列信息。
如何制作Genbank文件？
目前，我所知道的除了提取序列信息和序列能够用到Genbank文件，其他时候比如使用Mavue作比较基因组学分析的时候也能用到Genbank文件。加上有一次病毒所其他实验室的一个师兄问我怎么生成自己的gbk文件可能是他有其他用处吧。在Mauve分析中，因为比较基因组学在进行基因重排、共线性关系分析的时候需要知道序列上哪些位置是CDS或者基因序列，这样才有助于使用软件的算法进行重排分析、共线性分析等。当然，目前NCBI上的GenBank文件是有注释不全的。对于病毒基因组来说，例如噬菌体、或者其他研究并不太热门的病毒就没有完整的GenBank文件。因为缺少完整的注释，所以如果要使用这些缺少注释的文件进行分析那么结果肯定不好。为了得到好的结果就需要自己制作相对注释完整的Genbank文件。由于这种文件的第一部分基本信息只是包含提交文件的人、时间、地点和版本号等信息，这些基本信息是可要可不要的。我们需要的只是特征注释信息和基因序列信息。而特征注释信息可以用Prodigal软件预测基因得到。对于序列信息，我们需要将自己的FASTA序列制作为像前一部分序列信息一样的格式，即一行60个核苷酸序列，每10个核苷酸序列之间用分隔符分隔，每行第一列是每行第一个核苷酸在整个序列上的位置坐标。大体思路是将FASTA序列的核苷酸部分单独提取，合并为一行，然后进行字符串分割，对分割的字符串每60个核苷酸为一行，每行的每10个核苷酸用tab分隔，在最前面写上第一个核苷酸的坐标信息。
这部分我们需要结合perl语言和Shell语言制作。对于使用Prodigal输出注释信息部分的文件代码如下
对于处理序列信息的Perl语言代码和Shell代码如下
以上部分是直接在Shell中调用perl脚本，而perl脚本的代码如下
写在最后
持续工作的我和Соня不知不觉就到除夕了，不过工作中的我俩是最快乐的。这部分推文应该是我创立公众号以来写的涵盖知识点最广泛的推文了，因为这部分推文将我的R语言、Python语言、Perl语言和Linux的Shell编程都用上了，涉及的多语言编程广、深度相对以前较深。不过，这也是对我自己的鞭策，还要感谢我的Соня对我无私的帮助与支持！祝我的Соня新年快乐，猪年大吉，身体健康、多发高水平论文。早日博士毕业！！！！以后我还要继续提高我的多语言、多环境编程分析生物信息学数据的能力！
《Genbank文件的利用与制作》来自互联网公开内容，收录仅供学习使用，如侵权请联系删除。本文URL：https://www.ezixuan.com/1020641.html