ファイル内で最短と最長のシーケンスを得る

Question

私は、複数のgenbankのようなエントリを含むファイルで最短と最長のシーケンスを取得しようとしています。ファイルの例：

LOCUS  NM_182854    2912 bp mRNA linear PRI 20-APR-2016 
DEFINITION Homo sapiens mRNA. 
ACCESSION NM_182854 
SOURCE  Homo sapiens (human) 
    ORGANISM Homo sapiens 
      Eukaryota; Metazoa; Chordata; Craniata; Vertebrata; Euteleostomi; 
      Mammalia; Eutheria; Euarchontoglires; Primates; Haplorrhini; 
      Catarrhini; Hominidae; Homo. 

ORIGIN  
     1 gggcgatcag aagcaggtca cacagcctgt ttcctgtttt caaacgggga acttagaaag 
     61 tggcagcccc tcggcttgtc gccggagctg agaaccaaga gctcgaaggg gccatatgac 
     // 

LOCUS  NM_001323410   6992 bp mRNA linear PRI 20-APR-2016 
DEFINITION Homo sapiens mRNA. 
ACCESSION NM_001323410 
SOURCE  Homo sapiens (human) 
    ORGANISM Homo sapiens 
      Eukaryota; Metazoa; Chordata; Craniata; Vertebrata; Euteleostomi; 
      Mammalia; Eutheria; Euarchontoglires; Primates; Haplorrhini; 
      Catarrhini; Hominidae; Homo. 

ORIGIN  
     1 actacttccg gcttccccgc cccgccccgt ccccgggcgt ctccattttg gtctcaggtg 
     61 tggactcggc aagaaccagc gcaagaggga agcagagtta tagctacccc ggc 
     // 

私はこれまでのところ私のコード、受付番号を印刷するには、最短のシーケンスからの生物の種類と最も長いシーケンス

が欲しい：

#!/usr/bin/perl 

use strict; 
use warnings; 

print "enter file path\n"; 

while (my $line = <>){ 
    chomp $line; 
    my @record = ($line); 

    foreach my $file(@record){ 
    open(IN, "$file") or die "\n error opening file \n;/\n"; 

    $/="//"; 

    while (my $line = <IN>){ 
     my @gb_seq = split ("ORIGIN", $line); 
     my $definition = $gb_seq[0]; 
     my $sequence = $gb_seq[1]; 

     $definition =~ m/ORGANISM[\s\t]+(.+)[\n\s\t]+/; 
     my $organism = $1; 

     if ($definition =~ m/ACCESSION[\s\t]+(\D\D_\d\d\d\d\d\d(\d*))[\n\s\t]+/){ 
     my $accession = $1; 

      $sequence =~ s/\d//g; 
      $sequence =~ s/[\n\s\t]//g; 
      my $size = length($sequence); 
      my @sorted_keys = sort { $a <=> $b } keys my %size; 
      my $shortest = $sorted_keys[0]; 
      my $longest = $sorted_keys[-1]; 

      print "this is the shortest: $accession $organism size: $shortest\n"; 
      print "this is the longest: $accession $organism size: $longest\n"; 
    } 
    }}} 
    exit; 

を

私は長さをハッシュに入れて、最短と最長を得ることを考えましたが、何かが間違っています。これらのエラーが表示されます：

Use of uninitialized value $organism in concatenation (.) or string at test.pl line 39, <IN> chunk 1 
Use of uninitialized value $shortest in concatenation (.) or string at test.pl line 39, <IN> chunk 1. 
Use of uninitialized value $longest in concatenation (.) or string at test.pl line 40, <IN> chunk 1. 

私はどのような部分を変更しなければなりませんか？ありがとう

Answer 1

最長および最短シーケンスの2つのデータ（アクセスと生物）が必要であると述べています。つまり、ハッシュ値には2つの要素を格納する必要があります。レコードセパレータとして '//'を使用すると、各レコードの末尾に '//'が表示されます。だから、空白と数字を並べ替えると、最後に「//」が残っています。デバッガを使ってコードを実行したとき、私はこの長さが2であることを発見しました。

他の物事のカップル：あなたは$definitionを掘るとき、あなたは成功した試合を推定

readabillityのための空白を含めることができるように

1. regexsを使用して、、、/xを「拡張モード」を使う - よりよいですあなたの正規表現をテストし、マッチに割り当てる、ミスマッチで死ぬ
2. 長さをハッシュに格納するのではなく（シーケンス自体を失う）、シーケンスを保存して後で長さを計算することもできます。
3. それが最短と最長を計算し、resutsの印刷を行うにはすべてのものは、ループの外に移動する必要があるいくつかのライン
4. が含まれているように私は$chunkに変数$lineの名前を変更しました。その代わりに、単にハッシュへのエントリを作成する必要があります。前述のように、ハッシュ値は、アクセスと生物の2つの値を持つ配列である必要があります。
5. 1つのコマンドでシーケンスから数字を削除し、別のコマンドでシーケンスの空白を削除すると、それらの両方をtogeatherにすることができます。私たちがそれをしている間、レコードの最後に/を取り除くかもしれません。

上記の改造を考えればわかります。

use v5.14; 
use warnings; 

print "Enter file path: "; 
chomp(my $filename = <>); 
open(IN, $filename) or die "\n error opening file \n;/\n"; 

$/ = "//" ; 

my %organisms ; 
while (my $chunk = <IN>) { 
    next if $chunk =~ /^\s*\n\s*$/ ; 
    my ($definition , $sequence) = split "ORIGIN", $chunk ; 

    my $organism ; 
    $definition =~ m/ ORGANISM [\s\t]+ (.+) [\n\s\t]+ /x 
     ? $organism = $1 
     : die "Couldnt find ORGANISM line" ; 

    my $accession ; 
    $definition =~ m/ ACCESSION [\s\t]+ (\D\D _ \d{6} (\d*)) [\n\s\t]+ /x 
     ? $accession = $1 
     : die "Cant find ACCESSION line" ; 

    $sequence =~ s/[\d\n\s\t\/]//g; 
    $organisms{ $sequence } = [ $accession , $organism ] ; 
} 


my @sorted_keys = sort { length $a <=> length $b } keys %organisms ; 
my $shortest = $sorted_keys[0]; 
my $longest = $sorted_keys[-1]; 

say "this is the shortest: ", $organisms{$shortest}->[0], 
         ", ", $organisms{$shortest}->[1], 
        " size: ", length $shortest, "\n", 
       " sequence: ", $shortest ; 

say "this is the longest: ", $organisms{$longest}->[0], 
         ", ", $organisms{$longest}->[1], 
        " size: ", length $longest, "\n", 
       " sequence: ", $longest ; 

exit; 

データを実行すると、それが生成されます。

$ ./sequence.pl 
Enter file path: data.txt 
this is the shortest: NM_001323410, Homo sapiens size: 113 
sequence: actacttccggcttccccgccccgccccgtccccgggcgtctccattttggtctcaggtgtggactcggcaagaaccagcgcaagagggaagcagagttatagctaccccggc 
this is the longest: NM_182854, Homo sapiens size: 120 
sequence: gggcgatcagaagcaggtcacacagcctgtttcctgttttcaaacggggaacttagaaagtggcagcccctcggcttgtcgccggagctgagaaccaagagctcgaaggggccatatgac 

UPDATE 上記のコードの問題は、同じシーケンスが2つのチャンクに表示されている場合、そのデータは、ハッシュに上書きされて失われようとしているということです。以下は、問題を提起する配列の配列にデータを格納する更新されたバージョンです。

use v5.14; 
use warnings; 

print "Enter file path: "; 
chomp(my $filename = <>); 
open(IN, $filename) or die "\n error opening file \n;/\n"; 

$/ = "//" ; 

my @organisms ; 
while (my $chunk = <IN>) { 
    next if $chunk =~ /^\s*\n\s*$/ ; 
    my ($definition , $sequence) = split "ORIGIN", $chunk ; 

    my $organism ; 
    $definition =~ m/ ORGANISM [\s\t]+ (.+) [\n\s\t]+ /x 
     ? $organism = $1 
     : die "Couldnt find ORGANISM line" ; 

    my $accession ; 
    $definition =~ m/ ACCESSION [\s\t]+ (\D\D _ \d{6} (\d*)) [\n\s\t]+ /x 
     ? $accession = $1 
     : die "Cant find ACCESSION line" ; 

    $sequence =~ s/[\d\n\s\t\/]//g; 
    push @organisms, [$organism , $accession , $sequence] ; 
} 


my @sorted_organisms = sort { length $a->[2] <=> length $b->[2] } @organisms ; 

my ($organism , $accession , $sequence) = @{ $sorted_organisms[0] }; 
say "this is the shortest: $accession, $organism, size: ", 
    length $sequence, "\n", " sequence: ", $sequence ; 

($organism , $accession , $sequence) = @{ $sorted_organisms[-1] }; 
say "this is the longest: $accession, $organism, size: ", 
    length $sequence, "\n", " sequence: ", $sequence ; 

exit; 

Answer 2

我々は、彼らが所属するレコードを識別することができながら、極端な長さのエントリを検索する必要があります。それはまったく同じ出力を生成します。 //でレコードを読むことはまた素晴らしいアイデアです。しかし、各レコードは文字列であり、シーケンスを直接抜き出すことは、最初に行に分割するよりも困難です。したがって、必要なすべてのための明確なマーカーがあることを考えると、ラインごとに行こうとするかもしれません。

データ構造の選択は重要であり、目的によって異なります。で動作するように簡単になるように、ここで私は要素

%block = ('accession' => { 'type' => type, 'sequence' => sequence }, ...) 

「配列」（でこれを組織することではなく、することで、データを大幅に支援されるだろうが読み込まれたら実行するために、検索とハッシュに、データを整理します'アクセッション'）、それはそれを扱うことを非常に困難にするでしょう。私はこれがより多くのために使用されてしまうことがあり、スピードのわずかな損失は重要ではないと推測します。ここでの唯一の目的が、最適なパフォーマンスで特定の質問に答えることであれば、他のアプローチがより適しています。コメントはコードに従います。

use warnings; 
use strict; 
use feature qw(say); 

my $file = 'data_seqs.txt'; 
open my $fh, '<', $file or die "Can't open $file -- $!"; 

# Hash, helper variables, flag (inside a sequence?), sequence-end marker 
my (%block, $accession, $sequence); 
my $is_seq = 0; 
my $end_marker = qr(\s*//); # marks end of sequence: // 

while (my $line = <$fh>) 
{ 
    chomp($line); 
    next if $line =~ /^\s*$/;  # skip empty lines 

    if ($line =~ /$end_marker/) { # done with the sequence 
     $is_seq = 0; 
     $sequence = ''; 
     next; 
    } 

    if ($line =~ /^\s*ACCESSION\s*(\w+)/) { 
     $accession = $1; 
    } 
    elsif ($line =~ /^\s*ORGANISM\s*(.+)/) { 
     $block{$accession}{'type'} = $1; 
    } 
    elsif ($line =~ /^\s*ORIGIN/) { # start sequence on next line 
     $is_seq = 1; 
    } 
    elsif ($is_seq) {    # read (and add to) sequence 
     if ($line =~ /^\s*\d+\s*(.*)/) { 
      $block{$accession}{'sequence'} .= $1; 
     } 
     else { warn "Not sequence? Line: $line " } 
    } 
} 

# Identify keys for max and min lenght. Initialize with any keys 
my ($max, $min) = keys %block; 

foreach my $acc (keys %block) 
{ 
    my $current_len = length($block{$acc}{'sequence'}); 
    if ($current_len > length($block{$max}{'sequence'})) { 
     $max = $acc; 
    } 
    if ($current_len < length($block{$min}{'sequence'})) { 
     $min = $acc; 
    } 
} 

say "Maximum length sequence: ACCESSION: $max, ORGANISM: " . $block{$max}{'type'}; 
say "Minimum length sequence: ACCESSION: $min, ORGANISM: " . $block{$min}{'type'}; 

use Data::Dumper; 
print Dumper(\%block); 

この印刷物（ダンパのプリントアウトは省略）

 
Maximum length sequence: ACCESSION: NM_182854, ORGANISM Homo sapiens 
Minimum length sequence: ACCESSION: NM_001323410, ORGANISM Homo sapiens 

一つの一般的なアプローチから言って、最初のライブラリを使用し、その後、逆引き参照ハッシュを構築することです効率

の検索のコメントList::Utils、最大値と最小値を求め、それらがどこに属しているか調べます。このためには、検索ハッシュを構築する必要があり、ライブラリを2回使用しますが、上記のように手作業で検索すると、構造を1回通過するだけでなく、シンプルになります。もう1つの選択肢は、ハッシュのトップレベルキーをシーケンスにしてからmaxとminを直接見つけることです。しかし、このようなハッシュは作業するのがかなり難しくなります。

さらに別のアプローチは、おそらくアレイに基づいて、この特定の情報をより効率的に取得できる構造にデータを編成することです。

しかし、効率の向上は、利便性の大きな損失を正当化するようには見えません。速度が問題であることが判明した場合、これが考慮されるべきである。

複数のファイルを扱う必要がある場合は、ループをに変更し、コマンドラインで送信してください。それらのすべてからのすべての行が1行ずつ読み込まれ、コードは同じままになります。

私はいくつかの用語を誤解している可能性があります。私は "シーケンス"から空のスペースを削除しないで、最初の行の単語を "タイプ"のためだけに使用します。これらは調整しやすいです、私に知らせてください。