-1

我正在学习如何将 perl 用于基因组学应用程序。我正在尝试清理成对的末端读取(1 个正向,1 个反向)。这些存储在 2 个文件中,但行匹配。我遇到的麻烦是让相关子例程从第二个文件中读取(我得到的警告是针对未初始化的值)。

这些文件设置在 4 行块(fastq)中,其中第一行是运行 ID,第二行是序列,第三行是“+”,第四行保存第 2 行中序列的质量值。

当它仅应用于一个文件时,我对这段代码没有真正的问题,但我认为我误解了如何处理多个文件。

非常感谢任何指导!

我在这种情况下的警告是这样的:在 ./pairedendtrim.pl 第 137 行第 4 行的减法 (-) 中使用未初始化的值 $thisline。

#!/usr/bin/perl
#pairedendtrim.pl by AHU
use strict;
use warnings;

die "usage: readtrimmer.pl <file1> <file2> <nthreshold> " unless @ARGV == 3;
my $nthreshold = "$ARGV[2]";

open( my $fastq1, "<", "$ARGV[0]" );
open( my $fastq2, "<", "$ARGV[1]" );

my @forline;
my @revline;
while ( not eof $fastq2 and not eof $fastq1 ) {
    chomp $fastq1;
    chomp $fastq2;
    $forline[0] = <$fastq1>;
    $forline[1] = <$fastq1>;
    $forline[2] = <$fastq1>;
    $forline[3] = <$fastq1>;

    $revline[0] = <$fastq2>;
    $revline[1] = <$fastq2>;
    $revline[2] = <$fastq2>;
    $revline[3] = <$fastq2>;

    my $ncheckfor = removen( $forline[1] );

    my $ncheckrev = removen( $revline[1] );

    my $fortest = 0;
    if ( $ncheckfor =~ /ok/ ) { $fortest = 1 }

    my $revtest = 0;

    if ( $ncheckrev =~ /ok/ ) { $revtest = 1 }

    if ( $fortest == 1 and $revtest == 1 ) { print "READ 1 AND READ 2" }

    if ( $fortest == 1 and $revtest == 0 ) { print "Read 1 only" }

    if ( $fortest == 0 and $revtest == 1 ) { print "READ 2 only" }

}

sub removen {
    my ($thisline) = $_;
    my $ntotal = 0;
    for ( my $i = 0; $i < length($thisline) - 1; $i++ ) {
        my $pos = substr( $thisline, $i, 1 );
        #print "$pos\n";
        if ( $pos =~ /N/ ) { $ntotal++ }
    }
    my $nout;
    if ( $ntotal <= $nthreshold )    #threshold for N
    {
        $nout = "ok";
    } else {
        $nout = "bad";
    }
    return ($nout);
}
4

1 回答 1

0

子程序的参数在 中@_,而不是$_

sub removen {
    my ($thisline) = @_;

我还有一些其他的提示给你:

  1. use autodie;随时进行文件处理。
  2. 首先将值分配@ARGV给变量。这很快记录了持有的内容。
  3. 不要chomp一个文件句柄。这没有任何作用。而是应用于chomp读取返回的值。
  4. 不要使用字符串okbad作为布尔值。
  5. tr可用于计算字符在字符串中出现的次数。

以下是您的代码的清理版本:

#!/usr/bin/perl
#pairedendtrim.pl by AHU
use strict;
use warnings;
use autodie;

die "usage: readtrimmer.pl <file1> <file2> <nthreshold> " unless @ARGV == 3;

my ( $file1, $file2, $nthreshold ) = @ARGV;

open my $fh1, '<', $file1;
open my $fh2, '<', $file2;

while ( not eof $fh2 and not eof $fh1 ) {
    chomp( my @forline = map { scalar <$fh1> } ( 1 .. 4 ) );
    chomp( my @revline = map { scalar <$fh2> } ( 1 .. 4 ) );

    my $ncheckfor = removen( $forline[1] );
    my $ncheckrev = removen( $revline[1] );

    print "READ 1 AND READ 2" if $ncheckfor  and $ncheckrev;
    print "Read 1 only"       if $ncheckfor  and !$ncheckrev;
    print "READ 2 only"       if !$ncheckfor and $ncheckrev;
}

sub removen {
    my ($thisline) = @_;
    my $ntotal = $thisline =~ tr/N/N/;
    return $ntotal <= $nthreshold;    #threshold for N
}
于 2014-09-24T16:36:35.080 回答