主页 / unix / 问题 / 452602
Accepted
CN_229133
Asked: 2018-06-30 01:32:12 +0800 CST

计算最长的连续模式

  • 772

我有一个序列文件：

$猫文件
CACCGTTGCCAAACAATG
TTAGAAGCCTGTCAGCCT
CATTGCTCTCAGACCCAC
GATGTACGTCACATTAGA
ACACGGAATCTGCTTTTT
CAGAATTCCCAAAGATGG

我想计算最长的 C+T 段。我只能计算总 C+T，但我想要最长的延伸。

$猫文件 | awk '{ 打印 $0, gsub(/[cCtT]/,"",$1)}'
CACCGTTGCCAAACAATG 9
TTAGAAGCCTGTCAGCCT 10
CATTGCTCTCAGACCCAC 12
GATGTACGTCACATTAGA 8
ACACGGAATCTGCTTTTT 11
CAGAATTCCCAAAGATGG 7

预期结果将显示最长的 C+T 拉伸。

CACCGTTGCCAAACAATG 9 2
TTAGAAGCCTGTCAGCCT 10 3
CATTGCTCTCAGACCCAC 12 5
GATGTACGTCACATTAGA 8 2
ACACGGAATCTGCTTTTT 11 6
CAGAATTCCCAAAGATGG 7 5
awk pattern-matching

5 个回答

  1. Best Answer

    FWIW 这是一种在 perl 中执行此操作的方法，使用maxfromList::Util

    $ perl -MList::Util=max -lpe '$_ .= " " . max 0, map length, /[CT]+/gi' file
CACCGTTGCCAAACAATG 2
TTAGAAGCCTGTCAGCCT 3
CATTGCTCTCAGACCCAC 5
GATGTACGTCACATTAGA 2
ACACGGAATCTGCTTTTT 6
CAGAATTCCCAAAGATGG 5
    • 5 
  2. $ awk '{ split($0, a, "[^CTct]+"); m=0
         for (i in a) {
             len=length(a[i])
             if (len > m) m=len
         }
         print $0, m  }' file
CACCGTTGCCAAACAATG 2
TTAGAAGCCTGTCAGCCT 3
CATTGCTCTCAGACCCAC 5
GATGTACGTCACATTAGA 2
ACACGGAATCTGCTTTTT 6
CAGAATTCCCAAAGATGG 5

    该awk程序在运行任何不是大写或小写C或T. 然后它遍历分割结果的位，找到最长的位。然后它将原始行与找到的最大长度一起打印。

    由于 Roman 查看了各种解决方案的时间安排，因此这里有一个更快的解决方案：

    awk -F "[^CTct]+" '
    m = 0
    for (i = 1; i <= NF; ++i) {
        len = length($i)
        if (len > m) m = len
    }
    print m' file | paste file -

    它更快，因为它只分割线一次。第一个代码还尝试将输入行拆分为空格。

    mawk使用在 500000 行上显示 0.79 秒进行计时。第一个解决方案在相同的数据上使用了 1.69s，表明它可能是拆分操作花费的时间最多。

    • 4

  3. 更快的 GNUawk解决方案：

    awk -v FPAT='[ctCT]+' \
'{ 
     max_l = t_len = 0;
     for (i=1; i <= NF; i++) {
         len = length($i);
         if (len > max_l) max_l = len;
         t_len += len
     }
     print $0, t_len, max_l
 }' inputfile

    输出：

    CACCGTTGCCAAACAATG 9 2
TTAGAAGCCTGTCAGCCT 10 3
CATTGCTCTCAGACCCAC 12 5
GATGTACGTCACATTAGA 8 2
ACACGGAATCTGCTTTTT 11 6
CAGAATTCCCAAAGATGG 7 5

    时间性能对比（测试inputfile大概有120000行）：

    $ time awk -v FPAT='[ctCT]+' '{ max_l = t_len = 0; for (i=1; i <= NF; i++) { len = length($i); if (len > max_l) max_l = len; t_len += len } print $0, t_len, max_l }' inputfile > /dev/null

real    0m1.018s
user    0m0.948s
sys 0m0.012s

$ time awk '{ split($0, a, "[^CTct]+"); m=0; for (i in a) { len=length(a[i]); if (len > m) m=len } print $0, m }' inputfile > /dev/null

real    0m1.802s
user    0m1.688s
sys 0m0.028s

$ time perl -MList::Util=max -lpe '$_ = "$_ " . max map { length $_ } /[CT]*/gi' inputfile > /dev/null

real    0m1.216s
user    0m1.160s
sys 0m0.016s

$ time sed 'h;y/cCtT/xxxx/;x;H;s/./x/g;G; s/^\(x*\).*\n.*\1.*\n/\1 /; s/^x\{10\}/1/;s/$/:9876543210xxxxxxxxx/; s/^\(1*\)\(x*\) \(.*\):.*\(.\).\{9\}\2$/\3 \1\4/' inputfile > /dev/null

real    1m4.165s
user    1m2.784s
sys 0m0.352s
    • 3

  4. 使用sed（假设每行不超过 19 个字符），只是为了好玩并使用 RE 匹配的贪婪属性：

    sed '
  h;y/cCtT/xxxx/;x;H;s/./x/g;G
  s/^\(x*\).*\n.*\1.*\n/\1 /
  s/^x\{10\}/1/;s/$/:9876543210xxxxxxxxx/
  s/^\(1*\)\(x*\) \(.*\):.*\(.\).\{9\}\2$/\3 \1\4/'

    @Kusalananda 解决方案的一个变体：

    awk -F '[^cCtT]+' '
  {
    max = 0
    for (i = 1; i <= NF; i++)
      if ((l = length($i)) > max)
        max = l
    print $0, max
  }'
    • 3

  5. 也试试

    awk '
        {T0 = $0
         while (match (T0, /[CTct]+/))  {if (RLENGTH > MX) MX = RLENGTH
                                         T0 = substr (T0, RSTART+RLENGTH)
                                        }
         print $0, MX
        }
' file

    时间比其他awk提案快一点。

    • 0

相关问题