提取包含两种模式的线

Question

我有一个包含以下几行的文件：

>header1
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCGGGCCTCTTTTCCTGACGGCCGCCCCCACTGCCCCCACGACCGGCCCGTACAAC<pattern_2>
>header2
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_2>
>header3
<pattern_1>ATGGCCACCAACAACCAGAGCTCCC
>header4
GACCGGCACGTACAACCTCCAGGAAATCGTGCCCGGCAGCGTGTGGATGGAGAGGGACGTG
>header5
TGCCCCCACGACCGGCACGTACAAC<pattern_2>

我想提取包含 header 行的所有行。

我试过使用 grep，但它只提取序列行而不是 header 行。

grep <pattern_1> | grep <pattern_2> input.fasta > output.fasta

如何提取包含模式和 Linux 中的 header 的行？图案可以出现在线条的任何位置。不限于行首或行尾。

预期输出：

>header1
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCGGGCCTCTTTTCCTGACGGCCGCCCCCACTGCCCCCACGACCGGCCCGTACAAC<pattern_2>
>header2
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_2>

Answer 1

如果您希望 grep 打印匹配项周围的行，请使用 -B 标记表示匹配前的行，-A 表示匹配后的行，-C 表示匹配前后的行。

在你的情况下，grep -B 1 似乎可以完成这项工作。

Answer 2

$ grep -A 1 header[12] file
>header1
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCGGGCCTCTTTTCCTGACGGCCGCCCCCACTGCCCCCACGACCGGCCCGTACAAC<pattern_2>
>header2
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_2>

man grep:

   -A NUM, --after-context=NUM
          Print  NUM  lines  of  trailing  context  after  matching lines.
          Places  a  line  containing  a  group  separator  (--)   between
          contiguous  groups  of  matches.  With the -o or --only-matching
          option, this has no effect and a warning is given.

   -B NUM, --before-context=NUM
          Print NUM  lines  of  leading  context  before  matching  lines.
          Places   a  line  containing  a  group  separator  (--)  between
          contiguous groups of matches.  With the  -o  or  --only-matching
          option, this has no effect and a warning is given.

grep -B 1 pattern_[12] 也可以，但是样本数据中有几个 pattern_1，所以...这次不行。

Answer 3

您可以像这样使用 awk 轻松做到这一点：

awk '/^>/{h=[=10=];next}
     /<pattern_1>/&&/<pattern_2>/{print h;print}' input.fasta > output.fasta

这里还有一个 sed 解决方案，它也能产生所需的输出：

sed -n '/^>/{N;/<pattern_1>/{/<pattern_2>/p}}' input.fasta > output.fasta

如果很可能存在多行记录，可以这样用：

awk -v pat1='<pattern_1>' -v pat2='<pattern_2>' '
/^>/ {r=[=12=];p=0;next}
!p {r=r ORS [=12=];if(chk()){print r;p=1};next}
p

function chk(   tmp){
    tmp=gensub(/\n/,"","g",r)
    return (tmp~pat1&&tmp~pat2)
}' input.fasta > output.fasta

Answer 4

如果您的输入文件与 post 中描述的完全一致，那么您可以使用：

grep -B1 '^<pattern_1>.*<pattern_2>$' input 
>header1
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCGGGCCTCTTTTCCTGACGGCCGCCCCCACTGCCCCCACGACCGGCCCGTACAAC<pattern_2>
>header2
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_2>

其中 -B1 将在匹配行的顶部显示它之前的行。使用的正则表达式基于以下假设：您的 2 个模式在行首和行尾的顺序完全相同。如果不是这种情况：使用 '.*<pattern_1>.*<pattern_2>.*'。最后但并非最不重要的一点是，如果并不总是遵守 2 种模式的顺序，那么您可以使用：'^.*<pattern_1>.*<pattern_2>.*$\|^.*<pattern_2>.*<pattern_1>.*$'

在以下输入文件中：

cat input
>header1
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCGGGCCTCTTTTCCTGACGGCCGCCCCCACTGCCCCCACGACCGGCCCGTACAAC<pattern_2>
>header2
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_2>
>header2b
<pattern_2>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_1>
>header3
<pattern_1>ATGGCCACCAACAACCAGAGCTCCC
>header4
GACCGGCACGTACAACCTCCAGGAAATCGTGCCCGGCAGCGTGTGGATGGAGAGGGACGTG
>header5
TGCCCCCACGACCGGCACGTACAAC<pattern_2>

输出：

grep -B1 '^.*<pattern_1>.*<pattern_2>.*$\|^.*<pattern_2>.*<pattern_1>.*$' input 
>header1
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCGGGCCTCTTTTCCTGACGGCCGCCCCCACTGCCCCCACGACCGGCCCGTACAAC<pattern_2>
>header2
<pattern_1>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_2>
>header2b
<pattern_2>CGGCGGGCAGATGGCCACCAACAACCAGAGCTCCCTGGCCTGCAATCACTACTCGTGTTTTGCCACCACTGCCCCCACGACCGGCACGTACAAC<pattern_1>

Answer 5

您可能对 BioAwk 感兴趣，它是 awk 的改编版本，可以处理 fasta 文件

bioawk -c fastx -v seq1="pattern1" -v seq2="pattern2" \
       '($seq ~ seq1) && ($seq ~ seq2) { print ">"$name; print $seq }' file.fasta

如果想要开头是seq1，结尾是seq2，可以改成：

bioawk -c fastx -v seq1="pattern1" -v seq2="pattern2" \
       '($seq ~ "^"seq1) && ($seq ~ seq2"$") { print ">"$name; print $seq }' file.fasta

这对于处理 fasta 文件非常实用，因为序列通常分布在多行中。上面的代码很容易处理这个问题，因为变量 $seq 包含完整的序列。

如果你不想安装BioAwk，你可以使用下面的方法来处理你的FASTA文件。它将允许 multi-line 序列并执行以下操作：

一次读取一条记录（假设 header 中没有 >，第一个字符除外）
从记录中提取 header 并将其存储在 name 中（不是真的需要）
将完整序列合并为一个字符串，删除所有换行符和空格。这确保如果模式被拆分成多行，搜索 pattern1 或 pattern2 不会失败。
如果找到匹配项，则打印记录。

以下 awk 执行请求：

awk -v seq1="pattern1" -v seq2="pattern2" \
    'BEGIN{RS=">"; ORS=""; FS="\n"}
     { seq="";for(i=2;i<=NF;++i) seq=seq""$i; gsub(/[^a-zA-Z0-9]/,"",seq) }
     (seq ~ seq1 && seq ~ seq2){print ">" [=12=]}' file.fasta

如果记录 header 包含不在行首的其他 > 个字符，您必须采取稍微不同的方法（除非您使用 GNU awk）

awk -v seq1="pattern1" -v seq2="pattern2" \
    '/^>/ && (seq ~ seq1 && seq ~ seq2) {
         print name
         for(i=0;i<n;i++) print aseq[i]
     }
     /^>/ { seq=""; delete aseq; n=0; name=[=13=]; next }
     { aseq[n++] = [=13=]; seq=seq""[=13=]; sub(/[^a-zA-Z0-9]*$/,"",seq) }
     END { if (seq ~ seq1 && seq ~ seq2) {
              print name
              for(i=0;i<n;i++) print aseq[i]
            }
     }' file.fasta

注意： 我们在这里使用 sub 以防在 fasta 文件中引入意外字符（例如 spaces/tabs 或 CR (\r))

注意： BioAwk 基于 Brian Kernighan 的 awk，它记录在 "The AWK Programming Language", 作者：Al Aho、Brian Kernighan 和 Peter Weinberger （Addison-Wesley，1988 年，ISBN 0-201-07981-X）。我不确定这个版本是否与 POSIX.

兼容

提取包含两种模式的线

Extract lines containing two patterns

linux

awk

grep

fasta