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;
    }
}

发生的情况是，如果在一行数据的末尾和下一行的开头找到了主题，它将返回带有数据中换行符的主题。这种获取数据的方法以前效果很好。

样本结果：

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在脚本中的不同位置尝试了不同的方法，例如 or `s/\n//g。

4

1 回答 1

1

从您的描述中不清楚，但“在同一行也打印后半部分”听起来您的输出是覆盖在自身上的，因为它最后有一个回车符。

如果您在 Linux 系统上运行并且您只是chomp来自 Windows 的一行，则会发生这种情况。

您应该替换chomp为s/\s+\z//which 将删除所有尾随空格。而且因为回车和换行都算作“空白”，所以它将删除所有可能的终止字符。

顺便说一句，您误解了++运营商的目的。它还修改了它所应用的变量的内容，因此您只需要++$motif1_count等。您的代码按原样工作，因为运算符还返回递增变量的值，因此$motif1_count = ++$motif1_count首先递增变量，然后将其分配给自身。

此外，您\D在您的正则表达式中使用。你知道这匹配任何非数字字符吗？这似乎是一个非常模糊的分类有用。

于 2013-05-01T19:47:27.033 回答