22# -*- coding: utf-8 -*-
33
44import sys
5- < << << << HEAD
65import argparse
76import numpy as np
8- == == == =
97import matplotlib
108matplotlib .use ('Agg' )
119matplotlib .rcParams ['pdf.fonttype' ] = 42
1210matplotlib .rcParams ['svg.fonttype' ] = 'none'
13- > >> >> >> 4.0 .0
1411import matplotlib .pyplot as plt
1512from scipy import interpolate
1613from scipy .stats import poisson
@@ -405,6 +402,8 @@ def main(args=None):
405402 x = np .arange (total ).astype ('float' ) / total # normalize from 0 to 1
406403
407404 if args .plotFile is not None :
405+ plt .style .use ('ggplot' )
406+
408407 i = 0
409408 # matplotlib won't iterate through line styles by itself
410409 pyplot_line_styles = sum ([7 * ["-" ], 7 * ["--" ], 7 * ["-." ], 7 * [":" ]], [])
@@ -422,41 +421,39 @@ def main(args=None):
422421 plt .close ()
423422
424423 if args .outRawCounts is not None :
425- of = open (args .outRawCounts , "w" )
426- of .write ("#plotFingerprint --outRawCounts\n " )
427- of .write ("'" + "'\t '" .join (args .labels ) + "'\n " )
428- fmt = "\t " .join (np .repeat ('%d' , num_reads_per_bin .shape [1 ])) + "\n "
429- for row in num_reads_per_bin :
430- of .write (fmt % tuple (row ))
431- of .close ()
424+ with open (args .outRawCounts , "w" ) as of :
425+ of .write ("#plotFingerprint --outRawCounts\n " )
426+ of .write ("'" + "'\t '" .join (args .labels ) + "'\n " )
427+ fmt = "\t " .join (np .repeat ('%d' , num_reads_per_bin .shape [1 ])) + "\n "
428+ for row in num_reads_per_bin :
429+ of .write (fmt % tuple (row ))
432430
433431 if args .outQualityMetrics is not None :
434- of = open (args .outQualityMetrics , "w" )
435- of .write ("Sample\t AUC\t Synthetic AUC\t X-intercept\t Synthetic X-intercept\t Elbow Point\t Synthetic Elbow Point" )
436- if args .JSDsample :
437- of .write ("\t JS Distance\t Synthetic JS Distance\t % genome enriched\t diff. enrichment\t CHANCE divergence" )
438- else :
439- of .write ("\t Synthetic JS Distance" )
440- of .write ("\n " )
441- line = np .arange (num_reads_per_bin .shape [0 ]) / float (num_reads_per_bin .shape [0 ] - 1 )
442- for idx , reads in enumerate (num_reads_per_bin .T ):
443- counts = np .cumsum (np .sort (reads ))
444- counts = counts / float (counts [- 1 ])
445- AUC = np .sum (counts ) / float (len (counts ))
446- XInt = (np .argmax (counts > 0 ) + 1 ) / float (counts .shape [0 ])
447- elbow = (np .argmax (line - counts ) + 1 ) / float (counts .shape [0 ])
448- expected = getExpected (np .mean (reads )) # A tuple of expected (AUC, XInt, elbow)
449- of .write ("{0}\t {1}\t {2}\t {3}\t {4}\t {5}\t {6}" .format (args .labels [idx ], AUC , expected [0 ], XInt , expected [1 ], elbow , expected [2 ]))
432+ with open (args .outQualityMetrics , "w" ) as of :
433+ of .write ("Sample\t AUC\t Synthetic AUC\t X-intercept\t Synthetic X-intercept\t Elbow Point\t Synthetic Elbow Point" )
450434 if args .JSDsample :
451- JSD = getJSD (args , idx , num_reads_per_bin )
452- syntheticJSD = getSyntheticJSD (num_reads_per_bin [:, idx ])
453- CHANCE = getCHANCE (args , idx , num_reads_per_bin )
454- of .write ("\t {0}\t {1}\t {2}\t {3}\t {4}" .format (JSD , syntheticJSD , CHANCE [0 ], CHANCE [1 ], CHANCE [2 ]))
435+ of .write ("\t JS Distance\t Synthetic JS Distance\t % genome enriched\t diff. enrichment\t CHANCE divergence" )
455436 else :
456- syntheticJSD = getSyntheticJSD (num_reads_per_bin [:, idx ])
457- of .write ("\t {0}" .format (syntheticJSD ))
437+ of .write ("\t Synthetic JS Distance" )
458438 of .write ("\n " )
459- of .close ()
439+ line = np .arange (num_reads_per_bin .shape [0 ]) / float (num_reads_per_bin .shape [0 ] - 1 )
440+ for idx , reads in enumerate (num_reads_per_bin .T ):
441+ counts = np .cumsum (np .sort (reads ))
442+ counts = counts / float (counts [- 1 ])
443+ AUC = np .sum (counts ) / float (len (counts ))
444+ XInt = (np .argmax (counts > 0 ) + 1 ) / float (counts .shape [0 ])
445+ elbow = (np .argmax (line - counts ) + 1 ) / float (counts .shape [0 ])
446+ expected = getExpected (np .mean (reads )) # A tuple of expected (AUC, XInt, elbow)
447+ of .write ("{0}\t {1}\t {2}\t {3}\t {4}\t {5}\t {6}" .format (args .labels [idx ], AUC , expected [0 ], XInt , expected [1 ], elbow , expected [2 ]))
448+ if args .JSDsample :
449+ JSD = getJSD (args , idx , num_reads_per_bin )
450+ syntheticJSD = getSyntheticJSD (num_reads_per_bin [:, idx ])
451+ CHANCE = getCHANCE (args , idx , num_reads_per_bin )
452+ of .write ("\t {0}\t {1}\t {2}\t {3}\t {4}" .format (JSD , syntheticJSD , CHANCE [0 ], CHANCE [1 ], CHANCE [2 ]))
453+ else :
454+ syntheticJSD = getSyntheticJSD (num_reads_per_bin [:, idx ])
455+ of .write ("\t {0}" .format (syntheticJSD ))
456+ of .write ("\n " )
460457
461458
462459if __name__ == "__main__" :
0 commit comments