新聞中心
這篇文章主要介紹了fasta序列如何按指定格式輸出,具有一定借鑒價值,感興趣的朋友可以參考下,希望大家閱讀完這篇文章之后大有收獲,下面讓小編帶著大家一起了解一下。
十多年專注成都網(wǎng)站制作,成都企業(yè)網(wǎng)站建設,個人網(wǎng)站制作服務,為大家分享網(wǎng)站制作知識、方案,網(wǎng)站設計流程、步驟,成功服務上千家企業(yè)。為您提供網(wǎng)站建設,網(wǎng)站制作,網(wǎng)頁設計及定制高端網(wǎng)站建設服務,專注于成都企業(yè)網(wǎng)站建設,高端網(wǎng)頁制作,對混凝土攪拌罐車等多個方面,擁有豐富的網(wǎng)站維護經(jīng)驗。
前言:有時在處理fasta文件時,我們需要序列按照規(guī)定的格式排列。
很多人應該遇到過需要將序列排列到一行上,或者每行按照規(guī)定的bp數(shù)顯示。我也經(jīng)常遇到像60bp,70bp的不等長fasta序列共存于同一個fasta文件中的情況,為了避免不同長度對后面的處理造成影響,一般最好將格式統(tǒng)一。
fasta file format:
雖然是個小問題,但是卻有很多不同的方法來實現(xiàn)這些操作,那接下來還是以舉例說明,講解一些方法來實現(xiàn)上面講到的兩種格式排列。
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(換行?)數(shù)目:
我寫了一個biopython版本的,可以用它指定的參數(shù)nwrap完成上面的兩種操作,設置nwrap為0時即顯示到一行上。
wrap_xbp.py:
import argparse
from Bio import SeqIO
from Bio.SeqIO.FastaIO import FastaWriter
###usage description
describe=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 output
output_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_test1
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT
TTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATC
GCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCC
AACTAACC
>chr_test2
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT
TTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATC
GCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCC
AACTAACC
4、另外bbmap也有很快捷易用的reformat.sh來進行相同的操作。
按照50bp每行排列:
$ ~/tool/bbmap/reformat.sh in=test.fa out=test_out2.fa fastawrap=50
結果文件按照50bp每行排列:
$ cat test_out2.fa
>chr_test1
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT
TTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATC
GCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCC
AACTAACC
>chr_test2
GATCACAGGTCTATCACCCTATTAACCACTCACGGGAGCTCTCCATGCAT
TTGGTATTGCAGTATCTGTCTTTGATTCCTGCCTCATCCTATTATTTATC
GCACCTACGTTCAATATTTTTTAACAGTCACCCCCCAACTAACACATTCC
AACTAACC
當然也可以規(guī)劃到一行顯示,只需要設置大一些即可:fastawrap=50000000000
感謝你能夠認真閱讀完這篇文章,希望小編分享的“fasta序列如何按指定格式輸出”這篇文章對大家有幫助,同時也希望大家多多支持創(chuàng)新互聯(lián),關注創(chuàng)新互聯(lián)行業(yè)資訊頻道,更多相關知識等著你來學習!
新聞名稱:fasta序列如何按指定格式輸出
分享路徑:http://fisionsoft.com.cn/article/jihioh.html