Hello,

I am trying to do a BACI analysis on count data, and I am having trouble interpreting the output from multcomp::glht. I don't understand how the contrast's coefficients are related to effect size (if at all??).

I have 5 treatment conditions (one is a control), and I have counts from before the treatments were applied and after. Let's say that my model form is this ("Period" is the 'before' vs 'after' factor):

    m.pois <- glm(Y_count ~ Treatment + Period + Treatment:Period,
                       data = df.temp,
                       family = "poisson")

As in all BACI designs, I am interested in the interaction term, i.e., the differences of the differences. For example, I'd like to test whether TreatmentVR30 changed more than the Control did: (TreatmentVR30Later - TreatmentVR30Before) - (ControlLater - ControlBefore).

I have done the math, and I created all my planned contrasts, run them through the multcomp::glht, and I am struggling to interpret the output. As an example of my confusion, I ran the same contrast in two directions (A-B and B-A), which should give the same result (one positive, one negative):

    contr <- rbind(
      "VR30 vs Control" =     c(0, 0, 0, 0, 0, 0, -1, 0, 0, 1),
      "Control vs VR30" =     c(0, 0, 0, 0, 0, 0, 1, 0, 0, -1)  )
    m.pois.contr <- summary(glht(m.pois, contr))

which works perfectly, returning one positive and one negative estimate, as expected:

    Linear Hypotheses:
                         Estimate Std. Error z value Pr(>|z|)
    VR30 vs Control == 0   0.7354     0.5621   1.308    0.191
    Control vs VR30 == 0  -0.7354     0.5621  -1.308    0.191
    (Adjusted p values reported -- single-step method)

Understanding that the estimates are in log space (due to the link function of the poisson family in the glm), I back transformed using exp(coef(m.pois.contr) to get:

    VR30 vs Control Control vs VR30
          2.0862414       0.4793309

So, which is it? Did Control change more than VR30, or did VR30 change more than Control, and for both questions, by how much?

Clearly I am missing something here. I expect that this will be a simple fix, but surprisingly, I cannot find it anywhere online.

Thanks in advance to anyone who can help,

David Robichaud, Victoria, BC, Canada

______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Reply via email to