perl练习——FASTA格式文件中序列GC含量计算&perl数组排序如何获得下标或者键
2016-07-04 17:01
453 查看
一、关于程序:
FUN:计算FASTA文件中每条序列中G和C的含量百分比,输出最大值及其id
INPUT:FASTA格式文件
OUTPUT:最高含量的序列id及其含量(这是上面的结果)
二、编程思想及代码
当是注释行时(>……),获得序列 ID ,并跳过该次循环;当读到非注释行即序列行时,记录该行“G和C的含量”以及“序列的总含量”,这都可以利用perl上下文实现。(但是在这里有一些疑惑——当把14行@num换成$num会出现计算错误,知道的朋友欢迎留言)
三、技巧
神奇的perl,神奇的sort!!
对数组(或者哈希)排序获得下标的方式:
备注:贴一个感觉不错的代码(学习学习)
FUN:计算FASTA文件中每条序列中G和C的含量百分比,输出最大值及其id
INPUT:FASTA格式文件
>seq1 CGCCGAGCGCTTGACCTCCAGCAAGACGCCGTCTGGCACATGCAACGAGCTGTAGCAGAC >seq2 ATGCCTAGAACGTTCGAGACTTCTCGGGTGCGGTAGAATTAGCCATTCGACCGACTTCCA GCATCTGCGAGCCGCCTGTTGATTGCATCCGCCGGGGACGCAACAAGGCAAGGCCCTAAC
OUTPUT:最高含量的序列id及其含量(这是上面的结果)
seq1 63.333333%
二、编程思想及代码
当是注释行时(>……),获得序列 ID ,并跳过该次循环;当读到非注释行即序列行时,记录该行“G和C的含量”以及“序列的总含量”,这都可以利用perl上下文实现。(但是在这里有一些疑惑——当把14行@num换成$num会出现计算错误,知道的朋友欢迎留言)
use strict; my %GC_content; # id=>GC_content my %sequences; # id=>sequence my ($id, $sum); # id, 每个序列的字符个数 my @num; # 中间变量,用于存储单行中某字符的含量 while(my $seq = <>){ chomp($seq); if($seq =~ m/^>(.*)/){ $id = $1; next; } @num = ($seq =~ m/(G|C)/g); $GC_content{$id} += @num; @num = ($seq =~ m/(.)/g); $sequences{$id} += @num; } foreach(keys(%GC_content)){ $GC_content{$_} /= $sequences{$_}; } my @sort = sort{$GC_content{$b} <=> $GC_content{$a}} keys(%GC_content); printf("%s\n%.6f%\n", $sort[0], $GC_content{$sort[0]}*100);
三、技巧
神奇的perl,神奇的sort!!
对数组(或者哈希)排序获得下标的方式:
# 数字排序: my @arr = qw(2 3 41 2 34 ); my @result1 = sort{$a <=> $b} @arr; # 获得下标: my @result2 = sort{$arr[$a] <=> $arr[$b]} 0..$#arr; # 获得key: my %hash = ( one =>1, two =>5, tree=>9 ); my @result3 = sort{$hash{$a} <=> $hash{$b}} keys(%hash); print "数字排序:@result1\n获得下标:@result2\n获得key:@result3\n";
备注:贴一个感觉不错的代码(学习学习)
$/ = '>'; <>; # 读一次">"前的序列,以免下面代码出错 while (<>) { chomp; my ($id, @ary) = split '\n'; my $seq = join '', @ary; my $ratio = &GC_content($seq); if ($ratio > $highest) { $highest = $ratio; @result = ($id, $ratio); } } print join "\n", @result; sub GC_content { my ($seq) = @_; my $ratio = $seq =~ s/([CG])/$1/g / length($seq) * 100; return $ratio }
相关文章推荐
- 3.2 MEID和IMEI
- 提高Java 性能的5点 理解+实用
- C#中字符 '\0' 是a还是空格?
- 有用的鏈接
- 通过反射实现对象转JSON
- 3.1 4G FDD和TDD
- informix数据库的日志模式
- 多线程应用场景
- Machine Learning - 第6周(Advice for Applying Machine Learning、Machine Learning System Design)
- 物联网江湖 第三回 - 群鸦的盛宴 微软的阳谋
- 5个常用的Android自动化测试框架介绍
- [转]DataGridView 的右键菜单(ContextMenuStrip)
- linux 指令
- IOS 安全网络请求--HTTPS
- yii调用微信接口扫二维码
- HTML5 canvas画布(六)
- 非常好!!!Linux源代码阅读——内核引导【转】
- 消除performSelector:警告的方法
- linux下.run文件的安装与卸载
- 关于分布式系统的数据一致性问题(二)