用引用、哈希表和子类模拟限制酶的Perl程序

Question

我是Perl入门班的学生。我在寻找关于如何处理作业的建议。我的教授鼓励论坛。任务是：

编写了一个Perl程序，它将从命令行获取两个文件，一个酶文件和一个DNA文件。用限制酶读取文件，并将限制酶应用到DNA文件中。

输出的DNA片段将按照它们在dna文件中出现的顺序排列。输出文件的名称应该通过将限制酶的名称附加到DNA文件的名称，并在它们之间加上下划线来构造。

例如，如果酶是EcoRI，而DNA文件名为BC161026，则输出文件应该命名为BC161026_EcoRI。

我的方法是创建一个主程序和两个子程序，如下所示：

主语:不知道怎么把我的潜艇绑在一起？

子程序$DNA:获取DNA文件并删除任何新行以生成单个字符串。

子程序酶:从命令行读取和存储酶文件中的行，以将酶首字母缩略词与切割位置分开的方式解析该文件。将裁剪作为正则表达式存储在哈希表中，在哈希表中存储首字母缩略词的名称

关于酶文件格式的

注:酶文件遵循一种称为Staden的格式。示例：

AatI/AGG'CCT//

AatII/GACGT'C//

AbsI/CC'TCGAGG//

酶首字母缩略词由第一个斜杠之前的字符组成(AatI，在第一个示例中)。识别序列是在第一个斜杠和第二个斜杠之间的所有内容(在第一个示例中，同样是AGG‘’CCT)。切割点在识别序列中用撇号表示，酶中dna的标准缩写如下：

R=G或A=非A (C或G或T)等

除了一个主块的推荐之外，你有没有看到我遗漏的部分呢？你能推荐一些你认为在修补这个程序中有用的特定工具吗？

输入酶示例：TryII/RRR'TTT//

要读取的示例字符串：CCCCCCGGGTTTCCCCCCCCCCCCAAATTTCCCCCCCCCCCCAGATTTCCCCCCCCCCGAGTTTCCCCC

产出应是：

CCCCCCGGG

TTTCCCCCCCCCCCCAAA

TTTCCCCCCCCCCCCAGA

TTTCCCCCCCCCCGAG

TTTCCCCC

Answer 1

好吧，我知道我不应该只做你的家庭作业，但这个有一些有趣的技巧，所以我玩它。从中吸取教训，而不仅仅是抄袭。我没有很好的评论，所以如果有什么你不明白，请问。这里面有一些小魔术，如果你没有在课堂上讨论过，你的教授会知道的，所以一定要理解。

#!/usr/bin/env perl

use strict;
use warnings;

use Getopt::Long;

my ($enzyme_file, $dna_file);
my $write_output = 0;
my $verbose = 0;
my $help = 0;
GetOptions(
  'enzyme=s' => \$enzyme_file,
  'dna=s' => \$dna_file,
  'output' => \$write_output,
  'verbose' => \$verbose,
  'help' => \$help
);

$help = 1 unless ($dna_file && $enzyme_file);
help() if $help; # exits

# 'Main'
my $dna = getDNA($dna_file);
my %enzymes = %{ getEnzymes($enzyme_file) }; # A function cannot return a hash, so return a hashref and then store the referenced hash
foreach my $enzyme (keys %enzymes) {
  print "Applying enzyme " . $enzyme . " gives:\n";
  my $dna_holder = $dna;
  my ($precut, $postcut) = ($enzymes{$enzyme}{'precut'}, $enzymes{$enzyme}{'postcut'});

  my $R = qr/[GA]/;
  my $B = qr/[CGT]/;

  $precut =~ s/R/${R}/g;
  $precut =~ s/B/${B}/g;
  $postcut =~ s/R/${R}/g;
  $postcut =~ s/B/${B}/g;
  print "\tPre-Cut pattern: " . $precut . "\n" if $verbose;
  print "\tPost-Cut pattern: " . $postcut . "\n" if $verbose;

  #while(1){
  #  if ($dna_holder =~ s/(.*${precut})(${postcut}.*)/$1/ ) {
  #    print "\tFound section:" . $2 . "\n" if $verbose;
  #    print "\tRemaining DNA: " . $1 . "\n" if $verbose;
  #    unshift @{ $enzymes{$enzyme}{'cut_dna'} }, $2;
  #  } else {
  #    unshift @{ $enzymes{$enzyme}{'cut_dna'} }, $dna_holder;
  #    print "\tNo more cuts.\n" if $verbose;
  #    print "\t" . $_ . "\n" for @{ $enzymes{$enzyme}{'cut_dna'} };
  #    last;
  #  }
  #}
  unless ($dna_holder =~ s/(${precut})(${postcut})/$1'$2/g) {
    print "\tHas no effect on given strand\n" if $verbose;
  }
  @{ $enzymes{$enzyme}{'cut_dna'} } = split(/'/, $dna_holder);
  print "\t$_\n" for @{ $enzymes{$enzyme}{'cut_dna'} };

  writeOutput($dna_file, $enzyme, $enzymes{$enzyme}{'cut_dna'}) if $write_output; #Note that $enzymes{$enzyme}{'cut_dna'} is an arrayref already
  print "\n";
}

sub getDNA {
  my ($dna_file) = @_;

  open(my $dna_handle, '<', $dna_file) or die "Cannot open file $dna_file";
  my @dna_array = <$dna_handle>;
  chomp(@dna_array);

  my $dna = join('', @dna_array);

  print "Using DNA:\n" . $dna . "\n\n" if $verbose;
  return $dna;
}

sub getEnzymes {
  my ($enzyme_file) = @_;
  my %enzymes;

  open(my $enzyme_handle, '<', $enzyme_file) or die "Cannot open file $enzyme_file";
  while(<$enzyme_handle>) {
    chomp;
    if(m{([^/]*)/([^']*)'([^/]*)//}) {
      print "Found Enzyme " . $1 . ":\n\tPre-cut: " . $2 . "\n\tPost-cut: " . $3 . "\n" if $verbose;
      $enzymes{$1} = {
        precut => $2,
        postcut => $3,
        cut_dna => [] #Added to show the empty array that will hold the cut DNA sections
      };
    }
  }

  print "\n" if $verbose;
  return \%enzymes;
}

sub writeOutput {

  my ($dna_file, $enzyme, $cut_dna_ref) = @_;

  my $outfile = $dna_file . '_' . $enzyme;
  print "\tSaving data to $outfile\n" if $verbose; 
  open(my $outfile_handle, '>', $outfile) or die "Cannot open $outfile for writing";

  print $outfile_handle $_ . "\n" for @{ $cut_dna_ref };
}

sub help {

  my $filename = (split('/', $0))[-1];

  my $enzyme_text = <<'END';
AatI/AGG'CCT//
AatII/GACGT'C//
AbsI/CC'TCGAGG//
TryII/RRR'TTT//
Test/AAA'TTT//
END

  my $dna_text = <<'END';
CCCCCCGGGTTTCCCCCCC
CCCCCAAATTTCCCCCCCCCCCCAGATTTC
CCCCCCCCCGAGTTTCCCCC
END

  print <<END;
Usage: 
    $filename --enzyme (-e) <enzyme-filename> --dna (-d) <dna-filename> [options] (files may come in either order)
    $filename -h    (shows this help)

Options: 
    --verbose (-v)  print additional (debugging) information
    --output (-o)   output final data to files


Files:
The DNA file contains one DNA string which may be broken over many lines. E.G.:

$dna_text

The enzymes file constains enzyme definitions, one per line. E.G.:

$enzyme_text
END

exit;
}

编辑:显式地添加了cut_dna初始化，因为这是每个酶的最终结果保持器，所以我认为更清楚地看到它会更好。

编辑2:增加输出例程，呼叫，标志和相应的帮助。

编辑3:改变主例程，在去除循环的同时结合最好的canavanin方法。现在，它是一个全局替换，添加临时剪切标记(')，然后在剪切标记上拆分为数组。左旧方法作为注释，新方法是下面的5行。

编辑4:编写多个文件的附加测试用例。(下文)

my @names = ('cat','dog','sheep'); 
foreach my $name (@names) { #$name is lexical, ie dies after each loop
  open(my $handle, '>', $name); #open a lexical handle for the file, also dies each loop
  print $handle $name; #write to the handle
  #$handles closes automatically when it "goes out of scope"
}

Answer 2

注意，在酶中，当您在散列中存储一个酶时，酶的名称应该是关键，站点应该是值(因为在原则上，两种酶可以有相同的位置)。

在主例程中，您可以遍历散列；每个酶生成一个输出文件。最直接的方法是将站点转换为regex (通过其他regexs)并将其应用于DNA序列，但还有其他方法。(这可能是值得的，至少把它分成另外一艘潜艇。)

Answer 3

下面是我试图解决这个问题的方法(下面的代码)。

1)从参数中提取文件名，并创建相应的filehandles。

2)为以指定格式生成的输出文件创建了一个新的文件句柄

3)从第一个文件中提取“切点”。

( 4)第二个文件中的DNA序列被循环在步骤3中的切割点上。

#!/usr/bin/perl
use strict;
use warnings;
my $file_enzyme=$ARGV[0];
my $file_dna=$ARGV[1];

open DNASEQ, ">$file_dna"."_"."$file_enzyme";
open ENZYME, "$file_enzyme";
open DNA, "$file_dna";
while (<ENZYME>) {
 chomp;
  if( /'(.*)\/\//) { # Extracts the cut point between ' & // in the enzyme file
    my $pattern=$1;
    while (<DNA>) {
     chomp;
     #print $pattern;
     my @output=split/$pattern/,;
     print DNASEQ shift @output,"\n"; #first recognized sequence being output
     foreach (@output) {
        print DNASEQ "$pattern$_\n"; #prefixing the remaining sequences with the cut point pattern
     }
   }
 }
}
close DNA;
close ENZYME;
close DNASEQ;