Past exam of the mathematics course of the University of Cambridge / 2006 / iii / Paper 45 (Statistical and Population Genetics)

Treat $G_i$ as the entire phase-known genotype of individual $i$ at the genetic loci being studied, and let $\overline{F}$ denote the nonfounders. Assume unrelated pedigree founders have independent genotype priors, each nonfounder's genotype depends only on its parents' genotypes, and the observed phenotypes are conditionally independent given the genotypes. These are the standard assumptions behind the pedigree likelihood; shared environmental effects or related pedigree founders would require additional factors.

The directed ancestry graph is acyclic, even when there is inbreeding. Order individuals with parents preceding offspring and apply the chain rule for probabilities. The parental conditional independence assumption givesConditional phenotype independence gives $P(X_1,\dots ,X_n\mid G_1,\dots ,G_n)={\prod }_{i}P(X_i\mid G_i)$. Multiplying these two expressions and summing over every possible latent genotype configuration provesThis is the law of total probability, not a claim that relatives have independent unconditional genotypes.

For a binary disease phenotype, put $p_i(g)=P(X_i=\text{affected}\mid G_i=g)$. The penetrance function depends on the inheritance model, the disease alleles and genotype, and potentially age, sex, environment and other genetic modifiers. Phenocopies allow disease without the putative disease genotype, while incomplete penetrance allows a genetic carrier of that genotype to remain unaffected. Thus the phenotype factor is $p_i(g)$ for an affected individual and $1-p_i(g)$ for an unaffected one.

Under a fully penetrant model with no phenocopies the genotype determines disease status, so the factor is 1 for a compatible genotype and 0 for an incompatible genotype. For example, in fully penetrant recessive inheritance, $aa$ is affected and $AA,Aa$ are unaffected. Missing phenotype observations contribute the constant factor 1 rather than excluding genotypes.

Without genotyping error, the observed genetic marker genotype is a deterministic function $u(g)$ of the underlying phase-known genotype. In particular $u$ can forget chromosome phase and parental origin. ThereforeSeveral different phased haplotype pairs can receive factor 1 when the observed multilocus genotype is unphased. With laboratory error, missing alleles or uncertain calls, replace the indicator by an error or observation model. A wholly unobserved genetic marker contributes factor 1. Thus a 0-or-1 factor reflects perfect observation of a deterministic genetic marker state, not a general property of every genetic measurement.

The pedigree founder factor $P(G_i)$ depends on population allele and haplotype frequencies, ancestry, mating assumptions and any population linkage disequilibrium. At one autosomal genetic locus, Hardy-Weinberg equilibrium with allele frequencies $p,q$ gives genotype probabilities $p^2,2pq,q^2$. For a specified ordered pair of independent pedigree founder haplotypes $h_1,h_2$, the probability is $p_{h_1}{p}_{h_2}$; an unordered distinct pair has twice that probability. Pedigree founders assumed to be related require a joint prior rather than a product of unrelated-founder factors.

The nonfounder factor depends on the parental genotypes and phases, Mendelian segregation, and the recombination fractions between genetic loci. In a two-locus heterozygote of known phase, the two parental haplotypes have transmission probabilities $(1-\theta )/2$ each, and the two recombinant haplotypes $\theta /2$ each. Multiply the maternal and paternal gamete probabilities and sum any gamete combinations producing the child's genotype. Sex-specific genetic recombination, mutation or segregation distortion can be incorporated by changing these gamete probabilities. At one ordinary autosomal genetic locus, each parental copy is transmitted with probability $1/2$.

Call the parent genotype $g$ and child genotype $h$, each in $\{0,1\}$. DefineThere is one pedigree founder and one nonfounder, soExplicitly this is $a_0{p}_0{b}_0{t}_{0\mid 0}+a_0{p}_0{b}_1{t}_{1\mid 0}+a_1{p}_1{b}_0{t}_{0\mid 1}+a_1{p}_1{b}_1{t}_{1\mid 1}$. The genotype-transmission model can be arbitrary for the single-parent species; no two-parent Mendelian assumption is being introduced.

In the unrearranged likelihood there are four summands. Each is a product of four supplied probabilities and takes three multiplications. Summing the four results takes three additions. ThusThis count concerns direct evaluation without reusing repeated factors; zero-or-one phenotype values could of course allow extra special-case savings.

The Elston-Stewart algorithm first sums over the child:For each of the two parent states, the inner message requires two multiplications and one addition, totaling four multiplications and two additions. Multiplying each message by its two parent factors takes two further multiplications per parent state, totaling four. The final sum uses one addition. ThereforeMoving the sum inward saves four multiplications because the factors independent of the child genotype are evaluated only once per parent state.

Genetic linkage analysis follows co-segregation of genetic marker and disease alleles within pedigrees. It estimates or tests a recombination fraction, with independent assortment at $\theta =1/2$ as the ordinary null. It can locate a disease genetic locus even when different families have different genetic marker alleles on their disease chromosomes.

Genetic association analysis tests whether allele or genotype frequencies differ with disease status. It can compare unrelated cases and controls, or compare transmitted and untransmitted parental alleles. A noncausal genetic marker can be associated through population linkage disequilibrium with a causal genetic locus. Association therefore concerns population-frequency relationships, whereas linkage concerns inheritance; population stratification can confound the former without producing genuine within-family linkage.

Only heterozygous parents contribute an alternative-allele transmission contrast to the transmission disequilibrium test. Both parents are homozygous in rows 1, 13, 14 and 15, so rows 1, 13, 14 and 15 are uninformative. Their counts total 15 families. Every other row has at least one heterozygous parent. When both parents are heterozygous and the child is heterozygous, there is one transmission of each allele, regardless of the unresolved parent-specific phase; these transmissions contribute to both discordant cells.

Write $b$ for transmission of allele 1 with allele 2 untransmitted, and $c$ for the reverse. The informative transmissions give $b=20$, $c=30$. Under the fair-transmission null, conditional on $b+c=50$, $b$ is binomial with success probability $1/2$. Its standardized squared deviation isThe reference distribution is asymptotically chi-squared with one degree of freedom, giving a two-sided $p$-value about $0.157$. Allele 2 is overtransmitted and hence is the allele more associated with disease in this sample, although this result is not strong evidence against the null. The diagonal transmission cells contribute no discordant contrast.

Summing the transmission table by its rows and columns givesThus the marginal case-versus-pseudo-control odds ratio isThis is the requested approximate allelic effect. Keeping the parent-specific matching instead gives the discordant-pair odds ratio $c/b=30/20=1.5$, numerically close here. Neither calculation makes the within-family allele observations independent population samples.

For each family combine the two untransmitted parental alleles into one family pseudo-control genotype. A heterozygous child of two heterozygous parents has a heterozygous pseudo-control, whatever its parent-specific transmission phase. Combining the rows givesFor example, pseudo-control genotype $1/1$ receives contributions $1,2,5,5$ from the relevant rows; genotype $2/2$ receives $3,1,2,6$, and the remaining 25 are heterozygous. These counts also give the untransmitted allele totals 51 and 49.

The estimated allele frequencies are ${\hat{p}}_1=0.51$, ${\hat{p}}_2=0.49$. The Hardy-Weinberg equilibrium expected genotype counts areThey are almost identical to the observations, so these pseudo-control counts conform extremely closely to HWE proportions. The usual Pearson discrepancy is only about $8.0\times 10^{-6}$. This numerical agreement is not a general guarantee of HWE for family pseudo-controls, whose construction involves mating structure and case ascertainment in a genetic study.

Family comparisons condition on the parental genotypes, so cases are compared with the alleles they could have inherited. This protects a transmission test from ordinary population stratification, avoids having to recruit unrelated controls matched for ancestry, and can resolve haplotypes and identify Mendelian inconsistencies. The same pedigree may also provide linkage information.

The costs are recruiting and genotyping relatives, missing or unavailable parents, and fewer informative comparisons because homozygous parents contribute no transmission contrast. These problems are especially important for late-onset disease. Related cases or multiple siblings also require an analysis accounting for dependence, and genotyping error or transmission distortion can mimic overtransmission. Population-based case/control studies are usually easier to assemble and can offer greater power per genotyped person, but need careful control of ancestry and other sources of confounding. Family sampling does not remove every form of ascertainment in a genetic study or observation bias.

Fix one particular ancestral allele copy $q$ among the four copies in the founding couple. The pedigree founder carrying it transmits it independently to each of the two generation-2 siblings with probability $1/2$. ThereforeHere “one copy” means a specified copy, as required to follow “this ancestral copy” in the subsequent parts. If instead one asks whether the siblings share at least one of the four copies without specifying which, the probability is $1-(1/2{)}^2=3/4$: they fail to share either the paternal or maternal copy with probability $1/4$. These are different events.

Conditional on both generation-2 siblings carrying the specified copy $q$, each passes it to each of their three children with probability $1/2$. Their outside partners do not carry that ancestral copy. There are six independent meioses, one for each generation-3 individual forming the three cousin marriages. ConsequentlyThis calculation ignores phenotype ascertainment in a genetic study, as requested; conditioning on the later disease observations would force these transmissions rather than assign them their prior probabilities.

Each relevant cousin marriage is a cross between two genetic carriers of $q$. A child receives two copies of $q$ with probability $(1/2)(1/2)=1/4$, and fails to receive two with probability $3/4$. Conditional on the parental genotypes, distinct offspring transmissions are independent. There are eight specified double-copy children and eight specified remaining children, givingThere are no binomial coefficients because the subjects have been individually selected. If only the counts in the three sibships were specified, rather than their identities, the probability would instead be multiplied by $\left(\genfrac{}{}{0ex}{}{1}{1}\right)\left(\genfrac{}{}{0ex}{}{7}{4}\right)\left(\genfrac{}{}{0ex}{}{8}{3}\right)=1960$.

Let $E_q$ be the event that all eight selected subjects carry two copies of the same specified ancestral copy $q$, and that no other generation-4 subject is autozygous at this genetic locus. The selected subjects occur in every sibship, so this event forces both generation-2 transmissions and all six relevant generation-3 transmissions. Applying the preceding conditional probabilities,Once all six parents carry $q$, their other copies come from the two unrelated outside pedigree founders, one on each side of the pedigree. A child not receiving two $q$ copies therefore cannot become autozygous for a different founding copy. This verifies the “only” condition, rather than merely excluding homozygosity for $q$.

For the event $E$ without specifying which of the original four copies is shared, the four alternatives $E_q$ are disjoint: a selected child cannot simultaneously have both its copies descended from two different ancestral copies. HenceAll selected individuals are then autozygous for the same copy and mutually 2-IBD. The distinction from the $3/4$ probability in the alternative interpretation of part (a) is precisely this disjointness.

Let $D$ denote the specified disease pattern. With one introduced disease-allele copy, no new mutations or phenocopies, and fully penetrant recessive inheritance, every affected subject must carry two descendants of that one copy. Since every sibship has an affected child, all the upstream genetic carrier transmissions described above are forced. Conversely, unaffected siblings cannot carry two disease copies. Thus $D$ implies the inheritance event $E$ at the causal genetic locus, andThe rarity of the unconditional pattern in part (d) does not make its posterior probability small after observing the pattern that forces it. If the location of the one disease copy is initially uniform among the four founding copies, the likelihood of $D$ is the same for each, so its posterior location remains uniform; $P(E_q\mid D)=1/4$ for any named copy. If $q$ is already designated as the disease copy, its corresponding event has posterior probability 1.

Write $D$ for the disease observations and $M$ for the genetic marker observations. The LOD score for complete linkage isThe disease marginal cancels because the single-locus disease segregation law does not depend on the genetic marker–disease recombination fraction. The ten genetic markers form one nonrecombining haplotype; they do not supply ten independent copies of the pedigree's segregation evidence.

There is a distinction between identity by state and identity by descent. The reported observation that affected subjects all have $h/h$ establishes the former. An exact genetic marker likelihood requires pedigree founder haplotype frequencies, or the actual pedigree founder genotypes if conditioning on them. Here is an explicit expression using independent pedigree founder haplotypes, with frequency $p$ for the observed haplotype $h$, and linkage equilibrium between pedigree founder genetic marker and disease states. There are eight genetic marker copies in the four pedigree founders: the original couple and the two generation-2 outside spouses. For genetic marker inheritance configuration $v$, let $K(v)$ be the number of distinct pedigree founder copies ancestral to the sixteen chromosome copies in the eight selected subjects. Thenis the unlinked probability that all those copies have state $h$. At complete linkage the affected subjects all inherit the disease pedigree founder copy twice, and the probability that this copy has genetic marker state $h$ is $p$. Therefore

The denominator can be evaluated without an unspecified inheritance sum. Let $B_a=\left(\genfrac{}{}{0ex}{}{2}{a}\right)p^a(1-p{)}^{2-a}$ for $a=0,1,2$. If parental genetic marker genotypes contain $a,b$ copies of $h$, the child count distribution isFor a generation-2 couple defineFor a particular generation-3 parent with $j$ genetic marker copies, the chance it transmits $h$ to all its $r$ selected offspring is $(j/2{)}^r$, which explains each factor. Conditional on the original founding couple, its two generation-2 siblings and their independent outside spouses give independent left and right contributions. ThusExpanding givesIt satisfies $R(1)=1$, as it must for a monomorphic haplotype. This expression supplies a numerical LOD once $p$ is specified; the problem does not specify that frequency.

In the idealized limit where a matching rare haplotype identifies one ancestral copy, the leading coefficient is particularly simple. Dropping restrictions on the eight unaffected subjects, the probability that all selected subjects are autozygous for one common copy isAccordingly $R(p)/p\to 2^{-22}$ as $p\to 0$, and the affected-only IBD calculation gives $Z(0)=22{\log }_{10}2\simeq 6.623$.

If the genetic marker observations additionally establish the exact pattern $E$ of part (d), including that the eight unaffected subjects are not autozygous, its idealized IBD likelihood ratio is insteadThis latter value requires the additional negative genetic marker observations. It cannot be inferred merely from homozygosity in the affected subjects, and neither simplified IBD score is an exact identity-by-state likelihood for an unspecified common haplotype.

Under the original assumptions the extra affected child is impossible. A generation-3 woman can have at most one copy of the unique introduced disease allele, while an unrelated outside spouse has none. Their children cannot inherit two defective copies. Thus both linkage and no-linkage likelihoods of the enlarged data are zero, and the original single-entry model has no defined updated LOD score. At least one assumption must be relaxed: another introduced disease allele, a new mutation, a phenocopy, or the relationship of the outside spouse.

If the model is enlarged to allow an independent outside genetic carrier, the reported extra information contains disease status but no new genetic marker genotypes. In the simple extension with fixed outside heterozygous-carrier probability $c$, the old disease observations force both generation-2 family ancestors to be genetic carriers. A further generation-3 sister inherits the disease copy with probability $1/2$, independently of the existing offspring transmissions. Given that she and her spouse are genetic carriers, exactly one of their three children is affected with probabilityThe additional phenotype likelihood is therefore $27c/128$. It is the same at $\theta =0$ and $\theta =1/2$, so it cancels and the LOD remains unchanged in this phenotype-only, independent-carrier extension. If the identity of the affected child is specified, the factor is $9c/128$, which also cancels. No such conclusion applies automatically if new genetic marker data are supplied or the revised pedigree founder model couples this branch to the old genetic marker evidence.

In general an enlarged data set changes the score byafter choosing a coherent model with positive likelihood. The additional affected child alone does not justify assigning it the autozygosity evidence of an affected child of a cousin marriage.

Use standard neutral coalescent time units, so each pair of ancestral lineages merges at rate 1. With $j$ lineages, there are $\left(\genfrac{}{}{0ex}{}{j}{2}\right)$ possible pairs; consequently the successive waiting times are independent withDuring this epoch there are $j$ branches, each of length $T_j$. Hence the total branch length of a neutral coalescent is $L={\sum }_{j=2}^{n}j{T}_j$, and linearity of expectation givesThe root branch above the sample's common ancestor is excluded: its mutations would be shared by every chromosome and would not be segregating sites in the sample.

In the infinite sites mutation model, mutations occur as independent Poisson processes at rate $\theta /2$ per lineage per unit of the coalescent time used in part (a). Conditional on the genealogy, superposition over branches gives a Poisson distribution with mean $\theta L/2$. Since this depends only on $L$, conditioning on total length alone givesEach mutation occurs at a new site and below the common ancestor, so it contributes exactly one segregating site. The law of iterated expectation therefore gives, with $a_{n-1}={\sum }_{i=1}^{n-1}1/i$,Here $\theta$ is fixed when taking the expectation. If it is also random under a proper prior with finite mean, the unconditional expectation is $a_{n-1}E\theta$.

Assume $n\ge 2$, $k\ge 0$, a proper prior density $\pi$ supported on nonnegative mutation rates, and independence of the mutation-rate prior from the neutral genealogy. PutBy Bayes' theorem, the joint posterior isThe normalizing constant is the prior predictive probability $P(S=k)$ and must be positive. For $k=0$, interpret the factor $(\theta L/2{)}^0$ as 1 even at zero. If $n=1$, there are no segregating sites: $S=0$ surely, so that case supplies no information about $\theta$.

Use Bayesian rejection sampling for a segregating-site count. Independently propose $\theta \sim \pi$ and $T_j\sim \mathrm{Exp}({\lambda }_j)$ for $j=2,\dots ,n$, and calculate $L={\sum }_{j}j{T}_j$. Draw an independent uniform $U$ on $(0,1)$ and accept the whole proposed vector exactly whenRepeat after a rejection. The right side is a probability mass function value of a Poisson distribution and is therefore between zero and one. On acceptance, the probability density function of a proposed vector is its prior density multiplied by this acceptance probability, divided by the probability of acceptance. Part (c) shows that this is precisely $f(\theta ,T\mid S=k)$. Independent proposals and uniforms produce independent posterior draws. Equivalently, one could simulate a count from a Poisson distribution with mean $\theta L/2$ and retain the proposal when that count equals $k$; the uniform version avoids simulating that extra count.

The acceptance probability is the prior average of the Poisson likelihood:Thus the expected number of proposals per accepted draw is $1/Z_k$. A numerical rate depends on the specified prior, sample size and observed count; it is not determined by $k$ alone.

For a more explicit one-dimensional expression, condition on $\theta$. The probability generating function, integrating each independent exponential epoch, isLet $p_k(\theta )$ be its coefficient of $z^k$. Then $Z_k=\int \pi (\theta )p_k(\theta )d\theta$. This identifies the acceptance rate with a fully specified marginal likelihood average.

The rejection rule can use a tighter envelope than 1. For $k>0$, differentiate $\log (e^{-x}{x}^k/k!)=-x+k\log x-\log (k!)$. Its derivative $-1+k/x$ changes from positive to negative at $x=k$, soKeep the same joint prior proposal but accept with probability $w_k(\theta ,T)/M_k$. Multiplication of the acceptance function by this constant leaves the accepted probability density function unchanged. The improved rate isFor $k>0$, $M_k<1$, so this is a genuine improvement; for large $k$, Stirling's formula gives $M_k\sim (2\pi k{)}^{-1/2}$. This is the best global constant envelope for an unrestricted positive prior proposal, because $L$ has support $(0,\infty )$ and $\theta L/2$ can approach $k$. For $k=0$ there is no improvement from this constant scaling. Further gains can come from proposals concentrating near $\theta L/2=k$, or integrating out the times first, but such changes require the appropriate proposal-density ratio and envelope rather than an uncorrected alteration of the sampling law.

The accepted mutation parameters have marginal probability density functionEstimate this probability density function with a suitably smoothed histogram or a boundary-aware kernel density estimator from the accepted draws. On a parameter region where $\pi (\theta )>0$, divide the estimate by the known prior density and maximize:The unknown $Z_k$ is constant in $\theta$ and cancels. Maximizing the posterior density itself gives a posterior mode, not generally a maximum-likelihood estimate; the two coincide for a flat prior on a region containing the maximum. A prior excluding candidate maximizers cannot recover them through its posterior output. Boundary behavior also matters: for $k=0$, the likelihood is decreasing in $\theta$, so the unrestricted maximum is at $\theta =0$.

The exact likelihood provides a useful independent check. Put $W=L/2$. Its Laplace transform is ${\prod }_{i=1}^{n-1}i/(i+s)$, since $j{T}_j/2$ is exponential of rate $j-1$. The probability density functionhas that same transform: substitute $u=e^{-w}$ to obtain $(n-1){\int }_0^1{u}^s(1-u{)}^{n-2}du={\prod }_{i=1}^{n-1}i/(i+s)$. Thus, expanding its binomial factor and integrating each term,This is the segregating-site likelihood under the neutral coalescent. It can be maximized numerically to check the posterior-density estimate; the generating-function representation avoids cancellation in large alternating sums. For $n=2$, it reduces to ${\theta }^k/(1+\theta {)}^{k+1}$, whose maximum is at $\theta =k$ for $k>0$.

Alternatively, differentiating the integrated likelihood givesAn interior maximum therefore satisfies $\theta =k/E[W\mid S=k,\theta ]$. Conditional expectations of $W$ can be estimated by smoothing the accepted pairs $(\theta ,W)$ near each candidate parameter value, yielding another simulation-based score equation. The simple estimator $k/a_{n-1}$ is the unbiased Watterson estimator from part (b), and is not generally the likelihood maximizer.