首页 文章

重新复制FASTA，保留seq id

提问于
浏览 次
1

我需要格式化miRNA识别工具（miREAP）的文件 .

我有一个以下格式的fasta文件：

>seqID_1
CCCGGCCGTCGAGGC
>seqID_2
AGGGCACGCCTGCCTGGGCGTCACGC
>seqID_3
CCGCATCAGGTCTCCAAGGTGAACAGCCTCTGGTCGA
>seqID_4
CCGCATCAGGTCTCCAAGGTGAACAGCCTCTGGTCGA
>seqID_5
CCGCATCAGGTCTCCAAGGTGAACAGCCTCTGGTCGA
>seqID_6
AGGGCACGCCTGCCTGGGCGTCACGC

我想计算每个序列发生的次数，并将该数字附加到seqID行 . 每个序列的计数和引用序列的原始ID只需在文件中出现一次，如下所示：

>seqID_1 1
CCCGGCCGTCGAGGC
>seqID_2 2
AGGGCACGCCTGCCTGGGCGTCACGC
>seqID_3 3
CCGCATCAGGTCTCCAAGGTGAACAGCCTCTGGTCGA

Fastx_collapser几乎就像我想要的那样（http://hannonlab.cshl.edu/fastx_toolkit/index.html） . 但是，它不是维护seqID，而是返回：

>1 1
CCCGGCCGTCGAGGC
>2 2
AGGGCACGCCTGCCTGGGCGTCACGC
>3 3
CCGCATCAGGTCTCCAAGGTGAACAGCCTCTGGTCGA

这意味着我的序列，seqID和基因组作图位置之间的链接丢失了 . （每个seqID对应于我的fasta文件中的序列和单独的Bowtie2生成的.sam文件中的基因组映射点）

是否有一种简单的方法可以在命令行中执行所需的重复数据删除？

谢谢！

command-line formatting bioinformatics fasta

1 回答

  • 1

    线性化和排序/ uniq -c

    awk '/^>/ {if(N>0) printf("\n"); ++N; printf("%s ",$0);next;} {printf("%s",$0);} END { printf("\n");}'  input.fa  | \
  sort -t ' ' -k2,2 | uniq -f 1 -c |\
  awk '{printf("%s_%s\n%s\n",$2,$1,$3);}'

>seqID_2_2
AGGGCACGCCTGCCTGGGCGTCACGC
>seqID_1_1
CCCGGCCGTCGAGGC
>seqID_3_3
CCGCATCAGGTCTCCAAGGTGAACAGCCTCTGGTCGA
    回复于

相关问题