pima_report: pima_report.py comparison

comparison pima_report.py @ 14:95b1d1a9497d draft

Uploaded

author	greg
date	Fri, 17 Mar 2023 17:23:43 +0000
parents	f03c80bb22e9
children	02283aa193c3

comparison

equal deleted inserted replaced

-:f03c80bb22e9
+:95b1d1a9497d
 import sys
 from Bio import SeqIO
 from datetime import date
 from mdutils.mdutils import MdUtils
+# FIXME: TableOfContents doesn't work.
+# from mdutils.tools import TableOfContents
 CDC_ADVISORY = 'The analysis and report presented here should be treated as preliminary.  Please contact the CDC/BDRD with any results regarding _Bacillus anthracis_.'
 class PimaReport:
 def __init__(self, analysis_name=None, amr_deletions_file=None, amr_matrix_files=None, assembly_fasta_file=None,
 assembly_name=None, bedtools_version=None, blastn_version=None, circos_files=None,
 compute_sequence_length_file=None, contig_coverage_file=None, dbkey=None, dnadiff_snps_file=None,
 dnadiff_version=None, feature_bed_files=None, feature_png_files=None, flye_assembly_info_file=None,
-flye_version=None, genome_insertions_file=None, gzipped=None, illumina_fastq_file=None,
+flye_version=None, genome_insertions_file=None, gzipped=None, kraken2_report_file=None,
-kraken2_report_file=None, kraken2_version=None, minimap2_version=None, mutation_regions_bed_file=None,
+kraken2_version=None, minimap2_version=None, mutation_regions_bed_file=None, mutation_regions_tsv_files=None,
-mutation_regions_tsv_files=None, pima_css=None, plasmids_file=None, quast_report_file=None,
+ont_fastq_file=None, pima_css=None, plasmids_file=None, quast_report_file=None, reference_insertions_file=None,
-reference_insertions_file=None, samtools_version=None, varscan_version=None):
+samtools_version=None, varscan_version=None):
 self.ofh = open("process_log.txt", "w")
 self.ofh.write("amr_deletions_file: %s\n" % str(amr_deletions_file))
 self.ofh.write("amr_matrix_files: %s\n" % str(amr_matrix_files))
 self.ofh.write("analysis_name: %s\n" % str(analysis_name))
 self.ofh.write("feature_png_files: %s\n" % str(feature_png_files))
 self.ofh.write("flye_assembly_info_file: %s\n" % str(flye_assembly_info_file))
 self.ofh.write("flye_version: %s\n" % str(flye_version))
 self.ofh.write("gzipped: %s\n" % str(gzipped))
 self.ofh.write("genome_insertions_file: %s\n" % str(genome_insertions_file))
-self.ofh.write("illumina_fastq_file: %s\n" % str(illumina_fastq_file))
 self.ofh.write("kraken2_report_file: %s\n" % str(kraken2_report_file))
 self.ofh.write("kraken2_version: %s\n" % str(kraken2_version))
 self.ofh.write("minimap2_version: %s\n" % str(minimap2_version))
 self.ofh.write("mutation_regions_bed_file: %s\n" % str(mutation_regions_bed_file))
 self.ofh.write("mutation_regions_tsv_files: %s\n" % str(mutation_regions_tsv_files))
+self.ofh.write("ont_fastq_file: %s\n" % str(ont_fastq_file))
 self.ofh.write("pima_css: %s\n" % str(pima_css))
 self.ofh.write("plasmids_file: %s\n" % str(plasmids_file))
 self.ofh.write("quast_report_file: %s\n" % str(quast_report_file))
 self.ofh.write("reference_insertions_file: %s\n" % str(reference_insertions_file))
 self.ofh.write("samtools_version: %s\n" % str(samtools_version))
 self.flye_version = 'flye (version unknown)'
 else:
 self.flye_version = re.sub('_', '.', flye_version.rstrip(' _assembly info_'))
 self.gzipped = gzipped
 self.genome_insertions_file = genome_insertions_file
-self.illumina_fastq_file = illumina_fastq_file
 self.kraken2_report_file = kraken2_report_file
 if kraken2_version is None:
 self.kraken2_version = 'kraken2 (version unknown)'
 else:
 self.kraken2_version = re.sub('_', '.', kraken2_version.rstrip(' _report_'))
 self.minimap2_version = 'minimap2 (version unknown)'
 else:
 self.minimap2_version = re.sub('_', '.', minimap2_version)
 self.mutation_regions_bed_file = mutation_regions_bed_file
 self.mutation_regions_tsv_files = mutation_regions_tsv_files
-self.read_type = 'Illumina'
-self.ont_bases = None
-self.ont_n50 = None
-self.ont_read_count = None
 self.pima_css = pima_css
 self.plasmids_file = plasmids_file
 self.quast_report_file = quast_report_file
+self.read_type = 'ONT'
 self.reference_insertions_file = reference_insertions_file
 self.reference_insertions_file = reference_insertions_file
 if samtools_version is None:
 self.samtools_version = 'samtools (version unknown)'
 else:
 # Values
 self.assembly_size = 0
 self.contig_info = None
 self.did_medaka_ont_assembly = False
 self.feature_hits = pandas.Series(dtype='float64')
-self.illumina_length_mean = 0
+self.illumina_fastq_file = None
-self.illumina_read_count = 0
+self.illumina_length_mean = None
-self.illumina_bases = 0
+self.illumina_read_count = None
-# TODO: should the following 2 values  be passed as  parameters?
+self.illumina_bases = None
+self.ont_bases = None
+# TODO: should the following be passed as a parameter?
+self.ont_coverage_min = 30
+self.ont_fast5 = None
+self.ont_fastq_file = ont_fastq_file
+self.ont_n50 = None
+# TODO: should the following be passed as a parameter?
 self.ont_n50_min = 2500
-self.ont_coverage_min = 30
+self.ont_raw_fastq = self.analysis_name
+self.ont_read_count = None
 # Actions
 self.did_guppy_ont_fast5 = False
 self.did_qcat_ont_fastq = False
+self.info_ont_fastq(self.ont_fastq_file)
 self.info_illumina_fastq()
 self.load_contig_info()
 def run_command(self, command):
 self.ofh.write("\nXXXXXX In run_command, command:\n%s\n\n" % str(command))
 self.assembly_size += len(i.seq)
 self.assembly_size = self.format_kmg(self.assembly_size, decimals=1)
 def info_illumina_fastq(self):
 self.ofh.write("\nXXXXXX In info_illumina_fastq\n\n")
+if self.illumina_length_mean is None:
+return
 if self.gzipped:
 opener = 'gunzip -c'
 else:
 opener = 'cat'
 command = ' '.join([opener,
-self.illumina_fastq_file,
+self.ont_fastq_file,
 '| awk \'{getline;s += length($1);getline;getline;}END{print s/(NR/4)"\t"(NR/4)"\t"s}\''])
 output = self.run_command(command)
 self.ofh.write("output:\n%s\n" % str(output))
 self.ofh.write("re.split('\\t', self.run_command(command)[0]:\n%s\n" % str(re.split('\\t', self.run_command(command)[0])))
 values = []
 self.illumina_length_mean += values[0]
 self.ofh.write("values[1]:\n%s\n" % str(values[1]))
 self.illumina_read_count += int(values[1])
 self.ofh.write("values[2]:\n%s\n" % str(values[2]))
 self.illumina_bases += int(values[2])
-# The original PIMA code inserts self.illumina_fastq into
+# The original PIMA code inserts self.illumina_fastq_file into
 # a list for no apparent reason.  We don't do that here.
-# self.illumina_length_mean /= len(self.illumina_fastq)
+# self.illumina_length_mean /= len(self.illumina_fastq_file)
 self.illumina_length_mean /= 1
 self.illumina_bases = self.format_kmg(self.illumina_bases, decimals=1)
 def start_doc(self):
 header_text = 'Analysis of ' + self.analysis_name
 self.doc = MdUtils(file_name=self.report_md, title=header_text)
+def add_table_of_contents(self):
+self.doc.create_marker(text_marker="TableOfContents")
+self.doc.new_line()
+self.doc.new_line('<div style="page-break-after: always;"></div>')
+self.doc.new_line()
 def add_run_information(self):
 self.ofh.write("\nXXXXXX In add_run_information\n\n")
 self.doc.new_line()
 self.doc.new_header(1, 'Run information')
 "Category",
 "Information",
 "Date",
 date.today(),
 "ONT FAST5",
-"N/A",
+self.wordwrap_markdown(self.ont_fast5),
 "ONT FASTQ",
-"N/A",
+self.wordwrap_markdown(self.ont_raw_fastq),
 "Illumina FASTQ",
-self.wordwrap_markdown(self.analysis_name),
+self.wordwrap_markdown(self.illumina_fastq_file),
 "Assembly",
 self.wordwrap_markdown(self.assembly_name),
 "Reference",
 self.wordwrap_markdown(self.dbkey),
 ]
 self.doc.new_table(columns=2, rows=7, text=Table_list, text_align='left')
 self.doc.new_line()
+# FIXME: the following doesn't work.
+# self.add_table_of_contents()
 self.doc.new_line()
 def add_ont_library_information(self):
 self.ofh.write("\nXXXXXX In add_ont_library_information\n\n")
 if self.ont_n50 is None:
 "ONT reads",
 '{:,}'.format(self.ont_read_count),
 "ONT bases",
 '{:s}'.format(self.ont_bases),
 "Illumina FASTQ",
-self.wordwrap_markdown(self.illumina_fastq_file),
+self.wordwrap_markdown(self.ont_fastq_file),
 "Assembly",
 self.wordwrap_markdown(self.assembly_name),
 "Reference",
 self.wordwrap_markdown(self.dbkey),
 ]
 self.doc.new_table(columns=2, rows=7, text=Table_List, text_align='left')
 self.doc.new_line()
 def add_illumina_library_information(self):
 self.ofh.write("\nXXXXXX In add_illumina_library_information\n\n")
-if self.illumina_length_mean == 0:
+if self.illumina_length_mean is None:
 return
 self.doc.new_line()
 self.doc.new_header(2, 'Illumina library statistics')
 Table_List = [
 "Illumina Info.",
 '{if(NR == 1){sub(/+/, "", $1);s=$1}else{bp += $1;if(bp > s / 2 && f == 0){n50 = $1;f = 1}}}',
 'END{printf "%d\\t%d\\t%d\\n", n50, (NR - 1), s;exit}\''])
 result = list(re.split('\\t', self.run_command(command)[0]))
 if result[1] == '0':
 self.error_out('No ONT reads found')
-ont_n50, ont_read_count, ont_raw_bases = [int(i) for i in result]
+self.ont_n50, self.ont_read_count, ont_raw_bases = [int(i) for i in result]
 command = ' '.join([opener,
 fastq_file,
 '| awk \'{getline;print length($0);getline;getline;}\''])
 result = self.run_command(command)
 result = list(filter(lambda x: x != '', result))
-# TODO: the following are not yet used...
+self.ont_bases = self.format_kmg(ont_raw_bases, decimals=1)
-# ont_read_lengths = [int(i) for i in result]
+if self.ont_n50 <= self.ont_n50_min:
-# ont_bases = self.format_kmg(ont_raw_bases, decimals=1)
+warning = 'ONT N50 (%s) is less than the recommended minimum (%s)' % (str(self.ont_n50), str(self.ont_n50_min))
-if ont_n50 <= self.ont_n50_min:
-warning = 'ONT N50 (%s) is less than the recommended minimum (%s)' % (str(ont_n50), str(self.ont_n50_min))
 self.assembly_notes = self.assembly_notes.append(pandas.Series(warning))
 def wordwrap_markdown(self, string):
 if string:
 if len(string) < 35:
 if self.did_guppy_ont_fast5:
 methods += ['ONT reads were basecalled using guppy']
 if self.did_qcat_ont_fastq:
 methods += ['ONT reads were demultiplexed and trimmed using qcat']
 self.methods[self.basecalling_methods_title] = pandas.Series(methods)
-self.add_illumina_library_information()
+# self.add_illumina_library_information()
 self.add_contig_info()
 self.evaluate_assembly()
 self.add_assembly_information()
 if self.flye_assembly_info_file is not None:
 method = 'ONT reads were assembled using %s' % self.flye_version.rstrip('assembly info')
 warning = 'Flye reported {:d} open contigs ({:s}); assembly may be incomplete.'.format(open_contigs.shape[0], ', '.join(open_contig_ids))
 self.assembly_notes = self.assembly_notes.append(pandas.Series(warning))
 if self.did_medaka_ont_assembly:
 method = 'the genome assembly was polished using ont reads and medaka.'
 self.methods[self.assembly_methods_title] = self.methods[self.assembly_methods_title].append(pandas.series(method))
-self.info_ont_fastq(self.illumina_fastq_file)
 self.add_assembly_notes()
 def make_tex(self):
 self.doc.new_table_of_contents(table_title='detailed run information', depth=2, marker="tableofcontents")
 text = self.doc.file_data_text
 parser.add_argument('--feature_png_dir', action='store', dest='feature_png_dir', help='Directory of best feature hits png files')
 parser.add_argument('--flye_assembly_info_file', action='store', dest='flye_assembly_info_file', default=None, help='Flye assembly info tabular file')
 parser.add_argument('--flye_version', action='store', dest='flye_version', default=None, help='Flye version string')
 parser.add_argument('--genome_insertions_file', action='store', dest='genome_insertions_file', help='Genome insertions BED file')
 parser.add_argument('--gzipped', action='store_true', dest='gzipped', default=False, help='Input sample is gzipped')
-parser.add_argument('--illumina_fastq_file', action='store', dest='illumina_fastq_file', help='Input sample')
+parser.add_argument('--ont_fastq_file', action='store', dest='ont_fastq_file', help='Input sample')
 parser.add_argument('--kraken2_report_file', action='store', dest='kraken2_report_file', default=None, help='kraken2 report file')
 parser.add_argument('--kraken2_version', action='store', dest='kraken2_version', default=None, help='kraken2 version string')
 parser.add_argument('--minimap2_version', action='store', dest='minimap2_version', default=None, help='minimap2 version string')
 parser.add_argument('--mutation_regions_bed_file', action='store', dest='mutation_regions_bed_file', help='AMR mutation regions BRD file')
 parser.add_argument('--mutation_regions_dir', action='store', dest='mutation_regions_dir', help='Directory of mutation regions TSV files')
 feature_png_files,
 args.flye_assembly_info_file,
 args.flye_version,
 args.genome_insertions_file,
 args.gzipped,
-args.illumina_fastq_file,
 args.kraken2_report_file,
 args.kraken2_version,
 args.minimap2_version,
 args.mutation_regions_bed_file,
 mutation_regions_files,
+args.ont_fastq_file,
 args.pima_css,
 args.plasmids_file,
 args.quast_report_file,
 args.reference_insertions_file,
 args.samtools_version,

Mercurial > repos > greg > pima_report

comparison pima_report.py @ 14:95b1d1a9497d draft