Which test for this data set of chemical concentrations

Question

I have a very simple data set: a parameter (concentration of a chemical) was measured at day 0, day 1 and day 2 in three subjects (there is a control group as well, but here all the values are always 0). The results are as follows:

        S1    S2    S3
day 0    0     0     0
day 1    0     5    45
day 2   70   250    60

The zero hypothesis is that there is no increase in this parameter (in the untreated control group, the values are always 0 at any time point in any subject).

The alternative hypothesis is that there is a significant increase in the parameter with time.

As a biologist, I would say that no statistics is even necessary, the data are obvious: not a single observation of a non-zero concentration in the control group, and an observable increase in all of the treated subjects. However, the reviewer insists on having this thing called a "p-value".

At first, I considered using a mixed random model as follows. Treat the control group and the treatment group separately. Control group tells us only that the chemical in question is never observed without intervention, and that the H0 should be: no observation of increase in the chemical.

Now, fit a mixed random model to the treatment group and ask whether the slope is significantly greater than zero. In R:

d <- read.table( text=
  "day   c  S
    0    0  S1
    1    0  S1
    2   70  S1
    0    0  S2
    1    5  S2
    2  250  S2
    0    0  S3
    1   45  S3
    2   60  S3", header= TRUE )
require( nlme )
l.m <- lme( c ~ day, random= ~ 1 + day | S, data= d, 
            control= lmeControl( opt= 'optim' ) )
anova( l.m )

The result of this naive approach is as follows:

            numDF denDF  F-value p-value
(Intercept)     1     5 1.451292  0.2822
day             1     5 3.753617  0.1104

I don't get it. Intuitively, this seems all wrong.

I have driven myself in utmost confusion. There must be a simple way of handling these data! (my question is: how should I handle the above data?)

Answer 1

You write:

As a biologist, I would say that no statistics is even necessary, the data are obvious: not a single observation of a non-zero concentration in the control group, and an observable increase in all of the treated subjects.

As a statistician, I agree with you. Also because of the size of the values.

However, the reviewer insists on having this thing called a "p-value".

When I review articles, I am often asking for less emphasis on p values, but ... others are different.

Your data have two issues that make models and p values tricky: 1) Very small sample 2) Very nonlinear trend. Unfortunately, with only 3 time points, you can't estimate nonlinear trends very well. But you don't want to try to fit a linear trend - that won't work, the trend isn't linear.

Perhaps someone else has a better idea, but I think a randomization test may be the way out of this. Assuming your control group also has 3 subjects, you have 18 numbers. The basic idea of a randomization test is that you randomly assign those 18 numbers to treatment and control, calculate some statistic on that, and repeat a lot of times. Then you see how often the statistic is higher in random groups than in your data.

What statistic? That depends on what is of interest to you. What does theory say should happen to the concentration?

Answer 2

Since this appears to be count data (Poisson) and repeated measures I would use glmer in the lme4 package:

d <- read.table( text=
"day   c  S
0    0  S1
1    0  S1
2   70  S1
0    0  S2
1    5  S2
2  250  S2
0    0  S3
1   45  S3
2   60  S3", header= TRUE )

library(lme4)
mymodel<-(glmer(c ~ day + (1|S), data=d, family='poisson'))

summary(mymodel)

I believe these tests are asymptotic results requiring large sample size. Given that the p-value is so small it doesn't bother me too much, but you could look into bootstrapping methods with the lme4 package to get the p-value. Bootstrap methods also require you to get a good sample of the space of potential values so still might not be great. The best would be to use a Bayesian approach which doesn't make asymptotic assumptions. These are things you might pursue under review...

Congratulations on the great results!

Answer 3

You could forget the difference between days 1&2, and instead just look at binary variables defined as "does the measured concentration increase". (In the observed data, true for all 3 subjects in the treatment group and false for all subjects in the control group).

If the null hypothesis is indeed "there is no increase in the parameter", there is no need for p-values, the data is logically inconsistent with this hypothesis ($p=0$)! However, in this case, the control group would be useless. Presumably, we have the control group because we are worried that the increases in concentration might be caused by something completely else than your treatment. To address this, a better formulation of the null hypothesis would be "probability of increasing concentration under treatment $\leq$ probability of increasing concentration in the control". Then, perform some test on the difference of proportions (by, e.g., Fisher's exact test)