Ascertainment-Corrected Likelihoods for Selected Families

Why Selection Changes the Likelihood

Family studies are often assembled because one or more relatives have an affected status or an unusually high or low quantitative phenotype. The observed families are then not a random sample from the population of families: they are families that satisfied a selection rule. A likelihood that treats them as an unselected sample can attribute the unusual phenotypes that prompted recruitment to genetic variance, residual variance, or other model parameters.

The correction is to model the distribution of family phenotypes conditional on the selection event. The denominator is the probability that a family would have been selected under the same parameter values used in the numerator. It depends on the recruitment rule, the phenotype model, and any covariates or genotypes on which the analysis conditions.

Definition: Let \(Y\) be the phenotype vector for a family, let \(\theta\) denote the parameters of its phenotype model, and let \(S\) denote the event that the family meets the study’s ascertainment rule. For an observed family with phenotype vector \(y\) satisfying \(S\), its ascertainment-conditioned likelihood contribution is \(f(y\mid S;\theta)=f(y;\theta)/P_\theta(S)\). The density \(f(y;\theta)\) is the ordinary family likelihood contribution, and \(P_\theta(S)\) is the probability of selection under the model.

For continuous phenotypes, \(f(y;\theta)\) is a density, not the probability of observing the exact vector \(y\). The event probability in the denominator is different: it is obtained by integrating the model over all phenotype vectors that would have led to selection. For a family selected when at least one of its members is affected, for example, the event includes every phenotype configuration satisfying that rule.

If \(n\) independently ascertained families have vectors \(y_1,\ldots,y_n\) and selection events \(S_1,\ldots,S_n\), the conditional likelihood is

L_c(\theta)=\prod_{i=1}^{n}\frac{f(y_i;\theta)}{P_\theta(S_i)}, \qquad \ell_c(\theta)=\sum_{i=1}^{n}\left\{\log f(y_i;\theta)-\log P_\theta(S_i)\right\}.

The product requires that the families be independent under the sampling model, or that any dependence between families has been modeled. Different families may have different selection events, such as different probands or thresholds, so their denominators need not be identical.

Deriving the Conditional Likelihood

The formula follows directly from the conditional density. For a phenotype vector \(y\) inside the selected region, the joint density of observing \(y\) and selection is \(f(y;\theta)\), while the probability of selection is \(P_\theta(S)\). Dividing gives the density among selected families. Equivalently, if the selection region is \(R\), then \(P_\theta(S)=\int_R f(u;\theta)\,du\), and the conditional density integrates to one over \(R\).

\int_R \frac{f(y;\theta)}{P_\theta(S)}\,dy = \frac{\int_R f(y;\theta)\,dy}{P_\theta(S)} =1.

This normalization is also why the correction affects estimation. Differentiating the conditional log likelihood for one family gives

\frac{\partial \ell_c}{\partial\theta} = \frac{\partial \log f(y;\theta)}{\partial\theta} - \frac{\partial \log P_\theta(S)}{\partial\theta}.

The second term is generally not zero. In fact, under regularity conditions, it subtracts the mean ordinary score among selected families. To see this, write \(U(y;\theta)=\partial\log f(y;\theta)/\partial\theta\). Differentiating the selection probability gives

\frac{\partial\log P_\theta(S)}{\partial\theta} = \frac{\int_R f(y;\theta)U(y;\theta)\,dy}{P_\theta(S)} = E_\theta\{U(Y;\theta)\mid S\}.

Thus, the conditional score is the ordinary score minus its expectation under the selected-family distribution. The expected conditional score is zero at the true parameter, as a likelihood score should be. Omitting the denominator generally destroys that property for data sampled through selection.

Worked Example: One Phenotype Selected Above a Threshold

Suppose a hypothetical quantitative phenotype satisfies \(Y\sim N(0,\sigma^2)\), and a family enters the study only if \(Y>1\). A selected individual has observed value \(y=1.5\). Calculate the conditional density at this value when \(\sigma=1.5\).

The selection probability is the upper-tail probability at the standardized threshold: \(P_\sigma(Y>1)=P(Z>1/1.5)\), where \(Z\sim N(0,1)\). The numerator is the normal density at \(1.5\), with standard deviation \(1.5\):

f(1.5\mid Y>1;\sigma=1.5) = \frac{\frac{1}{1.5}\phi(1.5/1.5)} {P(Z>1/1.5)} \approx \frac{0.1613}{0.2525} =0.6389.

This is a conditional density, in inverse phenotype units; it is not the probability that \(Y\) equals \(1.5\). The threshold in the denominator must be the threshold used for recruitment. Using the unselected density \(0.1613\) alone would omit the fact that only observations above 1 could appear in this sample.

The Selection Rule Determines the Denominator

“Selected because the family is affected” is not a complete mathematical description. The analysis must specify who was screened and what event triggered recruitment. If a designated proband had to exceed a threshold, the denominator is the probability that this proband exceeds it. If the rule was that at least one of several relatives exceeded the threshold, the denominator is the probability of that family-level event. These probabilities are usually different.

For a family with phenotype vector \(Y\), some common forms are:

Proband rule: \(S=\{Y_1>c\}\), so the denominator is \(P_\theta(Y_1>c)\), even when the likelihood numerator includes phenotypes of other relatives.
At-least-one rule: \(S=\{\max_j Y_j>c\}\), so the denominator is a joint family probability. Dependence among relatives matters.
All-members rule: \(S=\{Y_j>c\text{ for every }j\}\), with a denominator that integrates over the region where every member meets the threshold.
Binary affected-member rule: if selection requires at least one affected relative, the denominator is the probability of that event under the binary phenotype model, or under a liability model if affected status is defined by crossing a liability threshold.

For a multivariate normal family model, an extreme-phenotype event is a region in the family phenotype space. The denominator is a univariate tail probability when selection uses only one observed proband phenotype; for a rule involving several relatives, it is generally a multivariate probability. Correlation induced by shared family effects must be retained when calculating that probability.

Worked Example: Selection Requires Both Relatives to Be Extreme

Consider two hypothetical independent relatives, each with \(Y_j\sim N(0,\sigma^2)\). A family is included only if both phenotypes exceed 1. The observed values are \(y_1=1.5\) and \(y_2=2\). Compare the conditional likelihood contributions under \(\sigma=1\) and \(\sigma=1.5\).

Independence makes the ordinary joint density the product of the two univariate densities, and it makes the selection probability the square of the single-relative upper-tail probability. For \(\sigma=1\), the uncorrected density is

\phi(1.5;0,1)\phi(2;0,1) \approx (0.1295)(0.05399)=0.006993.

The probability that both exceed 1 is \(P(Z>1)^2\approx(0.1587)^2=0.02517\). Therefore,

L_c(\sigma=1) = \frac{0.006993}{0.02517} \approx 0.2778.

For \(\sigma=1.5\), the uncorrected density is

\phi(1.5;0,1.5^2)\phi(2;0,1.5^2) \approx (0.1613)(0.1093)=0.01764.

Now the single-relative selection probability is \(P(Z>1/1.5)=P(Z>0.6667)\approx0.2525\), so the probability that both relatives exceed 1 is approximately \((0.2525)^2=0.06375\). Hence

L_c(\sigma=1.5) = \frac{0.01764}{0.06375} \approx 0.2767.

The uncorrected densities favor \(\sigma=1.5\), since \(0.01764>0.006993\). After conditioning on the rule that both relatives had to be extreme, the contributions are close: \(0.2778\) under \(\sigma=1\) and \(0.2767\) under \(\sigma=1.5\). These are likelihood contributions for two specified parameter values, not estimates by themselves; the example shows that selection can substantially change their comparison.

Worked Example: A Proband Rule in a Correlated Pair

Suppose two hypothetical standardized family phenotypes have a bivariate normal distribution with means zero, variances one, and correlation \(0.5\). Recruitment requires only that the designated proband satisfy \(X_1>1\). The observed pair is \((x_1,x_2)=(1.5,0.8)\). Find its conditional joint density given selection.

For a bivariate normal pair with correlation \(0.5\), \(X_2\mid X_1=1.5\) is normal with mean \(0.5(1.5)=0.75\) and variance \(1-0.5^2=0.75\). The joint density is the marginal density of the proband times the conditional density of the relative:

f(1.5,0.8) = \phi(1.5;0,1)\, \phi(0.8;0.75,\sqrt{0.75}) \approx (0.1295)(0.4599) \approx 0.05956.

The event is specifically \(X_1>1\), so its probability is \(P(Z>1)\approx0.1587\). Dividing the joint density by this probability gives

f(1.5,0.8\mid X_1>1) = \frac{0.05956}{0.1587} \approx 0.3754.

The denominator is not the probability that both relatives exceed 1: the second relative did not determine inclusion. If recruitment instead required at least one member to exceed 1, the denominator would be \(P(X_1>1\text{ or }X_2>1)\), which depends on their correlation and is a different event.

Using the Correction in Variance-Component and QTL Models

In a Gaussian family model, write the phenotype vector as \(Y\sim N(X\beta,V)\), where \(X\beta\) is the mean and \(V\) is the family covariance matrix. For example, a variance-component specification may be written

V=\sigma_A^2 K+\sigma_Q^2\Pi+\sigma_E^2 I.

Here \(K\) is the pedigree covariance multiplier for the polygenic component, \(\Pi\) is the IBD-sharing matrix for the QTL component, \(I\) is the identity matrix, and the \(\sigma^2\) terms are the corresponding variance components. The ordinary family density is the multivariate normal density with mean \(X\beta\) and covariance \(V\). Under phenotype-based selection, its contribution must be divided by \(P_{\beta,V}(S)\), calculated using that same mean and covariance model.

This matters for both estimation and testing. The selection probability can change with the mean parameters and with each variance component. For instance, increasing a variance may make extreme phenotypes more likely, increasing the denominator as well as changing the density at the observed vector. Maximizing only the numerator can therefore yield variance-component estimates that describe the selected families as though they had been sampled at random. The direction and size of the resulting distortion depend on the selection rule, model, and data; there is no universal direction of bias.

For a QTL test, the null and alternative must be compared using likelihoods that account for the same ascertainment design. Under a null model without a QTL component, the covariance excludes the QTL variance; under the alternative it includes that component. Each model uses its own parameter-dependent selection probability. The likelihood-ratio statistic is formed from the maximized conditional log likelihoods, not by correcting only the alternative or only the null. If the QTL variance is constrained to be nonnegative, it can lie on a parameter boundary under the null, so ordinary reference distributions may require care; ascertainment correction does not remove that issue.

When marker data leave IBD sharing uncertain, the treatment of that uncertainty must also match the conditioning scheme. If the analysis conditions on observed marker information, it should use the phenotype likelihood and selection probability conditional on that information, integrating over uncertain IBD consistently in both. In general, a ratio of mixtures is not the same as a mixture of ratios. The selected-data likelihood should therefore be derived for the actual observed and conditioned variables rather than assembled by applying a correction after an unrelated marginalization.

Assumptions and Common Pitfalls

Ascertainment correction is not a substitute for describing recruitment. It relies on a specified selection event and a phenotype model that can assign a probability to that event. It also assumes the modeled selection mechanism is the one that generated the sample; unmodeled dependence on family size, covariates, referral patterns, or other variables can leave residual selection bias. If the probability of selection is unknown or differs across families, that design information may need to be modeled explicitly.

Using the wrong event: dividing by the probability that any relative is extreme when only a designated proband triggered recruitment changes the sampling model. State the actual rule before computing the denominator.
Using an independent-relative denominator for correlated relatives: a product of marginal probabilities is valid only under the required independence assumptions. Family covariance generally makes a joint selection probability necessary.
Correcting only the observed proband: if the likelihood includes the full family phenotype vector, the numerator is the joint family density. The denominator still corresponds to the exact recruitment event, not automatically to a univariate truncation.
Leaving the denominator fixed while estimating parameters: \(P_\theta(S)\) generally changes with \(\theta\). Recalculate it as the model parameters change during likelihood optimization.
Treating selection as causal evidence: conditioning on a recruitment event corrects the sampling likelihood under its assumptions; it does not make an observed family association causal.

Key takeaway: An ascertainment-conditioned likelihood is the model density for the observed family divided by the model probability of the precise event that brought that family into the study. In variance-component and QTL analyses, both pieces must use the same phenotype model and conditioning choices.

Check Your Understanding

For each question, identify the selection event before deciding what belongs in the likelihood denominator.

A study recruits a family if its proband’s phenotype exceeds \(c\), then measures all relatives. Which phenotypes enter the likelihood numerator, and what event determines the denominator?
Why is \(f(y;\theta)\) not, by itself, the likelihood contribution for a family sampled only when \(Y>c\)?
Two correlated relatives are selected when at least one is affected. Why might multiplying their marginal non-affected probabilities give an incorrect selection probability?
In a QTL likelihood-ratio test, why must the ascertainment probability be evaluated under both the null and alternative models?
When may a selection probability remain unchanged during parameter estimation, and why is that not generally true for a phenotype-based threshold rule?