Score Tests for Genetic Models

Why Test a Genetic Effect Under the Null?

A genetic association analysis often asks whether a specified genetic parameter is zero. A Wald test answers this using an estimate and standard error from the fitted alternative model. A score test instead asks whether the likelihood is changing in the direction of the alternative when evaluated at the null model. It can therefore test a genetic effect without fitting the full model that includes that effect.

This can be useful when fitting a large alternative model is computationally demanding. The score test still requires fitting any nuisance parameters under the null, such as a baseline disease probability. It also depends on specifying a genetic model: an additive dosage effect and a carrier-versus-noncarrier effect are different hypotheses.

Definition: Let $\beta$ be the genetic parameter of interest and let $\eta$ denote nuisance parameters. A score test evaluates the score for $\beta$, the derivative of the log likelihood with respect to $\beta$, at $\beta=0$ and at the nuisance-parameter estimate under $H_0$. It scales that score using the information for $\beta$ after accounting for nuisance parameters. For a single parameter, the resulting statistic is approximately distributed as a chi-square random variable with one degree of freedom under $H_0$.

A Logistic Model for a Genetic Effect

Consider $n$ unrelated individuals with binary outcomes $Y_i$, where $Y_i=1$ indicates the outcome of interest and $Y_i=0$ indicates its absence. Let $G_i$ be a prespecified genotype coding. For additive dosage coding, $G_i$ is the number of copies of a specified allele, so $G_i\in\{0,1,2\}$. Suppose outcomes are conditionally independent given the genotypes and follow the logistic model

Y_i\mid G_i \sim \operatorname{Bernoulli}(\pi_i), \qquad \operatorname{logit}(\pi_i)=\alpha+\beta G_i.

Here $\alpha$ is the intercept and $\beta$ is the change in log odds of the outcome for a one-unit increase in the coded genotype. The null hypothesis of no genetic effect under this model is $H_0:\beta=0$; the usual two-sided alternative is $H_a:\beta\ne 0$. The odds ratio for a one-unit increase in $G_i$ is $\exp(\beta)$. The model describes association conditional on its assumptions; it does not by itself establish causation.

The assumptions for the calculations below are that the outcomes are independent conditional on the observed genotypes, each outcome follows the stated Bernoulli logistic model, the genotype coding is fixed before testing, and the null model has a finite fitted intercept. For the usual chi-square approximation, the sample must also be sufficiently informative for the likelihood to be locally well approximated by a quadratic. In particular, there must be genotype variation and a positive estimated null information.

Score, Information, and the Null Fit

For observed outcomes and genotypes, the log likelihood, apart from constants that do not involve the parameters, is

\ell(\alpha,\beta) = \sum_{i=1}^{n} \left[ Y_i(\alpha+\beta G_i) - \log\{1+\exp(\alpha+\beta G_i)\} \right].

The score for $\beta$ is the derivative of this log likelihood:

U_\beta(\alpha,\beta) = \frac{\partial\ell}{\partial\beta} = \sum_{i=1}^{n}G_i(Y_i-\pi_i).

Under $H_0$, all individuals have the same outcome probability, denoted by $p_0$, because the genotype term is zero. The null model is an intercept-only Bernoulli model. If there are both cases and non-cases, its maximum-likelihood estimate is $\hat p_0=\bar Y$, the observed proportion of outcomes equal to one. Thus the score for the genetic effect, evaluated at the null fit, is

U_\beta = \sum_{i=1}^{n}G_i(Y_i-\hat p_0).

The intercept is a nuisance parameter: it is not the target of the test, but it is estimated under the null. The information for $\beta$ must therefore account for the intercept. Write $I$ for the Fisher information matrix, the negative expected second-derivative matrix of the log likelihood. At the null fit, its entries are

I_{\alpha\alpha}=\hat p_0(1-\hat p_0)n,\qquad I_{\alpha\beta}=\hat p_0(1-\hat p_0)\sum_{i=1}^{n}G_i,\qquad I_{\beta\beta}=\hat p_0(1-\hat p_0)\sum_{i=1}^{n}G_i^2.

The information for the genetic effect after adjusting for the intercept is the Schur complement of $I_{\alpha\alpha}$:

I_{\beta\beta\cdot\alpha} = I_{\beta\beta} - I_{\beta\alpha}I_{\alpha\alpha}^{-1}I_{\alpha\beta} = \hat p_0(1-\hat p_0) \sum_{i=1}^{n}(G_i-\bar G)^2, \qquad \bar G=\frac{1}{n}\sum_{i=1}^{n}G_i.

This adjustment matters: a common shift in genotype coding is indistinguishable from a shift in the intercept, so the relevant information comes from variation in genotype around its mean. At the null maximum-likelihood estimate, the intercept score is zero, since $\sum_i(Y_i-\hat p_0)=0$. Consequently the adjusted genetic score is the score $U_\beta$ above.

Deriving and Applying the Score Statistic

Under regularity conditions, the genetic score at the null fit is approximately centered at zero under $H_0$, with variance given by the nuisance-adjusted information. Dividing by its standard deviation gives an approximately standard normal statistic. Squaring it gives the equivalent one-degree-of-freedom chi-square statistic.

Fit the null model.
Estimate the intercept under $\beta=0$; for the model here, $\hat p_0$ is the observed outcome proportion.

Evaluate the genetic score.
Calculate $U_\beta=\sum_iG_i(Y_i-\hat p_0)$ using the null fit, not an alternative-model estimate.

Calculate adjusted information.
Use $I_{\beta\beta\cdot\alpha}=\hat p_0(1-\hat p_0)\sum_i(G_i-\bar G)^2$.

Form the test statistic.
Compare the squared score divided by adjusted information with a chi-square distribution with one degree of freedom.

Theorem: For testing one interior parameter $\beta=0$, under the null model and standard large-sample regularity conditions, the score statistic $$ S=\frac{U_\beta^2}{I_{\beta\beta\cdot\alpha}} $$ has an approximate $\chi^2_1$ distribution under $H_0$. Its p-value is the probability, assuming $H_0$, of obtaining a $\chi^2_1$ statistic at least as large as the observed $S$. A small p-value is evidence against $H_0$, not the probability that $H_0$ is true.

Worked Example: Additive Dosage Score With Modest Evidence

Suppose hypothetical data include 60 unrelated individuals. Their genotype counts are 20 with $G=0$, 30 with $G=1$, and 10 with $G=2$. The corresponding numbers with the outcome are 4, 9, and 5. Assume the binary outcome model above and test $H_0:\beta=0$ against $H_a:\beta\ne0$ at $\alpha=0.05$.

There are $4+9+5=18$ outcomes, so the null estimate is $\hat p_0=18/60=0.3000$. The sum of genotype dosages is $30+2(10)=50$, and the dosage-weighted outcome total is $1(9)+2(5)=19$. Therefore,

U_\beta=19-(0.3000)(50)=4.000.

The mean dosage is $\bar G=50/60=0.8333$, while the sum of squared dosages is $30+4(10)=70$. Hence

\sum_i(G_i-\bar G)^2 = 70-\frac{50^2}{60} = 28.3333, \qquad I_{\beta\beta\cdot\alpha} = (0.3000)(0.7000)(28.3333) = 5.950.

The score statistic is $S=4.000^2/5.950=2.689$. Comparing with $\chi^2_1$ gives a two-sided p-value of approximately $0.101$. At the 0.05 level, these data do not provide sufficient evidence to reject the null of no additive genetic effect. This is not proof that the effect is zero.

Worked Example: A Stronger Additive-Effect Signal

Consider a second hypothetical sample of 100 unrelated individuals. The genotype counts are 40 with $G=0$, 40 with $G=1$, and 20 with $G=2$. The numbers with the outcome are 4, 12, and 12, respectively. Use the same logistic model and two-sided test.

The total number of outcomes is 28, so $\hat p_0=28/100=0.2800$. The sum of dosages is $40+2(20)=80$, and the dosage-weighted outcome total is $12+2(12)=36$. Thus

U_\beta=36-(0.2800)(80)=13.600.

The squared-dosage sum is $40+4(20)=120$. The centered sum of squares and adjusted information are

\sum_i(G_i-\bar G)^2 = 120-\frac{80^2}{100} = 56, \qquad I_{\beta\beta\cdot\alpha} = (0.2800)(0.7200)(56) = 11.2896.

Therefore $S=13.600^2/11.2896=16.384$, with an approximate $\chi^2_1$ p-value of $5.2\times10^{-5}$. The data provide strong evidence against no additive genetic effect under the stated model. The score test supplies evidence of association in the specified model; it does not establish a causal effect.

Worked Example: Testing a Dominant Carrier Coding

A genetic model need not use dosage. In hypothetical data on 80 unrelated individuals, suppose 30 are allele carriers and 50 are noncarriers. There are 12 outcomes among carriers and 10 among noncarriers. Set $G_i=1$ for a carrier and $G_i=0$ for a noncarrier, and test the coefficient of this carrier indicator in the logistic model.

There are 22 outcomes, giving $\hat p_0=22/80=0.2750$. Here $\sum_iG_i=30$, $\sum_iG_i^2=30$, and $\sum_iG_iY_i=12$. The score is

U_\beta=12-(0.2750)(30)=3.750.

Because the indicator is one for 30 people and zero for the remaining 50, $\sum_i(G_i-\bar G)^2=30-30^2/80=18.750$. It follows that

I_{\beta\beta\cdot\alpha} = (0.2750)(0.7250)(18.750) = 3.7383, \qquad S=\frac{3.750^2}{3.7383}=3.761.

The approximate two-sided p-value is $0.0524$. At $\alpha=0.05$, this result does not meet the chosen rejection threshold. The test concerns a carrier-versus-noncarrier effect; it is not the same test as an additive dosage effect. The model and genotype coding should therefore be chosen to represent the scientific hypothesis, rather than selected after inspecting which gives the smallest p-value.

What the Score Test Does—and Does Not—Require

The key computational feature is that the score and information are evaluated at the null fit. In the examples, this means estimating the outcome proportion under an intercept-only model, then calculating the genetic score and its adjusted information from the observed genotypes. The alternative model’s maximum-likelihood estimate of $\beta$ is not needed.

This advantage does not remove the need to specify a correct likelihood or to check whether the data support the approximation. If all individuals have the same genotype coding, then $\sum_i(G_i-\bar G)^2=0$, so the data contain no information for separating a genetic effect from the intercept. If the estimated null probability is 0 or 1, the intercept estimate is not finite in the ordinary logistic model, and the standard calculation above is not valid. Very sparse outcomes or genotype groups can also make the large-sample chi-square approximation unreliable.

The one-degree-of-freedom result applies to one tested genetic parameter. If several genetic parameters are tested jointly, the score vector and its information matrix must be used, and the degrees of freedom depend on the number of independently testable parameters. That is not the calculation shown here.

Common Mistake

Wrong reasoning: “A score test uses the null model, so I can use the raw genotype-weighted number of cases as the test statistic.”
Correction: The raw weighted count is not scaled for its variability and does not account for the fitted intercept. Calculate the score as $\sum_iG_i(Y_i-\hat p_0)$, then divide its square by the nuisance-adjusted information before comparing with the reference distribution.
Wrong reasoning: “A p-value above 0.05 shows there is no genetic effect.”
Correction: A p-value above the chosen threshold means the data do not provide sufficient evidence to reject the null at that threshold. It does not prove that the effect is absent.

Key takeaway: A score test for a genetic effect evaluates the genetic score and its nuisance-adjusted information at the fitted null model. For a single effect, $U_\beta^2/I_{\beta\beta\cdot\alpha}$ is compared with $\chi^2_1$; its validity depends on the likelihood, genotype coding, and large-sample conditions.

Check Your Understanding

Answer these questions using the definitions and calculations in this tutorial.

In the additive logistic model, what does $\beta=0$ mean for the relationship between genotype dosage and outcome odds?
Why does the adjusted information involve $\sum_i(G_i-\bar G)^2$, rather than only $\sum_iG_i^2$?
For a score statistic of 4.0 testing one parameter, which reference distribution is used under the null, and what are its degrees of freedom?
Why can a genetic effect not be tested from these data if every individual has the same genotype coding?
How does a dominant carrier coding define a different genetic hypothesis from additive dosage coding?

What You'll Learn