Created
September 14, 2022 21:42
-
-
Save StephenHwang/78ebf0125c08fb69d6d02e495d46ce5f to your computer and use it in GitHub Desktop.
Calculate N50 from list of contig sizes
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#!/usr/bin/env python3 | |
from numpy import cumsum, argmax | |
from sys import stdin | |
def readStdin(): | |
''' Read stdin line-by-line. ''' | |
for line in stdin: | |
line = line.rstrip() | |
if line == 'Exit': | |
break | |
yield line | |
def calculate_N50(contig_lengths): | |
''' Calculate contig N50. ''' | |
contig_lengths = sorted(contig_lengths, reverse=True) | |
csum = cumsum(contig_lengths) | |
N50_target_value = csum[-1] / 2 | |
N50_value = contig_lengths[argmax(csum >= N50_target_value)] | |
return N50_value | |
def main(): | |
contig_lengths = list(map(int, readStdin())) | |
print(calculate_N50(contig_lengths)) | |
main() |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment