2009年4月23日星期四

[PerlChina] Re: 编写提取染色体CDS的问题

以前写的,截取部分,你修改下估计可以用
sub parse_genbank {
    my $filename = shift;
    open GBFH, $filename;

    my ( $fflag, $features, $flag, $seqflag ) = ( 0, 0, 0, 0 );
    my ( $seq, $head, $id, $acc, $cds );
    my $data = {};

    foreach ( <GBFH> ) {
        if ( /^LOCUS/ && $flag == 0 ) {
            /^LOCUS\s*(.*\S)\s*$/;
            ( $id = ( $head = $1 ) ) =~ s/\s.*//;
            $features += 1;
            $flag = 1;
        }
        elsif ( /^DEFINITION\s*(.*)/ && $flag == 1 ) {
            $head .= " $1";
            $features += 2;
        }
        elsif ( /^ACCESSION/ && $flag == 1 ) {
            /^ACCESSION\s*(.*\S)\s*$/;
            $head .= " " . ( $acc = $1 );
            $features += 4;
        }
        elsif ( /^FEATURES/ && $flag == 1 ){
            $fflag = 1;
        }
        elsif ( /^ {5}CDS/ && $fflag == 1){
            /^ {5}CDS\s*(.*\S)\s*$/;
            $cds = $1;
        }
        elsif ( /^ORIGIN/ ) {
            $seqflag = 1;
            $features += 8;
        }
        elsif (/^\/\//) {       # end of genbank file
            $features += 16;
        }
        elsif ( $seqflag == 0 ) {
            ;
        }
        elsif ( $seqflag == 1 ) {
            $seq .= $_;
        }
        if ( $features == 31 ) {
            $seq =~ tr/a-zA-Z//cd;
            $seq =~ tr/A-Z/a-z/;
            my @pos = split /\.+/, $cds;
            my ($begin,$end) = @pos;
            $data->{_sequence}  = $seq;
            $data->{_header}    = $head;
            $data->{_id}        = $id;
            $data->{_accession} = $acc;
            $data->{_cds}       = substr($seq,$begin-1,$end-$begin+1);
            return \$data;
        }
    }
    return;
}


2009/4/23 guo <sczz520@126.com>
我在windows上下载安装了activeperl,想编写给程序把人类染色体序列上的编码区也就是CDS给提取出来。有朋友说用bioperl模块
来编程效率更高。可是我不会用。自己试着使用ppm来安装bioperl也没有成功。
现在主要是想请教一下,有没有人会bioperl,或者activperl的,能不能帮我写一个程序来把人类染色体序列中的CDS提取出来。
染色体序列我已经从NCBI上下载下来了。就是不知道如何编写程序。请教高手。
即使程序付费,我也愿意。
希望广大perl友帮忙。



--~--~---------~--~----~------------~-------~--~----~
您收到此信息是由于您订阅了 Google 论坛"PerlChina Mongers 讨论组"论坛。
 要在此论坛发帖,请发电子邮件到 perlchina@googlegroups.com
 要退订此论坛,请发邮件至 perlchina+unsubscribe@googlegroups.com
 更多选项,请通过 http://groups.google.com/group/perlchina?hl=zh-CN 访问该论坛

-~----------~----~----~----~------~----~------~--~---

没有评论: