主页 / unix / 问题 / 563901
Accepted
Ricardo Guerreiro
Asked: 2020-01-25 07:56:27 +0800 CST

Unix折叠命令行为异常

  • 772

所以我有这个fasta（生物学）文件，看起来像这样：

>m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATCAAATAATTTCATTCGCGTCAGAACACCTGGTTTCACGACC
ATAAATAATTTACCAGTGAATCGAGGCTCAATTATAGATCCTCGGACGCGAGTTCTCGGTTGACGAGTGG
GATTCGAATTATTTTTCACCGAAAATTTTAGTCGACGAGTTCAGATAAATTTGTTCGGGATAAAATCATC
TGAGTAGGTCGGGCTTCTGAATTTCGTATTCTTGCGAGCAATGAATTTTAAATAATCATCGGACATACCA
ATTTTTGGAACAATAATGTTCCGAACATCCCGAAAATATAGGAAGAGCCCGGATAGATAAAAATAAACAC

每行最长为 70 个字符。通常，如果我想将其格式化为最多 50 个字符，我使用：

fold -50 input.fasta > output.fasta # 也试过 -b 和 -w args

但不知何故，这是行不通的。该文件看起来与我见过的许多其他文件完全相同。输出现在如下所示：

>m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATCAAATAATTTCATTCGCGTCA
GAACACCTGGTTTCACGACC
ATAAATAATTTACCAGTGAATCGAGGCTCAATTATAGATCCTCGGACGCG
AGTTCTCGGTTGACGAGTGG
GATTCGAATTATTTTTCACCGAAAATTTTAGTCGACGAGTTCAGATAAAT
TTGTTCGGGATAAAATCATC
TGAGTAGGTCGGGCTTCTGAATTTCGTATTCTTGCGAGCAATGAATTTTA
AATAATCATCGGACATACCA
ATTTTTGGAACAATAATGTTCCGAACATCCCGAAAATATAGGAAGAGCCC

它剪切了突出的 20 个字符并将它们正确放置在下面，但是它没有加入下一行并将其剪切到最多 50 个字符上。

我回到以前创建的 fasta 文件， fold 命令仍然正常工作。如果我复制新文件的一部分并将其粘贴到另一个文件中，问题仍然存在。

我认为可能存在我不知道的编码问题。任何人都可以帮忙吗？

干杯，

编辑：很好的答案，谢谢！

bioinformatics fold

3 个回答

  1. Best Answer

    您的问题与文件的编码无关。该fold实用程序非常原始，可以以特定长度断开线条，但它不会连接线条。

    您可能还需要小心保留 fasta 标题行（即，不要折叠这些）。

    awk -v W=50 '
    /^>/ { if (seq != "") print seq; print; seq = ""; next }
    {
        seq = seq $1
        while (length(seq) > W) {
            print substr(seq, 1,W)
            seq = substr(seq, 1+W)
        }
    }
    END { if (seq != "") print seq }' file.fa

    此awk命令会将您的序列重新格式化为 50 个字符，而标题行保持不变。宽度 50 可以通过W变量进行调整，并且可以设置为任何正整数。

    代码中的第一个块处理标题行，并将输出前一个序列中累积的序列位（如果有任何剩余要输出），然后将标题行未经修改地传递到输出。

    第二个块累积一行序列，如果它足够长，可能会以适当的块输出累积的序列。

    最后一个块 ( END) 在到达输入末尾时输出任何剩余的序列。

    在包含两个序列副本的文件上运行此命令将产生

    >m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATCAAATAATTTCATTCGCGTCA
GAACACCTGGTTTCACGACCATAAATAATTTACCAGTGAATCGAGGCTCA
ATTATAGATCCTCGGACGCGAGTTCTCGGTTGACGAGTGGGATTCGAATT
ATTTTTCACCGAAAATTTTAGTCGACGAGTTCAGATAAATTTGTTCGGGA
TAAAATCATCTGAGTAGGTCGGGCTTCTGAATTTCGTATTCTTGCGAGCA
ATGAATTTTAAATAATCATCGGACATACCAATTTTTGGAACAATAATGTT
CCGAACATCCCGAAAATATAGGAAGAGCCCGGATAGATAAAAATAAACAC
>m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATCAAATAATTTCATTCGCGTCA
GAACACCTGGTTTCACGACCATAAATAATTTACCAGTGAATCGAGGCTCA
ATTATAGATCCTCGGACGCGAGTTCTCGGTTGACGAGTGGGATTCGAATT
ATTTTTCACCGAAAATTTTAGTCGACGAGTTCAGATAAATTTGTTCGGGA
TAAAATCATCTGAGTAGGTCGGGCTTCTGAATTTCGTATTCTTGCGAGCA
ATGAATTTTAAATAATCATCGGACATACCAATTTTTGGAACAATAATGTT
CCGAACATCCCGAAAATATAGGAAGAGCCCGGATAGATAAAAATAAACAC

    改为W30 给出

    >m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATC
AAATAATTTCATTCGCGTCAGAACACCTGG
TTTCACGACCATAAATAATTTACCAGTGAA
TCGAGGCTCAATTATAGATCCTCGGACGCG
AGTTCTCGGTTGACGAGTGGGATTCGAATT
ATTTTTCACCGAAAATTTTAGTCGACGAGT
TCAGATAAATTTGTTCGGGATAAAATCATC
TGAGTAGGTCGGGCTTCTGAATTTCGTATT
CTTGCGAGCAATGAATTTTAAATAATCATC
GGACATACCAATTTTTGGAACAATAATGTT
CCGAACATCCCGAAAATATAGGAAGAGCCC
GGATAGATAAAAATAAACAC
>m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATC
AAATAATTTCATTCGCGTCAGAACACCTGG
TTTCACGACCATAAATAATTTACCAGTGAA
TCGAGGCTCAATTATAGATCCTCGGACGCG
AGTTCTCGGTTGACGAGTGGGATTCGAATT
ATTTTTCACCGAAAATTTTAGTCGACGAGT
TCAGATAAATTTGTTCGGGATAAAATCATC
TGAGTAGGTCGGGCTTCTGAATTTCGTATT
CTTGCGAGCAATGAATTTTAAATAATCATC
GGACATACCAATTTTTGGAACAATAATGTT
CCGAACATCCCGAAAATATAGGAAGAGCCC
GGATAGATAAAAATAAACAC

    您可能还对CSHL的 FASTX-Toolkit感兴趣。我自己从未使用过它，但它似乎包含一个“FASTA 格式化程序（更改 FASTA 文件中序列行的宽度）”。这些工具的最新版本来自 2014 年（相当旧），因此您可能希望自己从源代码编译它们，而不是使用提供的预编译二进制文件之一，除非您的特定 Unix 发行版提供了一个包（检查您的包存储库）。

    • 5

  2. 尝试这个：

    <input.fasta tr -d '\n'|fold -w 50 >output.fasta

    这用于tr删除现有的行尾，然后将生成的单行格式化为多行，每行最大长度为 50。

    要将第一行保持在当前长度，而不是将其与以下行合并，这应该可以工作（并以行尾结束输出）：

    awk '{if (NR==1) {print $0 gensub(/ /, " ", "g", sprintf("%*s", 50-length($0), ""))} else print $0}' input.fasta|tr -d '\n'|sed '$s/$/\n/'|fold -w 50|awk '{$1=$1};1' >output.fasta
    • 2

  3. 这就是fold工作原理。你以前从未播过它，因为你以前没有这种长度的线。折叠分别发生在每一行。因此，如果线的长度不是您要折叠到的大小的精确倍数，您将获得这种输出。例如：

    $ perl -le 'for (0..2){print "A" x 12}' 
AAAAAAAAAAAA
AAAAAAAAAAAA
AAAAAAAAAAAA
$ perl -le 'for (0..2){print "A" x 12}' | fold -w 6
AAAAAA
AAAAAA
AAAAAA
AAAAAA
AAAAAA
AAAAAA
$ perl -le 'for (0..2){print "A" x 12}' | fold -w 7
AAAAAAA
AAAAA
AAAAAAA
AAAAA
AAAAAAA
AAAAA

    现在，这实际上不是问题。这仍然是一个有效的 fasta 文件，但它不是很漂亮。作为一种解决方法，您可以使用我之前发布的脚本FastaToTbl并执行以下操作：TblToFasta

    $ FastaToTbl input.fasta | TblToFasta 
>m64093_191209_130050/133911/ccs_64 
TTCAGGCTGTGTTCCATTTGATTTAAAATCAAATAATTTCATTCGCGTCAGAACACCTGG
TTTCACGACCATAAATAATTTACCAGTGAATCGAGGCTCAATTATAGATCCTCGGACGCG
AGTTCTCGGTTGACGAGTGGGATTCGAATTATTTTTCACCGAAAATTTTAGTCGACGAGT
TCAGATAAATTTGTTCGGGATAAAATCATCTGAGTAGGTCGGGCTTCTGAATTTCGTATT
CTTGCGAGCAATGAATTTTAAATAATCATCGGACATACCAATTTTTGGAACAATAATGTT
CCGAACATCCCGAAAATATAGGAAGAGCCCGGATAGATAAAAATAAACAC

    该TblToFasta脚本将确保输出是标准的每行 60bp。如果你真的需要 50，你可以这样做（假设 GNU sed）：

    $ FastaToTbl input.fasta | sed -E 's/^/>/;s/\t/\n/ ' | sed -E 's/(.{50})/\1\n/g'
>m64093_191209_130050/133911/ccs_64
TTCAGGCTGTGTTCCATTTGATTTAAAATCAAATAATTTCATTCGCGTCA
GAACACCTGGTTTCACGACCATAAATAATTTACCAGTGAATCGAGGCTCA
ATTATAGATCCTCGGACGCGAGTTCTCGGTTGACGAGTGGGATTCGAATT
ATTTTTCACCGAAAATTTTAGTCGACGAGTTCAGATAAATTTGTTCGGGA
TAAAATCATCTGAGTAGGTCGGGCTTCTGAATTTCGTATTCTTGCGAGCA
ATGAATTTTAAATAATCATCGGACATACCAATTTTTGGAACAATAATGTT
CCGAACATCCCGAAAATATAGGAAGAGCCCGGATAGATAAAAATAAACAC
    • 2

相关问题