专栏首页生物信息学fasta序列按指定格式输出

fasta序列按指定格式输出

前言:有时在处理fasta文件时,我们需要序列按照规定的格式排列。

很多人应该遇到过需要将序列排列到一行上,或者每行按照规定的bp数显示。我也经常遇到像60bp,70bp的不等长fasta序列共存于同一个fasta文件中的情况,为了避免不同长度对后面的处理造成影响,一般最好将格式统一。

fasta file format:

虽然是个小问题,但是却有很多不同的方法来实现这些操作,那接下来还是以举例说明,讲解一些方法来实现上面讲到的两种格式排列。

1、这里我使用全长158bp,60bp每行显示,最后一行38bp排列的两条fasta序列组成的fasta文件来举例。

test.fa:

$ cat test.fa          >chr_test1GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC>chr_test2GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC

2、首先是使用awk排列到一行:

$ awk '/^>/ { if(NR>1) print "";  printf("%s\n",$0); next; } { printf("%s",$0);}  END {printf("\n");}' test.fa >chr_test1GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC>chr_test2GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC

3、另外biopython处理fasta、fastq文件也很方便,也有相应的解决办法。

biopython中默认是按照60bp每行输出的,如果去查查它的帮助文档,可以查到FastaWriter可以在写出文件中指定fasta序列的wrap(换行?)数目:

我写了一个biopython版本的,可以用它指定的参数nwrap完成上面的两种操作,设置nwrap为0时即显示到一行上。

wrap_xbp.py:

import argparsefrom Bio import SeqIOfrom Bio.SeqIO.FastaIO import FastaWriter
###usage descriptiondescribe=argparse.ArgumentParser(description="Make Fasta Sequence in a Single Line or Wrap N bp One Line")describe.add_argument("-nwrap",help="n base per line;default=0 means seq in one line",default=0,type=int)describe.add_argument("orgf",help="Original fasta")#原始fasta文件describe.add_argument("optf",help="Output fasta")#修改格式后的输出文件args=describe.parse_args()
###handle to output and FastaWriter to make normalized outputoutput_fasta = open(args.optf,"w")#打开文件句柄用于写出文件writer = FastaWriter(output_fasta,wrap=args.nwrap)#设置写出格式writer.write_file(SeqIO.parse(args.orgf,"fasta"))#读取原始文件并按照要求格式写出output_fasta.close()#关闭文件句柄

运行得到50bp每行的输出文件test_50wrap.fa

$ python3 wrap_xbp.py -nwrap 50 test.fa test_50wrap.fa
$ cat test_50wrap.fa>chr_test1GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC>chr_test2GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC

4、另外bbmap也有很快捷易用的reformat.sh来进行相同的操作。

按照50bp每行排列:

$ ~/tool/bbmap/reformat.sh in=test.fa out=test_out2.fa fastawrap=50

结果文件按照50bp每行排列:

$ cat test_out2.fa                            >chr_test1GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC>chr_test2GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCATTTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATCGCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCCAACTAACC

当然也可以规划到一行显示,只需要设置大一些即可:fastawrap=50000000000

本文分享自微信公众号 - 生物信息学(swxxx1),作者:saber

原文出处及转载信息见文内详细说明,如有侵权,请联系 yunjia_community@tencent.com 删除。

原始发表时间:2019-12-24

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • linux-两个文件求交集、并集、差集

    (用sort将a.txt b.txt文件进行排序,uniq使得两个文件中的行唯一,使用-d输出两个文件中次数大于1的内容,即是得到交集)

    阿凡亮
  • 植物miRNA的鉴定原理

    由上面的介绍我们知道miRNA是由可以折叠成茎环结构的的RNA产生的,由此我们可以想到只要预测基因组序列中可以折叠成茎环结构的区域就可以对miRNA进行预测了,...

    阿凡亮
  • 在线功能富集分析与qPCR引物设计

    之前写的RNA-seq数据差异表达分析一文中提到,筛选得到差异表达基因list后,需要进一步分析这些基因参与了哪些功能,因此要进行后续的一些分析,比如功能富集分...

    阿凡亮
  • Android.mk文件中添加第三方jar文件的方法

    下面给大家介绍Android.mk文件中添加第三方jar文件的方法,具体内容详情如下所示:

    砸漏
  • eclipse中将java文件打成jar包

    软件开发的最后一步就是软件的打包与发布,这也是很重要的一步。这几天在Eclipse中做了一个小软件,准备将其打成可运行的jar包进行发布,在网上搜了好多关于在E...

    闵开慧
  • Java工程中添加依赖jar包不起作用问题总结

    Java工程中添加依赖jar包不起作用问题总结 此次总结两种方式的依赖问题 1 在Eclipse中添加依赖jar包不起作用问题     这种方式可能是Eclip...

    闵开慧
  • 斯坦福大学NLP组Python深度学习自然语言处理工具Stanza试用

    众所周知,斯坦福大学自然语言处理组出品了一系列NLP工具包,但是大多数都是用Java写得,对于Python用户不是很友好。几年前我曾基于斯坦福Java工具包和N...

    汤贤
  • Maven 优势

    如果项目非常庞大,就不适合使用 package 来划分模块,最好是每个模块对应一个工程,利于分工协作,而借助于 maven 就可以将项目拆分成多个工程

    happyJared
  • 1、在eclipse中导入Java的jar包方法---JDBC【图文说明】

    Eclipse环境下jar包导入 在Eclipse环境下编写Java程序,常常会借用到各种jar包。如:连接数据库时,导入jar包是必须的。导入方法如下: 1....

    YGingko
  • SpringBoot Jar包瘦身 - 跟大文件说再见!

    SpringBoot部署起来配置非常少,如果服务器部署在公司内网,上传速度还行,但是如果部署在公网(阿里云等云服务器上),部署起来实在头疼、就是 编译出来的 J...

    yunlgonn

扫码关注云+社区

领取腾讯云代金券