Handle cases with no actionable mutations (resolves #214)

resolves #214

ProTECT now fails gracefully if Transgene receives an empty file,
 Transgene finds no actionable mutations, or merge_mhc_peptide_calls
produces no peptides. This is done by propagating a None value.

Author: davidstew

src/protect/binding_prediction/common.py

@@ -58,6 +58,10 @@ def spawn_antigen_predictors(job, transgened_files, phlat_files, univ_options, m
                             +- 'normal': fsID
     :rtype: tuple(dict, dict)
     """
+    if transgened_files is None:
+        job.fileStore.logToMaster('spawn_antigen_predictors received no peptides from Transgene for %s. Skipping.'
+                                  % univ_options['patient'])
+        return None
     work_dir = os.getcwd()
     mhci_options, mhcii_options = mhc_options
     pept_files = {
@@ -498,6 +502,10 @@ def merge_mhc_peptide_calls(job, antigen_predictions, transgened_files, univ_opt
     :rtype: dict
     """
     job.fileStore.logToMaster('Merging MHC calls')
+    if antigen_predictions is None and transgened_files is None:
+        job.fileStore.logToMaster('merge_mhc_peptide_calls received no peptides from Transgene (and no pMHC binding '
+                                  'predictions) for %s. Skipping.' % univ_options['patient'])
+        return None
     work_dir = os.getcwd()
     pept_files = {
         '10_mer.faa': transgened_files['transgened_tumor_10_mer_peptides.faa'],
@@ -568,7 +576,9 @@ def merge_mhc_peptide_calls(job, antigen_predictions, transgened_files, univ_opt
                 print_mhc_peptide(pred, peptides, pepmap, mhcii_resfile,
                                   netmhc=mhcii_preds[key]['predictor'] == 'netMHCIIpan')
     if not(mhci_called or mhcii_called):
-        raise RuntimeError('No peptides available for ranking')
+        job.fileStore.logToMaster('No peptides available for ranking after running merge_mhc_peptide_calls on %s.'
+                                  % univ_options['patient'])
+        return None
     output_files = defaultdict()
     for mhc_file in [mhci_resfile.name, mhcii_resfile.name]:
         output_files[os.path.split(mhc_file)[1]] = job.fileStore.writeGlobalFile(mhc_file)

src/protect/mutation_translation.py

@@ -90,7 +90,8 @@ def run_transgene(job, snpeffed_file, rna_bam, univ_options, transgene_options,
                   '--pep_lens', '9,10,15',
                   '--cores', str(transgene_options['n']),
                   '--genome', input_files['genome.fa'],
-                  '--annotation', input_files['annotation.gtf']]
+                  '--annotation', input_files['annotation.gtf'],
+                  '--log_file', '/data/transgene.log']
 
     if snpeffed_file is not None:
         parameters.extend(['--snpeff', input_files['snpeffed_muts.vcf']])
@@ -111,13 +112,19 @@ def run_transgene(job, snpeffed_file, rna_bam, univ_options, transgene_options,
         parameters += ['--transcripts', fusion_files['transcripts.fa'],
                        '--fusions', fusion_files['fusion_calls']]
 
-    docker_call(tool='transgene',
-                tool_parameters=parameters,
-                work_dir=work_dir,
-                dockerhub=univ_options['dockerhub'],
-                tool_version=transgene_options['version'])
+    try:
+        docker_call(tool='transgene',
+                    tool_parameters=parameters,
+                    work_dir=work_dir,
+                    dockerhub=univ_options['dockerhub'],
+                    tool_version=transgene_options['version'])
+    finally:
+        logfile = os.path.join(os.getcwd(), 'transgene.log')
+        export_results(job, job.fileStore.writeGlobalFile(logfile), logfile, univ_options,
+                       subfolder='mutations/transgened')
 
     output_files = defaultdict()
+    peptides_not_found = False
     for peplen in ['9', '10', '15']:
         for tissue_type in ['tumor', 'normal']:
             pepfile = '_'.join(['transgened', tissue_type, peplen, 'mer_peptides.faa'])
@@ -128,9 +135,13 @@ def run_transgene(job, snpeffed_file, rna_bam, univ_options, transgene_options,
                 if tissue_type == 'tumor':
                     os.rename(os.path.join(work_dir, old_pepfile + '.map'),
                               os.path.join(work_dir, pepfile + '.map'))
-
+            if not os.path.exists(pepfile):
+                peptides_not_found = True
+                break
             output_files[pepfile] = job.fileStore.writeGlobalFile(os.path.join(work_dir, pepfile))
             export_results(job, output_files[pepfile], pepfile, univ_options, subfolder='peptides')
+        if peptides_not_found:
+            break
         mapfile = '_'.join(['transgened_tumor', peplen, 'mer_peptides.faa.map'])
         output_files[mapfile] = job.fileStore.writeGlobalFile(os.path.join(work_dir, mapfile))
         export_results(job, output_files[mapfile], mapfile, univ_options, subfolder='peptides')
@@ -144,6 +155,10 @@ def run_transgene(job, snpeffed_file, rna_bam, univ_options, transgene_options,
         os.rename('transgened_transgened.bedpe', 'fusions.bedpe')
         export_results(job, job.fileStore.writeGlobalFile('fusions.bedpe'), 'fusions.bedpe',
                        univ_options, subfolder='mutations/transgened')
-
-    job.fileStore.logToMaster('Ran transgene on %s successfully' % univ_options['patient'])
-    return output_files
+    if peptides_not_found:
+        job.fileStore.logToMaster('Transgene failed to find any peptides for %s.'
+                                  % univ_options['patient'])
+        return None
+    else:
+        job.fileStore.logToMaster('Ran transgene on %s successfully' % univ_options['patient'])
+        return output_files

src/protect/rankboost.py

@@ -36,6 +36,10 @@ def wrap_rankboost(job, rsem_files, merged_mhc_calls, transgene_out, univ_option
                 +- 'mhci_rankboost_detailed_results.txt': fsID
     :rtype: dict
     """
+    if merged_mhc_calls is None:
+        job.fileStore.logToMaster('Rankboost was provided no peptides from sample %s for ranking. Skipping.'
+                                  % univ_options['patient'])
+        return None
     rankboost = job.addChildJobFn(boost_ranks, rsem_files['rsem.isoforms.results'],
                                   merged_mhc_calls, transgene_out, univ_options, rankboost_options)
 

src/protect/test/ci/test_protect.py

@@ -142,10 +142,10 @@ def _test_ran_successfully(self, expected_dirs):
                                              ['mutations.vcf']),
                                 'transgened_f': ('/mnt/ephemeral/done/TEST/mutations/transgened',
                                                  [],
-                                                 ['fusions.bedpe', 'mutations.vcf']),
+                                                 ['fusions.bedpe', 'mutations.vcf', 'transgene.log']),
                                 'transgened_v': ('/mnt/ephemeral/done/TEST/mutations/transgened',
                                                  [],
-                                                 ['mutations.vcf']),
+                                                 ['mutations.vcf', 'transgene.log']),
                             },
                             'peptides': ('/mnt/ephemeral/done/TEST/peptides',
                                          [],