Haplotype Inference

Pipeline step: compute_haploid_genotypes (Step 5)

Module: grid.utils.hi_inference

Diploid copy number estimates from Step 4 (denoted IRR — individual repeat ratio) cannot be decomposed into per-haplotype contributions without phase information. GRiD infers haplotype-specific copy numbers iteratively by identifying haplotype-matched neighbors, following Hujoel et al. (2026).

Two neighbor sources are supported via the method config parameter:

  • IBS (method: "ibs") — neighbors from computeIBSpbwt PBWT matching output. Each neighbor is weighted equally.

  • IBD (method: "ibd") — neighbors from iLASH IBD segment output. Set weighted: True to apply a Lorentzian distance+match weight to each segment.

For more information on the neighbor file formats, see Computing IBS/IBD Neighbors.

Note

This step is locus-agnostic — it operates on diploid copy number estimates and a haplotype neighbors file, with no assumptions about which locus produced those estimates. It was developed and validated on LPA KIV-2, but can be applied to any VNTR locus for which phased neighbor data are available.

Notation

For \(N\) individuals, define:

  • \(\mathrm{IRR}_i\) — diploid copy number estimate for individual \(i\)

  • \(h_{i,k}\) — haplotype copy number for individual \(i\), haplotype \(k \in \{1, 2\}\)

  • \(\mathcal{H}_{i,k}\) — set of IBD-matched haplotype neighbors for haplotype \(k\) of individual \(i\)

The constraint \(h_{i,1} + h_{i,2} = \mathrm{IRR}_i\) must hold at every iteration.

Initialization

Each haplotype is initialized to half the diploid estimate:

\[h_{i,1}^{(0)} = h_{i,2}^{(0)} = \frac{\mathrm{IRR}_i}{2}\]

Individuals with fewer than \(\mathrm{MIN\_NBR}\) IBD neighbors on either haplotype are excluded from iterative phasing and receive mean-imputed values.

Iterative Update

At each iteration \(t\), a weighted mean haplotype copy number across matched neighbors is computed for each haplotype:

\[\bar{h}^{(t)}_{\mathcal{H}_{i,k}} = \frac{\displaystyle\sum_{(j,\, k') \in \mathcal{H}_{i,k}} w_{j,k'} \cdot h^{(t)}_{j,\, k'}} {\displaystyle\sum_{(j,\, k') \in \mathcal{H}_{i,k}} w_{j,k'}}\]

where \((j, k')\) denotes haplotype \(k'\) of neighbor \(j\) matched to haplotype \(k\) of individual \(i\).

IBS method: all weights \(w_{j,k'} = 1\).

IBD method (unweighted): all weights \(w_{j,k'} = 1\).

IBD method (weighted): the weight accounts for segment proximity to the target region and IBD match quality:

\[w_{j,k'} = \frac{\lambda}{\delta_{j,k'} + \lambda} \cdot m_{j,k'}\]

where \(\delta_{j,k'}\) is the physical distance (bp) from the IBD segment to the target region (0 if overlapping), \(m_{j,k'}\) is the iLASH match score, and \(\lambda\) is the weight_scale parameter (default 1,000,000 bp).

Each haplotype is updated proportionally, preserving the diploid constraint:

\[h^{(t+1)}_{i,k} = \mathrm{IRR}_i \cdot \frac{\bar{h}^{(t)}_{\mathcal{H}_{i,k}}} {\bar{h}^{(t)}_{\mathcal{H}_{i,1}} + \bar{h}^{(t)}_{\mathcal{H}_{i,2}}}\]

This redistributes the diploid total between haplotypes in proportion to their neighbors’ copy numbers. In practice the algorithm runs for a fixed number of iterations (default: \(T = 100\)).

Mean Imputation for Unphased Samples

For individuals with insufficient IBD neighbors, haplotype copy numbers are mean-imputed using the cohort mean on each haplotype:

\[h^{\mathrm{imp}}_{i,k} = \bar{h}^{(T)}_{\mathcal{H}_{i,k}}\]

where the mean is taken over any available neighbors (or the cohort mean if no neighbors exist).

Output

The output file contains one row per individual with six columns:

Column

Description

ID

Sample identifier

IRRs

Original diploid copy number estimate

hap1phased

Inferred haplotype 1 copy number (iterative phasing)

hap2phased

Inferred haplotype 2 copy number (iterative phasing)

hap1imp

Mean-imputed haplotype 1 copy number (neighbor mean)

hap2imp

Mean-imputed haplotype 2 copy number (neighbor mean)