gnomad_qc.v5.variant_qc.isolation_forest

Script to run the GATK isolation forest variant QC model on an AS-annotated sites VCF.

SNPs and INDELs are modeled separately (GATK --mode): for each variant type the workflow runs ExtractVariantAnnotations -> TrainVariantAnnotationsModel (PYTHON_IFOREST) -> ScoreVariantAnnotations (scattered over intervals). The two per-mode scored VCFs are merged into a single result HT in the --load-iforest step.

usage: gnomad_qc.v5.variant_qc.isolation_forest.py [-h] [-o]
                                                   [--environment {rwb,batch}]
                                                   --model-id MODEL_ID
                                                   [--test]
                                                   [--test-chrom TEST_CHROM [TEST_CHROM ...]]
                                                   [--test-on-v4]
                                                   [--export-v4-test-vcf]
                                                   [--export-only]
                                                   --batch-billing-project
                                                   BATCH_BILLING_PROJECT
                                                   --gcp-billing-project
                                                   GCP_BILLING_PROJECT
                                                   [--batch-suffix BATCH_SUFFIX]
                                                   [--app-name APP_NAME]
                                                   [--driver-cores DRIVER_CORES]
                                                   [--driver-memory DRIVER_MEMORY]
                                                   [--worker-cores WORKER_CORES]
                                                   [--worker-memory WORKER_MEMORY]
                                                   [--gatk-image GATK_IMAGE]
                                                   [--scatter-count SCATTER_COUNT]
                                                   --out-vcf-name OUT_VCF_NAME
                                                   [--hyperparameters-json HYPERPARAMETERS_JSON]
                                                   [--transmitted-singletons]
                                                   [--sibling-singletons]
                                                   [--adj] [--load-iforest]
                                                   [--load-only]
                                                   [--n-partitions N_PARTITIONS]
                                                   [--header-path HEADER_PATH]
                                                   [--array-elements-required]

Named Arguments

-o, --overwrite

Whether to overwrite existing outputs (with –load-only, only overwrite the result HT).

Default: False

--environment

Possible choices: rwb, batch

Compute environment.

Default: “batch”

--model-id

Model ID for the isolation forest run. Must start with ‘if’ (e.g. if_test1).

testing options

--test

Filter to –test-chrom and use a small scatter count for testing.

Default: False

--test-chrom

Contig(s) to score under –test/–test-on-v4. Must match the contigs in the test info VCF.

Default: [‘chr21’]

--test-on-v4

Use gnomAD v4 sites/true-positive VCFs as inputs (v5 data unavailable).

Default: False

--export-v4-test-vcf

Export the v4 split info HT filtered to –test-chrom as the –test-on-v4 sites VCF. Requires –test-on-v4.

Default: False

--export-only

Stop after –export-v4-test-vcf instead of continuing to the Batch.

Default: False

batch/GCP configuration

--batch-billing-project

Hail Batch billing project.

--gcp-billing-project

GCP billing project for requester-pays buckets.

--batch-suffix

String to append to the Batch name.

Default: “”

--app-name

Job name for the batch/QoB backend.

--driver-cores

Number of driver cores. The load step has been run with 4-8.

--driver-memory

Driver memory. Must be ‘highmem’ for –load-iforest/–load-only.

--worker-cores

Number of worker cores.

--worker-memory

Worker memory.

GATK scoring

--gatk-image

GATK docker image.

Default: “us.gcr.io/broad-gatk/gatk:4.6.1.0”

--scatter-count

Number of intervals to scatter scoring across.

Default: 100

--out-vcf-name

Base name for scored VCF shards.

--hyperparameters-json

Optional GATK model hyperparameters JSON path.

training/calibration set

--transmitted-singletons

Include transmitted singletons as a training/calibration set.

Default: False

--sibling-singletons

Include sibling singletons as a training/calibration set.

Default: False

--adj

Use adj genotypes for the true-positive singletons VCF.

Default: False

load results HT

--load-iforest

Merge the scored VCFs into a variant QC result HT.

Default: False

--load-only

Skip the GATK Batch and only run the load step.

Default: False

--n-partitions

Number of partitions for the result HT.

Default: 5000

--header-path

Optional header file to use when importing the scored VCFs.

--array-elements-required

Set array_elements_required=True when importing the scored VCFs.

Default: False

Module Functions

gnomad_qc.v5.variant_qc.isolation_forest.DEFAULT_GATK_IMAGE

GATK image.

gnomad_qc.v5.variant_qc.isolation_forest.extract_variant_annotations_job(b, ...)

Run GATK ExtractVariantAnnotations for one variant type mode.

gnomad_qc.v5.variant_qc.isolation_forest.train_variant_annotations_model_job(b, ...)

Run GATK TrainVariantAnnotationsModel (PYTHON_IFOREST) for one variant type mode.

gnomad_qc.v5.variant_qc.isolation_forest.split_intervals_job(b, ...)

Split the calling contigs into scatter_count intervals, excluding telomeres/centromeres.

gnomad_qc.v5.variant_qc.isolation_forest.score_variant_annotations_job(b, ...)

Run GATK ScoreVariantAnnotations for one variant type mode over one interval.

gnomad_qc.v5.variant_qc.isolation_forest.check_run_config(...)

Guard reuse of outputs under run_prefix against a changed run configuration.

gnomad_qc.v5.variant_qc.isolation_forest.isolation_forest_workflow(b, ...)

Add the per-mode extract/train/score jobs (SNP and INDEL) to the Batch.

gnomad_qc.v5.variant_qc.isolation_forest.export_exclusion_intervals(...)

Export telomere/centromere intervals to a GATK-readable chr:start-end file.

gnomad_qc.v5.variant_qc.isolation_forest.export_v4_test_sites_vcf(...)

Export the v4 split info HT, filtered to contigs, as a GATK-readable sites VCF.

gnomad_qc.v5.variant_qc.isolation_forest.run_batch(b, ...)

Run a Batch to completion and raise if it did not succeed.

gnomad_qc.v5.variant_qc.isolation_forest.merge_iforest_result(...)

Merge the per-mode scored VCFs into a single variant QC result HT.

gnomad_qc.v5.variant_qc.isolation_forest.reconcile_scored_sites(...)

Raise if any input variant in the scored region is missing from the output or unscored.

gnomad_qc.v5.variant_qc.isolation_forest.main(args)

Run the isolation forest variant QC workflow.

gnomad_qc.v5.variant_qc.isolation_forest.get_script_argument_parser()

Get script argument parser.

Script to run the GATK isolation forest variant QC model on an AS-annotated sites VCF.

SNPs and INDELs are modeled separately (GATK --mode): for each variant type the workflow runs ExtractVariantAnnotations -> TrainVariantAnnotationsModel (PYTHON_IFOREST) -> ScoreVariantAnnotations (scattered over intervals). The two per-mode scored VCFs are merged into a single result HT in the --load-iforest step.

gnomad_qc.v5.variant_qc.isolation_forest.DEFAULT_GATK_IMAGE = 'us.gcr.io/broad-gatk/gatk:4.6.1.0'

GATK image. 4.6.1.0 ships the conda env (scikit-learn, dill) needed for PYTHON_IFOREST.

gnomad_qc.v5.variant_qc.isolation_forest.extract_variant_annotations_job(b, mode, sites_only_vcf, features, resource_args, calling_intervals_args, exclude_intervals, out_root, gatk_image, gcp_billing_project)[source]

Run GATK ExtractVariantAnnotations for one variant type mode.

Parameters:
  • b (Batch) – Batch to add the job to.

  • mode (str) – Variant type mode, ‘SNP’ or ‘INDEL’.

  • sites_only_vcf (str) – AS-annotated sites-only input VCF.

  • features (List[str]) – Features to extract for this mode.

  • resource_args (str) – GATK --resource args for the labeled sets.

  • calling_intervals_args (str) – GATK -L args restricting the extracted sites.

  • exclude_intervals (str) – Telomere/centromere intervals excluded from the training/calibration set (-XL); these sites are not scored either.

  • out_root (str) – Output prefix (files written as {out_root}.*).

  • gatk_image (str) – GATK docker image.

  • gcp_billing_project (str) – GCP billing project for requester-pays buckets.

Return type:

Job

Returns:

Job with output ResourceGroup extract.

gnomad_qc.v5.variant_qc.isolation_forest.train_variant_annotations_model_job(b, mode, annotations_hdf5, out_root, gatk_image, gcp_billing_project, hyperparameters_json=None)[source]

Run GATK TrainVariantAnnotationsModel (PYTHON_IFOREST) for one variant type mode.

Parameters:
  • b (Batch) – Batch to add the job to.

  • mode (str) – Variant type mode, ‘SNP’ or ‘INDEL’.

  • annotations_hdf5 (ResourceFile) – Extracted annotations HDF5 from the extract step.

  • out_root (str) – Output prefix (files written as {out_root}.{mode}.*).

  • gatk_image (str) – GATK docker image.

  • gcp_billing_project (str) – GCP billing project for requester-pays buckets.

  • hyperparameters_json (Optional[str]) – Optional model hyperparameters JSON path.

Return type:

Job

Returns:

Job with output ResourceGroup model.

gnomad_qc.v5.variant_qc.isolation_forest.split_intervals_job(b, calling_intervals_args, exclude_intervals, scatter_count, gatk_image, gcp_billing_project)[source]

Split the calling contigs into scatter_count intervals, excluding telomeres/centromeres.

Parameters:
  • b (Batch) – Batch to add the job to.

  • calling_intervals_args (str) – GATK -L args covering the calling contigs.

  • exclude_intervals (str) – Path to an intervals file to exclude (-XL).

  • scatter_count (int) – Number of interval shards to produce.

  • gatk_image (str) – GATK docker image.

  • gcp_billing_project (str) – GCP billing project for requester-pays buckets.

Return type:

Job

Returns:

Job with output ResourceGroup intervals.

gnomad_qc.v5.variant_qc.isolation_forest.score_variant_annotations_job(b, mode, sites_only_vcf, features, extracted_vcf, model, resource_args, interval, idx, out_root, gatk_image, gcp_billing_project)[source]

Run GATK ScoreVariantAnnotations for one variant type mode over one interval.

Parameters:
  • b (Batch) – Batch to add the job to.

  • mode (str) – Variant type mode, ‘SNP’ or ‘INDEL’.

  • sites_only_vcf (str) – AS-annotated sites-only input VCF.

  • features (List[str]) – Features for this mode (must match the trained model).

  • extracted_vcf (ResourceFile) – Extracted training/calibration VCF from the extract step.

  • model (Union[ResourceGroup, str]) – --model-prefix: the trained model ResourceGroup from the train job, or a prefix string built from a read-in group (see isolation_forest_workflow).

  • resource_args (str) – GATK --resource args for the labeled sets.

  • interval (ResourceFile) – Interval file to restrict scoring to.

  • idx (int) – Shard index, used to identify the job in the Batch.

  • out_root (str) – Output prefix (files written as {out_root}.*).

  • gatk_image (str) – GATK docker image.

  • gcp_billing_project (str) – GCP billing project for requester-pays buckets.

Return type:

Job

Returns:

Job with output ResourceGroup output_score.

Note

GATK traversal returns every record overlapping the shard interval, and the default output mode (ANYWHERE) writes all of them, so a deletion straddling a shard boundary is emitted by both neighbors. --variant-output-filtering STARTS_IN writes a record only if its start is in the interval, making emission exactly-once across shards. It also makes -XL start-based, which is the definition reconcile_scored_sites uses.

gnomad_qc.v5.variant_qc.isolation_forest.check_run_config(run_prefix, config, overwrite, write)[source]

Guard reuse of outputs under run_prefix against a changed run configuration.

Extract/model/score outputs are reused by path alone, so a rerun with different contigs, scatter count, or inputs would silently reuse region- or split-mismatched files. Compare against the config recorded by the first run and raise on any difference unless overwrite.

Parameters:
  • run_prefix (str) – GCS prefix for this run’s GATK outputs.

  • config (Dict[str, Any]) – Current run configuration.

  • overwrite (bool) – Whether outputs are being overwritten. Skips the comparison only in the Batch step; –load-only reuses whatever is under run_prefix (its –overwrite only replaces the result HT), so it always compares.

  • write (bool) – Whether to record config (the Batch step; the load step only compares).

Return type:

None

Returns:

None.

gnomad_qc.v5.variant_qc.isolation_forest.isolation_forest_workflow(b, sites_only_vcf, run_prefix, out_vcf_name, singletons_vcf, calling_intervals_args, exclude_intervals, scatter_count, gatk_image, gcp_billing_project, hyperparameters_json=None, overwrite=False)[source]

Add the per-mode extract/train/score jobs (SNP and INDEL) to the Batch.

Extract and train outputs are reused if they already exist (unless overwrite).

Parameters:
  • b (Batch) – Batch to add jobs to.

  • sites_only_vcf (str) – AS-annotated sites-only input VCF.

  • run_prefix (str) – GCS prefix for this run’s GATK outputs.

  • out_vcf_name (str) – Base name for scored VCF shards.

  • singletons_vcf (Optional[str]) – Optional true-positive singletons VCF.

  • calling_intervals_args (str) – GATK -L args for the calling contigs.

  • exclude_intervals (str) – Path to the telomere/centromere exclusion intervals file.

  • scatter_count (int) – Number of interval shards for scoring.

  • gatk_image (str) – GATK docker image.

  • gcp_billing_project (str) – GCP billing project for requester-pays buckets.

  • hyperparameters_json (Optional[str]) – Optional model hyperparameters JSON path.

  • overwrite (bool) – Whether to rerun extract/train/score even if outputs exist.

Return type:

None

Returns:

None.

gnomad_qc.v5.variant_qc.isolation_forest.export_exclusion_intervals(out_path)[source]

Export telomere/centromere intervals to a GATK-readable chr:start-end file.

Written to a stable path (not a temp file) so the same file used by GATK -XL can be read back during reconciliation to define the excluded region exactly.

telomeres_and_centromeres is a BED import, so its intervals are half-open [start, end), while GATK reads chr:start-end as 1-based inclusive on both ends. Emit end - 1 so the file covers exactly the interval.

Parameters:

out_path (str) – Destination path for the intervals file.

Return type:

str

Returns:

out_path.

gnomad_qc.v5.variant_qc.isolation_forest.export_v4_test_sites_vcf(contigs, out_path)[source]

Export the v4 split info HT, filtered to contigs, as a GATK-readable sites VCF.

v4’s --export-split-info-vcf output is misnamed: it is exported from the unsplit info HT, so it cannot stand in for the split v5 info VCF. Read the split info HT directly instead.

Parameters:
  • contigs (List[str]) – Contigs to keep.

  • out_path (str) – Destination VCF path.

Return type:

None

Returns:

None.

gnomad_qc.v5.variant_qc.isolation_forest.run_batch(b, description)[source]

Run a Batch to completion and raise if it did not succeed.

Batch.run waits but only prints the final state, so a failed batch would otherwise fall through to the next step.

Parameters:
  • b (Batch) – Batch to run.

  • description (str) – Batch description used in the error message.

Return type:

None

Returns:

None.

gnomad_qc.v5.variant_qc.isolation_forest.merge_iforest_result(run_prefix, out_vcf_name, model_id, scatter_count, n_partitions, header_path, array_elements_required, transmitted_singletons, sibling_singletons, adj)[source]

Merge the per-mode scored VCFs into a single variant QC result HT.

The SNP-pass VCF contributes SNP alleles and the INDEL-pass VCF contributes all non-SNP alleles, so the union covers every scored allele exactly once.

Only the current run’s scatter_count shards are read (by explicit index), so a rerun with a smaller scatter count does not pick up stale shards from a prior run.

Parameters:
  • run_prefix (str) – GCS prefix for this run’s GATK outputs.

  • out_vcf_name (str) – Base name used for scored VCF shards.

  • model_id (str) – Isolation forest model ID (starts with if_).

  • scatter_count (int) – Number of score shards per mode (must match the scoring run).

  • n_partitions (int) – Number of partitions for the imported HTs.

  • header_path (Optional[str]) – Optional VCF header file for import.

  • array_elements_required (bool) – Value passed to hl.import_vcf.

  • transmitted_singletons (bool) – Whether transmitted singletons were a training set.

  • sibling_singletons (bool) – Whether sibling singletons were a training set.

  • adj (bool) – Whether the singletons VCF used adj genotypes.

Return type:

Table

Returns:

Merged variant QC result HT.

gnomad_qc.v5.variant_qc.isolation_forest.reconcile_scored_sites(sites_only_vcf, contigs, exclude_intervals_path, scored_ht, array_elements_required=False)[source]

Raise if any input variant in the scored region is missing from the output or unscored.

GATK writes every input record within its interval to the scored VCF, so any variant absent from the merged output was dropped (e.g. a truncated shard). GATK classifies every allele as either SNP or INDEL, so every variant should also carry a SCORE.

Input and output are both split, so variants are compared on the full key.

Parameters:
  • sites_only_vcf (str) – Input sites VCF that was scored.

  • contigs (List[str]) – Calling contigs passed to GATK (-L).

  • exclude_intervals_path (str) – Telomere/centromere exclusion file passed to GATK (-XL).

  • scored_ht (Table) – Merged scored result HT.

  • array_elements_required (bool) – Value passed to hl.import_vcf for the input.

Return type:

None

Returns:

None.

gnomad_qc.v5.variant_qc.isolation_forest.main(args)[source]

Run the isolation forest variant QC workflow.