1

私は少し乱雑なコードで試したばかりの新しい perl です。

猫input1.txt

  ##gff-version 2
  ##source-version geneious 5.6.4
  Xm_ABL1 Geneious  CDS 1   168 .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 169 334 .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 335 628 .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 629 901 .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 902 985 .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 986 1165    .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 1166    1350    .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  CDS 1351    1504    .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
  Xm_ABL1 Geneious  BLAST Hit   169 334 .   +   .   
  Xm_ABL1 Geneious  extracted region    1   168 .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="351297 -> 351464"
  Xm_ABL1 Geneious  extracted region    169 334 .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="371785 -> 371950"
  Xm_ABL1 Geneious  extracted region    335 628 .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="372554 -> 372847"
  Xm_ABL1 Geneious  extracted region    629 901 .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="374760 -> 375032"
  Xm_ABL1 Geneious  extracted region    902 985 .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="375230 -> 375313"
  Xm_ABL1 Geneious  extracted region    986 1165    .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="375992 -> 376171"
  Xm_ABL1 Geneious  extracted region    1166    1350    .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="376575 -> 376759"
  Xm_ABL1 Geneious  extracted region    1351    1504    .   +   .   Name=Extracted region from gi|371443098|gb|JH556762.1|;Extracted interval="376914 -> 377067"

入力ファイルに (->) 進む矢印が含まれている場合。 if($array[7]=~/.*interval=\"\d+ -> \d+\"$/gm){ $array[5]= "+"; }

猫出力1.txt

gi_371443098_gb_JH556762.1  gene    351297  377067  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 351297  351464  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 371785  371950  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 372554  372847  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 374760  375032  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 375230  375313  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 375992  376171  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 376575  376759  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 376914  377067  .   +   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
###

cat output1.txt 入力ファイルに (<-) 逆矢印が含まれている場合。if($array[7]=~/.*interval=\"\d+ <- \d+\"$/gm){ $array[5]="-"; }

gi_371443098_gb_JH556762.1  gene    351297  377067  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 351297  351464  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 371785  371950  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 372554  372847  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 374760  375032  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 375230  375313  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 375992  376171  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 376575  376759  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
gi_371443098_gb_JH556762.1  CDS 376914  377067  .   -   .   Name=Xm_ABL1;created by=User;modified by=User;ID=w0IVHutPuN4H4FVDCg4sFVRaJjQ.1340919460469.4
###

私は初心者なので、少し厄介なコードで試しました。

#usr/bin/perl;
use strict;
open(FH,"$ARGV[0]");

while(<FH>){
chomp $_; 
my @array=split("\t");
my $key="$array[2]-$array[0]-$array[1]-$array[2]-$array[3]"; 
if($array[1] eq "CDS"){
$cds_cnt{$key}++;
$cds{$key}="$array[4]\t$array[5]\t$array[6]\t$array[7]";
    }
if($array[1] eq "extracted region"){ 
    (my $pos1,my $pos2)=($array[7]=~/.*interval=\"(\d+) -> (\d+)\"$/gm);
$extract_cnt{$key}++;
$extract{$key}="$pos1\t$pos2";
            }
}

 foreach $i (  sort {$a<=>$b} keys %cds){ 
 my $a=$i; #print "$i\n";
 $a=~s/CDS/extracted region/g;
 if($cds_cnt{$i} == $extract_cnt{$a}){
 #print "$i\t$cds{$i}\n$a\t$extract{$a}\n"; 
 my @array=split /\-/,$i;
 my @pos=split "\t",$extract{$a};
 print "$array[1]\t$array[2]\t$pos[0]\t$pos[1]\t$cds{$i}\n";
   }
  }
 print "###";

アップデート

コードで変更する必要があるもの

1.抽出された領域の行から値を取得するには (つまり、array[7]=/gi|371443098|gb|JH556762.1|/)、任意の値にすることができます。アンダースコアを追加し (つまり、gi_371443098_gb_JH556762.1)、出力します。示されているように、output1.txt の array[0] にあります。

2. 印刷中に最初の行として新しい行を追加 (gi_371443098_gb_JH556762.1 遺伝子)、列 3 で CDS の開始値 (つまり 351297) を取得し、列 4 で CDS の終了値 (つまり 377067) を取得し、最初の行のように印刷します。 output1.txt に表示

3. If /extracted region/ block-all rows for.egExtracted interval="351297 -> 351464" (つまり、前方矢印) 配列 [5] を "+" 記号として出力に遺伝子ヘッダーを含めて出力します。if egExtracted interval="351297 <- 351464" (逆矢印) 出力に遺伝子ヘッダーを含む "-" 記号として array[5] を出力します。

4

1 に答える 1

2

達成しようとしているのは、CDSというラベルの付いた行の詳細を、抽出された領域というラベルの付いた一致する行とマージし、マージされた結果を、いくつかの最小値と最大値に基づいて、名前。あれは正しいですか?

$ array [0](Xm_ABL1 Geneious)と$ array [2](169、335など)と呼ぶものは、それらを結合するのに十分であると仮定しますが、それはあなたの例ではあまり明確ではありません。

最初の質問は単なる正規表現です。これは一般的なコツを持っていると思います。問題は、データのキャプチャをどのように行ったかにあると思います。

2番目に要求することを行うには、最初のパスでhi値とlo値をキャプチャし、それらを保存します。

私は完全なソリューションを書くことを計画していませんでしたが、ここにあります...

use strict;
use warnings;

my $metadata = {}; # hashref to store CDS info in..
my $group = {}; # hashref to store summary/detail in..
my $arrow = { "->" => '+', "<-" => '-' }; # decode arrow to pos/neg

open(FH,"$ARGV[0]");
while(<FH>){
    chomp;
    next if /^#/;
    my @array=split("\t");
    my $key = join(":", $array[0], $array[2]);
    if ($array[1] =~ /CDS/){
        $metadata->{$key} = $array[7];
    }
    if ($array[1] =~ /extracted region/){
        #assert CDS already processed..
        die "No CDS record for $key!\n" unless $metadata->{$key};
        (my $label = $array[7]) =~ s/.*region from (.*)\|;.*/$1/;
        $label =~ s/\|/_/g;
        $group->{$label} ||= { #seed summary if not exists
                pos1 => 1e10,
                pos2 => 0,
                metadata => $metadata->{$key},
                sequences => [],
        };
        (my $pos1, my $arr, my $pos2) = ($array[7]=~/.*interval=\"(\d+) (<?->?) (\d+)\"$/gm);
        # capture hi/lo values for group
        $group->{$label}->{pos1} = $pos1 if $pos1 < $group->{$label}->{pos1};
        $group->{$label}->{pos2} = $pos2 if $pos2 > $group->{$label}->{pos2};
        # push this sequence onto the group's array
        push(@{ $group->{$label}->{sequences} }, [ $pos1, $pos2, $arrow->{$arr} ]);
    }
}
for my $gene (sort keys %{ $group }){
    #write out header
    printf "%s\t%s\t%d\t%d\t.\t%s\t.\t%s\n",
        $gene, 'gene',
        $group->{$gene}->{pos1}, $group->{$gene}->{pos2},
        $group->{$gene}->{sequences}->[0]->[2],
        $group->{$gene}->{metadata};
    foreach my $sequence ( @{ $group->{$gene}->{sequences} } ){
        # write out details
        printf "%s\t%s\t%d\t%d\t.\t%s\t.\t%s\n",
            $gene, 'CDS',
            $sequence->[0], $sequence->[1], $sequence->[2],
            $group->{$gene}->{metadata};
    }
}
print "###\n";

私はそれが意味をなすのに十分にコメントされていることを望みます。このように書かれているので、6か月が経過した後に戻って変更する必要がある場合は、保守がはるかに簡単になります

アップデート

4番目のコメントに続いてこのコードを変更しました。$ array [7] regexpブロックは、$ array [5]の値をキャプチャし、それをシーケンスarrayrefに格納するようになりました。

更新#2

  • ->デコードおよび<-必要なシンボルへの$arrowhashrefの使用に注意してください。(6行目)
  • 遺伝子ヘッダーは、最初のシーケンスの3番目のフィールドの値に基づいて+または-を示します。遺伝子のすべての配列が同じ方向を持っているという仮定があります。(最後から11行)

私たちは今、Q&Aハイウェイからフリーソフトウェア開発局に迷い込んだと思います。私が書いたのは複雑なコードではなく、コメントと構造があります。それはあなたがそれの論理を理解するようになる時です。そして、私の答えに賛成票を投じてください。

于 2012-07-12T02:11:11.217 に答える