A Simple Example (Exponential)

say \(X_1, ..., X_n\) iid with

\[f(x)=\beta \exp(-\beta x), x>0\]

and we want to test

\[H_0:\beta=\beta_0\text{ vs }H_a:\beta \ne \beta_0\]

Wald test

we can often find such a test by first finding the method of moments estimator:

\(E[X]=1/\beta\), so \(\hat{\beta}_1=1/\bar{x}\)

so if \(1/\bar{x}\) is close to \(\beta_0\) we should accept the null hypothesis, otherwise we should reject it. This is of course equivalent to \(\bar{x}\) is close to \(1/\beta_0\). Now \(var(X)=1/\beta^2\) and so

\[\sqrt{n}\frac{\bar{X}-1/\beta_0}{1/\beta_0} =\sqrt{n}(\beta_0\bar{X}-1)\sim N(0,1)\]

and we reject \(H_0\) if

\[ \sqrt{n}|\beta_0\bar{x}-1|>z_{\alpha/2} \]

LRT test (5.4.1)

\[ \begin{aligned} &f(\pmb{x}\vert\beta) = \prod_{i=1}^n \beta e^{-\beta x_i} = \beta^n e^{-\beta \sum x_i}\\ &l(\beta\vert\pmb{x}) = n\log \beta-\beta \sum x_i\\ & \frac{d l(\beta\vert\pmb{x})}{d \beta} = \frac{n}{\beta} -\sum x_i=0 \\ &\hat{\beta}=1/\bar{x} \end{aligned} \] for the likelihood ratio test statistic we find

\[ \begin{aligned} &\lambda(\pmb{x}) = \frac{\beta_0^n e^{-\beta_0 \sum x_i}}{(1/\bar{x})^n e^{-(1/\bar{x}) \sum x_i}} = \left(\beta_0\bar{x}\right)^ne^{n(1-\beta_0\bar{x})}\\ &-2\log \lambda(\pmb{x}) = -2n\left[\log \left(\beta_0\bar{x}\right)+(1-\beta_0\bar{x})\right]\\ & = \\ \end{aligned} \]

and we reject \(H_0\) if

\[-2n\left[ \log (\beta_0 \bar{x}) +1 - \beta_0 \bar{x}) \right] > \text{qchisq}(1-\alpha, 1)\]

How about using the LRT without the chisquare approximation? As always we have

“reject \(H_0\) if \(\lambda(\pmb{x})\) is small”

iff

“reject \(H_0\) if \(-2\log\lambda(\pmb{x})\) is large”"

iff

“reject \(H_0\) if \(\bar{x}\) is small or large”"

this can be seen by drawing the graph of the lrt as a function of \(\bar{x}\):

n <- 100; beta0 <- 1
fun <- function(xbar) 
  (-2)*n*(log(xbar*beta0)+1-beta0*xbar)
ggcurve(fun=fun, A=0.1, B=5)

so we reject \(H_0\) if \(\bar{x} < cv_1\) or \(\bar{x}>cv_2\).

Now under the null hypothesis \(X_i\sim Exp(\beta_0)\), so \(\sum_{i=1}^n X_i\sim Gamma(n,1/\beta_0\)) and we have

\[\alpha/2=P(\bar{X}<cv_1)=P(\sum_{i=1}^n X_i<n\times cv_1)=\text{pgamma}(n\times cv_1,n, \beta_0)\] and so \(cv_1=\text{qgamma}(\alpha/2,n, \beta_0)/n\). And also clearly we have \(cv_2=\text{qgamma}(1-\alpha/2,n, \beta_0)/n\).

This solution only works because we know the distribution of a sum of exponential variables (gamma), in another example we might be stuck with using the chisquare approximation.

Which test is better?

for that we need to find the power of the tests:

• Wald test:

\[ \begin{aligned} &Pow_A(\beta_1) = P \left( \sqrt{n}\vert\beta_0\bar{X}-1\vert>z_{\alpha/2}\vert\beta_1\right) = \\ &1-P \left(-z_{\alpha/2}< \sqrt{n}(\beta_0\bar{X}-1)<z_{\alpha/2}\vert\beta_1\right)= \\ &1-P \left(-z_{\alpha/2}+\sqrt{n}< \sqrt{n}\beta_0\bar{X})<z_{\alpha/2}\beta_0\vert\beta_1\right)= \\ &1-P \left([-z_{\alpha/2}+\sqrt{n}]\frac{\beta_0}{\beta_1}< \sqrt{n}\beta_1\bar{X})<[z_{\alpha/2}+\sqrt{n}]\frac{\beta_0}{\beta_1}\vert\beta_1\right)= \\ &1-P \left([-z_{\alpha/2}+\sqrt{n}]\frac{\beta_0}{\beta_1}-\sqrt{n}< \sqrt{n}(\beta_1\bar{X}-1)<[z_{\alpha/2}+\sqrt{n}]\frac{\beta_0}{\beta_1}-\sqrt{n}\vert\beta_1\right)= \\ &1-\left(\Phi([z_{\alpha/2}+\sqrt{n}]\frac{\beta_0}{\beta_1}-\sqrt{n})-\Phi([-z_{\alpha/2}+\sqrt{n}]\frac{\beta_0}{\beta_1}-\sqrt{n})\right) \end{aligned} \]

• likelihood ratio test:

\[ \begin{aligned} &Pow_B(\beta_1) = P(cv_1<\bar{X}<cv_2\vert\beta_1) = \\ &P(n\cdot cv_1<\sum X_i<n\cdot cv_2\vert\beta_1) = \\ &pgamma(n\cdot cv_2;n;\beta_1)-pgamma(n\cdot cv_1;n;\beta_1) \end{aligned} \] here are the power graphs:

alpha <- 0.05
beta1 <- seq(0.7, 1.5, length = 250)
Pow1 = 1-(pnorm((qnorm(1-alpha/2) + sqrt(n)) * beta1/beta0 - 
          sqrt(n)) - pnorm((-qnorm(1 - alpha/2) + sqrt(n)) * 
          beta1/beta0 - sqrt(n)))
cv <- qgamma(c(alpha/2, 1 - alpha/2), n, beta0)/n
Pow2 <- 1-(pgamma(n*cv[2], n, beta1) - 
             pgamma(n*cv[1], n, beta1))
df <- data.frame(beta=c(beta1, beta1),
          Power=c(Pow1, Pow2),
          Test=rep(c("Wald", "Lrt"), each=250))
ggplot(data=df, aes(beta, Power, color=Test)) +
  geom_line() 

and so we see that it depends on the true value of \(\beta\). Of course the big advantage of the LRT test here is that we were able to turn it into an exact test and so don’t have to worry about the sample size.