# Multiple testing and logistic regression

I want to perform a number of univariate regressions with different symptoms (e.g. fever, cough, sneezing) as the response variable and one categorical variable (which is always the same each time) as an explanatory variable - age group (0-4, 5-14, 15-64, 65+). In total there are 18 symptoms so hence, I want to do 18 regressions and want to correct for multiple testing using Holm's method.

For each regression I get 4 p-values - one overall age group p-value and 3 p-values comparing individual age groups to a reference age group (15-64). I am unsure about how to calculate confidence intervals and adjusted p-values using Holm's method in this scenario. Does anybody know how to do this?

The R rms package supports general simultaneous confidence intervals. The analyst can specify a series of comparisons, and simultaneous adjustment will be done over all the current set of comparisons. For more information see the rms contrast.rms function, along with Predict and other functions. Here's an example using ols. Substitute lrm for ols to do binary or proportional odds logistic regression.

require(rms)
set.seed(1)
n <- 800
treat <- factor(sample(c('drug','placebo'), n,TRUE))
sex   <- factor(sample(c('female','male'),  n,TRUE))
age   <- rnorm(n, 50, 10)
y     <- .05*age + (sex=='female')*(treat=='drug')*.05*abs(age-50) + rnorm(n)
f     <- ols(y ~ rcs(age,4)*treat*sex)

# show separate estimates by treatment and sex

plot(Predict(f, age, treat, sex='female'))
plot(Predict(f, age, treat, sex='male'))
ages  <- seq(35,65,by=5); sexes <- c('female','male')
w     <- contrast(f, list(treat='drug',    age=ages, sex=sexes),
list(treat='placebo', age=ages, sex=sexes),
conf.type='simultaneous'))
xYplot(Cbind(Contrast, Lower, Upper) ~ age | sex, data=w,
ylab='Drug - Placebo')
xYplot(Cbind(Contrast, Lower, Upper) ~ age, groups=sex, data=w,
ylab='Drug - Placebo', method='alt bars')


For confidence interval, I think you can just use the logistic regression outputs. For adjust p-values, refer to this wiki link: http://en.wikipedia.org/wiki/Holm%E2%80%93Bonferroni_method

To my understanding, you are first test an overall age group, then testing a nested hypothesis. I would recommend you test overall group effects first, then go deeper to test individual groups.

• No, 95% confidence intervals are constructed such that, with a single data sample, the actual point is within the interval with probability $0.05$. It does not account for multiple testing. The probability of finding a significant effect for at least one measure (when there is none, assuming measures are independent) would be $1-0.95^{18}\approx 0.6>>0.05$. – Nameless Jun 25 '13 at 14:56
• @Nameless I'm pretty sure that your interpretation of a 95%-confidence interval is wrong. In a frequentist framework, the value is fixed and the endpoints of the confidence interval are random variables. For a single sample, the confidence interval either includes the true value or not. To say a single confidence interval includes the true value with a probability is not applicable. See here and here. – COOLSerdash Jun 25 '13 at 17:47
• Yes, please replace "probability" with "percent of cases" above. Which, as it happens, is the same in a frequentist framework.. – Nameless Jun 25 '13 at 19:36

So lets imagine there are only 2 symptoms and I got the following results:

### Fever: ###

Overall p-value
Wald test:
----------

Chi-squared test:
X2 = 16.8, df = 3, P(> X2) = 0.00078

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept)                   1.5849     0.1103  14.366  < 2e-16 ***
relevel(agegp, "15-64")0-4    0.6323     0.3361   1.881  0.05994 .
relevel(agegp, "15-64")5-14   0.7273     0.2182   3.333  0.00086 ***
relevel(agegp, "15-64")65+   -1.0253     0.6364  -1.611  0.10718

### Cough ###

Overall p-value
Wald test:
----------

Chi-squared test:
X2 = 8.5, df = 3, P(> X2) = 0.036

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept)                   2.2693     0.1416  16.024  < 2e-16 ***
relevel(agegp, "15-64")0-4   -0.3234     0.3189  -1.014  0.31051
relevel(agegp, "15-64")5-14  -0.5850     0.2050  -2.854  0.00432 **
relevel(agegp, "15-64")65+   -0.7652     0.7945  -0.963  0.33544


What would be my adjusted p-values (using Holm's method) in this scenario of having only 2 symptoms and how would I calculate confidence intervals?