[R] mixed models

James Henson jfhenson1 at gmail.com
Fri May 27 23:25:48 CEST 2016 

Greetings David,
I am new to R and neglected to check vigorously for missing values.
Apologize for posting without checking and finding the one NA.

I appreciate your help.
Thanks.
James F. Henson

On Fri, May 27, 2016 at 1:49 PM, David Winsemius <dwinsemius at comcast.net> wrote:
>
>> On May 27, 2016, at 10:07 AM, James Henson <jfhenson1 at gmail.com> wrote:
>>
>> Greetings Jeff,
>> You are correct that the unequal number of levels is not the problem.
>> I revised the data frame so that the number of levels was equal and
>> the same error message occurred.  The code is below, and the
>> Eboni2.txt file is attached. This problem baffles me.  I appreciate
>> any help.
>> Best regards,
>> James
>> Eboni2 <- read.csv("Eboni2.csv", header = TRUE)
>>
>> library("nlme")
>>
>> str(Eboni2)
>>
>> head(Eboni2)
>>
>> model1 <- lme(preDawn ~ Irrigation, random=~season_order|treeNo, data=Eboni2)
>
> I downloaded the attached file to your first posting that was called a "csv" file but it was tab-separated (as could be clearly seen with the str output, so would only load properly with read.delim rather than read.csv. Running then with the lme call, it produced this message
>
>> model1 <- lme(preDawn ~ Irrigation, random=~season_order|treeNo, data=Eboni2)
> Error in na.fail.default(list(season_order = c(5L, 5L, 5L, 5L, 5L, 5L,  :
>   missing values in object
>
> And looking at the str result made it clear that there were many NA's in the file.
>
>> head(Eboni2)
>   number Location   Season season_order Month  treeID treeNo preDawn midday
> 1      1      UCC November            5   Nov UCCLO 1     60     1.4    1.3
> 2      2      UCC November            5   Nov UCCLO 2     72     1.2    1.3
> 3      3      UCC November            5   Nov UCCLO 3     78     1.1    1.2
> 4      4      UCC November            5   Nov UCCLO 4     79     1.1    2.1
> 5      5      UCC November            5   Nov UCCLO 5     80     1.4    1.3
> 6      6      UCC November            5   Nov UCCLO 6     81     0.6    1.8
>   Irrigation Pnet        Gs        E      WUE d15N   d13C Nper  Cper include2
> 1          N    9 0.2907004 3.766207 2.389672   NA     NA   NA    NA       no
> 2          N   11 0.3262582 3.120574 3.524993   NA     NA   NA    NA       no
> 3          N    8 0.2870957 1.693821 4.723050 3.00 -27.44 2.12 52.12      yes
> 4          N   10 0.2475180 1.839343 5.436724 3.61 -29.50 1.42 51.97      yes
> 5          N   13 0.3009228 3.082278 4.217660   NA     NA   NA    NA       no
> 6          N   17 0.3487337 2.534550 6.707304 2.79 -30.50 1.49 49.94      yes
>
> And even more importantly, there was one NA in your outcome variable:
>> sum( is.na(Eboni2$Irrigation))
> [1] 0
>> sum( is.na(Eboni2$preDawn))
> [1] 1
>
> So after restricting to complete.cases, I then formed the hypothesis that you reversed the order of the variables in the formula for the random parameter:
>
>> table(Eboni2$season_order)
>
>  1  2  3  4  5
> 83 83 83 83 83
>> length( Eboni2$treeNo)
> [1] 415
>
> So it seemed unreasonable to have a "grouping" on variable with only one item per group.
>
>> model1 <- lme(preDawn ~ Irrigation, random=~treeNo|season_order, data=Eboni2[ complete.cases( Eboni2[ , c('preDawn','Irrigation','season_order','treeNo')]), ] )
>> model1
> Linear mixed-effects model fit by REML
>   Data: Eboni2[complete.cases(Eboni2[, c("preDawn", "Irrigation", "season_order",      "treeNo")]), ]
>   Log-restricted-likelihood: -183.4708
>   Fixed: preDawn ~ Irrigation
> (Intercept) IrrigationY
>  1.04520145 -0.06037706
>
> Random effects:
>  Formula: ~treeNo | season_order
>  Structure: General positive-definite, Log-Cholesky parametrization
>             StdDev      Corr
> (Intercept) 0.140239324 (Intr)
> treeNo      0.003766019 -0.725
> Residual    0.365678898
>
> Number of Observations: 414
> Number of Groups: 5
>
> (Warning, I'm not a frequent user of this package or any of the mixed effects packages.)
>
>
>
> Just to correct some misinformation that appeared earlier: You can attach "csv" or "tsv" files as long as you name them with an .txt extension so the mail clients and servers consider them to be MIME-text.
>
> --
> David.
>
>
>> On Wed, May 25, 2016 at 6:23 PM, Jeff Newmiller
>> <jdnewmil at dcn.davis.ca.us> wrote:
>>> Please keep the mailing list in the loop by using reply-all.
>>>
>>> I don't think there is a requirement that the number of levels is equal, but
>>> there may be problems if you don't have the minimum number of records
>>> corresponding to each combination of levels specified in your model.
>>>
>>> You can change the csv extension to txt and attach for the mailing list. Or,
>>> better yet, you can use the dput function to embed the data directly in your
>>> sample code.
>>>
>>> Also, please learn to post plain text email to avoid corruption of R code by
>>> the HTML formatting.
>>> --
>>> Sent from my phone. Please excuse my brevity.
>>>
>>> On May 25, 2016 2:26:54 PM PDT, James Henson <jfhenson1 at gmail.com> wrote:
>>>>
>>>> Good afternoon Jeff,
>>>> The sample sizes for levels of the factor "Irrigation" are not equal. If
>>>> 'nlme' requires equal sample sizes this may be the problem. The same data
>>>> frame runs in 'lme4' without a problem.
>>>>
>>>> Best regards,
>>>> James
>>>>
>>>>
>>>> On Wed, May 25, 2016 at 3:41 PM, James Henson <jfhenson1 at gmail.com> wrote:
>>>>>
>>>>> Good afternoon Jeff,
>>>>>
>>>>> When working with this data frame, I just open the .csv file in R Studio.
>>>>> But, we should not send .csv file to R_help.  What should I send?
>>>>>
>>>>> Best regards,
>>>>> James
>>>>>
>>>>> On Wed, May 25, 2016 at 2:52 PM, Jeff Newmiller
>>>>> <jdnewmil at dcn.davis.ca.us> wrote:
>>>>>>
>>>>>> You forgot to show the commands to us that you used to read the data in
>>>>>> with (your example is not "reproducible"). This step can make all the
>>>>>> difference in the world as to whether your analysis commands will work or
>>>>>> not.
>>>>>> --
>>>>>> Sent from my phone. Please excuse my brevity.
>>>>>>
>>>>>> On May 25, 2016 11:59:06 AM PDT, James Henson <jfhenson1 at gmail.com>
>>>>>> wrote:
>>>>>>>
>>>>>>> Greetings R community,
>>>>>>>
>>>>>>> My aim is to analyze a mixed-effects model with temporal
>>>>>>> pseudo-replication
>>>>>>> (repeated measures on the same experimental unit) using ‘nlme’.
>>>>>>> However,
>>>>>>> my code returns the error message “Error in na.fail.default’, even
>>>>>>> though
>>>>>>> the data frame does not contain missing values. My code is below, and
>>>>>>> the
>>>>>>> data file is attached as ‘Eboni2.txt.
>>>>>>>
>>>>>>> library("nlme")
>>>>>>>
>>>>>>> str(Eboni2)
>>>>>>>
>>>>>>> head(Eboni2)
>>>>>>>
>>>>>>> model1 <- lme(preDawn ~ Irrigation, random=~season_order|treeNo,
>>>>>>> data=Eboni2)
>>>>>>>
>>>>>>> I am genuinely confused.  Hope someone can help.
>>>>>>>
>>>>>>> Best regards,
>>>>>>>
>>>>>>> James F. Henson
>>>>>>>
>>>>>>> ________________________________
>>>>>>>
>>>>>>> R-help at 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.
>>>>>
>>>>>
>>>>
>>>
>>
>> ______________________________________________
>> R-help at 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.
>
> David Winsemius
> Alameda, CA, USA
>

More information about the R-help mailing list