【Linux】fasta标准格式转换

接到一个非标准格式保存的fasta序列格式的文本,看起来像是直接从NCBI copy下来的。内含多条序列,后续需要对序列进行对位排列,而手工改不现实。有不想写perl脚本,故寻求awk这一宝藏命令。

eg. 序列input.fa

>AF371111
        1 atggatgccg acaagattgt attcaaagtc aataatcaag tggtttcctt gaagcccgag
       61 atcattgtag accagtatga gtacaggtac ccgcgtatca aagacctgaa ggaacccagt
      121 ataaccttag ggaaggctcc tgacctaaac aaagcataca agtctgtttt gtcgggcatg
      181 aatgctgcca aacttgatcc agatgatgtg tgctcctatt tggcagctgc aatgcagttt
      241 ttcgagggat cctgtcctga ggactggacc agctacggga tcctgattgc acgaaaagga
      301 gacaagatca ctccagattc tcttgtggat ataaaacgta ctgatgtaga gggaagttgg
      361 gccctgacag gaggaatgga gttaacgaga gaccccactg tttccgagca tgcatcttta
      421 gtcggtcttc tcttgagtct gtatcggttg agcaaaatat cggggcaaaa cactggcaac
      481 tacaagacaa acattgcaga taggatagag cagattttcg agacagcccc tttcgtcaag
>DQ281111
        1 acagtcgaca atggatgccg acaagattgt gttcaaagtc aataatcagg tggtctcttt
       61 gaagcctgag attatcgtgg atcaatatga gtacaagtac cctgccatca aggatttgaa
      121 aaagccttgt atcaccctag ggaaagcccc cgacttgaac aaagcataca aatcagtttt
      181 atcaggcatg aatgccgcca aacttgatcc ggatgatgta tgctcctact tggcagcagc
      241 aatgcagttc tttgagggga catgtccgga agactggacc agctatggaa tcctgattgc
      301 acgaaaagga gataggatca ccccaaactc tctagtggag ataaagcgta ctgatgtaga
      361 agggaattgg gctctgacag gaggcatgga attgacaagg gaccccactg tctctgaaca
      421 tgcatcttta gtcggtcttc tcctgagtct gtacaggttg agcaaaatat caggacagaa
      481 cactggtaac tataagacaa acattgcaga taggatagag cagattttcg agacagcacc
      541 ttttgttaag atcgtggaac accataccct aatgacaact cacaagatgt gtgctaattg

序列存在多个问题:

    1. 一条序列分很多行;
    1. 序列开头用数字标识,且数字前是空格还是tab不明;
    1. 每行序列中,每10bp用空格分隔;
    1. 来源文件是非生信操作人员,熟悉Windows,故换行符很可能是\r\n;
    1. 碱基是使用的小写字母,大写字母为常态。

基于以上几点的考虑,一个功能一个功能的实现:

awk '{if($0~/^>/) ID=$0 ; else seq[ID]=seq[ID]$0;}END{for (i in seq) {gsub(" ","",seq[i]);gsub("[0-9]","",seq[i]);print i"\n"toupper(seq[i])}}' 264_N_complete.txt

    1. 使用vim -b 264_N_complete.txt打开文件,输入:%s/^M//将\r换行符替换掉,^M为ctrl+v ctrl+m建的组合,不是直接键盘打出来的。
    1. {if($0~/^>/) ID=$0 ; else seq[ID]=seq[ID]$0;} 将多行序列转换成单行,并将序列ID和Seq存如hash;
    1. gsub(" ","",seq[i]) 将空格替换成空;
    1. gsub("[0-9]","",seq[i])将数字转化为空
    1. toupper(seq[i] 将小写字母转换成大写字母
    1. print i"\n"toupper(seq[i]) 打印标准格式的fasta序列文件。序列后没有加“\n”是因为最后一行$0是带着\n的。若没有的话就再加个“\n”。

至此,序列转化为:

>AF371111
ATGGATGCCGACAAGATTGTATTCAAAGTCAATAATCAAGTGGTTTCCTTGAAGCCCGAGATCATTGTAGACCAGTATGAGTACAGGTACCCGCGTATCAAAGACCTGAAGGAACCCAGTATAACCTTAGGGAAGGCTCCTGACCTAAACAAAGCATACAAGTCTGTTTTGTCGGGCATGAATGCTGCCAAACTTGATCCAGATGATGTGTGCTCCTATTTGGCAGCTGCAATGCAGTTTTTCGAGGGATCCTGTCCTGAGGACTGGACCAGCTACGGGATCCTGATTGCACGAAAAGGAGACAAGATCACTCCAGATTCTCTTGTGGATATAAAACGTACTGATGTAGAGGGAAGTTGGGCCCTGACAGGAGGAATGGAGTTAACGAGAGACCCCACTGTTTCCGAGCATGCATCTTTAGTCGGTCTTCTCTTGAGTCTGTATCGGT
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
平台声明:文章内容(如有图片或视频亦包括在内)由作者上传并发布,文章内容仅代表作者本人观点,简书系信息发布平台,仅提供信息存储服务。

推荐阅读更多精彩内容