Bio Immunoinformatics Neoantigen Prediction

Version Compatibility

Reference examples tested with: Ensembl VEP 111+, pVACtools 4.1+, MHCflurry 2.1+, VAtools 5+, pandas 2.2+

Before using code patterns, verify installed versions match. If versions differ:

• Python: pip show <package> then help(module.function) to check signatures
• CLI: <tool> --version then <tool> --help to confirm flags

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Notes specific to this skill: pVACtools is now at 7.x; positional CLI args are stable across 4.x-7.x but defaults and the supported-algorithm list change between releases — always run pvacseq run --help against the installed build. pVACseq requires the Wildtype and Frameshift VEP plugins; the Downstream plugin was replaced by Frameshift in pVACtools 2.0, so 4.x+ pipelines must NOT use Downstream. A local IEDB install (--iedb-install-directory) is strongly preferred over the rate-limited public API for patient data.

Neoantigen Prediction

"Find neoantigens from my tumor mutations" -> Translate somatic variants into mutant peptides, predict patient-HLA presentation, and rank by tumor-specific quality for vaccine/biomarker use.

• CLI: pvacseq run on a VEP-annotated, expression/readcount-annotated somatic VCF + patient HLA (pVACtools)
• CLI: pvacfuse (fusions via AGFusion/Arriba), pvacbind (arbitrary peptides), pvacview (manual re-tiering)
• Python: VAtools annotation, LOHHLA/CCF integration, aggregate-report parsing

The Single Most Important Modern Insight -- binding is the easy part; PPV lives downstream

The visible surface of the field — NetMHCpan, MHCflurry, the IC50 column everyone sorts on — is the binding step, and binding is the one step the field has genuinely cracked. The positive predictive value of a binding-only neoantigen pipeline is single-digit percent: of peptides confidently called strong binders, the large majority are never presented, and of those presented, the large majority never elicit a T-cell response (TESLA; Wells 2020). This is structural, not a bad IC50 cutoff — each step of the presentation-and-recognition cascade multiplies a low conditional probability. The corrective: spend the analysis on the filters and features that govern the predicted->presented->immunogenic attrition (clonality/CCF, HLA LOH, expression, agretopicity, foreignness, processing, validation tiers) and treat the choice of binding algorithm as a near-afterthought with sane defaults. TESLA's five features that actually separated immunogenic peptides from binders: HLA binding affinity, source-gene expression ("tumor abundance"), peptide-HLA binding stability, hydrophobicity, and the two recognition features — agretopicity and foreignness.

The pVACtools Suite

Upstream Chain (every link can silently poison the output)

Decision Tree by Scenario

Run VEP, Then pVACseq

Goal: Produce the VEP annotation pVACseq actually consumes, then call neoantigens.

Approach: Run VEP with the Wildtype + Frameshift plugins and a protein FASTA; annotate expression and read counts with VAtools; supply a phased proximal-variants VCF; then pvacseq run with the patient HLA and sane filters.

pvacseq install_vep_plugin $VEP_PLUGINS          # installs Wildtype + Frameshift
vep --input_file somatic.vcf --output_file somatic.vep.vcf --format vcf --vcf \
    --symbol --terms SO --tsl --hgvs --fasta GRCh38.fa --offline --cache --dir_cache $VEP_CACHE \
    --plugin Frameshift --plugin Wildtype --pick

vcf-expression-annotator somatic.vep.vcf kallisto.tsv custom transcript -s TUMOR \
    --id-column target_id --expression-column tpm -o somatic.vep.expn.vcf

pvacseq run somatic.vep.expn.vcf TUMOR \
    "HLA-A*02:01,HLA-A*24:02,HLA-B*07:02,HLA-B*44:02,HLA-C*07:02,DRB1*01:01" \
    all_class_i pvac_out/ \
    -e1 8,9,10,11 --iedb-install-directory $IEDB \
    --phased-proximal-variants-vcf phased.vcf.gz \
    --normal-vaf 0.02 --tdna-vaf 0.25 --trna-vaf 0.25 --expn-val 1.0 -t 8

Key flags: -e1/-e2 epitope lengths; -b/--binding-threshold (default 500 nM); --percentile-threshold (recommend 2); -m/--top-score-metric median|lowest; --allele-specific-binding-thresholds (preferred over flat 500 nM); --net-chop-method/--netmhc-stab (processing + stability features).

Compute Agretopicity (DAI) Correctly

Goal: Quantify how much more foreign the mutant looks than its wild-type counterpart.

Approach: Agretopicity (the fitness-model amplitude; Łuksza 2017) is the WT/MT binding ratio. A high value means the mutant binds while the WT does not — the surface is new to the immune system, so reactive T cells were not deleted in the thymus. The original differential agretopicity index (DAI; Duan 2014) is the difference form; both forms share the traps below. Requires the matched WT peptide (the Wildtype plugin), so pVACbind cannot compute it.

import pandas as pd

def add_agretopicity(df, wt='Median WT IC50 Score', mt='Median MT IC50 Score'):
    '''Agretopicity (amplitude) = IC50_WT / IC50_MT (ratio > 1 = mutant binds better -> favorable).
    Anchor-position mutations inflate DAI without changing the TCR-facing surface, so
    pair DAI with anchor evaluation rather than trusting it alone.'''
    out = df.copy()
    out['agretopicity'] = out[wt] / out[mt]
    out['dai_favorable'] = out['agretopicity'] > 1
    return out

Drop Candidates on Lost HLA Alleles (LOHHLA)

Goal: Remove neoantigens predicted to be presented by an HLA allele the tumor has deleted.

Approach: HLA LOH is an immune-escape mechanism in ~40% of NSCLC (McGranahan 2017) and is invisible to binding/expression/clonality filters. Run LOHHLA (or a subclonal-sensitive equivalent like DASH) with the HLA type and tumor purity/ploidy, then filter the aggregate report. This step sits outside pVACtools and errors silently if skipped.

def drop_lost_allele_candidates(df, lost_alleles, allele_col='HLA Allele'):
    '''lost_alleles: set of alleles called as LOH-lost by LOHHLA. A peptide assigned
    to a lost allele is not weakly presented - it is not presented at all.'''
    return df[~df[allele_col].isin(set(lost_alleles))].copy()

Per-Method Failure Modes

HLA LOH silent invalidation

Trigger: ranking candidates without running LOHHLA. Mechanism: tumor deletes the haplotype that would present its neoantigens; upstream signals all look fine. Symptom: beautiful candidates on an absent allele. Fix: mandatory separate LOHHLA step; drop lost-allele candidates.

Subclonal mis-tiering from raw VAF

Trigger: using VAF as clonality without purity/CN correction. Mechanism: clonality needs cancer cell fraction (CCF = f(VAF, purity, local CN)). Symptom: clonal mutation in low-purity sample read as subclonal (and vice versa in amplified regions). Fix: estimate purity (ASCAT/Sequenza/PURPLE) and CCF (PyClone) before tiering.