Dear R-listers,

I'm trying to run a mixed effect model using the lmer() function and have
run into some issues in interpreting the p-values generated by
pvals.fnc(). The design is a between-subjects design, with two fixed
effects (condition & block; each with two levels), and one random effect
(subject). Additionally, I have a set of weights that I want to include.

When looking at the pvals.fnc() output,there appears to be a large
discrepancy between the pMCMC values and the t-statistic p-values. Whereas
one of the main effects and the interaction are far from significant
judging by the pMCMC values, they are highly significant when looking at
the t-statistic p-values (e.g. Condition: pMCMC = 0.2294; Pr(>|t|) = 0.0000
& Condition*Block: pMCMC = 0.3296; Pr(>|t|) = 0.0000) . I have read that
the t-statistic based p-values are less conservative, but the difference
between these two values seems really extreme.

Below some code that simulates the model and the data. The original data
set has two precise characteristics that might influence the results, so I
tried to simulate those characteristics in the mock data. That is: 1)
there's fewer observations in block A than in block B; and 2) the weights
for observations in block A generally are lower than those for block B.

Running this code reproduces the original observation of conflicting pMCMC
and p-T-test values. However, when excluding the weights argument from the
lmer model, these values seem to converge, suggesting that the weights
specification might be underlying these problems.

In short, my question is whether anyone knows why these values diverge and
what I could do to address this issue.

Many thanks in advance!

Tom

block <- as.factor(c(rep('a', times = 20), rep('b', times = 200)))
condition <- as.factor(c(rep(c('x', 'y'), each = 10), rep(c('x','y'), each
= 100)))
contrasts(block) <- c(-0.5, 0.5)
contrasts(condition) <- c(-0.5, 0.5)

subject <- c(rep(1:4, each = 5), rep(1:4, each = 50))

intercept <- 100
block.me <- 20
condition.me <- 30
err <- rnorm(length(block), sd = 20)
weights <- c(rep(1, times = 20), rep(10, times = 200))

y <- intercept + ifelse(block == 'a', block.me, 0) + ifelse(condition ==
'x', condition.me, 0) +
    ifelse(block == 'a' & condition == 'x', 30, 0) + (subject * 10) + err


fm.1 <- lmer(y ~ block * condition + (1 | as.factor(subject)),
             weights = weights, REML = FALSE)
fm.1.mcmc <- pvals.fnc(fm.1, addPlot=F)

	[[alternative HTML version deleted]]


