Multiple Tests, an Introduction

[This article was first published on Freakonometrics » R-english, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Last week, a student asked me about multiple tests. More precisely, she ran an experience over – say – 20 weeks, with the same cohort of – say – 100 patients. An we observe some{i,t}


(here, I just generate some fake data). I can visualize some trajectories, over the 20 weeks,

for(i in 1:length(cl)){
lines(1:20,X[i,],type="b",col=cl[i])  }


She wanted to compare the averages, over the experiments.

for(i in 1:length(cl)){
lines(1:20,X[i,],type="b",col=cl[i])  }

The question is what do you want to test ? Is it,jin{1,cdots,%20T},%20mu_i=mu_j


Because that one might be not very informative (if we reject, we don’t really know why). Here, I did consider completely independent observations, over individuals, and over time. Nevertheless, if we perform 190 tests on pairs, we should not be suprised – even if is valid – that we reject it many times !

> for(i in 1:(nb-1)){ for(j in (i+1):nb){
+ P=c(P,t.test(X[,i],X[,j])$p.value) }}
> mean(P>.05)
[1] 0.9368421

i.e. for 12 pairs, we reject the assumption of identical mean. Even if that assumption is true, here. That’s the problem with multiple tests.

Actually, observe that over those 190 tests, the are uniformly distributed over the unit interval. Thus, for (almost) 95% of the tests, the exceeds 5%.

Now, what if we relax the assumption of independence ? Here, we need to be more specific. A first idea might be to assume that{i,t}=mu_t%20+%20varepsilon_{i,t}

for some trend We can actually assume that the trend depends on the individuals, and write, for instance,{i,t}=alpha_i%20+%20beta_i%20t%20+%20varepsilon_{i,t}

Standard regression techniques can be considered, under those models. We’ll probably discuss that point before the end of this year.

My point was the following : if the assume that{X}_i=(X_{i,1},cdots,X_{i,T})

consist of independent observations, where all{i,t}‘s are assumed to be independent, if we use pairwise comparison tests for the means, in 95% of the tests, the exceeds 5%, and in 90% of the tests, the exceeds 10%. That the interpreation of the property that are uniformly distributed on the unit interval. Now, what if we relax the independence assumption. What if obsevations were – somehow – correlated.

A first idea could be to consider a linear trend,


Here, we still compute pairwise comparison tests

> for(i in 1:(nb-1)){ for(j in (i+1):nb){
+ P=c(P,t.test(X[,i],X[,j])$p.value) }}

and we counts the number of pairs where the exceeds some given values


If we look at the proportions, as a function of the slope, we have the following graph

In the middle, we have the independent case, e.g. 97.5% of the the exceeds 2.5% – the upper curve. But, as the slope increase (or decrease), the proportion decreases. With a 0.2 slope, 90% of the the exceeds 2.5%.

One can also assume that{X}_i=(X_{i,1},cdots,X_{i,T})

is some time series. An autoregressive one, for instance,{i,t}=varphi%20X_{i,t-1}+varepsilon_{i,t}

Here we use to following code to generate that vector

for(i in 1:(n-2)){

The graph below is the evolution of the proportion of exceeding some threshold, as a function of the autoregressive coefficient.

The value on the left is the independent case (and are uniformly distributed). The stronger the autocorrelation, the larger the proportion of exceeding any value.

Finally, consider the case where components of vector{X}_i=(X_{i,1},cdots,X_{i,T})

are exchangeable.


In that case, we have the following graph, with the proportion of exceeding some value, as the function of the correlation of (any) pair{i,s},X_{i,t}).

This is the idea that is used in any correction-type method, with multiple tests, for instance Bonferroni correction. When a run tests, instead of comparing the with, we use a correction on, But as we can see here, it will depend on the kind of dependence we assume there might be.

Next time, I will discuss my student’s data, see what could be said.

To leave a comment for the author, please follow the link and comment on their blog: Freakonometrics » R-english. offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Never miss an update!
Subscribe to R-bloggers to receive
e-mails with the latest R posts.
(You will not see this message again.)

Click here to close (This popup will not appear again)