Replicate post-hoc test performed with glht

41 Views Asked by At

I´m performing an ancova followed by post-hoc Dunnet´s test and arrive at a set of p-values. I´m using glht in R to perform the Dunnet´s test.

I would like to replicate the Dunnet´s test, and understand which underlying values are used for deriving at the p-value.

Data consists of several groups (Treatment) for which I have a measurement (value) and a continuous co-variate (covariate).

head(data)
  Treatment  value covariate
1         A  95.76      58.0
2   Control  95.16      59.8
3         B 143.92      41.7
4         C 119.80      51.4
5         D 100.32      43.5
6         D 107.44      51.3
 
model=aov(value~Treatment+covariate, data=data)
postHocs <- glht(model, linfct = mcp(Treatment = "Dunnet"))
summary(postHocs)

     Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Dunnett Contrasts


Fit: aov(formula = value ~ Treatment + covariate, data = data)

Linear Hypotheses:
                 Estimate Std. Error t value Pr(>|t|)  
A - Control == 0  -1.7318     6.8485  -0.253   0.9986  
B - Control == 0  -3.0630     8.7941  -0.348   0.9941  
C - Control == 0  11.7043     8.3848   1.396   0.4483  
D - Control == 0  -0.2675     8.4407  -0.032   1.0000  
E - Control == 0 -44.1218    17.5675  -2.512   0.0529 .
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)

It does not seem to be the fitted values which the Dunnet´s test uses as the following does not yield the same results:

fitted_values=model$fitted.values
treatments=data$Treatment
DescTools::DunnettTest(fitted_values, treatments)

  Dunnett's test for comparing several treatments with a control :  
    95% family-wise confidence level

$Control
                 diff      lwr.ci    upr.ci   pval    
A-Control  -0.7236364 -18.3798733 16.932601 1.0000    
B-Control   4.3463636 -12.9381214 21.630849 0.9486    
C-Control  18.2763636   0.9918786 35.560849 0.0343 *  
D-Control   6.1781818 -11.4780551 23.834419 0.8332    
E-Control -23.3156364 -41.4078954 -5.223377 0.0066 ** 

---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Very best, M

###A reproducible example:

> df=data.frame(Treatment=rep(c("A","B"), 10),
+               covariate=rnorm(10, 20),
+               value=rnorm(10,50))
> df
   Treatment covariate    value
1          A  19.15914 51.22763
2          B  21.38436 49.19822
3          A  18.74451 48.91961
4          B  20.07014 49.84247
5          A  21.71144 48.92824
6          B  19.39709 49.86101
7          A  19.52783 49.40269
8          B  19.36463 47.81603
9          A  19.71423 50.24082
10         B  20.13811 49.74064
11         A  19.15914 51.22763
12         B  21.38436 49.19822
13         A  18.74451 48.91961
14         B  20.07014 49.84247
15         A  21.71144 48.92824
16         B  19.39709 49.86101
17         A  19.52783 49.40269
18         B  19.36463 47.81603
19         A  19.71423 50.24082
20         B  20.13811 49.74064
> 
> #Ancova
> library(multcomp)
> model=aov(value~Treatment+covariate, data=df)
> posthocs=glht(model, linfct=mcp(Treatment="Dunnet"))
> summary(posthocs)

     Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Dunnett Contrasts


Fit: aov(formula = value ~ Treatment + covariate, data = df)

Linear Hypotheses:
           Estimate Std. Error t value Pr(>|t|)
B - A == 0  -0.4135     0.4054   -1.02    0.322
(Adjusted p values reported -- single-step method)

> 
> #Trying to reprocude stats from glht
> model=aov(value~Treatment+covariate, data=df)
> fitted_values=model$fitted.values
> treatments=df$Treatment
> 
> DescTools::DunnettTest(fitted_values, treatments)

  Dunnett's test for comparing several treatments with a control :  
    95% family-wise confidence level

$A
          diff    lwr.ci     upr.ci    pval    
B-A -0.4521208 -0.565957 -0.3382846 1.3e-07 ***

---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
0

There are 0 best solutions below