n50.pl

# Differences

This shows you the differences between two versions of the page.

 — n50.pl [2019/11/12 15:16] (current)hyjeong created 2019/11/12 15:16 hyjeong created 2019/11/12 15:16 hyjeong created Line 1: Line 1: + #​!/​usr/​bin/​perl + # + # Source: + #    http://​genomics-array.blogspot.com/​2011/​02/​calculating-n50-of-contig-assembly-file.html + #    Modified by Haeyoung Jeong + ## Read Fasta File and compute length ### + # + my \$length; + my \$totalLength;​ + my @arr; + my \$num; + + + while(<>​){ + ​chomp;​ + ​if(/​^>​([^\s]+)/​){ + print " \$length\n"​ unless \$num == 0; + ​\$num++;​ + print "SEQ: \$1 "; + push (@arr, \$length); + ​\$totalLength += \$length; + ​\$length=0;​ + } else { + s/\s//g; + s/\t//g; + \$length += length(\$_);​ + } + } + push (@arr, \$length); # for the last contig + print " \$length\n";​ + \$totalLength += \$length; + + close(FH);​ + + my @sort = sort {\$b <=> \$a} @arr; + my \$n50; + print "Total \$num contigs (\$totalLength)\n";​ + foreach my \$val(@sort){ + ​\$n50+=\$val;​ + if(\$n50 >= \$totalLength/​2){ + print "N50 contig length is \$n50 and N50 value is: \$val\n";​ + ​last;​ + } + } + print "​Average contig length is ", \$totalLength / \$num, "​\n";​ + print "Max. contig length is ", \$sort, "​\n";​ 