Merge pull request #1015 from griffithlab/issue_1006

During reference proteome similarity step, handle cases where a variant's full peptide sequence is not in the fasta file
griffithlab · Aug 7, 2023 · f93fbd9 · f93fbd9
2 parents 86ecee5 + e91c29c
commit f93fbd9
Show file tree

Hide file tree

Showing 2 changed files with 5 additions and 2 deletions.
diff --git a/pvactools/lib/calculate_reference_proteome_similarity.py b/pvactools/lib/calculate_reference_proteome_similarity.py
@@ -279,7 +279,8 @@ def _get_full_peptide(self, line, mt_records_dict, wt_records_dict):
             (parsed_aa_change, pos, wt_aa, mt_aa) = index_to_aggregate_report_aa_change(aa_change, variant_type)
             if line['Best Transcript'] == transcript and line['AA Change'] == parsed_aa_change:
                 return (mt_records_dict[record_id], wt_records_dict[record_id], variant_type, mt_aa, wt_aa)
-        raise Exception("Unable to find full_peptide for variant {}".format(line['ID']))
+        print("Unable to find full_peptide for variant {}".format(line['ID']))
+        return (None, None, variant_type, mt_aa, wt_aa)
 
     def _get_peptide(self, line, mt_records_dict, wt_records_dict):
         ## Get epitope, peptide and full_peptide
@@ -293,6 +294,8 @@ def _get_peptide(self, line, mt_records_dict, wt_records_dict):
             if self._input_tsv_type(line) == 'aggregated':
                 epitope = line['Best Peptide']
                 (full_peptide, wt_peptide, variant_type, mt_amino_acids, wt_amino_acids) = self._get_full_peptide(line, mt_records_dict, wt_records_dict)
+                if full_peptide is None:
+                    return None, None
                 if variant_type != 'FS':
                     if line['Pos'] == 'NA':
                         mt_pos = None

diff --git a/pvactools/tools/pvacseq/generate_protein_fasta.py b/pvactools/tools/pvacseq/generate_protein_fasta.py
@@ -90,7 +90,7 @@ def convert_vcf(input_vcf, temp_dir, sample_name, phased_proximal_variants_vcf,
         convert_params['proximal_variants_vcf'] = phased_proximal_variants_vcf
         proximal_variants_tsv = os.path.join(temp_dir, 'proximal_variants.tsv')
         convert_params['proximal_variants_tsv'] = proximal_variants_tsv
-        convert_params['flanking_bases'] = flanking_sequence_length * 4
+        convert_params['flanking_bases'] = (flanking_sequence_length + 1)* 4
     else:
         proximal_variants_tsv = None
     if pass_only: