0

我有这样的数据，有 6 列

LINES XY1  XY2  XY3  XY4  XY5
P1    Z/Z  T/T  -/-  T/T  T/T
P2    A/A  A/A  G/G  Z/Z  T/T
1     G/G  T/T  G/G  T/T  G/G
2     T/T  A/A  C/C  C/C  T/T
3     T/T  G/G  T/T  G/G  T/T
4     A/A  C/C  A/A  A/A  A/A
5     A/A  A/A  T/T  T/T  A/A

  1. 我想找出哪些列（从XY1到XY5）对于行P1和P2. 相等意味着P1和P2包含相同的字母（等位基因）或它们中的任何一个是Z/Z或-/-。

  2. 我将比较从行1到跨列到5的列。如果它们与输出匹配，则应包含else 。我继续这个我的程序遇到第二组和行。P2XY1XY510P1P2

  3. 我将计算行的总行数，1以仅包括在和5之间不同的列。P1P2

  4. 我将通过将总和除以 和 之间不同的列数来1计算行的百分比。5P1P2

我期待这样的输出

LINES XY1  XY2  XY3  XY4  XY5       
P1    eq   nq   eq   eq   eq   SUM  %
P2                                  1   
1     0    0    1    0    0    0    0
2     0    1    0    0    1    1    100
3     0    0    0    0    1    0    0
4     1    0    0    0    0    0    0
5     1    1    0    0    0    1    100

我有超过 5,000 行的数据，目前我正在使用不同的公式在 Excel 2010 中工作，但这占用了我很多精力。

我想做这个 Perl，但我是 Perl 的新手。我已成功将文件读取到屏幕上。

这是我写到文件读取部分的代码。

#!/usr/bin/perl
use strict;
use warnings;

use Text::CSV;

my $file = 'csv.csv';

my $csv = Text::CSV->new();

open(CSV, "<", $file) or die $!;

while (<CSV>) {
  if ($csv->parse($_)) {
    my @columns = $csv->fields();
    print "@columns\n";
  }
  else {
    my $err = $csv->error_input;
    print "Failed to parse line: $err";
  }
}

close CSV;
4

1 回答 1

2

该程序似乎可以满足您的需求。它期望输入文件的路径作为命令行上的参数。生成的输出被发送到STDOUT.

很多人不相信我，但是用它来解析制表符分隔的数据是错误的。Text::CSV它需要设置分隔符，并禁用引用和转义选项，如果做得对，它最终的行为与split /\t/. 您没有说您的数据是制表符分隔的还是只是空格分隔的，但由于您似乎没有包含空格的数据，所以我假设后者。使用 just 可以非常简单地解析此类数据split。

它通过检查输入的每一行来查看它是P1一行、P2一行、一个编号行还是其他任何东西（假设是标题行）。

  • 当遇到标题行时，printf后续行的格式由标题字段的间距得出，并将该行复制到输出

  • 当P1遇到一行时，数据只是简单地保存在@p1

  • 当P2遇到一行时，将保存数据@p2以与后续数据进行比较。适当的P1和P2行被打印到输出，并@unequal计算一个数组，其中包含数据与数据不匹配的列的P1索引P2

  • 将编号的数据行与保存在 中的数据逐列进行比较@p2，并相应1地0插入相应的输出列中。的值$sum是通过将数组中列出的输出列的值相加来计算的@unequal。通过除以 中的条目数计算百分比@unequal，并将数据打印到输出。

笔记

目前还不清楚“第二组P1和P2”如何。启动，因此此代码可能无法正确处理它。您也没有说您希望如何显示小数百分比，因此此代码仅打印33.3333333333333.

此外，您没有说出1输出行末尾的P2表示什么，所以我只是按字面意思复制了它。

 

use strict;
use warnings;

sub compare_alleles {
  return 1 if grep {$_ eq '-/-' or $_ eq 'Z/Z' } @_;
  return $_[0] eq $_[1] ? 1 : 0;
}

my $format;
my (@p1, @p2);
my @unequal;

while (<>) {

  unless (/^(P?\d)/) {
    my @widths;
    push @widths, $+[1] - $-[1] while /(\S+\s*)/g;
    pop @widths;
    push @widths, $widths[-1], $widths[-1];
    $format = join '', map("%-${_}s", @widths, ''), "\n";
    print;
    next;
  }

  my @fields = split;

  if ($fields[0] eq 'P1') {
    @p1 = @fields;
  }
  elsif ($fields[0] eq 'P2') {
    @p2 = @fields;
    printf $format, 'P1', map (compare_alleles($p1[$_], $p2[$_]) ? 'eq' : 'nq', 1..5), 'SUM', '%';
    printf $format, 'P2', map('', 1..5), '', '1';
    @unequal = grep { not compare_alleles($p1[$_], $p2[$_]) } 1..5;
  }
  else {
    my @columns = ($fields[0], map { $fields[$_] eq $p2[$_] ? 1 : 0 } 1..5);
    my $sum = 0;
    $sum += $_ for @columns[@unequal];
    my $percent = $sum == 0 ? 0 : $sum * 100 / @unequal;
    printf $format, @columns, $sum, $percent;
  }
}

输出

LINES XY1  XY2  XY3  XY4  XY5
P1    eq   nq   eq   eq   eq   SUM  %
P2                                  1
1     0    0    1    0    0    0    0
2     0    1    0    0    1    1    100
3     0    0    0    0    1    0    0
4     1    0    0    0    0    0    0
5     1    1    0    0    0    1    100
于 2013-09-21T16:31:24.673 回答