Categorized | Biocompute

快速提取序列的Perl脚本

Posted on 04 一月 2010 by 柳城 ,阅读 777

是这么回事的,整理一下一些有用的提问。感谢大家的提问。

网友sun:

我从miRbase数据库下载到了,matrue.fa,即里面有测出来的物种的成熟miRNA,里面有human 我想用perl把human的提出来另外存到一个文件.这是其中人的一条的
>hsa-miR-96 MIMAT0000095 Homo sapiens miR-96
UUUGGCACUAGCACAUUUUUGCU
怎么把所有的都提出来.
ps:正则表达式不是这个$s=~/^>hsa/期待你的回信.

有接触过miRbase数据库的,都知道miRNA的命名Human是用hsa开头,Mouse是用mmu开头,Rat是用rno开头。

解释一下

忘了大伙对Fasta格式可能有不明白的,我来稍微解释一下上面的话。比如有一个文件:matrue.fa。里面存放的内容的格式是这样子的:

>hsa-miR-96
UUUGGCACUAGCACAUUUUUGCU
>mmu-miR-96
AAAGGCACUAGCACUCGAUCGA
>rno-miR-85
UAGUCUGCAUGUUGCACAGUGUA

上面只是举一个为例子,数据都是批量的。想要取出以hsa命名开头的数据,取出来的结果是:

>hsa-miR-96
UUUGGCACUAGCACAUUUUUGCU
……

 

这样子就有了下面的Perl脚本:

提取Fasta序列的Perl脚本

#!/usr/bin/perl
$/ = '>';
print '>';

while(<>){
    print if(/hsa-/);
}

这样就实现了从Fasta格式的文件中提取序列的目的,如果你再稍微改改,类似GenBank格式的也是可以提取的。还有好多,只要有唯一分隔符的都可以,不一一列出。

提醒一下:

$/ = '>';

强大之处是用了这句,自行慢慢体会。

转载请注明 : 来源于 快速提取序列的Perl脚本 | 柳城

赞助商

23条评论 于 “快速提取序列的Perl脚本”

  1. 卢松松 卢松松 Says:

    嘿嘿,没看懂,
    新年快乐

    [回复]

    A.shun
    A.shun 回复:

    坐一下老卢的大腿,也没看懂

    [回复]

  2. evlos evlos Says:

    汗,看得云里雾里 ~
    对了,博主应该有从php转到perl的经历吧 ~
    能否简述下过程 ~ 小邪有点心动 ~ [呲牙]

    [回复]

    柳城
    柳城 回复:

    你也想学Perl?
    我是要用到Perl才学的啊. [呲牙]

    [回复]

  3. Blinux Blinux Says:

    佩服数据分析 [呲牙]

    [回复]

  4. niche niche Says:

    问问楼主有什么好方法可以进行批量数据的PAML分析?比bioperl好的方法,呵呵 [握手]

    [回复]

    柳城
    柳城 回复:

    PAML不是建进化树的么? 不是太能理解你的意思~

    [回复]

  5. gahoo gahoo Says:

    这句是什么意思?
    $/ = '>';

    [回复]

    柳城
    柳城 回复:

    你可以自己测试一下啊。你找一个fasta格式的序列。

    #!/usr/bin/perl
    while(<>){
    print;
    }

    第二种是:

    #!/usr/bin/perl
    $/='>';
    while(<>){
    print;
    }

    试了你就知道了呗~~ [ok]

    [回复]

  6. niche niche Says:

    非常感谢您的回复,PAML是算核苷酸替换率的,http://www.bioperl.org/wiki/HOWTO:PAML
    非常喜欢你的博客,很多内容都很有用。

    [回复]

  7. 娱乐博客 娱乐博客 Says:

    好好加油吧,更上一层楼。

    [回复]

  8. 先看看 先看看 Says:

    真是太专业了,我完全看不懂

    [回复]

  9. 土狼妹妹 土狼妹妹 Says:

    好专业,学习了

    [回复]

  10. 淘宝皇冠店铺大全 淘宝皇冠店铺大全 Says:

    太专业了,我看不了啊

    [回复]

  11. 电脑资讯网 电脑资讯网 Says:

    [偷笑] [偷笑] [偷笑] [偷笑] 专业

    [回复]

  12. cricketlyz cricketlyz Says:

    您好,请问一下,怎么把提取出来的序列保存在一个新的文件里?谢谢了。

    [回复]

    柳城
    柳城 回复:

    大于号 > 输出啊

    [回复]

  13. motianfeng motianfeng Says:

    想请问一下各位大虾,我现在有一个达6G的txt数据,也是用'>'分开的,采用上面的方法提取每条染色体下面的数据,结果生成了几千个只有'>'的文件,不知是咋回事?另外,如果想要达到我的目的,有没有什么不太耗内存的方法?多谢啦先!~~

    [回复]

    柳城
    柳城 回复:

    具体问题具体分析嘛。可到yunbio.com提问哦

    [回复]

  14. tempo8 tempo8 Says:

    类似于http://liucheng.name/1205/
    #Usage:perl fasta2.pl in_fasta out_file
    in_fasta,out_file写在哪行?如何读文件以及如何输出结果?

    [回复]

  15. tempo8 tempo8 Says:

    楼主还有一个问题:
    如何批量定点提取序列,即CDS-Genome作BLAST比对之后,要提取ATG上游2000bp的序列!
    如:
    正链/正链(plus/plus)
    >subject1
    Query 1 ATCG
    Subject 5186 tagc
    这里有个特征是Query 1下一行对应的Subject后的数字(5186)减2000,即subject1中提取3186到5186

    此外,如何提取CDS下游2000bp的序列?
    若(plus/plus),在下一个">"到上一个">"之间在所有Subject这一行取最大数,直接加2000
    若(plus/minus)在Query这一行取最大数,再看Subject这一行对应的数值,然后减2000

    望不吝赐教!谢谢

    [回复]

    柳城
    柳城 回复:

    请到yunbio.com提问。。

    [回复]

Leave a Reply

广告招租

[强] [握手] [可爱] [ok] [呲牙] :) [偷笑] [流泪] [疑问] [亲亲] [擦汗] [得意] [衰] [可怜] [抱拳] [坏笑] more »

无觅相关文章插件,快速提升流量