5

私はperlを初めて使用し、rtfファイルに保存されているDNA配列に対する基本的な文字列操作を実行したいと思います。

基本的に、私のファイルは次のようになります（ファイルはFASTA形式です）：

>LM1
AAGTCTGACGGAGCAACGCCGCGTGTATGAAGAAGGTTTTCGGATCGTAA
AGTACTGTCCGTTAGAGAAGAACAAGGATAAGAGTAACTGCTTGTCCCTT
GACGGTATCTAACCAGAAAGCCACGGCTAACTACGTGCCAGCAGCCGCGG
TAATACGTAGGTGGCAAGCGTTGTCCGGATTTATTGGGCGTAAAGCGCGC
GCAGGCGGTCTTTTAAGTCTGATGTGAAAGCCCCCGGCTTAACCGGGGAG
GGTCATTGGAAACTGGAAGACTGGAGTGCAGAAGAGGAGAGTGGAATTCC
ACGTGTAGCGGTGAAATGCGTAGATATGTGGAGGAACACCAGTGGCGAAG
GCGACTCTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCA
AACAGGATTAGATACCCTGGTAGTCCACGCCGT

私がやりたいのは、ファイルを読み込んでヘッダー（ヘッダーは> LM1）を印刷し、次のDNA配列GTGCCAGCAGCCGCと一致させてから、前のDNA配列を印刷することです。
したがって、私の出力は次のようになります。

>LM1 
AAGTCTGACGGAGCAACGCCGCGTGTATGAAGAAGGTTTTCGGATCGTAA
AGTACTGTCCGTTAGAGAAGAACAAGGATAAGAGTAACTGCTTGTCCCTT
GACGGTATCTAACCAGAAAGCCACGGCTAACTAC

私は次のプログラムを書きました：

#!/usr/bin/perl

use strict; use warnings;

open(FASTA, "<seq_V3_V6_130227.rtf") or die "The file could not be found.\n";

while(<FASTA>) {
    chomp($_);
    if ($_ =~  m/^>/ ) {
        my $header = $_;
        print "$header\n";
    }

    my $dna = <FASTA>;
    if ($dna =~ /(.*?)GTGCCAGCAGCCGC/) {
        print "$dna";
    }

}
close(FASTA);

問題は、私のプログラムがファイルを1行ずつ読み取り、受信する出力が次のようになることです。

>LM1
GACGGTATCTAACCAGAAAGCCACGGCTAACTAC

基本的に、DNAシーケンス全体を$ dna変数に割り当てる方法がわかりません。また、最終的には、DNAシーケンスを1行ずつ読み取らないようにする方法もわかりません。また、次の警告が表示されます。stacked.pl行14、行1113のパターン一致（m //）で初期化されていない値$dnaを使用しています。

誰かが私にもっと良いコードを書くのを手伝ってくれるか、正しい方向に私を向けることができれば、それは大いにありがたいです。

4

4 に答える 4

3

pos関数の使用：

use strict;
use warnings;

my $dna = "";
my $seq = "GTGCCAGCAGCCGC";
while (<DATA>) {
  if (/^>/) {
    print;
  } else {
    if (/^[AGCT]/) {
      $dna .= $_;
    }
  }

}

if ($dna =~ /$seq/g) {
  print substr($dna, 0, pos($dna) - length($seq)), "\n";
}

__DATA__
>LM1

AAGTCTGACGGAGCAACGCCGCGTGTATGAAGAAGGTTTTCGGATCGTAA
AGTACTGTCCGTTAGAGAAGAACAAGGATAAGAGTAACTGCTTGTCCCTT
GACGGTATCTAACCAGAAAGCCACGGCTAACTACGTGCCAGCAGCCGCGG
TAATACGTAGGTGGCAAGCGTTGTCCGGATTTATTGGGCGTAAAGCGCGC
GCAGGCGGTCTTTTAAGTCTGATGTGAAAGCCCCCGGCTTAACCGGGGAG
GGTCATTGGAAACTGGAAGACTGGAGTGCAGAAGAGGAGAGTGGAATTCC
ACGTGTAGCGGTGAAATGCGTAGATATGTGGAGGAACACCAGTGGCGAAG
GCGACTCTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCA
AACAGGATTAGATACCCTGGTAGTCCACGCCGT

次のように、複数のエントリを持つファイルを処理できます。

while (<DATA>) {
  if (/^>/) {
    if ($dna =~ /$seq/g) {
      print substr($dna, 0, pos($dna) - length($seq)), "\n";
      $dna = ""; 
    }   
    print;
  } elsif (/^[AGCT]/) {
    $dna .= $_; 
  }   
}

if ($dna && $dna =~ /$seq/g) {
  print substr($dna, 0, pos($dna) - length($seq)), "\n";
}
于 2013-03-04T21:08:06.407 に答える
2

whileステートメントは、ファイルの終わりまで読み取ります。つまり、ループが繰り返されるたびに、$_がの次の行になり<FASTA>ます。だから$dna = <FASTA>、あなたが思っていることをしていません。おそらくあなたが望む以上に読んでいます。

while(<FASTA>) { #Reads a line here
  chomp($_);
  if ($_ =~  m/^>/ ) {
    my $header = $_;
    print "$header\n";
  }
  $dna = <FASTA> # reads another line here - Causes skips over every other line
}

ここで、シーケンスをに読み込む必要があります$dna。elseステートメントを使用してwhileループを更新できます。したがって、見出しの場合は印刷し、そうでない場合はに追加し$dnaます。

while(<FASTA>) {
  chomp($_);
  if ($_ =~  m/^>/ ) {
    # It is a header line, so print it
    my $header = $_;
    print "$header\n";
  } else {
    # if it is not a header line, add to your dna sequence.
    $dna .= $_;
  }
}

ループの後、正規表現を実行できます。

注：このソリューションは、fastaファイルにシーケンスが1つしかないことを前提としています。複数ある場合、$dna変数にはすべてのシーケンスが1つになります。

編集：複数のシーケンスを処理する簡単な方法を追加

my $dna = "";
while(<FASTA>) {
  chomp($_);
  if ($_ =~  m/^>/ ) {

    # Does $dna match the regex?
    if ($dna =~ /(.*?)GTGCCAGCAGCCGC/) {
      print "$1\n";
    }

    # Reset the sequence
    $dna = "";

    # It is a header line, so print it
    my $header = $_;
    print "$header\n";

  } else {
    # if it is not a header line, add to your dna sequence.
    $dna .= $_;
  }
}

# Check the last sequence
if ($dna =~ /(.*?)GTGCCAGCAGCCGC/) {
  print "$1\n";
}
于 2013-03-04T21:10:52.633 に答える
2

私はBioSeqIO（およびBioPerlディストリビューションのBioSeqのtruncメソッド）を使用して解決策を考え出しました。正規表現を使用するのではなく、インデックスを使用してサブシーケンスを検索しました。

このソリューションは、サブシーケンスが見つからなかった場合、またはサブシーケンスが最初の位置で始まった場合（したがって、前の文字がない場合）、id （行は>で始まります）を出力しません。

#!/usr/bin/perl
use strict;
use warnings;
use Bio::SeqIO;

my $in  = Bio::SeqIO->new( -file   => "fasta_junk.fasta" ,
                           -format => 'fasta');

my $out = Bio::SeqIO->new( -file   => '>test.dat',
                           -format => 'fasta');

my $lookup = 'GTGCCAGCAGCCGC';

while ( my $seq = $in->next_seq() ) {
    my $pos = index $seq->seq, $lookup;

    # if $pos != -1, ($lookup not found),
    # or $pos != 0, (found $lookup at first position, thus
    #   no preceding characters).
    if ($pos > 0) {
        my $trunc = $seq->trunc(1,$pos);
        $out->write_seq($trunc);
    }
}

__END__
*** fasta_junk.fasta
>LM1
AAGTCTGACGGAGCAACGCCGCGTGTATGAAGAAGGTTTTCGGATCGTAA
AGTACTGTCCGTTAGAGAAGAACAAGGATAAGAGTAACTGCTTGTCCCTT
GACGGTATCTAACCAGAAAGCCACGGCTAACTACGTGCCAGCAGCCGCGG
TAATACGTAGGTGGCAAGCGTTGTCCGGATTTATTGGGCGTAAAGCGCGC
GCAGGCGGTCTTTTAAGTCTGATGTGAAAGCCCCCGGCTTAACCGGGGAG
GGTCATTGGAAACTGGAAGACTGGAGTGCAGAAGAGGAGAGTGGAATTCC
ACGTGTAGCGGTGAAATGCGTAGATATGTGGAGGAACACCAGTGGCGAAG
GCGACTCTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCA
AACAGGATTAGATACCCTGGTAGTCCACGCCGT

*** contents of test.dat
>LM1
AAGTCTGACGGAGCAACGCCGCGTGTATGAAGAAGGTTTTCGGATCGTAAAGTACTGTCC
GTTAGAGAAGAACAAGGATAAGAGTAACTGCTTGTCCCTTGACGGTATCTAACCAGAAAG
CCACGGCTAACTAC
于 2013-03-05T01:57:18.683 に答える
0

ファイル全体をメモリに読み込み、正規表現を探します

while(<FASTA>) {
    chomp($_);
    if ($_ =~  m/^>/ ) {
        my $header = $_;
        print "$header\n";
    } else {
    $dna .= $_;
    }
}
if ($dna =~ /(.*?)GTGCCAGCAGCCGC/) {
    print $1;
}
于 2013-03-04T21:05:23.303 に答える