美文网首页
计算全基因组GC content

计算全基因组GC content

作者: 余绕 | 来源:发表于2020-07-16 20:46 被阅读0次

    计算每条染色体上5000bp 一个bin的GC百分比。

    代码如下:

    open FA,"$ARGV[0]";
    
    $/=">";
    
    while(<FA>){
        chomp;
        
        my($ID,$seq)=split /\n/,$_,2;
        $seq=~s/\n//g;
        $seq=uc($seq);
        #print "$ID\n$seq\n";
        $hash{$ID}=$seq;
    }
    
    open OU,">$ARGV[1]";
    foreach (keys %hash){
        $key=$_;
        $seq1=$hash{$key};
        #print  OU "$key\n";
        $len=length($seq1);
        #print OU "$len\n";
        
        $num=int($len/5000);
        #print OU "--- $num --- \n";
        $num=$num+1 if ($len/5000 > $num);
        #print OU "---- $num ---- \n";
        $start=0;
        $end=0;
        
        for ($i=1;$i<= $num;$i+=1){
            #print "$seq1\n";
        
            $seq2=substr($seq1,$start,5000);
            #print OU "$start\n$seq2\n";
            $start+=5000;
            $G=($seq2=~s/G/G/g);
            $C=($seq2=~s/C/C/g);
            $GC=($G+$C)/5000;
            #$HASH{"$key.$i"}=$GC; #可以创建hash,用于后面输出进行排序!
            print OU "$key\_$i\t$GC\n";
            
        }
        
    }
    
    #foreach (sort keys %HASH ){
    #   
    #   print OU "$_\t$HASH{$_}\n";
    #   
    #}
    
    
    

    运行code:perl GC.pl IRGSP-1.0_genome.fasta OUT.GC
    结果展示:


    image.png

    相关文章

      网友评论

          本文标题:计算全基因组GC content

          本文链接:https://www.haomeiwen.com/subject/ecimhktx.html