【问题标题】:Adding columns to a file based on existing columns根据现有列向文件添加列
【发布时间】:2013-01-18 22:36:51
【问题描述】:

我正在尝试修改这样设置的文件:

chr start   ref alt 
chr1    18884   C   CAAAA
chr1    135419  TATACA  T
chr1    332045  T   TTG
chr1    453838  T   TAC
chr1    567652  T   TG
chr1    602541  TTTA    T
chr1    614937  C   CTCTCTG
chr1    654889  C   CA
chr1    736800  AC  A

我想修改它: 如果列“ref”是一个字符串 >1(即第 2 行),那么我生成 2 个新列,其中:

第一个新列 = 起始坐标-1 第二个新列=起始坐标+(参考字符串的长度)+1

因此,第 2 行的输出如下所示:

chr1 135419 TATACA T 135418 135426

或: 如果 "ref" 中的字符串长度 = 1 且列 "alt"=长度>1 的字符串(即第 1 行)则

第一个新列 = 起始坐标 第二个新列 = 起始坐标+2

所以,第 1 行的输出将是:

chr1 18884 C CAAAA 18884 18886

我在 awk 中尝试过,但没有成功 我的 perl 不存在,但这是最好的方法吗?或者也许在 R 中?

【问题讨论】:

  • 如果我正确阅读了规范,你的两个 if 会合二为一(如 1 + 1 = 2),除非你确实有可能缺少 alt 的情况。如果是这样(alt 可能会丢失),那么您可能有一个固定的记录字段输入并且缺少规范 - 因为下面的每个解决方案都会在空白处拆分。

标签: perl awk calculated-columns tsv


【解决方案1】:

Perl 解决方案。请注意,您的规范没有提及如果两个字符串的长度均为 1 时该怎么做。

#!/usr/bin/perl
use warnings;
use strict;
use feature qw(say);

#use Data::Dumper;
<DATA>; # Skip the header;
while (<DATA>) {
    my ($chr, $start, $ref, $alt) = split;
    my @cols;
    if (1 < length $ref) {
          @cols = ( $start - 1, $start + 1 + length $ref);
    } elsif (1 < length $alt) {
        @cols = ($start, $start + 2);
    } else {
        warn "Don't know what to do at $.\n";
    }
    say join "\t", $chr, $start, $ref, $alt, @cols;
}


__DATA__
chr start   ref alt
chr1    18884   C   CAAAA
chr1    135419  TATACA  T
chr1    332045  T   TTG
chr1    453838  T   TAC
chr1    567652  T   TG
chr1    602541  TTTA    T
chr1    614937  C   CTCTCTG
chr1    654889  C   CA
chr1    736800  AC  A

【讨论】:

    【解决方案2】:

    这是使用awk 的一种方式。运行如下:

    awk -f script.awk file | column -t
    

    script.awk的内容:

    NR==1 {
        next
    }
    
    length($3)>1 && length($4)==1 {
        print $0, $2-1, $2+length($3)+1
        next
    }
    
    length($3)==1 && length($4)>1 {
        print $0, $2, $2+2
        next
    }1
    

    结果:

    chr1  18884   C       CAAAA    18884   18886
    chr1  135419  TATACA  T        135418  135426
    chr1  332045  T       TTG      332045  332047
    chr1  453838  T       TAC      453838  453840
    chr1  567652  T       TG       567652  567654
    chr1  602541  TTTA    T        602540  602546
    chr1  614937  C       CTCTCTG  614937  614939
    chr1  654889  C       CA       654889  654891
    chr1  736800  AC      A        736799  736803
    

    或者,这里是单行:

    awk 'NR==1 { next } length($3)>1 && length($4)==1 { print $0, $2-1, $2+length($3)+1; next } length($3)==1 && length($4)>1 { print $0, $2, $2+2; next }1' filem | column -t
    

    代码应该是不言自明的。脚本末尾的1 仅启用每行的默认打印(即“1”返回真)。 HTH。

    【讨论】:

      【解决方案3】:

      在 perl 中执行它是微不足道的(但在 awk 中也是如此):

      #!/usr/bin/perl
      while (<>) {
        chmop;
        my ($chr,$start,$ref,$alt)=split(/\s+/,$_);
        if (len($ref) > 1) {
      print STDOUT
        "$chr\t$start\t$ref\t$alt\t",
          $start+len($ref)+1,"\n";
        } elsif (len($ref)==1) {
      print STDOUT
        "$chr\t$start\t$ref\t$alt\t",
          $start+2,"\n";
        } else {
      print STDERR "ERROR: ???\n"; #actually impossible
        }
      }
      

      将它粘贴到文件 morecols.pl 中,chmod +x morecols.pl,运行 morecols.pl。 (请注意,此代码/说明中有很多假设)。我感觉您的实际问题更多是编程/文本处理,而不是工具或语言。如果是这样,这段代码只是权宜之计......

      干杯。

      【讨论】:

        猜你喜欢
        • 2018-08-13
        • 2020-11-07
        • 1970-01-01
        • 1970-01-01
        • 2015-11-19
        • 1970-01-01
        • 2022-01-22
        • 1970-01-01
        • 2019-09-29
        相关资源
        最近更新 更多