3. Detect HRD

Homologous-recombination deficiency leaves marks in three places: Signature 3 in the point mutations, loss of heterozygosity, telomeric allelic imbalance and large-scale state transitions in the copy-number profile, and microhomology at the breakpoints of deletions. SigMA2's HRD model reads these marks and predicts, per tumor, whether its mutations and copy number resemble those of tumors with biallelic BRCA1/2 inactivation: HRD, HRD-low or HRP. It needs the tumor's SNVs and indels and uses copy-number segments when you give them; what a call means states what a call supports and what it does not establish. A Signature 3 call answers another question, whether the point mutations carry Signature 3, and is reported beside the HRD call, not in its place:

Get the example files

The MMRD/POLE step an HRD run starts with needs the repeat-locus files (sigma2 data fetch --all, or from the public SigMA (R) repository); without them the run skips that step without a message. The deletion features need only the genome.

curl -LO https://raw.githubusercontent.com/parklab/SigMA/master/inst/extdata/examples/test_mutations_50sample.maf
curl -LO https://sigscape.org/sigma2/examples/msk_breast_50_segments.seg

Run the HRD model

import sigma2 as sg2

hrd = sg2.run(
    input="test_mutations_50sample.maf", ref_genome="hg19",
    tumor_type="breast", seq_type="msk", soi="HRD",
    segments="msk_breast_50_segments.seg",
)
hrd.df[[
    "total_snvs",
    "sig3_score",
    "fga",
    "lst_lr",
    "mmej",
    "nhej",
    "hrd_variant",
    "hrd_score",
    "hrd_call",
    "MMRD_category",
]].head(3)

                   total_snvs  sig3_score       fga  lst_lr  mmej  nhej hrd_variant  hrd_score hrd_call MMRD_category
P-0001484-T02-IM5        <NA>         NaN  0.006771     0.0   NaN   NaN         tcn   0.001089      HRP           NaN
P-0001939-T02-IM5          13    0.219434  0.277629     0.0   0.0   0.0         tcn   0.010756      HRP          MMRP
P-0003375-T01-IM5          10    0.104040  0.334849     0.0   0.0   0.0         tcn   0.003901      HRP          MMRP

Provide the three inputs

Three sources, each optional when given another way:

The model variant follows what is present:

hrd0 = sg2.run(
    input="test_mutations_50sample.maf",
    ref_genome="hg19",
    tumor_type="breast",
    seq_type="msk",
    soi="HRD",
    verbose=1,
)

HRD: no copy-number input (segments= or gis=); scoring with the model trained without copy number.

The gis variant reads the three GIS components; a GIS sum alone (gis=) falls back to a weaker variant with a warning. On msk, segments are read on the model's IMPACT gene grid, which is hg19, and hg38 segments are refused. The refusal suggests passing a GIS sum as gis=, but the msk model's gis variant reads the three components, so a sum alone sends the run to the variant without copy number. To keep copy number, lift the segments to hg19, compute the three components on the grid with sigma2.hrd.genomic_instability(lifted, genome="hg19", grid=sg2.hrd_model(seq_type="msk").gis_grid), write hrd_loh, tai, lst and gis into the table sg2.hrd_features() returns, and score it with sg2.score_hrd() (Optional: run in two steps); on hg19 input this route gives the same scores as passing the segments.

Read the output

An HRD run keeps every sample in out.df. Below the SNV floor a tumor has no Signature 3 features and no MMRD/POLE category, and the output shows its deletion counts as NA, although the model read them. It still receives a call: with segments, from its copy number and deletion counts; without segments, from its deletion counts alone, or from no feature of the tumor at all, when its HRP is the model's prior for an empty row. Report such tumors as not assessable and leave them out of an HRD prevalence; in the example, 12 of the 50 panels are below the floor and all 12 read HRP. segments= restricts nothing either: a tumor without segments is scored without copy number, with its own hrd_variant.

hrd.df["hrd_call"].value_counts()

hrd_call
HRP        41
HRD-low     4
HRD         4
MMRD        1

Cohorts. An HRD run costs about a second per tumor on one worker, most of it the MMRD/POLE pre-step. Pass mxw="max" to use every core:

sg2.run(
    input="cohort.maf",
    ref_genome="hg19",
    tumor_type="breast",
    seq_type="tcga_mc3",
    soi="HRD",
    segments="ascat.tsv",
    mxw="max",
)

Note HRD-low means the two models disagree: one of the pan-cancer and the tumor-type model scores the tumor above 0.5 and the other does not. The result is indeterminate, and an orthogonal assay settles it; the standard columns count it as a loose positive (pass_mva True, confidence low), so read hrd_call. A Signature 3 call with HRP is a discordant case: the spectrum carries Signature 3 while the model, reading the deletions and copy number as well, does not place the tumor with the BRCA1/2-deficient ones; look at the spectrum and the copy number. MMRD and POLE tumors keep a score, but their call is replaced: a hypermutated spectrum's Signature 3 features mean something else, so the pre-step (Detect MMRD or POLE) runs first, and its label becomes hrd_call; ignore hrd_score for them. The pre-step needs the indels and the repeat-locus files. Without the files it is skipped with no warning, the table has no MMRD_category column, and a hypermutated tumor gets an HRD call like any other: check for that column.

The spectrum of the panel the model scores highest, drawn by the package:

from sigma2 import plot as sp

scored = hrd.df[hrd.df["total_snvs"] >= 5]
top = scored["hrd_score"].idxmax()
row = hrd.df.loc[top]
sp.sbs(
    hrd.df.loc[[top]].iloc[:, :96],
    title=f"{top}: {row['hrd_call']}, hrd_score {row['hrd_score']:.2f}",
    path="top_hrd.svg",
)

Check that the tumor type is scored

An HRD model is shipped for msk, fo, op, seqcap and tcga_mc3. Within each, a tumor type has its own model when it had more than ten positives in the first training step, its BRCA1/2-deficient tumors plus three neighbors of each, so an own model can rest on as few as three BRCA1/2-deficient tumors; the others share the pooled one, whose call rests on other tumor types' data. The bundled classifiers page gives, per tumor type and platform, the model and its number of BRCA1/2-deficient training tumors. hrd_group names which scored your sample; the bundled classifiers page lists own and pooled per platform.

Some tumor types the models saw in training have no Signature 3 signature set and cannot be scored; run(soi="HRD") refuses them with the reason:

Tumor types absent from the HRD training (bone_other, ewing, liposarc, lms, medullo, myeloma, oligo, osteo, panc_en, piloastro) are refused too.

m = sg2.hrd_model("msk")
sorted(m.covered_tumor_types()), sorted(m.unscorable_tumor_types())

hrd_model is called with its defaults. All parameters of hrd_model

HRD models are cross-validated and shipped without a release gate or a weak tier. Their caveats: a pooled type's call rests on other types' data; the exome-derived platforms (fo, op, seqcap, tcga_mc3) share one training cohort; without copy number the weaker nogis variant runs; the msk gene grid is hg19. Reference: what each HRD model was trained on gives the cohorts and their truth.

Choose between SBS3 and HRD

Report hrd_call with a Signature 3 call beside it, and name the classifier: in an HRD run SBS3_pass_mva is SigMA v1's, the one the model reads; a soi="SBS3" run gives SigMA2's for the cell. Neither is a clinical determination on its own.

Optional: run in two steps

For a feature table you can inspect or score more than once:

feats = sg2.hrd_features(
    "test_mutations_50sample.maf",
    tumor_type="breast",
    seq_type="msk",
    segments="msk_breast_50_segments.seg",
    ref_genome="hg19",
)
scores = sg2.score_hrd(feats, seq_type="msk", tumor_type="breast")

sg2.genomic_instability(segments, genome="hg19") computes HRD-LOH, TAI and LST from allele-specific segments on their own; sg2.copy_number_instability(segments) the log-ratio features. HRD on exomes runs both on ASCAT segments.

Optional: run from the command line

sigma2 hrd score --input test_mutations_50sample.maf --tumor-type breast \
    --seq-type msk --segments msk_breast_50_segments.seg -o hrd.tsv
sigma2 hrd features --input test_mutations_50sample.maf --tumor-type breast \
    --seq-type msk --segments msk_breast_50_segments.seg -o features.tsv
sigma2 hrd models

Command line lists every option.

Optional: train your own HRD model

sg2.tune_hrd(features, labels, seq_type=...) trains an HRD model from an hrd_features() table of MMR-proficient tumors with biallelic BRCA1/2 labels and a per-row tumor type, for groups with their own assay and their own truth. Its validation is cross-validation only: there is no release gate and no experimental tier for HRD models. The API reference lists the options.

Reference: what each HRD model was trained on

The models follow the two-step procedure of the preprint of Gulhan et al. (2025), which describes the TCGA exome model; the panel models apply it to other cohorts, with the truth each cohort allows. The positives are tumors with biallelic BRCA1/2 inactivation, plus three nearest wild-type neighbors of each; mismatch-repair deficient and POLE tumors were left out.

Four consequences for reading a call. On msk, an HRP call is partly learned from germline carriers labeled negative, and the model's probabilities are expressed at the class balance of the TCGA rows, whose truth is complete; its training also left the BRCA1, BRCA2 and PALB2 indels out of the deletion counts, which scoring keeps. On fo, op and seqcap, copy number from the panel itself resolves fewer breakpoints than the SNP-array segments the model learned from, so its scar counts run lower than the training's; report the call without copy number (nogis) beside it. And where SigMA v1 has no Signature 3 classifier for the tumor type and platform, the model never reads a Signature 3 score: on msk this holds for every lung, skin, head and neck and cervix tumor, and most MSK-IMPACT training panels of the other types had no score either, below SigMA v1's floor. Last, an HRD run computes its Signature 3 features with SigMA2's port of SigMA v1's classifier (sg2.hrd_model("msk").sig3 names it); where the port passes a feature as missing (Check the legacy grid), those features differ from what SigMA v1 in R gives, so compare an HRD run's SBS3_* columns with a SigMA2 legacy run, not with R.