[R-sig-ME] glmer formula with partially-crossed fixed effects

[Dan McCloy]

You may want to consider Helmert coding applied to a single new factor made
from A and B, with 3 levels: A1, A2B1, A2B2. Make sure that A1 ends up as
the highest value in the factor.

foo <- factor(1:3, labels=c("A2B1", "A2B2", "A1"))
contrasts(foo) <- contr.helmert
contrasts(foo) <- contrasts(foo) / rep(2:3, each=3)

The last line is to get "sum to one" factor weights.  See ?contr.helmert
for more info.
-- dan

Daniel McCloy
http://dan.mccloy.info/
Postdoctoral Research Fellow
Institute for Learning and Brain Sciences
University of Washington





On Wed, Jul 29, 2015 at 10:05 AM, Becky Gilbert <beckyannegilbert at gmail.com>
wrote:

> Dear List,
>
> I'm trying to create a logistic linear mixed effect model with 3 fixed
> factors, 2 of which are partially crossed (I think this is the right
> term?).  Factor A has 2 levels, factor B has 2 levels that only occur
> within one level of A, and factor C has 2 levels.  There are also 2 crossed
> random factors, subject and item.  I've provided an example data set with
> the same structure at the end of this email.
>
> My question is, how do I code the partially crossed fixed effects (A and B)
> in order to get all the fixed effects of interest?
>
> Here's what I'd like to get from the model:
> 1. main effect of A
> 2. effect of B (which occurs within one level of A)
> 3. main effect of C
> 4. A:C interaction
> 5. B:C interaction (which occurs within one level of A)
>
> So far I have tried the following:
>
> 1. Code factor B as two levels, which are only valid within level 1 of A.
> Within level 2 of A, factor B = NA. This gives me the following
> combinations:
> A   B
> 1   1
> 1   2
> 2   NA
>
> I get an error using this formula (random intercepts only, for
> convenience):
> glmer(dv ~ 1 + A + B + C + A:C + B:C +
>     (1|subject) +
>     (1|item),
>     data = mydata,
>     family = binomial)
>
> Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
>   contrasts can be applied only to factors with 2 or more levels
>
> 2. Use dummy coding for factor B, where A1B1 is coded as 1 and everything
> else is 2, and then look at the effect of B within A1.  Recoding the B
> factor in this way gives me:
> A  BnoNA
> 1   1
> 1   2
> 2   2
>
> glmer(dv ~ 1 + A + BnoNA + C + A:C + BnoNA:C +
>                 (1|subject) +
>                 (1|item),
>               data = mydata,
>               family = binomial)
>
> This model converges without errors/warnings, but I don't know how to get
> the effect of B within A1.  Also, it seems like I don't want to model the
> main effects and interactions of 'BnoNA' (i.e. the effect of A1B1), but I
> don't know how else to write the formula...
>
> 3. Remove factor A from the model, recode B with 3 levels (i.e. 1 = A1B1, 2
> = A1B2, 3 = A2), then use contrasts of B levels to get the main effect of A
> (i.e. (B = 1 and B = 2) vs B = 3).
> A  B   B3levels
> 1  1    1
> 1  2    2
> 2  NA 3
>
> glmer(dv ~ 1 + B3levels + C + B3levels:C +
>                 (1|subject) +
>                 (1|item),
>               data = mydata,
>               family = binomial)
>
> This model converges without errors/warnings, but I don't know how to add
> the contrasts for (B1 and B2) vs B3, or the interaction between this
> contrast coding (i.e. main effect of A) and C.
>
> 4. Run 2 separate models. The first model uses the full data set to get the
> effects of A, C and the interaction.  The second model uses a subset of the
> data (only level 1 of A) to get the effect of B within A1, and the B:C
> interaction.  My concern with this solution is that the random effects
> structure would differ between the two models, because the random effects
> for subjects and items would be based on data from both levels of A in the
> first model and only A1 in the second model.
>
> Apologies for the novice question, messy data simulation code and possible
> confusion of terms!  I searched the list archives and I googled extensively
> but couldn't find an answer.
>
> Thanks in advance.
> Becky
>
>
> ### Here is code to produce some example data of the same structure.
>
> A <- c(rep(rep(c(1,1,1,1,0,0), each = 13), 13),
>        rep(rep(c(0,1,1,1,1,0), each = 13), 13),
>        rep(rep(c(0,0,1,1,1,1), each = 13), 13),
>        rep(rep(c(1,0,0,1,1,1), each = 13), 13),
>        rep(rep(c(1,1,0,0,1,1), each = 13), 13),
>        rep(rep(c(1,1,1,0,0,1), each = 13), 13))
> B <- c(rep(rep(c(1,1,0,0,NA,NA), each = 13), 13),
>        rep(rep(c(NA,1,1,0,0,NA), each = 13), 13),
>        rep(rep(c(NA,NA,1,1,0,0), each = 13), 13),
>        rep(rep(c(0,NA,NA,1,1,0), each = 13), 13),
>        rep(rep(c(0,0,NA,NA,1,1), each = 13), 13),
>        rep(rep(c(1,0,0,NA,NA,1), each = 13), 13))
> C <- c(rep(rep(c(1,0,1,0,1,0), each = 13), 13),
>        rep(rep(c(0,1,0,1,0,1), each = 13), 13),
>        rep(rep(c(1,0,1,0,1,0), each = 13), 13),
>        rep(rep(c(0,1,0,1,0,1), each = 13), 13),
>        rep(rep(c(1,0,1,0,1,0), each = 13), 13),
>        rep(rep(c(0,1,0,1,0,1), each = 13), 13))
> mydata <- data.frame(A = factor(A),
>                      B = factor(B),
>                      C = factor(C),
>                      subject = factor(rep(1:78, each = 78)),
>                      item = factor(rep(1:78, 78)),
>                      dv = factor(sample(0:1, 78*78, replace = TRUE)))
>
> # recode B factor - change NAs into 0s
> mydata$BnoNA <- mydata$B
> mydata$BnoNA[which(is.na(mydata$BnoNA))] <- 0
>
> # recode B factor - change NAs into a 3rd level
> mydata$B3levels <- as.numeric(mydata$B)
> mydata$B3levels[which(is.na(mydata$B3levels))] <- 3
> mydata$B3levels <- as.factor(mydata$B3levels)
>
> > str(mydata)
> 'data.frame': 6084 obs. of  8 variables:
>  $ A       : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
>  $ B       : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
>  $ C       : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
>  $ subject : Factor w/ 78 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 1 1
> ...
>  $ item    : Factor w/ 78 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10
> ...
>  $ dv      : Factor w/ 2 levels "0","1": 2 2 1 1 1 2 1 2 2 1 ...
>  $ BnoNA   : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
>  $ B3levels: Factor w/ 3 levels "1","2","3": 2 2 2 2 2 2 2 2 2 2 ...
>
> #Because B is only partially crossed with A, there should be no valid
> entries of B in one of the levels of A, as confirmed by:
> > table(mydata$A, mydata$B)
>        0    1
>   0    0    0
>   1 2028 2028
> > table(mydata$A,mydata$BnoNA)
>        0    1
>   0 2028    0
>   1 2028 2028
> > table(mydata$A,mydata$B3levels)
>
>        1    2    3
>   0    0    0 2028
>   1 2028 2028    0
>
> ### My R environment
> R version 3.0.2 (2013-09-25)
> Platform: x86_64-w64-mingw32/x64 (64-bit)
> lme4_1.1-7
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>

	[[alternative HTML version deleted]]