Skip to content

Update generate protein fasta and create peptide ordering form commands to use transcripts.fa #1453

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
susannasiebert wants to merge 5 commits into
base: 8.0.0
Choose a base branch
from
Open

Update generate protein fasta and create peptide ordering form commands to use transcripts.fa #1453

Show file tree
Hide file tree
Changes from 3 commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
7d4d190
Add a standalone command to generate the transcripts fasta
susannasiebert Aug 17, 2026
418f32b
Update generate protein fasta commands to take in a transcripts.fa fi…
susannasiebert Aug 18, 2026
dc4725f
Update create peptide ordering form commands to take in a transcripts…
susannasiebert Aug 18, 2026
b92ced1
Add documentation for the generate transcripts fasta command and upda…
susannasiebert Aug 19, 2026
267efc5
Apply batched suggestions from code review
susannasiebert Aug 31, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
The table of contents is too big for display.
Diff view
Diff view
  •  
  •  
  •  
157 changes: 20 additions & 137 deletions pvactools/lib/create_peptide_ordering_form.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,10 +6,11 @@
from pvactools.lib.generate_protein_fasta import PvacseqGenerateProteinFasta, PvacspliceGenerateProteinFasta, PvacfuseGenerateProteinFasta
from pvactools.lib.generate_reviews_files import main as run_generate_reviews_files
from pvactools.lib.color_peptides51mer import main as run_color_peptides
from pvactools.lib.run_argument_utils import downstream_sequence_length, aggregate_report_evaluations, pvacsplice_anchors
from pvactools.lib.run_argument_utils import aggregate_report_evaluations

class CreatePeptideOrderingForm:
def __init__(self, **kwargs):
self.transcripts_fasta = kwargs['transcripts_fasta']
self.flanking_sequence_length = kwargs['flanking_sequence_length']
self.classI_aggregated_tsv = kwargs['classI_aggregated_tsv']
self.classII_aggregated_tsv = kwargs['classII_aggregated_tsv']
Expand All @@ -23,12 +24,8 @@ def __init__(self, **kwargs):
sys.exit(f"Error: {self.output_path} must specify a directory.")
else:
os.makedirs(self.output_path)
self.phased_proximal_variants_vcf = kwargs.pop('phased_proximal_variants_vcf', None)
self.input_vcf = kwargs.pop('input_vcf', None)
self.external_vcf = kwargs.pop('external_vcf', None)
self.pass_only = kwargs.pop('pass_only', False)
self.biotypes = kwargs.pop('biotypes', ['protein_coding'])
self.allow_incomplete_transcripts = kwargs.pop('allow_incomplete_transcripts', False)
self.downstream_sequence_length = kwargs.pop('downstream_sequence_length', 1000)
self.aggregate_report_evaluation = kwargs.pop('aggregate_report_evaluation', ['Accept'])
self.classI_IC50 = kwargs.pop('classI_IC50', 1000.0)
self.classI_percent = kwargs.pop('classI_percent', 2.0)
Expand All @@ -47,42 +44,12 @@ def parser(cls, tool):
description="Generate peptide ordering files (FASTA, annotated ordering Excel spreadsheet, and review template Excel spreadsheet) to streamline preparation of peptides for synthesis and review.",
formatter_class=argparse.ArgumentDefaultsHelpFormatter
)
if tool == 'pvacseq':
parser.add_argument(
"input_vcf",
help="A VEP-annotated single- or multi-sample VCF containing genotype, transcript, "
+"Wildtype protein sequence, and Frameshift protein sequence information. "
+"The VCF may be gzipped (requires tabix index). This VCF will be used to extract "
+"peptide sequences for processable variants with 25 flanking amino acids on either "
+"side of the mutation. These sequences will be included in the peptide ordering spreadsheet."
)
elif tool == 'pvacsplice':
parser.add_argument(
"input_file",
help="RegTools junctions output TSV file"
)
parser.add_argument(
"annotated_vcf",
help="A VEP-annotated single- or multi-sample VCF containing genotype and transcript information."
+ "The VCF may be gzipped (requires tabix index)."
)
parser.add_argument(
"ref_fasta",
help="A reference FASTA file. Note: this input should be the same as the RegTools vcf input."
)
parser.add_argument(
"gtf_file",
help="A reference GTF file. Note: this input should be the same as the RegTools gtf input."
)
elif tool == 'pvacfuse':
parser.add_argument(
"input",
help="An AGFusion output directory or Arriba fusion.tsv output file."
)
parser.add_argument(
"ref_fasta",
help="A reference CDS FASTA file. Note: this input should match the build and Ensembl version used to create the fusion annotations."
)

parser.add_argument(
"transcripts_fasta",
help="A pVACsplice transcripts.fa file with transcript protein sequences of splicing events and matching wildtypes. "
Comment thread
susannasiebert marked this conversation as resolved.
Outdated
Comment thread
susannasiebert marked this conversation as resolved.
Outdated
+ f"This file can be found in the top-level output directory of your {tool.replace('vac', 'VAC')} run can be generated using the `{tool} generate_transcripts_fasta` command."
Comment thread
susannasiebert marked this conversation as resolved.
Outdated
)
parser.add_argument(
"flanking_sequence_length",
help="Number of amino acids to add on each side of the mutation when creating the FASTA.",
Expand Down Expand Up @@ -113,64 +80,21 @@ def parser(cls, tool):
)
if tool == 'pvacseq':
parser.add_argument(
"-p", "--phased-proximal-variants-vcf",
help="A VCF with phased proximal variant information to incorporate into the predicted fasta sequences "
+"generated from the input_vcf. Must be gzipped and tabix indexed."
"--input-vcf",
help='The original VEP-annotated VCF used in your pVACseq run to check variants against. '
+'Any variant with a PASS filter or no other filter applied in the input VCF and external VCF '
+'will be marked as called in the '
+'"Variant Called in External VCF" column of the updated aggregated report '
+'"<sample_name>.Annotated.Neoantigen_Candidates.xlsx"'
)
parser.add_argument(
'--external-vcf',
help='A VCF file from an external provider to check variants against. Any variant '
+'with a PASS filter or no other filter applied will be marked as called in the '
+'with a PASS filter or no other filter applied in the input VCF and external VCF '
+'will be marked as called in the '
+'"Variant Called in External VCF" column of the updated aggregated report '
+'"<sample_name>.Annotated.Neoantigen_Candidates.xlsx"'
)
if tool in ['pvacseq', 'pvacsplice']:
parser.add_argument(
'--pass-only',
help="Only process VCF entries with a PASS status.",
default=False,
action='store_true',
)
parser.add_argument(
"--biotypes", type=lambda s:[a for a in s.split(',')],
help="A list of biotypes to use for pre-filtering transcripts when generating peptide sequences from "
+"the input_vcf.",
default=['protein_coding']
)
parser.add_argument(
"--allow-incomplete-transcripts",
help="By default, transcripts annotated with incomplete CDS (i.e., 'cds_start_NF' or 'cds_end_NF' flags in the VEP CSQ field) "
+ "are excluded from analysis, as they often produce invalid protein sequences. "
+ "Use this flag to allow candidates from such transcripts. Only peptides that do not contain 'X' will be included. "
+ "These candidates will be deprioritized relative to those from transcripts without incomplete CDS flags.",
default=False,
action='store_true'
)
parser.add_argument(
"-d", "--downstream-sequence-length",
default="1000",
help="Cap to limit the downstream sequence length for frameshifts when creating the fasta file. "
+ "Use 'full' to include the full downstream sequence.",
type=downstream_sequence_length()
)
if tool == 'pvacsplice':
parser.add_argument(
"-j", "--junction-score", type=int,
help="Junction Coverage Cutoff. Only sites above this read depth cutoff will be considered.",
default=10
)
parser.add_argument(
"-v", "--variant-distance", type=int,
help="Regulatory variants can lie inside or outside of splicing junction."
+ "Maximum distance window (upstream and downstream) for a variant outside the junction.",
default=100
)
parser.add_argument(
"--anchor-types", type=pvacsplice_anchors(),
help="The anchor types of junctions to use. Multiple anchors can be specified using a comma-separated list."
+ "Choices: A, D, NDA, DA, N",
default=['A', 'D', 'NDA'],
)
parser.add_argument(
"--aggregate-report-evaluation",
help="Only include variants where the Evaluation column in the classI_aggregated_tsv matches this evaluation. "
Expand Down Expand Up @@ -247,19 +171,9 @@ def execute(self):
os.remove(self.combined_fasta_output_file)

class PvacseqCreatePeptideOrderingForm(CreatePeptideOrderingForm):
def __init__(self, **kwargs):
self.input_vcf = kwargs['input_vcf']
super().__init__(**kwargs)

def create_fastas(self):
params = {
'input_vcf': self.input_vcf,
'sample_name': self.sample_name,
'pass_only': self.pass_only,
'phased_proximal_variants_vcf': self.phased_proximal_variants_vcf,
'biotypes': self.biotypes,
'allow_incomplete_transcripts': self.allow_incomplete_transcripts,
'downstream_sequence_length': self.downstream_sequence_length,
'transcripts_fasta': self.transcripts_fasta,
'flanking_sequence_length': self.flanking_sequence_length,
'mutant_only': True,
'aggregate_report_evaluation': self.aggregate_report_evaluation,
Expand All @@ -278,31 +192,9 @@ def create_fastas(self):
os.remove(generator.manufacturability_file)

class PvacspliceCreatePeptideOrderingForm(CreatePeptideOrderingForm):
def __init__(self, **kwargs):
self.input_file = kwargs['input_file']
self.annotated_vcf = kwargs['annotated_vcf']
self.ref_fasta = kwargs['ref_fasta']
self.gtf_file = kwargs['gtf_file']
self.junction_score = kwargs.pop('junction_score', 10)
self.variant_distance = kwargs.pop('variant_distance', 100)
self.anchor_types = kwargs.pop('anchor_types', ['A', 'D', 'NDA'])
self.input_vcf = None
super().__init__(**kwargs)

def create_fastas(self):
params = {
'input_file': self.input_file,
'annotated_vcf': self.annotated_vcf,
'ref_fasta': self.ref_fasta,
'gtf_file': self.gtf_file,
'junction_score': self.junction_score,
'variant_distance': self.variant_distance,
'anchor_types': self.anchor_types,
'sample_name': self.sample_name,
'pass_only': self.pass_only,
'biotypes': self.biotypes,
'allow_incomplete_transcripts': self.allow_incomplete_transcripts,
'downstream_sequence_length': self.downstream_sequence_length,
'transcripts_fasta': self.transcripts_fasta,
'flanking_sequence_length': self.flanking_sequence_length,
'mutant_only': True,
'aggregate_report_evaluation': self.aggregate_report_evaluation,
Expand All @@ -321,18 +213,9 @@ def create_fastas(self):
os.remove(generator.manufacturability_file)

class PvacfuseCreatePeptideOrderingForm(CreatePeptideOrderingForm):
def __init__(self, **kwargs):
self.input = kwargs['input']
self.ref_fasta = kwargs['ref_fasta']
self.input_vcf = None
super().__init__(**kwargs)

def create_fastas(self):
params = {
'input': self.input,
'ref_fasta': self.ref_fasta,
'sample_name': self.sample_name,
'downstream_sequence_length': self.downstream_sequence_length,
'transcripts_fasta': self.transcripts_fasta,
'flanking_sequence_length': self.flanking_sequence_length,
'mutant_only': True,
'aggregate_report_evaluation': self.aggregate_report_evaluation,
Expand Down
104 changes: 7 additions & 97 deletions pvactools/lib/generate_protein_fasta.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,26 +13,18 @@

class GenerateProteinFasta:
def __init__(self, **kwargs):
self.sample_name = kwargs.pop('sample_name', 'tmp')
if self.sample_name is None:
self.sample_name = 'tmp'
self.downstream_sequence_length = kwargs.pop('downstream_sequence_length', 1000)
self.pass_only = kwargs.pop('pass_only', False)
self.biotypes = kwargs.pop('biotypes', ['protein_coding'])
self.allow_incomplete_transcripts = kwargs.pop('allow_incomplete_transcripts', False)
self.transcripts_fasta = kwargs['transcripts_fasta']
self.flanking_sequence_length = kwargs.pop('flanking_sequence_length')
self.mutant_only = kwargs.pop('mutant_only', False)
self.aggregate_report_evaluation = kwargs.pop('aggregate_report_evaluation', ['Accept'])
self.temp_dir = tempfile.mkdtemp()
self.fasta_file_path = kwargs.pop('fasta_file_path', os.path.join(self.temp_dir, f"{self.sample_name}.transcripts.fa"))
self.trimmed_fasta_file_path = kwargs.pop('trimmed_fasta_file_path', os.path.join(self.temp_dir, f"{self.sample_name}.transcripts.trimmed.fa"))
self.tsv_filtered_fasta_file_path = os.path.join(self.temp_dir, f"{self.sample_name}.transcripts.filtered.fa")
self.mutant_only_fasta_file_path = os.path.join(self.temp_dir, f"{self.sample_name}.transcripts.filtered.mt_only.fa")
self.trimmed_fasta_file_path = kwargs.pop('trimmed_fasta_file_path', os.path.join(self.temp_dir, f"tmp.transcripts.trimmed.fa"))
self.tsv_filtered_fasta_file_path = os.path.join(self.temp_dir, f"tmp.transcripts.filtered.fa")
self.mutant_only_fasta_file_path = os.path.join(self.temp_dir, f"tmp.transcripts.filtered.mt_only.fa")
self.input_tsv = kwargs.pop('input_tsv', None)
self.output_file = kwargs.pop('output_file', None)

def execute(self):
self.generate_fasta()
self.trim_sequences()
self.tsv_filter_sequences()
if self.mutant_only:
Expand All @@ -46,9 +38,6 @@ def execute(self):
CalculateManufacturability(self.output_file, self.manufacturability_file, 'fasta').execute()
print("Completed")

def generate_fasta(self):
raise Exception("Implement in child class")

def trim_sequences(self):
raise Exception("Implement in child class")

Expand Down Expand Up @@ -111,32 +100,11 @@ def filter_wildtype_sequences(self):
SeqIO.write(output_records, self.mutant_only_fasta_file_path, "fasta")

class PvacseqGenerateProteinFasta(GenerateProteinFasta):
def __init__(self, **kwargs):
self.input_vcf = kwargs.pop('input_vcf', None)
self.phased_proximal_variants_vcf = kwargs.pop('phased_proximal_variants_vcf', None)
super().__init__(**kwargs)

def generate_fasta(self):
from pvactools.lib.variant_to_kmer_pipeline import VariantToKmerPipeline
params = {
'output_dir' : self.temp_dir,
'input_file' : self.input_vcf,
'sample_name' : self.sample_name,
'pass_only' : self.pass_only,
'proximal_variants_vcf' : self.phased_proximal_variants_vcf,
'biotypes' : self.biotypes,
'allow_incomplete_transcripts': self.allow_incomplete_transcripts,
'downstream_sequence_length' : self.downstream_sequence_length,
'flanking_bases' : self.flanking_sequence_length,
}
pipeline = VariantToKmerPipeline(**params)
pipeline.generate_fasta()

def trim_sequences(self):
print("Trimming Variant Peptide FASTA")
records = {}
keys = set()
for record in SeqIO.parse(self.fasta_file_path, "fasta"):
for record in SeqIO.parse(self.transcripts_fasta, "fasta"):
records[record.id] = str(record.seq)
keys.add(record.id.split('.', 1)[1])

Expand Down Expand Up @@ -195,51 +163,10 @@ def trim_sequences(self):
print("Completed")

class PvacspliceGenerateProteinFasta(GenerateProteinFasta):
def __init__(self, **kwargs):
self.input_file = kwargs.pop('input_file', None)
self.annotated_vcf = kwargs.pop('annotated_vcf', None)
self.ref_fasta = kwargs.pop('ref_fasta', None)
self.gtf_file = kwargs.pop('gtf_file', None)
self.junction_score = kwargs.pop('junction_score', 10)
self.variant_distance = kwargs.pop('variant_distance', 100)
self.anchor_types = kwargs.pop('anchor_types', ['A', 'D', 'NDA'])
super().__init__(**kwargs)

def generate_fasta(self):
from pvactools.lib.junction_to_kmer_pipeline import JunctionToKmerPipeline
junction_arguments = {
'input_file_type' : 'junctions',
'junctions_dir' : self.temp_dir,
'input_file' : self.input_file,
'gtf_file' : self.gtf_file,
'save_gtf' : False,
'sample_name' : self.sample_name,
'ref_fasta' : self.ref_fasta,
'annotated_vcf' : self.annotated_vcf,
'pass_only' : self.pass_only,
'biotypes' : self.biotypes,
'allow_incomplete_transcripts' : self.allow_incomplete_transcripts,
'junction_score' : self.junction_score,
'variant_distance' : self.variant_distance,
'anchor_types' : self.anchor_types,
'downstream_sequence_length' : self.downstream_sequence_length,
'normal_sample_name' : None,
'keep_tmp_files' : False,
'class_i_epitope_length' : [],
'class_ii_epitope_length' : [],
'class_i_hla' : [],
'class_ii_hla' : [],
}

pipeline = JunctionToKmerPipeline(**junction_arguments)
pipeline.generate_fasta()

return pipeline.create_file_path('fasta')

def trim_sequences(self):
records = {}
keys = set()
for record in SeqIO.parse(self.fasta_file_path, "fasta"):
for record in SeqIO.parse(self.transcripts_fasta, "fasta"):
records[record.id] = record.seq
keys.add(record.id.split('.', 1)[1])

Expand All @@ -263,28 +190,11 @@ def trim_sequences(self):
SeqIO.write(output_records, self.trimmed_fasta_file_path, "fasta")

class PvacfuseGenerateProteinFasta(GenerateProteinFasta):
def __init__(self, **kwargs):
self.input = kwargs.pop('input', None)
self.ref_fasta = kwargs.pop('ref_fasta', None)
super().__init__(**kwargs)

def generate_fasta(self):
from pvactools.lib.fusion_to_kmer_pipeline import FusionToKmerPipeline
params = {
'input_file': self.input,
'output_dir': self.temp_dir,
'sample_name': self.sample_name,
'transcript_fasta': self.ref_fasta,
'downstream_sequence_length': self.downstream_sequence_length,
}
pipeline = FusionToKmerPipeline(**params)
pipeline.generate_fasta()

def trim_sequences(self):
print("Trimming Variant Peptide FASTA")
records = {}
keys = set()
for record in SeqIO.parse(self.fasta_file_path, "fasta"):
for record in SeqIO.parse(self.transcripts_fasta, "fasta"):
records[record.id] = str(record.seq)
keys.add(record.id.split('.', 1)[1])

Expand Down
Loading
Loading