0

Ensembl FASTA ファイルからタンパク質モチーフを見つけようとしています。シーケンス ID やシーケンス自体の取得など、スクリプトの大部分は完了しましたが、おかしな結果が返ってきました。

#!/usr/bin/perl
use strict;
use warnings;
use autodie;

my $motif1 = qr/(HE(\D)(\D)H(\D{18})E)/x;
my $motif2 = qr/(AMEN)/x;
my $input;
my $output;
my $count_total     = 0;
my $count_processed = 0;
my $total_run       = 0;
my $id;
my $seq;
my $motif1_count    = 0;
my $motif2_count    = 0;
my $motifboth_count = 0;

############################################################################################################################
# FILEHANDLING - INPUT/OUTPUT
# User input prompting and handling
print "**********************************************************\n";
print "Question 3\n";
print "**********************************************************\n";

#opens the user input file previously assigned to varible to new variable or kills script.
open my $fh, '<', "chr2.txt" || die "Error! Cannot open file:$!\n";

#Opens and creates output file previously assigned to variable to new variable or kills script
#open(RESULTS, '>', $output)||die "Error! Cannot create output file:$!\n";

# FILE and DATA PROCESSING
############################################################################################################################

while (<$fh>) {

    if (/^>(\S+)/) {
        $count_total = ++$count_total;    # Plus one to count
        find_motifs($id, $seq) if $seq;   # Passing to subroutine
        $id = substr($1, 0, 15);          # Taking only the first 16 characters for the id
        $seq = '';
    }
    else {
        chomp;
        $seq .= $_;
    }
}

print "Total proteins: $count_total \n";
print "Proteins with both motifs: $motifboth_count \n";
print "Proteins with motif 1: $motif1_count \n";
print "Proteins with motif 2: $motif2_count \n";

exit;

######################################################################################################################################
# SUBROUTINES
#
# Takes passed variables from special array
# Finds the position of motif within seq
# Checks for motif 1 presence and if found, checks for motif 2. If not found, prints motif 1 results
# If no motif 1, checks for motif 2

sub find_motifs {
    my ($id, $seq) = @_;
    if ($seq =~ $motif1) {
        my $motif_position = index $seq, $1;
        my $motif = $1;
        if ($seq =~ $motif2) {
            $motif1_count    = ++$motif1_count;
            $motif2_count    = ++$motif2_count;
            $motifboth_count = ++$motifboth_count;
            print "$id, $motif_position, \n$motif \n";
        }
        else {
            $motif1_count = ++$motif1_count;
            print "$id, $motif_position,\n $motif\n\n";
        }
    }
    elsif ($seq =~ $motif2) {
        $motif2_count = ++$motif2_count;
    }
}

何が起こっているかというと、モチーフがデータの 1 行の終わりと次の行の先頭にある場合、データ内の改行を含むモチーフが返されます。データを丸呑みするこの方法は、以前はうまく機能していました。

サンプル結果:

ENSG00000119013, 6,  HEHGHHKMELPDYRQWKIEGTPLE (CORRECT!)

ENSG00000142327, 123,  HEVAHSWFGNAVTNATWEEMWLSE (CORRECT!) 

ENSG00000151694, 410, **AECAPNEFGAEHDPDGL**

これが問題です。モチーフは一致しますが、前半の改行を返し、後半も同じ行に出力します (これは、より大きな問題の症状です - 改行を取り除く!)

Total proteins: 13653  
Proteins with both motifs: 1  
Proteins with motif 1: 12  
Proteins with motif 2: 22

@seq =~ s/\r//gまたは `s/\n//gなどのさまざまな方法を、スクリプト内のさまざまな場所で試しました。

4

1 に答える 1

1

あなたの説明からは明らかではありませんが、「後半も同じ行に出力します」は、最後に改行文字があるため、出力がそれ自体にオーバーレイされているように聞こえます。

これは、Linux システムで実行しchompていて、Windows からの行だけを使用している場合に発生します。

末尾の空白をすべて削除する which にchomp置き換える必要があります。s/\s+\z//また、キャリッジ リターンとラインフィードの両方が「空白」としてカウントされるため、すべての可能な終了文字が削除されます。

ところで、あなたは++演算子の目的を誤解しています。また、適用される変数の内容を変更するため、必要なものはすべて変更されます++$motif1_count。演算子はインクリメントされた変数の値も返すため、コードはそのまま機能し、$motif1_count = ++$motif1_count最初に変数をインクリメントしてからそれ自体に割り当てます。

また、\D正規表現で使用します。これが数字以外の文字と一致することを知っていますか? 非常に曖昧な分類が役立つようです。

于 2013-05-01T19:47:27.033 に答える