[R] lmer( ) parameter estimates changing depending on dummy-coded reference level
Bert Gunter
gunter.berton at gene.com
Wed Aug 7 19:14:40 CEST 2013
I think it's fair to say that this should really be posted on the
mixed-model specific (especially using lme4) list r-sig-mixed-models
and not here.
Cheers,
Bert
On Wed, Aug 7, 2013 at 9:01 AM, Alan Mishler <amishler at umd.edu> wrote:
> Hi everyone,
>
> I'm writing with a question about lmer( ) in the lme4 package. I've
> searched around for
> answers and done quite a bit of experimentation with toy data sets to
> figure out my issue, and I haven't been able to resolve it.
>
> I'm running linear mixed effects models on a large, sparse dataset in
> which I'm regressing reaction time (a continuous variable) on several
> categorical factors: Block (Block1/Block2/Block3), Group
> (monolingual/bilingual), and Type (target/nontarget).
>
> As a way of examining simple effects, I am dummy-coding specific
> factors, setting each level of a given factor as the reference level
> in turn. For example, I generate three models with each of the three
> levels of Block coded as the reference level, without changing the
> codings of the other factors:
>
> ## Model with Block 1 as reference level
> contrasts(nback.low$Group) <- c(1, 0) # monoling ref
> contrasts(nback.low$Type) <- c(1, -1)
> contrasts(nback.low$Block) <- matrix(c(0, 1, 0, 0, 0, 1), ncol=2) #B1
> ref, B2: 1, B3: 2
> glmerNL.RS.SI_RTgxb1 <- glmer(WinRTs~(Group*Block*Type) +
> (1+Block+Type|Subject) + (1+Group+Block|Item),data=nback.low)
>
> ## Model with Block 2 as reference level
> contrasts(nback.low$Group) <- c(1, 0) # monoling ref
> contrasts(nback.low$Type) <- c(1, -1)
> contrasts(nback.low$Block) <- matrix(c(1, 0, 0, 0, 0, 1), ncol=2) #B2
> ref, B1: 1, B3: 2
> glmerNL.RS.SI_RTgxb2 <- glmer(WinRTs~(Group*Block*Type) +
> (1+Block+Type|Subject) + (1+Group+Block|Item),data=nback.low)
>
> ## Model with Block 3 as reference level
> contrasts(nback.low$Group) <- c(1, -1) # monoling 'ref'
> contrasts(nback.low$Type) <- c(1, 0) # target ref
> contrasts(nback.low$Block) <- matrix(c(1, 0, 0, 0, 1, 0), ncol=2) #B3
> ref, B1: 1, B2: 2
> glmerNL.RS.SI_RTbxt3 <- glmer(WinRTs~(Group*Block*Type) +
> (1+Block+Type|Subject) + (1+Group+Block|Item),data=nback.low)
> summary(glmerNL.RS.SI_RTbxt3)
>
> The issue I'm having is that contrasts that I believe should be
> identical are not. Below are summaries of the three
> models. You can see that the estimate of the fixed effect of Block1
> (the contrast between Block 1 and Block 2) is -117.98 in the first
> model and 118.478 in the second model. To my understanding, they
> should be identical except for the sign. Similar discrepancies can be
> seen in the other Block contrasts.
>
> There are two subjects who have no data at Block 1, so I removed them
> and re-ran the models, but the same issue occurred. Separately, I
> removed the random effects for Item, without removing those two
> subjects, and when I did that the discrepancies disappeared. I have a
> feeling this means that my models are too complex for my data, but I'm
> not sure what I should look at to (dis)confirm this hunch or how
> exactly to proceed if that is the case. (As an example of the
> sparseness of the data, items are repeated across subjects, but each
> subject has only one data point per item, or zero data points per item
> for trials where they didn't respond correctly. However, I didn't get
> any warnings about model convergence, or any warnings at
> all.)
>
> Any clues as to why I'm getting this results would be
> very much appreciated.
>
> Thanks in advance,
>
> Alan Mishler
> Research Assistant
> University of Maryland
> --
> ## Model 1 output: Block 1 as reference level ##
>> glmerNL.RS.SI_RTgxb1
> Linear mixed model fit by REML
> Formula: WinRTs ~ (Group * Block * Type) + (1 + Block + Type |
> Subject) + (1 + Group + Block | Item)
> Data: nback.low
> AIC BIC logLik deviance REMLdev
> 181965 182210 -90949 181990 181899
>
> Random effects:
> Groups Name Variance Std.Dev. Corr
> Item (Intercept) 6134.64 78.324
> Group1 2336.31 48.335 -1.000
> Block1 1759.56 41.947 1.000 -1.000
> Block2 846.06 29.087 0.771 -0.771 0.771
> Subject (Intercept) 132462.22 363.954
> Block1 6011.83 77.536 -0.034
> Block2 10883.42 104.324 -0.127 0.915
> Type1 13653.97 116.850 -0.194 0.260 0.109
> Residual 97048.54 311.526
> Number of obs: 12640, groups: Item, 288; Subject, 52
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 1016.73 73.39 13.855
> Group1 12.39 101.67 0.122
> Block1 -117.98 18.97 -6.220
> Block2 -175.96 23.60 -7.455
> Type1 -136.36 24.96 -5.463
> Group1:Block1 46.43 26.05 1.782
> Group1:Block2 76.70 32.53 2.358
> Group1:Type1 20.55 34.15 0.602
> Block1:Type1 62.16 10.47 5.934
> Block2:Type1 96.66 10.46 9.243
> Group1:Block1:Type1 -18.39 14.05 -1.309
> Group1:Block2:Type1 -39.94 14.13 -2.826
>
> Correlation of Fixed Effects:
> (Intr) Group1 Block1 Block2 Type1 Gr1:B1 Gr1:B2 Gr1:T1
> Bl1:T1 Bl2:T1 G1:B1:
> Group1 -0.721
> Block1 -0.065 0.049
> Block2 -0.146 0.106 0.815
> Type1 -0.186 0.134 0.219 0.106
> Grop1:Blck1 0.053 -0.074 -0.715 -0.587 -0.160
> Grop1:Blck2 0.108 -0.150 -0.586 -0.721 -0.077 0.818
> Group1:Typ1 0.136 -0.189 -0.160 -0.078 -0.721 0.226 0.110
> Block1:Typ1 0.012 -0.009 -0.086 -0.038 -0.163 0.063 0.028 0.131
> Block2:Typ1 0.013 -0.009 -0.051 -0.070 -0.189 0.037 0.051 0.144
> 0.533
> Grp1:Bl1:T1 -0.009 0.013 0.064 0.028 0.153 -0.090 -0.043 -0.216
> -0.702 -0.371
> Grp1:Bl2:T1 -0.010 0.014 0.038 0.052 0.156 -0.055 -0.072 -0.218
> -0.371 -0.717 0.527
>
> ## Model 2 output: Block 2 as reference level ##
>> glmerNL.RS.SI_RTgxb2 [Block 2 as reference level]
> Linear mixed model fit by REML
> Formula: WinRTs ~ (Group * Block * Type) + (1 + Block + Type |
> Subject) + (1 + Group + Block | Item)
> Data: nback.low
> AIC BIC logLik deviance REMLdev
> 181931 182177 -90933 181957 181865
>
> Random effects:
> Groups Name Variance Std.Dev. Corr
> Item (Intercept) 14663.50 121.093
> Group1 2400.09 48.991 -1.000
> Block1 6561.73 81.004 -0.582 0.582
> Block2 590.58 24.302 -0.677 0.677 -0.205
> Subject (Intercept) 136638.85 369.647
> Block1 5868.05 76.603 -0.172
> Block2 2128.70 46.138 -0.146 -0.388
> Type1 13788.61 117.425 -0.140 -0.259 -0.188
> Residual 95743.04 309.424
> Number of obs: 12640, groups: Item, 288; Subject, 52
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 898.925 74.610 12.048
> Group1 58.779 103.098 0.570
> Block1 118.478 19.267 6.149
> Block2 -58.059 13.680 -4.244
> Type1 -74.323 25.553 -2.909
> Group1:Block1 -46.515 25.799 -1.803
> Group1:Block2 30.232 18.734 1.614
> Group1:Type1 2.166 34.141 0.063
> Block1:Type1 -64.291 11.246 -5.717
> Block2:Type1 33.892 10.045 3.374
> Group1:Block1:Type1 19.392 13.995 1.386
> Group1:Block2:Type1 -21.160 13.634 -1.552
>
> Correlation of Fixed Effects:
> (Intr) Group1 Block1 Block2 Type1 Gr1:B1 Gr1:B2 Gr1:T1
> Bl1:T1 Bl2:T1 G1:B1:
> Group1 -0.720
> Block1 -0.183 0.126
> Block2 -0.152 0.108 -0.033
> Type1 -0.132 0.096 -0.174 -0.095
> Grop1:Blck1 0.126 -0.176 -0.699 0.019 0.130
> Grop1:Blck2 0.106 -0.147 0.021 -0.722 0.070 -0.030
> Group1:Typ1 0.099 -0.137 0.130 0.071 -0.714 -0.188 -0.102
> Block1:Typ1 0.009 -0.007 -0.079 -0.049 -0.239 0.059 0.036 0.148
> Block2:Typ1 0.010 -0.008 -0.036 -0.107 -0.219 0.027 0.079 0.152
> 0.416
> Grp1:Bl1:T1 -0.008 0.010 0.063 0.040 0.137 -0.092 -0.050 -0.194
> -0.655 -0.345
> Grp1:Bl2:T1 -0.008 0.010 0.026 0.079 0.141 -0.036 -0.101 -0.200
> -0.315 -0.722 0.482
>
> ## Model 3 output: Block 3 as reference level ##
>> glmerNL.RS.SI_RTgxb3 [Block 3 as reference level]
> Linear mixed model fit by REML
> Formula: WinRTs ~ (Group * Block * Type) + (1 + Block + Type |
> Subject) + (1 + Group + Block | Item)
> Data: nback.low
> AIC BIC logLik deviance REMLdev
> 181932 182177 -90933 181958 181866
>
> Random effects:
> Groups Name Variance Std.Dev. Corr
> Item (Intercept) 11296.69 106.286
> Group1 2419.00 49.183 -1.000
> Block1 7418.92 86.133 -0.522 0.522
> Block2 494.56 22.239 0.596 -0.596 0.373
> Subject (Intercept) 133759.20 365.731
> Block1 10713.51 103.506 -0.153
> Block2 2123.58 46.082 0.021 0.733
> Type1 13764.21 117.321 -0.165 -0.107 0.189
> Residual 95762.39 309.455
> Number of obs: 12640, groups: Item, 288; Subject, 52
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 840.86 73.76 11.399
> Group1 89.07 102.02 0.873
> Block1 176.58 23.92 7.381
> Block2 58.07 13.66 4.251
> Type1 -40.40 25.30 -1.597
> Group1:Block1 -76.71 32.29 -2.376
> Group1:Block2 -30.30 18.72 -1.618
> Group1:Type1 -18.88 34.08 -0.554
> Block1:Type1 -98.00 11.46 -8.555
> Block2:Type1 -33.93 10.03 -3.382
> Group1:Block1:Type1 40.12 14.06 2.853
> Group1:Block2:Type1 21.06 13.63 1.545
>
> Correlation of Fixed Effects:
> (Intr) Group1 Block1 Block2 Type1 Gr1:B1 Gr1:B2 Gr1:T1
> Bl1:T1 Bl2:T1 G1:B1:
> Group1 -0.720
> Block1 -0.170 0.119
> Block2 -0.031 0.024 0.595
> Type1 -0.156 0.113 -0.073 0.138
> Grop1:Blck1 0.119 -0.164 -0.706 -0.434 0.054
> Grop1:Blck2 0.026 -0.034 -0.429 -0.723 -0.101 0.603
> Group1:Typ1 0.116 -0.161 0.054 -0.102 -0.717 -0.079 0.142
> Block1:Typ1 0.009 -0.006 -0.063 -0.045 -0.231 0.047 0.034 0.148
> Block2:Typ1 0.009 -0.007 -0.032 -0.107 -0.177 0.024 0.079 0.138
> 0.461
> Grp1:Bl1:T1 -0.007 0.009 0.052 0.037 0.141 -0.073 -0.048 -0.195
> -0.651 -0.357
> Grp1:Bl2:T1 -0.007 0.009 0.025 0.079 0.144 -0.030 -0.101 -0.200
> -0.324 -0.723 0.490
>
> ______________________________________________
> R-help at r-project.org mailing list
> 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.
--
Bert Gunter
Genentech Nonclinical Biostatistics
Internal Contact Info:
Phone: 467-7374
Website:
http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm
More information about the R-help
mailing list