6
$\begingroup$

I received this question from a colleague which left me with a desire to find out a bit more about it.

If you have two groups A and B with means--for example, 7.9 and 5.8--is it possible that the confidence interval for the difference between the groups can be such that it doesn't encompass the current difference? In this case the 95% CI is 0.9-1.79 when the difference is 2.1.

(Edit: For clarification this was from an M.Sc. statistics question where 150 angina pasients were divided into two groups (A and B) and treated with differing combinations of physical activity and medication. Cholesterol levels were measured after the intervention: Group A had an average of 7.9 +/- 4 mmol/l and group B 5.8 +/- 2.9 mmol/l. The CI stated is for the difference between the group averages)

I believe this is possible but I would like to be able to give him a concrete example (preferably with some data that he can process in R) to show him where this can occur.

Cheers!

$\endgroup$
8
$\begingroup$

It is possible for a confidence interval of a mean not to include the sample mean.

It is not part of the definition of a CI that it must always cover the sample mean. Thus one may, in theory, construct a CI procedure that never covers the sample mean. But most people would consider that a bad procedure.

Let us therefore examine CI procedures that not only have been seriously proposed, but have been studied and found to be good. One of them is the "generalized lognormal confidence interval" procedure explained and studied (via simulation) by Ulf Olsson in the Journal of Statistics Education (Volume 13, Number 1 (2005), http://www.amstat.org/publications/jse/v13n1/olsson.html). This procedure is a reasonable one to use when the (natural) logarithms of the $n$ data are assumed to be independent and identically distributed with a Normal distribution.

Recall that when the population mean log is $\mu$ and the population standard deviation of the logarithms is $\sigma$, then the population mean is $\exp(\mu + \sigma^2/2)$. (This relationship uses the lognormal assumption.) We will obtain confidence limits for $\mu+\sigma^2/2$; exponentiating them will give confidence limits for the population mean. The procedure is based on a "generalized confidence interval" known to produce good confidence intervals for complicated combinations of parameters like this one.

Calculate the mean logarithm $\bar y$ and the sample variance of the logarithms $s^2$ from the data. A symmetric confidence interval of size $\alpha$ for $\mu + \sigma^2/2$ is found by identifying the middle $100 - 100\alpha\%$ of the distribution of

$$T_{2} = \bar y - Z \sqrt{A^2/n} + A^2/2$$

where $Z$ and $A^2$ are independent variates, $Z$ has a standard Normal distribution,

$$A^2 = \frac{s^2}{U^2 / (n-1)},$$

and $U^2$ has a chi-squared distribution with $n-1$ degrees of freedom. Because the distribution of $T_2$ is difficult to work with analytically, we may estimate it through simulation. When exponentiated, the $\alpha/2$ and $1-\alpha/2$ quantiles of $T_2$ are the lower and upper confidence limits for $\exp(\mu+\sigma^2/2)$.

Olsson's work indicates that once $n \ge 20$ or so, this procedure tends to achieve its nominal characteristics for $\alpha=0.05$. That is, about $2.5\%$ of the time it is less than $\mu+\sigma^2/2$ and $2.5\%$ of the time it is greater than $\mu+\sigma^2/2$.

One does not have to look hard to find datasets that (a) appear to meet the assumptions of this test yet for which (b) the confidence interval does not include the sample mean. Here is one with $n=50$:

#0.08 0.14 0.21 0.25 0.28 0.3 0.35 0.37 0.39 0.41 0.46 0.51 0.55 0.55 0.66 0.66 0.69 0.71 0.74 0.74 0.77 0.81 0.85 1.04 1.09 1.1 1.17 1.18 1.19 1.25 1.29 1.38 1.54 1.62 1.62 1.68 1.74 1.87 2.11 2.29 2.37 2.42 2.93 2.99 4.8 5.12 5.94 7.09 11.26 120

The generalized CI of $(1.56, 3.95)$ does not include the sample mean of $4.03$.

(This was computed using ten million simulated values of the distribution of $T_2$, so it should be pretty accurate. Twenty independent simulations using just a million simulated values never produced an upper limit larger than $4.02$, still below the sample mean.)

Although the last few data values ($11.26, 120$) may look like outliers, their logarithms are not. Here is a histogram of their logs:

Figure

OK, that final value of $\log(120)$ looks a wee bit high. But the (very) powerful Shapiro-Wilk test does not strongly reject the Normal hypothesis ($p = 0.012$). This provides some insight: lognormal (and other heavy-tailed) distributions frequently produce unusually large values by their very nature. Such values can strongly influence the sample mean but should have less influence on estimates of the underlying distributional properties. We shouldn't find anything paradoxical about this.

(Although this example concerns a single group, it could be generalized to compare the difference of means between two groups, at some cost in complexity of the calculations. Nothing really changes, though: we may think a CI must include the sample mean only when we become so accustomed to using Normal-theory calculations that we come to believe, through sheer repetition, that all CIs must share their properties.)


The following R code will reproduce these calculations and allow you to explore the properties of the generalized lognormal confidence interval. In particular, if you are wondering whether the failure of this CI sometimes to cover the sample mean might be due to an error in coding (a possibility that always worries me!), or if it isn't really a CI for the population mean in the first place, you can reproduce a part of Olsson's work by simulating the coverage of this CI, as in

set.seed(17)
x.mean <- exp(mu + sigma^2/2)                 # The true lognormal (population) mean
sim <- replicate(1e3, {ci.lognormal(exp(rnorm(n, mu, sigma)))})
mean(sim[1, ] <= x.mean & sim[2, ] >= x.mean) # Fraction of times covering the true mean

The output of $0.949$ shows that this nominal $95\%$ interval has covered the true mean $94.9\%$ of the time, which is excellent. (I chose this particular CI procedure specifically because it is so good.) By contrast, you could check how often this interval covers (say) the geometric mean:

mean(sim[1, ] <= exp(mu) & sim[2, ] >= exp(mu)) # Fraction of times covering the GM

The output of $0.899$ confirms it is not a $95\%$ confidence interval for the geometric mean.

Here is the full code (which you will need to compile before running the preceding lines).

#
# Generalized confidence intervals for lognormal means.
#
ci.generalized <- function(y, alpha=0.05, n.boot=1e4) {
  n <- length(y); m <- mean(y); s2 <- var(y)
  z <- rnorm(n.boot)
  u2 <- rchisq(n.boot, n-1)
  a2 <- s2 / (u2 / (n-1))
  sim <- m - z * sqrt(a2 / n) + a2 / 2         # Simulated distribution of T^2
  return(quantile(sim, c(alpha/2, 1-alpha/2))) # CI for mu + sigma^2/2
}
ci.lognormal <- function(x, ...) exp(ci.generalized(log(x), ...))
#
# Experiment with simulated data.
#
n <-50
mu <- 0                    # Mean log (the actual value is irrelevant)
sigma <- 1/5               # SD of logs (affects the shape of `x`)
set.seed(968)              # Reproduces the example in the text
y <- rnorm(n, mu, sigma)
y <- (y - mean(y)) / sd(y) # (Make sure the logs start out looking very Normal)
x <- round(exp(y), 2)      # Model a limited-precision data collection process
x[1] <- 120                # Tweak the data to give them a large sample mean
#
# Display the data and look at the CI of the mean.
#
hist(log(x))
ci <- ci.lognormal(x, alpha=0.05, n.boot=1e7) # Takes a few seconds
c(p.value=shapiro.test(log(x))$p.value, ci.lognormal(x), sample.mean=mean(x))
$\endgroup$
  • 1
    $\begingroup$ Please read the paper I referenced, @John: This is an honest-to-goodness CI for the arithmetic mean, not the geometric mean. $\endgroup$ – whuber Feb 20 '15 at 16:49
  • 1
    $\begingroup$ @John, no: it is not an interval for the GM. That would be of the form $\exp(\bar y \pm t_{\alpha/2, n-1} \sqrt{s^2/n})$, which it manifestly is not. I entreat you to look at the reference, or read the formulas in my post, or study the code therein, before you comment further. $\endgroup$ – whuber Feb 20 '15 at 17:03
  • 1
    $\begingroup$ @John I have inserted some code to help you study the characteristics of this generalized confidence interval procedure. $\endgroup$ – whuber Feb 20 '15 at 17:23
  • 1
    $\begingroup$ @John The problem lies not with ci.generalized but with the distribution itself. With such an enormous SD of the logs (sigma in the code), we can no longer trust a simulation to give useful guidance. In practice, lognormal data have log SDs around $0.1$ to $2$ or so, but any value above $1$ can be difficult to handle unless you do enormous simulations and take care to handle floating point overflows. Your example has a mean of $\exp(20+20^2/2)\approx 3.5\times 10^{95}$ and an SD that (in IEEE double precision) is infinite! (It's actually around $10^{356}$.) $\endgroup$ – whuber Feb 20 '15 at 18:05
  • 2
    $\begingroup$ This a great answer. Thanks for taking the time. I should thank @John as the discussion between the two of you was really insightful. Thanks again gents. $\endgroup$ – Hantabaru Feb 20 '15 at 20:39
0
$\begingroup$

Presumably, the confidence interval you are looking at is of the form $\bar{X}_1 - \bar{X}_2 \pm M$ where $M$ is some margin of error measure.This of course includes $\bar{X}_1 - \bar{X}_2$.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.