Bio Applied Hla Typing

HLA Typing and Antigen Presentation

When to Use

• Typing HLA class I (-A/-B/-C) or class II (-DRB1/-DQ/-DP) alleles from WGS, WES, or RNA-seq reads
• Predicting peptide-MHC binding affinity for a peptide list or tumor somatic mutations (neoantigen prioritization, pVACseq-style pipelines)
• Detecting tumor HLA loss of heterozygosity (LOH) as an immune-escape mechanism
• Screening a patient's HLA genotype for known drug-hypersensitivity or disease associations (e.g. HLA-B57:01/abacavir, HLA-B15:02/carbamazepine)
• Summarizing HLA allele frequencies / homozygosity across a cohort

Version Compatibility

OptiType 1.3.5, arcasHLA 0.6.0, HLA-LA 1.0.3, NetMHCpan 4.1 / NetMHCIIpan 4.3, pVACtools ≥4.0, Python ≥3.10, pandas ≥2.0, numpy ≥1.26, matplotlib ≥3.8.

Prerequisites

pip install pandas numpy matplotlib. OptiType, arcasHLA, HLA-LA, and NetMHCpan are separate CLI tools (install via conda/bioconda or Docker/Singularity; NetMHCpan requires a license from DTU Health Tech). Prior concepts: bio-variant-calling-vcf-basics (somatic VCF for neoantigen input), bio-read-alignment-bwa-alignment (BAM inputs for typing).

HLA Typing from NGS Reads

Goal: Call 4-digit HLA-A/B/C (and optionally class II) genotypes from a BAM or FASTQ. Approach: Use a read-based typer against the IMGT/HLA reference; RNA-seq samples use arcasHLA, DNA (WGS/WES) samples use OptiType or HLA-LA.

# OptiType: DNA or RNA reads -> 4-digit class I genotype (ILP over HLA-mapped reads)
OptiTypePipeline.py -i tumor_R1.fastq tumor_R2.fastq \
    --dna -v -o optitype_out/
# Output: optitype_out/<timestamp>/<sample>_result.tsv
#   A*02:01,A*03:01,B*07:02,B*44:02,C*05:01,C*07:02

# arcasHLA: RNA-seq BAM -> class I + II genotype, good when only RNA-seq exists
arcasHLA extract --unmapped -o hla_reads/ tumor_rna.bam
arcasHLA genotype hla_reads/sample.extracted.1.fq.gz \
                   hla_reads/sample.extracted.2.fq.gz \
                   -g A,B,C,DPB1,DQB1,DQA1,DRB1 \
                   -o hla_typing/
# Output: hla_typing/sample.genotype.json

Tumor HLA LOH is called by comparing tumor vs matched-normal genotypes (or allele-specific copy number with a tool like LOHHLA), not by typing the tumor alone — apparent homozygosity in tumor-only RNA-seq can also come from allele-specific expression, not true LOH.

Cohort HLA Summary and Neoantigen Ranking

Goal: Summarize HLA allele frequencies/homozygosity across a cohort, and rank candidate neoantigens by predicted binding affinity. Approach: Parse per-sample typing calls into a DataFrame for cohort-level stats; parse NetMHCpan -BA output into a ranked table filtered by %Rank_EL/IC50 thresholds.

import numpy as np
import pandas as pd


def simulate_hla_cohort(n_patients=20, seed=42):
    """Build a toy HLA-A/B/C genotype table for n_patients (stand-in for
    parsing real OptiType/arcasHLA result.tsv files into one DataFrame).
    Allele pools/frequencies are approximate European population values.
    """
    rng = np.random.default_rng(seed)
    hla_a_pool = ['A*02:01', 'A*01:01', 'A*03:01', 'A*24:02', 'A*11:01',
                  'A*29:02', 'A*23:01', 'A*26:01', 'A*31:01', 'A*32:01']
    hla_b_pool = ['B*07:02', 'B*08:01', 'B*44:02', 'B*44:03', 'B*35:01',
                  'B*51:01', 'B*40:01', 'B*15:01', 'B*18:01', 'B*57:01']
    hla_c_pool = ['C*07:01', 'C*07:02', 'C*03:04', 'C*05:01', 'C*04:01',
                  'C*06:02', 'C*01:02', 'C*02:02', 'C*08:02', 'C*16:01']
    a_probs = np.array([0.283, 0.161, 0.143, 0.098, 0.072,
                         0.054, 0.041, 0.038, 0.031, 0.026])
    a_probs /= a_probs.sum()
    b_probs = np.array([0.126, 0.099, 0.093, 0.068, 0.061,
                         0.058, 0.048, 0.043, 0.039, 0.036])
    b_probs /= b_probs.sum()

    rows = []
    for i in range(n_patients):
        rows.append({
            'Patient': f'P{i+1:02d}',
            'HLA-A1': rng.choice(hla_a_pool, p=a_probs),
            'HLA-A2': rng.choice(hla_a_pool, p=a_probs),
            'HLA-B1': rng.choice(hla_b_pool, p=b_probs),
            'HLA-B2': rng.choice(hla_b_pool, p=b_probs),
            'HLA-C1': rng.choice(hla_c_pool),
            'HLA-C2': rng.choice(hla_c_pool),
        })
    df = pd.DataFrame(rows)
    df['A_homozygous'] = df['HLA-A1'] == df['HLA-A2']
    df['B_homozygous'] = df['HLA-B1'] == df['HLA-B2']
    return df


def rank_neoantigens(mutations, ic50_nM, strong_cutoff=50, weak_cutoff=500):
    """Rank candidate neoantigens by NetMHCpan-style binding affinity.

    mutations : list[str] e.g. ['KRAS_G12D', 'TP53_R175H', ...]
    ic50_nM   : list[float] predicted IC50 (nM) per mutation, from the
                Affinity(nM) column of `netMHCpan ... -BA` output.
    Returns a DataFrame sorted best-binder-first with a Binding_Level call.
    In production, use NetMHCpan's own %Rank_EL column rather than deriving
    a rank from IC50 (the two scales are not linearly related).
    """
    df = pd.DataFrame({'mutation': mutations, 'IC50_nM': ic50_nM})
    df['Binding_Level'] = pd.cut(
        df['IC50_nM'], bins=[0, strong_cutoff, weak_cutoff, np.inf],
        labels=['Strong Binder', 'Weak Binder', 'Non-binder'],
    )
    return df.sort_values('IC50_nM').reset_index(drop=True)


hla_df = simulate_hla_cohort()
candidate_neoantigens = rank_neoantigens(
    mutations=['KRAS_G12D', 'TP53_R175H', 'EGFR_L858R', 'PIK3CA_H1047R',
               'BRAF_V600E', 'NRAS_Q61K', 'IDH1_R132H', 'CTNNB1_S45F'],
    ic50_nM=[12, 340, 28, 890, 45, 210, 1500, 60],
)

print(hla_df[['Patient', 'HLA-A1', 'HLA-A2', 'HLA-B1', 'HLA-B2']].head(5).to_string(index=False))
print(f"\nA homozygous: {hla_df['A_homozygous'].sum()}/{len(hla_df)}")
print(f"\n{candidate_neoantigens}")

Visualizing HLA and Neoantigen Results

Goal: Produce one figure combining a neoantigen binding waterfall, cohort allele distribution, and a pharmacogenomic frequency check. Approach: Three matplotlib panels driven by the DataFrames above.

import matplotlib.pyplot as plt
import matplotlib.patches as mpatches


def plot_hla_neoantigen_summary(hla_df, candidate_neoantigens, out_png=None):
    """Three-panel HLA/neoantigen summary figure.

    Panel 1: neoantigens ranked by predicted binding (lower IC50 = tighter).
    Panel 2: HLA-A allele distribution across the cohort.
    Panel 3: HLA-B*57:01 population frequency (abacavir screen relevance).
    """
    colors = {'Strong Binder': 'firebrick', 'Weak Binder': 'orange', 'Non-binder': 'lightgray'}
    fig, axes = plt.subplots(1, 3, figsize=(16, 5))

    sorted_df = candidate_neoantigens.sort_values('IC50_nM')
    bar_colors = [colors[b] for b in sorted_df['Binding_Level']]
    axes[0].barh(sorted_df['mutation'], -np.log10(sorted_df['IC50_nM'] + 1),
                 color=bar_colors, edgecolor='black', linewidth=0.5)
    axes[0].axvline(-np.log10(51), color='red', linestyle='--')
    axes[0].axvline(-np.log10(501), color='orange', linestyle='--')
    axes[0].set_xlabel('-log10(IC50 nM)')
    axes[0].set_title('Neoantigen MHC-I Binding')
    axes[0].legend(handles=[mpatches.Patch(color=v, label=k) for k, v in colors.items()],
                   fontsize=8, loc='lower right')

    all_a = pd.concat([hla_df['HLA-A1'], hla_df['HLA-A2']]).value_counts()
    axes[1].bar(all_a.index, all_a.values, color='steelblue', edgecolor='black', linewidth=0.5)
    axes[1].set_xticklabels(all_a.index, rotation=45, ha='right', fontsize=8)
    axes[1].set_title(f'HLA-A Allele Distribution (n={len(hla_df)})')

    populations = ['European', 'East Asian', 'South Asian', 'African', 'Latin American']
    b5701_freq = [3.6, 0.2, 1.8, 0.1, 1.2]  # % carriers, approx. from AFND
    axes[2].bar(populations, b5701_freq, color='mediumpurple', edgecolor='black', linewidth=0.5)
    axes[2].axhline(1.0, color='red', linestyle='--', label='1% threshold')
    axes[2].set_xticklabels(populations, rotation=30, ha='right', fontsize=9)
    axes[2].set_ylabel('HLA-B*57:01 carriers (%)')
    axes[2].set_title('Abacavir Pharmacogenomic Screen')
    axes[2].legend(fontsize=9)

    plt.tight_layout()
    if out_png:
        plt.savefig(out_png, dpi=120, bbox_inches='tight')
    return fig

HLA Biology Reference

• Class I (HLA-A/-B/-C): all nucleated cells; presents 8–11 aa intracellular peptides (proteasome → TAP → ER loading) to CD8+ T cells.
• Class II (HLA-DR/-DQ/-DP): professional APCs only; presents 13–25 aa extracellular/endosomal peptides to CD4+ T cells.
• Nomenclature: A*02:01 = gene A, field 1 (02) = allele group/serotype, field 2 (01) = protein sequence; 6/8-digit fields add synonymous/non-coding variants. HLA-A*02:01 is the most common HLA-A allele in Europeans (~28%).
• Disease/PGx associations: HLA-B57:01 → abacavir hypersensitivity (OR >1000, mandatory pre-screen); HLA-B15:02 → carbamazepine SJS/TEN (~OR 80); HLA-DQ2/DQ8 → celiac disease; HLA-A*02:01 → improved melanoma immunotherapy response via neoantigen presentation.
• Tumor HLA LOH: ~40% of NSCLC tumors show HLA LOH as an immune-escape route (McGranahan et al. 2017, TRACERx).

Pitfalls

• Reference bias in typing: standard linear reference genomes collapse HLA diversity in the hyper-polymorphic exons 2–3; always type against an HLA-specific reference (IMGT/HLA), not the primary GRCh38 HLA region.
• Coverage at the HLA locus: <30x at HLA-A/B/C in WES/WGS causes false homozygosity calls; check locus-level depth before trusting a typing result.
• RNA-seq apparent LOH: allele-specific expression can make a tumor look HLA-homozygous in RNA-seq even when both alleles are present in DNA — confirm LOH with matched-normal DNA typing or a copy-number-aware tool (e.g. LOHHLA), not RNA-seq alone.
• %Rank_EL vs %Rank_BA: NetMHCpan reports both an eluted-ligand rank and a binding-affinity rank; they use different score distributions — don't apply a %Rank_EL < 0.5% cutoff to a %Rank_BA column or vice versa.
• Class II is less reliable: NetMHCIIpan has fewer training epitopes and an open-ended binding groove, so predictions are noisier than class I — treat class II neoantigen calls as lower confidence.
• Nomenclature mismatches: typing tools, NetMHCpan (HLA-A02:01, no *), and IMGT/HLA (A*02:01) each expect slightly different allele string formats — normalize before joining tables.

See Also

• bio-clinical-databases-hla-typing
• bio-immunoinformatics-mhc-binding-prediction
• bio-immunoinformatics-neoantigen-prediction
• bio-workflows-neoantigen-pipeline