[R-sig-ME] Unidentifiable model in lmer

Thierry Onkelinx thierry.onkelinx at inbo.be
Fri Sep 18 17:15:56 CEST 2015


Yes, this works.

library(lme4)
result <- lmer(y ~ A * B + (1|ID), data = dataset)
summary(result) # works fine
library(lmerTest)
summary(result) # works fine
result <- lmer(y ~ A * B + (1|ID), data = dataset)
summary(result) # yields "Model is not identifiable..."

This might be due to the fact that ID is nested in A. I recommend that
you contact the authors of lmerTest about this problem.

Best regards,
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature
and Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium

To call in the statistician after the experiment is done may be no
more than asking him to perform a post-mortem examination: he may be
able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does
not ensure that a reasonable answer can be extracted from a given body
of data. ~ John Tukey


2015-09-18 17:08 GMT+02:00 Takahiro Fushimi <taka.6765 op gmail.com>:
> Sorry about it. Does this work?
>
> dataset <- data.frame(y = c(60625.19, 51088.38, 67283.03, 56550.54, NA,
> NA, 58720.44,63939.01, 66298.32, 62662.15, 53604.41, 88338.55, 62818.78,
> 68696.71,57433.98, 59105.99, 44801.97, 68655.48, 62645.16, 69409.15,
> 43568.70,55693.21, 57322.63, 51095.55, 57660.41, 64128.63, 55216.42,
> 59945.17,67284.31, 55401.13, 60471.59, NA, 63633.07, 65939.14,NA, 69488.38,
> 43950.50, 52782.99, 55674.94, 70130.25, 72560.25, 69297.95, 53188.98,
> 75687.53, 68242.23, 60696.15, 65087.04, 59377.71, 55055.45, 66673.40,
> 65933.37, 62636.56, 49606.80, 58668.13, 56694.13, NA, NA, 61928.95,63208.84,
> 83786.71, 72727.53, 78873.45, 47482.32, 65967.84, 53228.25, 65699.59,
> 61792.81, 44816.72, 71048.06, 59671.84, 89513.43, 56972.30, 45439.25,
> 51357.83, NA, 60146.86, 56845.88, 63181.13, 63704.80, 49140.22, 54538.14,
> 49411.32, 66059.73, 73147.55, 55665.96, 66821.27, 66935.12, 54350.63,
> 48116.25, 67664.66, 64278.54, 64555.03, 62463.60, 55831.53, 57392.61,
> 75727.00, 64839.50, 51009.78, 65274.59, 63339.91, 62276.49, 66159.42,
> 58333.76, 60275.09) , A = as.factor(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
> 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2,
> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2,
> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2)) , B = as.factor(c(1, 1, 1,
> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 ,1 ,1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0,
> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
> 0)), ID = as.factor(c(72,  73,  74,  75,  76, 78,  79,  86,  88,  89,  90,
> 91,  92,  93,  94,  95,  96,  121, 122, 123, 124, 77,  80,  81, 82,  83,
> 84,  98,  116, 117, 118, 119, 120, 125, 128, 100, 101, 102, 103, 104, 105,
> 106, 107, 108, 109, 110, 112, 113, 114, 115, 126, 127, 72, 73, 74, 75, 76,
> 78, 79, 86, 88, 89, 90, 91, 92, 93, 94, 95, 96, 121, 122, 123, 124, 77, 80,
> 81, 82, 83, 84, 98, 116, 117, 118, 119, 120, 125, 128, 100, 101, 102, 103,
> 104, 105, 106, 107, 108, 109, 110, 112, 113, 114, 115, 126, 127)))
>
> Best regards,
> Takahiro
>
>
> On 9/18/15 10:41 AM, Thierry Onkelinx wrote:
>>
>> Dear Takahiro,
>>
>> Please send the output of dput(dataset). That's a lot easier to
>> copy-paste into an R session.
>>
>> Best regards,
>> ir. Thierry Onkelinx
>> Instituut voor natuur- en bosonderzoek / Research Institute for Nature
>> and Forest
>> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
>> Kliniekstraat 25
>> 1070 Anderlecht
>> Belgium
>>
>> To call in the statistician after the experiment is done may be no
>> more than asking him to perform a post-mortem examination: he may be
>> able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
>> The plural of anecdote is not data. ~ Roger Brinner
>> The combination of some data and an aching desire for an answer does
>> not ensure that a reasonable answer can be extracted from a given body
>> of data. ~ John Tukey
>>
>>
>> 2015-09-18 16:40 GMT+02:00 Takahiro Fushimi <taka.6765 op gmail.com>:
>>>
>>> The data set is as follows:
>>>>
>>>> dataset
>>>
>>>             y A B  ID
>>> 1   60625.19 0 1  72
>>> 2   51088.38 0 1  73
>>> 3   67283.03 0 1  74
>>> 4   56550.54 0 1  75
>>> 5         NA 0 1  76
>>> 6         NA 0 1  78
>>> 7   58720.44 0 1  79
>>> 8   63939.01 0 1  86
>>> 9   66298.32 0 1  88
>>> 10  62662.15 0 1  89
>>> 11  53604.41 0 1  90
>>> 12  88338.55 0 1  91
>>> 13  62818.78 0 1  92
>>> 14  68696.71 0 1  93
>>> 15  57433.98 0 1  94
>>> 16  59105.99 0 1  95
>>> 17  44801.97 0 1  96
>>> 18  68655.48 0 1 121
>>> 19  62645.16 0 1 122
>>> 20  69409.15 0 1 123
>>> 21  43568.70 0 1 124
>>> 22  55693.21 1 1  77
>>> 23  57322.63 1 1  80
>>> 24  51095.55 1 1  81
>>> 25  57660.41 1 1  82
>>> 26  64128.63 1 1  83
>>> 27  55216.42 1 1  84
>>> 28  59945.17 1 1  98
>>> 29  67284.31 1 1 116
>>> 30  55401.13 1 1 117
>>> 31  60471.59 1 1 118
>>> 32        NA 1 1 119
>>> 33  63633.07 1 1 120
>>> 34  65939.14 1 1 125
>>> 35        NA 1 1 128
>>> 36  69488.38 2 1 100
>>> 37  43950.50 2 1 101
>>> 38  52782.99 2 1 102
>>> 39  55674.94 2 1 103
>>> 40  70130.25 2 1 104
>>> 41  72560.25 2 1 105
>>> 42  69297.95 2 1 106
>>> 43  53188.98 2 1 107
>>> 44  75687.53 2 1 108
>>> 45  68242.23 2 1 109
>>> 46  60696.15 2 1 110
>>> 47  65087.04 2 1 112
>>> 48  59377.71 2 1 113
>>> 49  55055.45 2 1 114
>>> 50  66673.40 2 1 115
>>> 51  65933.37 2 1 126
>>> 52  62636.56 2 1 127
>>> 53  49606.80 0 0  72
>>> 54  58668.13 0 0  73
>>> 55  56694.13 0 0  74
>>> 56        NA 0 0  75
>>> 57        NA 0 0  76
>>> 58  61928.95 0 0  78
>>> 59  63208.84 0 0  79
>>> 60  83786.71 0 0  86
>>> 61  72727.53 0 0  88
>>> 62  78873.45 0 0  89
>>> 63  47482.32 0 0  90
>>> 64  65967.84 0 0  91
>>> 65  53228.25 0 0  92
>>> 66  65699.59 0 0  93
>>> 67  61792.81 0 0  94
>>> 68  44816.72 0 0  95
>>> 69  71048.06 0 0  96
>>> 70  59671.84 0 0 121
>>> 71  89513.43 0 0 122
>>> 72  56972.30 0 0 123
>>> 73  45439.25 0 0 124
>>> 74  51357.83 1 0  77
>>> 75        NA 1 0  80
>>> 76  60146.86 1 0  81
>>> 77  56845.88 1 0  82
>>> 78  63181.13 1 0  83
>>> 79  63704.80 1 0  84
>>> 80  49140.22 1 0  98
>>> 81  54538.14 1 0 116
>>> 82  49411.32 1 0 117
>>> 83  66059.73 1 0 118
>>> 84  73147.55 1 0 119
>>> 85  55665.96 1 0 120
>>> 86  66821.27 1 0 125
>>> 87  66935.12 1 0 128
>>> 88  54350.63 2 0 100
>>> 89  48116.25 2 0 101
>>> 90  67664.66 2 0 102
>>> 91  64278.54 2 0 103
>>> 92  64555.03 2 0 104
>>> 93  62463.60 2 0 105
>>> 94  55831.53 2 0 106
>>> 95  57392.61 2 0 107
>>> 96  75727.00 2 0 108
>>> 97  64839.50 2 0 109
>>> 98  51009.78 2 0 110
>>> 99  65274.59 2 0 112
>>> 100 63339.91 2 0 113
>>> 101 62276.49 2 0 114
>>> 102 66159.42 2 0 115
>>> 103 58333.76 2 0 126
>>> 104 60275.09 2 0 127
>>>
>>> Best regards,
>>> Takahiro
>>>
>>>
>>> On 9/18/15 3:23 AM, Thierry Onkelinx wrote:
>>>>
>>>> The formula shouldn't be the problem. model.matrix() is clever enough
>>>> to handle this. Have a look at the examples below.
>>>>
>>>> model.matrix(~ A + B + A * B, data= data.frame(A = 1, B = 1))
>>>> model.matrix(~ A + B + A * B + B + A:B + B:A, data= data.frame(A = 1, B
>>>> =
>>>> 1))
>>>>
>>>> I suggest to give use a reproducible example of the problem so that we
>>>> can invesigate what when wrong.
>>>> ir. Thierry Onkelinx
>>>> Instituut voor natuur- en bosonderzoek / Research Institute for Nature
>>>> and Forest
>>>> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
>>>> Kliniekstraat 25
>>>> 1070 Anderlecht
>>>> Belgium
>>>>
>>>> To call in the statistician after the experiment is done may be no
>>>> more than asking him to perform a post-mortem examination: he may be
>>>> able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
>>>> The plural of anecdote is not data. ~ Roger Brinner
>>>> The combination of some data and an aching desire for an answer does
>>>> not ensure that a reasonable answer can be extracted from a given body
>>>> of data. ~ John Tukey
>>>>
>>>>
>>>> 2015-09-18 7:23 GMT+02:00 Ché Lucero <chelucero op uchicago.edu>:
>>>>>
>>>>> I don't know if this is causing the problem you're seeing, but A*B
>>>>> expands
>>>>> to A + B + A:B, so your model right now is R ~ A + B + A + B + A:B +
>>>>> (1|S).
>>>>>
>>>>> On Thu, Sep 17, 2015 at 5:53 PM Takahiro Fushimi <taka.6765 op gmail.com>
>>>>> wrote:
>>>>>
>>>>>> Hi everyone,
>>>>>>
>>>>>> I have been working on a linear mixed effect model by using lmer()
>>>>>> function, but I got an error saying that the fitted model is not
>>>>>> identifiable.
>>>>>>
>>>>>> The data set includes the following variables:
>>>>>> y = a numeric variable
>>>>>> factorA = a 3-level categorical variable
>>>>>> factorB = a 2-level categorical variable
>>>>>> subjectID = subject id number. 2 measurements of y for each subject
>>>>>>
>>>>>> R code and output are as follows:
>>>>>>    > (result <- lmer(y ~ factorA + factorB + factorA*factorB +
>>>>>> (1|subjectID)))
>>>>>> Linear mixed model fit by REML ['merModLmerTest']
>>>>>> Formula: y ~ factorA + factorB + factorA * factorB + (1 | subjectID)
>>>>>> REML criterion at convergence: 1928.966
>>>>>> Random effects:
>>>>>>     Groups    Name        Std.Dev.
>>>>>>     subjectID (Intercept) 4711
>>>>>>     Residual              7688
>>>>>> Number of obs: 97, groups:  subjectID, 51
>>>>>> Fixed Effects:
>>>>>>          (Intercept)           factorA1           factorA2 factorB1
>>>>>> factorA1:factorB1  factorA2:factorB1
>>>>>>                62411              -2700              -1124
>>>>>> -1037               1279               2482
>>>>>>    > summary(result)
>>>>>> Model is not identifiable...
>>>>>> summary from lme4 is returned
>>>>>> some computational error has occurred in lmerTest
>>>>>> Linear mixed model fit by REML ['lmerMod']
>>>>>> Formula: y ~ factorA + factorB + factorA * factorB + (1 | subjectID)
>>>>>>
>>>>>> REML criterion at convergence: 1929
>>>>>>
>>>>>> Scaled residuals:
>>>>>>         Min       1Q   Median       3Q      Max
>>>>>> -1.93423 -0.62611  0.01837  0.48887  2.73380
>>>>>>
>>>>>> Random effects:
>>>>>>     Groups    Name        Variance Std.Dev.
>>>>>>     subjectID (Intercept) 22194074 4711
>>>>>>     Residual              59108495 7688
>>>>>> Number of obs: 97, groups:  subjectID, 51
>>>>>>
>>>>>> Fixed effects:
>>>>>>                      Estimate Std. Error t value
>>>>>> (Intercept)          62411       2065  30.227
>>>>>> factorA1             -2700       3238  -0.834
>>>>>> factorA2             -1124       3008  -0.374
>>>>>> factorB1             -1037       2512  -0.413
>>>>>> factorA1:factorB1     1279       4013   0.319
>>>>>> factorA2:factorB1     2482       3642   0.681
>>>>>>
>>>>>> Correlation of Fixed Effects:
>>>>>>                (Intr) fctrA1 fctrA2 fctrB1 fA1:B1
>>>>>> factorA1    -0.638
>>>>>> factorA2    -0.687  0.438
>>>>>> factorB1    -0.608  0.388  0.418
>>>>>> fctrA1:fcB1  0.381 -0.601 -0.262 -0.626
>>>>>> fctrA2:fcB1  0.420 -0.268 -0.606 -0.690  0.432
>>>>>>    >
>>>>>>
>>>>>> Could anyone give me some idea of why this unidentifiability problem
>>>>>> happens and how to fix it?
>>>>>> Any help would be appreciated.
>>>>>>
>>>>>> Best regards,
>>>>>> Takahiro
>>>>>>
>>>>>> _______________________________________________
>>>>>> R-sig-mixed-models op r-project.org mailing list
>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>>
>>>>>           [[alternative HTML version deleted]]
>>>>>
>>>>> _______________________________________________
>>>>> R-sig-mixed-models op r-project.org mailing list
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>>
>>> --
>>> Takahiro Fushimi
>>> Columbia University in the City of New York
>>> Master of Arts in Statistics
>>>
>
> --
> Takahiro Fushimi
> Columbia University in the City of New York
> Master of Arts in Statistics
>



More information about the R-sig-mixed-models mailing list