#2
fall_bernana2019-11-29 17:50
回复 楼主 h17785687955
|
如何用perl计算密码子不同位点的GC含量???急用!!!!!
我是想计算每一条CDS密码子三个位点的GC含量,用下面这个脚本但没有数据输出,我的文本是fasta格式的,文件名是M.fasta,fasta结构的文本是以下形式:
>a
ATG GAA GAA GAG AAC GCG CAA CTC CGT ACC GAA TTG
>b
ATG GTT TTG TCA CAA GAA AAT GTT GAT CAC AAT GAT
“ $li=~s/^\w+://; # 去除每条序列前的序列编号和冒号”怎么修改?
脚本:
#!/usr/bin/perl
sub gc {
my $x;
if($_[0] eq "G" or $_[0] eq "C" {
$x=1;
}
$x;
}
my $input=$ARGV[0];
my ($li,$i,$len,$s,$tm,@t,@g);
open TMP, $M.fasta|| die;
while($li=<TMP> {
my @t;
my @g;
chomp $li;
#print "$li\n";
$li=~s/^\w+://; # 去除每条序列前的序列编号和冒号
$li=~s/\s//g;
#print "$li\n";
$len=length $li;
$i=0;
while($i<$len) {
$s=substr($li,$i,1);
$y=$i%3; # 判断是三个位点中的哪个
$t[$y]++; # 该位点的总碱基数目
$g[$y]+=gc($s); # 该位点的GC碱基数目
$i++;
}
$tm=$g[0]/$t[0]*100;
printf "The GC content in site1 is %.3f\%\n",$tm;
$tm=$g[1]/$t[1]*100;
printf "The GC content in site2 is %.3f\%\n",$tm;
$tm=$g[2]/$t[2]*100;
printf "The GC content in site3 is %.3f\%\n",$tm;
}
close TMP;