如何解析文件,创建记录并对记录执行操作,包括术语频率和距离计算

Koa*_*ala 6 perl hash bioinformatics

我是一个介绍Perl课程的学生,正在寻找建议和反馈我的方法来编写一个小的(但棘手的)程序来分析有关原子的数据.我的教授鼓励论坛.我没有使用Perl subs或模块(包括Bioperl),因此请将响应限制在适当的"初学者级别",以便我能够理解并从您的建议和/或代码中学习(也请限制"魔术").

该计划的要求如下:

  1. 从命令行读取一个文件(包含有关Atoms的数据)并创建一个原子记录数组(每个换行符一个记录/原子).对于程序需要存储的每条记录:

    •原子的序列号(第7 - 11栏)
    •它所属的氨基酸的三个字母的名称(第18 - 20栏)
    •原子的三个坐标(x,y,z)(第31-54栏)
    •原子的单字母或双字母元素名称(例如C,O,N,Na)(cols 77-78)

  2. 提示三个命令之一:freq,length,density d(d是某个数字):

    •freq - 文件中每种原子的数量(例如氮,钠等)将显示如下:N:918 S:23
    •长度 - 坐标之间的距离
    •密度d(其中d是数字) - 程序将提示输入文件的名称以保存计算,并将包含该原子与每个其他原子之间的距离.如果该距离小于或等于数字d,则增加原子数的计数在该距离内,除非该计数为零到文件.输出将看起来像:
    1:5
    2:3
    3:6
    ...(非常大的文件),并将在它完成时关闭.

我正在寻找下面代码中我写的(并且需要写)的反馈.我特别感谢有关如何编写我的潜艇的任何反馈.我在底部包含了示例输入数据.

我看到的程序结构和功能描述:

$^W = 1; # turn on warnings
use strict; # behave!

my @fields;
my @recs;

while ( <DATA> ) {
 chomp;
 @fields = split(/\s+/);
 push @recs, makeRecord(@fields);
}

for (my $i = 0; $i < @recs; $i++) {
 printRec( $recs[$i] );
}
    my %command_table = (
 freq => \&freq,
 length => \&length,
 density => \&density,
 help => \&help, 
 quit => \&quit
 );

print "Enter a command: ";
while ( <STDIN> ) {
 chomp; 
 my @line = split( /\s+/);
 my $command = shift @line;
 if ($command !~ /^freq$|^density$|length|^help$|^quit$/ ) {
    print "Command must be: freq, length, density or quit\n";
    }
  else {
    $command_table{$command}->();
    }
 print "Enter a command: ";
 }

sub makeRecord 
    # Read the entire line and make records from the lines that contain the 
    # word ATOM or HETATM in the first column. Not sure how to do this:
{
 my %record = 
 (
 serialnumber => shift,
 aminoacid => shift,
 coordinates => shift,
 element  => [ @_ ]
 );
 return\%record;
}

sub freq
    # take an array of atom records, return a hash whose keys are 
    # distinct atom names and whose values are the frequences of
    # these atoms in the array.  

sub length
    # take an array of atom records and return the max distance 
    # between all pairs of atoms in that array. My instructor
    # advised this would be constructed as a for loop inside a for loop. 

sub density
    # take an array of atom records and a number d and will return a
    # hash whose keys are atom serial numbers and whose values are 
    # the number of atoms within that distance from the atom with that
    # serial number. 

sub help
{
    print "To use this program, type either\n",
          "freq\n",
          "length\n",
          "density followed by a number, d,\n",
          "help\n",
          "quit\n";
}

sub quit
{
 exit 0;
}

# truncating for testing purposes. Actual data is aprox. 100 columns 
# and starts with ATOM or HETATM.
__DATA__
ATOM   4743  CG  GLN A 704      19.896  32.017  54.717  1.00 66.44           C  
ATOM   4744  CD  GLN A 704      19.589  30.757  55.525  1.00 73.28           C  
ATOM   4745  OE1 GLN A 704      18.801  29.892  55.098  1.00 75.91           O 
Run Code Online (Sandbox Code Playgroud)

FMc*_*FMc 5

看起来你的Perl技能正在很好地推进 - 使用引用和复杂的数据结构.以下是一些一般建议的提示和部分.

  • 使用use warnings而不是启用警告$^W = 1.前者是自我记录的,并且具有封闭块本地的优势,而不是全局设置.

  • 使用命名良好的变量,这将有助于记录程序的行为,而不是依赖于Perl的特殊功能$_.例如:

    while (my $input_record = <DATA>){
    }
    
    Run Code Online (Sandbox Code Playgroud)
  • 在用户输入方案中,无限循环提供了一种避免重复指令的方法,例如"输入命令".见下文.

  • 您的正则表达式可以简化,以避免重复锚点的需要.见下文.

  • 作为一般规则,肯定测试比否定测试更容易理解.请参阅if-else下面的修改结构.

  • 将程序的每个部分都包含在自己的子例程中.由于一系列原因,这是一个很好的通用做法,所以我只是开始习惯.

  • 相关的良好做法是尽量减少全局变量的使用.作为练习,您可以尝试编写程序,使其根本不使用全局变量.相反,任何所需的信息都将在子例程之间传递.对于小程序,人们不一定需要对避免使用全局变量来保持僵硬,但是记住理想并不是一个坏主意.

  • length子例程指定一个不同的名称.该名称已被内置length函数使用.

  • 关于你的问题makeRecord,一种方法是忽略内部的过滤问题makeRecord.相反,makeRecord可以包括一个额外的哈希字段,过滤逻辑将驻留在其他地方.例如:

    my $record = makeRecord(@fields);
    push @recs, $record if $record->{type} =~ /^(ATOM|HETATM)$/;
    
    Run Code Online (Sandbox Code Playgroud)

以上一些要点的说明:

use strict;
use warnings;

run();

sub run {
    my $atom_data = load_atom_data();
    print_records($atom_data);
    interact_with_user($atom_data);
}

...

sub interact_with_user {
    my $atom_data = shift;
    my %command_table = (...);

    while (1){
        print "Enter a command: ";
        chomp(my $reply = <STDIN>);

        my ($command, @line) = split /\s+/, $reply;

        if ( $command =~ /^(freq|density|length|help|quit)$/ ) {
            # Run the command.
        }
        else {
            # Print usage message for user.
        }
    }
}

...
Run Code Online (Sandbox Code Playgroud)