parameters: (contig_seqs, read_seqs, work_dir, num_iters, num_threads, read_platform,
             read_type, output_progress)
return_statement: return prev_assembly, stats_file
docstring: High-level polisher interface
docstring_summary: High-level polisher interface
docstring_tokens: High - level polisher interface
function: def polish(contig_seqs, read_seqs, work_dir, num_iters, num_threads, read_platform,
             read_type, output_progress):
      """
      High-level polisher interface
      """
      logger_state = logger.disabled
      if not output_progress:
          logger.disabled = True
  
      subs_matrix = os.path.join(cfg.vals["pkg_root"],
                                 cfg.vals["err_modes"][read_platform]["subs_matrix"])
      hopo_matrix = os.path.join(cfg.vals["pkg_root"],
                                 cfg.vals["err_modes"][read_platform]["hopo_matrix"])
      use_hopo = cfg.vals["err_modes"][read_platform]["hopo_enabled"]
      use_hopo = use_hopo and (read_type == "raw")
      stats_file = os.path.join(work_dir, "contigs_stats.txt")
  
      bam_input = read_seqs[0].endswith("bam")
  
      prev_assembly = contig_seqs
      contig_lengths = None
      coverage_stats = None
      for i in range(num_iters):
          logger.info("Polishing genome (%d/%d)", i + 1, num_iters)
  
          ####
          if not bam_input:
              logger.info("Running minimap2")
              alignment_file = os.path.join(work_dir, "minimap_{0}.bam".format(i + 1))
              make_alignment(prev_assembly, read_seqs, num_threads,
                             work_dir, read_platform, alignment_file,
                             reference_mode=True, sam_output=True)
          else:
              logger.info("Polishing with provided bam")
              alignment_file = read_seqs[0]
  
          #####
          logger.info("Separating alignment into bubbles")
          contigs_info = get_contigs_info(prev_assembly)
          bubbles_file = os.path.join(work_dir,
                                      "bubbles_{0}.fasta".format(i + 1))
          coverage_stats, mean_aln_error = \
              make_bubbles(alignment_file, contigs_info, prev_assembly,
                           read_platform, num_threads,
                           bubbles_file)
  
          logger.info("Alignment error rate: %f", mean_aln_error)
          consensus_out = os.path.join(work_dir, "consensus_{0}.fasta".format(i + 1))
          polished_file = os.path.join(work_dir, "polished_{0}.fasta".format(i + 1))
          if os.path.getsize(bubbles_file) == 0:
              logger.info("No reads were aligned during polishing")
              if not output_progress:
                  logger.disabled = logger_state
              open(stats_file, "w").write("#seq_name\tlength\tcoverage\n")
              open(polished_file, "w")
              return polished_file, stats_file
  
          #####
          logger.info("Correcting bubbles")
          _run_polish_bin(bubbles_file, subs_matrix, hopo_matrix,
                          consensus_out, num_threads, output_progress, use_hopo)
          polished_fasta, polished_lengths = _compose_sequence(consensus_out)
          fp.write_fasta_dict(polished_fasta, polished_file)
  
          #Cleanup
          os.remove(bubbles_file)
          os.remove(consensus_out)
          if not bam_input:
              os.remove(alignment_file)
  
          contig_lengths = polished_lengths
          prev_assembly = polished_file
  
      with open(stats_file, "w") as f:
          f.write("#seq_name\tlength\tcoverage\n")
          for ctg_id in contig_lengths:
              f.write("{0}\t{1}\t{2}\n".format(ctg_id,
                      contig_lengths[ctg_id], coverage_stats[ctg_id]))
  
      if not output_progress:
          logger.disabled = logger_state
  
      return prev_assembly, stats_file
function_tokens: def polish ( contig_seqs , read_seqs , work_dir , num_iters , num_threads , read_platform , read_type , output_progress ) : logger_state = logger . disabled if not output_progress : logger . disabled = True subs_matrix = os . path . join ( cfg . vals [ "pkg_root" ] , cfg . vals [ "err_modes" ] [ read_platform ] [ "subs_matrix" ] ) hopo_matrix = os . path . join ( cfg . vals [ "pkg_root" ] , cfg . vals [ "err_modes" ] [ read_platform ] [ "hopo_matrix" ] ) use_hopo = cfg . vals [ "err_modes" ] [ read_platform ] [ "hopo_enabled" ] use_hopo = use_hopo and ( read_type == "raw" ) stats_file = os . path . join ( work_dir , "contigs_stats.txt" ) bam_input = read_seqs [ 0 ] . endswith ( "bam" ) prev_assembly = contig_seqs contig_lengths = None coverage_stats = None for i in range ( num_iters ) : logger . info ( "Polishing genome (%d/%d)" , i + 1 , num_iters ) #### if not bam_input : logger . info ( "Running minimap2" ) alignment_file = os . path . join ( work_dir , "minimap_{0}.bam" . format ( i + 1 ) ) make_alignment ( prev_assembly , read_seqs , num_threads , work_dir , read_platform , alignment_file , reference_mode = True , sam_output = True ) else : logger . info ( "Polishing with provided bam" ) alignment_file = read_seqs [ 0 ] ##### logger . info ( "Separating alignment into bubbles" ) contigs_info = get_contigs_info ( prev_assembly ) bubbles_file = os . path . join ( work_dir , "bubbles_{0}.fasta" . format ( i + 1 ) ) coverage_stats , mean_aln_error = make_bubbles ( alignment_file , contigs_info , prev_assembly , read_platform , num_threads , bubbles_file ) logger . info ( "Alignment error rate: %f" , mean_aln_error ) consensus_out = os . path . join ( work_dir , "consensus_{0}.fasta" . format ( i + 1 ) ) polished_file = os . path . join ( work_dir , "polished_{0}.fasta" . format ( i + 1 ) ) if os . path . getsize ( bubbles_file ) == 0 : logger . info ( "No reads were aligned during polishing" ) if not output_progress : logger . disabled = logger_state open ( stats_file , "w" ) . write ( "#seq_name\tlength\tcoverage\n" ) open ( polished_file , "w" ) return polished_file , stats_file ##### logger . info ( "Correcting bubbles" ) _run_polish_bin ( bubbles_file , subs_matrix , hopo_matrix , consensus_out , num_threads , output_progress , use_hopo ) polished_fasta , polished_lengths = _compose_sequence ( consensus_out ) fp . write_fasta_dict ( polished_fasta , polished_file ) #Cleanup os . remove ( bubbles_file ) os . remove ( consensus_out ) if not bam_input : os . remove ( alignment_file ) contig_lengths = polished_lengths prev_assembly = polished_file with open ( stats_file , "w" ) as f : f . write ( "#seq_name\tlength\tcoverage\n" ) for ctg_id in contig_lengths : f . write ( "{0}\t{1}\t{2}\n" . format ( ctg_id , contig_lengths [ ctg_id ] , coverage_stats [ ctg_id ] ) ) if not output_progress : logger . disabled = logger_state return prev_assembly , stats_file