2014-01-13 45 views
0

請幫助改進以下代碼。我無法在單行中打印序列。想輸出印有四行的核苷酸頻率的四個字符之一。提前致謝。 enter code here使用Perl計算fasta文件中的核苷酸頻率

#!/usr/bin/perl 
use strict; 
use warnings; 
my $A;  
my $T; 
my $G; 
my $C; 
my $fileIN; 
my $fileOUT; 

my $seq ; 
open ($fileIN, "basecount.nfasta") or die "can't open file "; 
open ($fileOUT, ">basecount.out") or die "can't open file "; 

while (<$fileIN>) 
{ 

      if ($_ =~/^>/) #ignore header line 
      {next;} 

      else 
        { 
        $seq = $_; #copy the all line with only nucleotide characters ATGC 
        } 
      $seq =~ s/\n//g; #create one single line containing all ATGC characters 

      print "$seq\n"; # verify previous step 

      my @dna = split ("",$seq); #create an array to include each nucleotide as array element 

      foreach my $element (@dna) 

      { 
      if ($element =~/A/) # match nucleotide pattern and countstrong text 
          { 
          $A++; 
          } 
      if ($element =~/T/) 
          { 
          $T++; 
          } 
      if ($element =~/G/) 
          { 
          $G++; 
          } 
      if ($element =~/C/) 
          { 
          $C++; 
          } 

      } 

      print $fileOUT "A=$A\n"; 
      print $fileOUT "T=$T\n"; 
      print $fileOUT "G=$G\n"; 
      print $fileOUT "C=$C\n"; 
} 

close ($fileIN); 
close ($fileOUT); 
+0

究竟是什麼問題?如果您想要將所有行的頻率合併到一起,只需將循環外的打印語句移動到循環外 –

+0

,謝謝將印刷語句移出循環工作。我擔心的是,雖然我已經替換了\ n字符,但爲什麼仍然看到核苷酸在多行中而不是在一行中? – harsh

+0

因爲您正在逐行讀取文件while語句。如果你想讀一個標量的所有行,你需要改變它。例如,使用file :: slurp模塊(有更多的方法只需在谷歌搜索它們) –

回答

1

首先,我會使用一些快捷方式。它更容易閱讀:

use strict; 
use warnings; 
use feature 'say'; 
my $A; 
my $T; 
my $G; 
my $C; 
my $fileIN; 
my $fileOUT; 

open $fileIN, '<',"basecount.nfasta" or die "can't open file basecount.nfasta for reading"; 
open $fileOUT, '>','basecount.out' or die "can't open file basecount.out for writing"; 

while (my $seq = <$fileIN>) { 

    next if $seq =~ /^>/; 
    $seq =~ s/\n//g; 
    say $seq; 

    my @dna = split //, $seq; 

    foreach my $element (@dna) { 
    $A++ if $element =~ m/A/; 
    $T++ if $element =~ m/T/; 
    $G++ if $element =~ m/G/; 
    $C++ if $element =~ m/C/; 
    } 

    say $fileOUT "A=$A"; 
    say $fileOUT "T=$T"; 
    say $fileOUT "G=$G"; 
    say $fileOUT "C=$C"; 
} 

close $fileIN; 
close $fileOUT; 

使用3語句打開也建議(以及良好的死亡警告以及)。

編輯: 我在這裏使用use feature 'say'因爲所有的打印結束了一個換行符。 sayprint完全相同,只是在最後添加換行符。

+0

感謝您的改進;特別是對'說';它的幫助。 – harsh

相關問題