Heterozygosity Outlier Diagnostics

Why Diagnose Heterozygosity Outliers?

Observed heterozygosity summarizes how often a sample has a heterozygous genotype across a chosen set of loci. Comparing this summary across samples can identify values that deserve quality-control follow-up. The goal is to detect unusual observations reliably, not to decide from one number that a sample must be removed. The earlier material on inbreeding coefficients from heterozygosity gives a related sample-level summary; here the emphasis is on finding candidates using robust summaries and plots.

A useful comparison requires a clearly defined marker set and consistent genotype handling. Different locus sets, missing calls, or cohort subgroups can make otherwise valid sample summaries difficult to compare. A candidate flag is therefore a prompt to inspect the data and its context, not a diagnosis of its cause.

Define the Sample Summary Before Comparing Samples

Definition: Let $g_{ij}$ denote the genotype call for sample $i$ at locus $j$, and let $J_i$ be the set of eligible loci with a nonmissing call for that sample. For diploid biallelic loci, define $I(g_{ij}\text{ is heterozygous})=1$ for a heterozygous call and 0 otherwise. The sample’s observed heterozygosity is $$ H_i=\frac{\sum_{j\in J_i} I(g_{ij}\text{ is heterozygous})}{|J_i|}. $$ It is a proportion, so it has no physical units. Its denominator is the number of eligible called loci for sample $i$.

The eligible loci should be specified before interpreting the comparison. For example, an analysis might use a consistent panel of autosomal, biallelic SNPs that passed variant-level quality filters. The same locus criteria should be applied to every sample. Missing genotypes are not heterozygous calls and do not enter the numerator; they also do not enter the denominator. A sample with few called loci may consequently have a less stable summary and should be checked alongside its sample call rate.

The calculation presumes that the calls being counted are meaningful for this purpose: genotype coding is interpreted consistently, missing values are correctly encoded, and loci with unreliable calls have been handled appropriately. It also presumes that the summaries being compared refer to a suitably comparable locus panel. These conditions do not make all samples biologically identical, and a difference in $H_i$ alone does not establish why samples differ.

Robust Summaries and an Outlier Score

The mean and standard deviation can be pulled toward extreme observations, which can make them less useful for identifying those same observations. A common robust alternative uses the median and the median absolute deviation (MAD). For sample values $H_1,\ldots,H_n$, define $$ \widetilde{H}=\operatorname{median}(H_1,\ldots,H_n), \qquad \operatorname{MAD}=\operatorname{median}_i |H_i-\widetilde{H}|. $$

To put the MAD on a scale comparable to a standard deviation under a normal distribution, multiply it by approximately 1.4826. A robust outlier score is then

R_i=\frac{H_i-\widetilde{H}}{1.4826\,\operatorname{MAD}}.

Positive scores indicate values above the cohort median; negative scores indicate values below it. A commonly used screening rule flags $|R_i|>3.5$ as a candidate outlier. This is a heuristic, not a universal biological or statistical boundary: it is not a hypothesis test, and it does not assign a p-value.

Justification for the scale factor: If $Z\sim N(0,1)$, its median is 0. The population median absolute deviation $m$ satisfies $P(|Z|\leq m)=0.5$. By symmetry, $2\Phi(m)-1=0.5$, so $\Phi(m)=0.75$ and $m\approx0.67449$. Thus dividing MAD by 0.67449, or multiplying by approximately 1.4826, estimates the standard deviation when the data are normal. This calibration explains the scale factor; it does not require treating the sample heterozygosities as normally distributed or the 3.5 rule as a formal test.

Worked Example: Flagging Extreme Values With the MAD

Suppose ten hypothetical samples have heterozygosities calculated from the same eligible marker panel, with adequate call rates for this illustration: $0.298, 0.302, 0.299, 0.301, 0.300, 0.303, 0.297, 0.301, 0.327, 0.273$.

Find the median. In sorted order the values are $0.273, 0.297, 0.298, 0.299, 0.300, 0.301, 0.301, 0.302, 0.303, 0.327$. The median is the average of the fifth and sixth values: $(0.300+0.301)/2=0.3005$.

Find the MAD. The absolute deviations from $0.3005$, sorted, are $0.0005, 0.0005, 0.0005, 0.0015, 0.0015, 0.0025, 0.0025, 0.0035, 0.0265, 0.0275$. Therefore $\operatorname{MAD}=(0.0015+0.0025)/2=0.0020$, and the scaled MAD is $1.4826(0.0020)=0.002965$.

Score the extremes. For the sample with $H_i=0.327$, $R_i=(0.327-0.3005)/0.002965=8.94$. For the sample with $H_i=0.273$, $R_i=(0.273-0.3005)/0.002965=-9.27$. Both exceed the heuristic cutoff in absolute value. The corresponding screening limits are $0.3005\pm3.5(0.002965)$, or approximately $0.2901$ and $0.3109$. These two samples are candidates for follow-up; the calculation alone does not establish why their values differ.

Plots and Follow-Up Checks

A robust score is most informative when paired with a plot. A histogram or density plot shows the overall shape and whether the sample values form separate clusters. A boxplot gives a compact view of the center and spread, while a dot plot can make individual candidates easy to identify. Where relevant, color or panel the display by batch or another known grouping, and compare heterozygosity with sample call rate. A candidate that is also unusually low in call rate deserves a check of its call completeness and the loci contributing to its summary.

A follow-up should verify that the result is not an artifact of calculation or data handling. Check that the marker panel, genotype coding, and missing-value encoding were consistent; confirm that the candidate has enough called loci; and inspect whether the pattern persists on a suitable shared set of loci. Also examine whether the cohort contains known subgroups or batches with different distributions. If the marker set or sample set changes after QC, recalculate the summaries and document the inputs used.

Fix the comparison set.
Specify eligible loci and consistent genotype and missing-call handling before comparing samples.

Summarize and flag.
Calculate $H_i$, the cohort median, MAD, and robust scores. Treat a threshold crossing as a candidate flag, not a final decision.

Inspect the pattern.
Use plots and check call rate, number of called loci, marker-panel consistency, and relevant cohort groupings.

Record the decision.
Document what was checked and whether the candidate is retained, investigated further, or handled under the project’s QC criteria.

Worked Example: Checking a Candidate With Unequal Call Rates

In a hypothetical panel of 10,000 eligible loci, sample A has 9,500 called genotypes, 2,850 of them heterozygous. Sample B has 8,000 called genotypes, 2,600 of them heterozygous. Their call rates and observed heterozygosities are $$ \text{Call rate}_A=\frac{9500}{10000}=0.950,\qquad H_A=\frac{2850}{9500}=0.3000, $$ $$ \text{Call rate}_B=\frac{8000}{10000}=0.800,\qquad H_B=\frac{2600}{8000}=0.3250. $$

Sample B has both a lower call rate and a higher observed heterozygosity. That combination makes it important to check the calls and the locus set rather than compare the two proportions without qualification. Suppose a check finds 7,800 loci called in both samples, with 2,340 heterozygous calls for A and 2,496 for B on those shared loci. Then $2340/7800=0.3000$ for A and $2496/7800=0.3200$ for B. The difference remains on this shared set, so it is not explained simply by the different denominators. It still warrants checking the underlying calls and cohort context; this calculation by itself does not resolve the cause.

The example also shows why call rate belongs in the diagnostic review. A shared-locus comparison can help assess comparability, but it is not a substitute for examining whether the sample’s calls and the chosen marker panel are suitable for the analysis.

Worked Example: A Pooled Flag That Needs a Group Check

Suppose eight hypothetical samples from batch A have heterozygosities $0.298, 0.299, 0.300, 0.300, 0.301, 0.301, 0.302, 0.303$. Two samples from batch B have values $0.340$ and $0.341$. Assume these values were calculated using the same marker criteria.

Across all ten samples, the median is $0.301$. The sorted absolute deviations are $0, 0, 0.001, 0.001, 0.001, 0.002, 0.002, 0.003, 0.039, 0.040$, so the MAD is $(0.001+0.002)/2=0.0015$. The scaled MAD is $1.4826(0.0015)=0.002224$. The scores for the batch B samples are $(0.340-0.301)/0.002224=17.54$ and $(0.341-0.301)/0.002224=17.99$. Both are far beyond the pooled screening cutoff.

The plot should make the grouping visible: the two batch B observations form a distinct cluster rather than an isolated point within batch A. The pooled calculation accurately signals unusual values relative to the combined set, but does not establish that the two samples are erroneous. With only two observations, there is also too little information to estimate a stable within-batch MAD for batch B. The appropriate diagnostic response is to investigate the grouping and data quality, not to treat the pooled flag as a removal rule.

Alternatives, Limitations, and Common Mistakes

Approach	What it does	Diagnostic caution
Mean and standard deviation	Flags values by distance from the mean in standard-deviation units.	Extreme observations can shift the mean and inflate the standard deviation, masking their own extremeness.
Median and scaled MAD	Measures distance from the median relative to a robust spread estimate.	May be uninformative when the MAD is zero, and a pooled summary can flag a whole subgroup.
IQR-based fences	Uses the first and third quartiles and the interquartile range to mark unusually distant values.	Fences are also screening rules, not formal tests; results can differ from MAD-based flags.

No threshold works equally well for every cohort. A small number of samples makes the center and spread estimates unstable. A highly skewed distribution, a mixture of groups, or many tied values can also make a single pooled robust score difficult to interpret. If the MAD is zero, $R_i$ is undefined; do not divide by zero or add an arbitrary small constant. Inspect the observed values and consider a clearly stated alternative, such as IQR-based summaries, alongside the plot and data checks. Any alternative still needs an appropriate comparison set.

Common mistake: “A sample has $|R_i|>3.5$, so it is a bad sample and should be removed.” The score only says that the observed value is distant from the cohort’s robust center under a chosen screening rule. It is not a probability, a test result, or proof of a technical error. Check call rate, locus comparability, plots, and relevant groupings before making a QC decision.

Key takeaway: Heterozygosity outlier diagnostics are a screening step. Define a comparable marker set, calculate robust summaries, inspect plots and related QC information, and treat each flag as a candidate for follow-up rather than an automatic exclusion.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions below.

For 1,200 called loci, a hypothetical sample has 372 heterozygous calls. What is its observed heterozygosity?
A sample’s robust score is $R_i=-4.1$ under the 3.5 screening rule. What does this indicate, and what does it not establish?
Why is a zero MAD a problem for the robust score, and what should be done instead of adding an arbitrary constant?
A pooled analysis flags several samples that all belong to one known batch. Name two checks that should precede any exclusion decision.

Solutions. (1) $H_i=372/1200=0.3100$, or 31.00% of the sample’s called eligible loci. (2) The sample is a candidate with a value below the cohort median, unusually distant under this screening rule. The score is not a p-value and does not show that the sample is erroneous or should be removed. (3) The score divides by $1.4826\,\operatorname{MAD}$, so a zero MAD makes the denominator zero and the score undefined. Inspect the tied values and distribution, then report an appropriate alternative such as IQR-based summaries; do not conceal the problem with an arbitrary constant. (4) Check sample call rates and the number and identity of loci contributing to the summaries; also inspect plots by batch and verify marker and genotype handling. These checks help determine whether the pooled flag reflects data quality or a group pattern that needs investigation.