[R-sig-ME] Repeated measures - general question
Ben Bolker
bbolker @ending from gm@il@com
Tue Oct 30 18:15:22 CET 2018
Don't know about matching SAS, but nlme::lme() with
correlation=corCAR1() or glmmTMB::glmmTMB() with an ou() correlation
structure should both work for accounting for autocorrelation with
unevenly spaced data. Maybe you could let us know how you're analyzing
the data in SAS?
On 2018-10-30 11:45 a.m., Nik Tuzov wrote:
> What is the best package to analyze RM data in R, given that time is on continuous scale and the time pointsare not equally spaced. Ideally, that package should match SAS results.
> Regards,Nik
>
> From: "r-sig-mixed-models-request using r-project.org" <r-sig-mixed-models-request using r-project.org>
> To: r-sig-mixed-models using r-project.org
> Sent: Monday, October 29, 2018 6:58 PM
> Subject: R-sig-mixed-models Digest, Vol 142, Issue 28
>
> Send R-sig-mixed-models mailing list submissions to
> r-sig-mixed-models using r-project.org
>
> To subscribe or unsubscribe via the World Wide Web, visit
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> or, via email, send a message with subject or body 'help' to
> r-sig-mixed-models-request using r-project.org
>
> You can reach the person managing the list at
> r-sig-mixed-models-owner using r-project.org
>
> When replying, please edit your Subject line so it is more specific
> than "Re: Contents of R-sig-mixed-models digest..."
>
>
> Today's Topics:
>
> 1. Re: Simulations in Fisher's Exact test (Alexandre Courtiol)
> 2. Re: glmmTMB (Ben Bolker)
> 3. Re: Random intercept model- unbalanced cluster (Yashree Mehta)
> 4. Re: Random intercept model- unbalanced cluster (Ben Bolker)
> 5. zero-inflated glmmTMB with poly() - confidence band (John Wilson)
>
> ----------------------------------------------------------------------
>
> Message: 1
> Date: Mon, 29 Oct 2018 12:11:24 +0100
> From: Alexandre Courtiol <alexandre.courtiol using gmail.com>
> To: r-sig-mixed-models using r-project.org
> Subject: Re: [R-sig-ME] Simulations in Fisher's Exact test
> Message-ID:
> <CAERMt4dwmh1rLe_gdLBw8_K-cutyFVFphwHBL85tQCNtHswH1w using mail.gmail.com>
> Content-Type: text/plain; charset="utf-8"
>
>>
>> Message: 2
>> Date: Sun, 28 Oct 2018 12:18:44 -0400
>> From: Ben Bolker <bbolker using gmail.com>
>> To: r-sig-mixed-models using r-project.org
>> Subject: Re: [R-sig-ME] Simulations in Fisher's Exact test
>> Message-ID: <b9fd6b4b-e9e6-d368-4148-d9762544e169 using gmail.com>
>> Content-Type: text/plain; charset="utf-8"
>>
>>
>> This is indeed the wrong list; r-help using r-project.org, or StackOverflow,
>> might be more appropriate. I am guessing this is an assignment for a
>> class? If so, it might be more useful to get clarification from your
>> instructor or teaching assistant (or a colleague in your class). The
>> help page for ?fisher.test says:
>>
>> simulate.p.value: a logical indicating whether to compute p-values by
>> Monte Carlo simulation, **in larger than 2 by 2 tables**.
>>
>> Emphasis (**) added. Since you're using a 2x2 table, I'm guessing that
>> simulate.p.value has no effect ... R probably should warn you, but oh
>> well ..
>>
>>
>> On 2018-10-28 7:47 a.m., Adeela Munawar wrote:
>>> hi all,
>>> Probably I am posted in wrong mailing list. I am getting a problem in
>>> applying Fisher's exact test. I am applying Fisher's exact test as
>>>
>>> ntable<- array(data = c(3, 1, 8,11), dim = c(2,2))
>>> fisher.test(ntable)
>>>
>>> now, I have to repeat the same 10000 times and have to report p-values.
>>> Using the arguments simulate.p.value in the command is producing the same
>>> results
>>> test<-fisher.test(ntable,workspace=20000,simulate.p.value=T,B=10000)
>>>
>>> what changes I have to made in my model.
>>>
>>> regards
>>> Adeela
>>>
>>> [[alternative HTML version deleted]]
>>>
>>> _______________________________________________
>>> R-sig-mixed-models using r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>
>> Adeela, with simulate.p.value = FALSE, the Fisher exact test will attempt
> to create all possible contingency tables (keeping marginal sum constant)
> to compute the p-value of the test. With simulate.p.value = TRUE (when it
> has an effect -> see Ben's comment), it will only sample the space of
> possible contingency tables. If the number of possible contingency tables
> is not too large, there is not need to use simulate.p.value and if it is
> lower than the number of table you simulate, then you should obtain nearly
> the same results anyhow. In other words, I don't really understand what you
> are trying to achieve but a simple call to fisher.test(ntable) should do
> the job.
> ++
> Alex
>
> [[alternative HTML version deleted]]
>
>
>
>
> ------------------------------
>
> Message: 2
> Date: Mon, 29 Oct 2018 08:48:05 -0400
> From: Ben Bolker <bbolker using gmail.com>
> To: mairafatoretto using gmail.com
> Cc: R SIG Mixed Models <r-sig-mixed-models using r-project.org>
> Subject: Re: [R-sig-ME] glmmTMB
> Message-ID:
> <CABghstQFdj44Awgfhj66ZTSEAarx1xTGOVZk=XzhJ5B0vd_uBg using mail.gmail.com>
> Content-Type: text/plain; charset="utf-8"
>
> glmmTMB implements REML by treating the fixed effect parameters as
> 'random variables', i.e. turning on the Laplace approximation
> machinery; see http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#reml-for-glmms
> for starting points in the literature.
>
> https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/R/glmmTMB.R#L191
> On Mon, Oct 29, 2018 at 6:26 AM Maira Fatoretto
> <mairafatoretto using gmail.com> wrote:
>>
>> Hello,
>>
>> I am using glmmTMb with binomial family. I chose REML=TRUE because I have
>> three random effects and I would like to know which is the restriction used
>> in this estimation in the package.
>> The REML = TRUE have better estimation than ML in my case.
>>
>>
>>
>> Thank you.
>> Best regards.
>>
>>
>>
>> --
>> *Maíra Blumer Fatoretto*
>> Estatística - Universidade Estadual de Campinas- UNICAMP
>> Mestra em Ciências(Estatística e Experimentação Agronômica) - ESALQ/USP -
>> Piracicaba - SP
>> Doutoranda em Estatística e Experimentação Agronômica
>> Telefone:(19) 988309481
>> E-mail: mairafatoretto using gmail.com
>>
>> [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-mixed-models using r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
>
>
> ------------------------------
>
> Message: 3
> Date: Mon, 29 Oct 2018 22:13:06 +0100
> From: Yashree Mehta <yashree19 using gmail.com>
> To: r-sig-mixed-models using r-project.org
> Subject: Re: [R-sig-ME] Random intercept model- unbalanced cluster
> Message-ID:
> <CAOE=hq+3u9pAr5Xp-BMWeQ4JNE-ZVpm6YYS=wEmUXhGVihvG8w using mail.gmail.com>
> Content-Type: text/plain; charset="utf-8"
>
> Or is there an alternative method of modeling this subset of households who
> only own one plot?
>
> thank you,
>
> Regards,
> Yashree
>
> On Thu, Oct 25, 2018 at 6:00 PM Yashree Mehta <yashree19 using gmail.com> wrote:
>
>> Hi,
>>
>> I am working with a random intercept model on a cluster dataset (Repeated
>> measurements of plots per household). I have the usual "X" vector
>> of covariates and one id variable which will make up the random
>> intercept. For example,
>>
>> Response variable: Production of maize
>> Covariates: Size, input quantities, soil fertility dummies etc..
>> ID variable: Household_ID
>>
>> However, about 40% of the households own one plot. The number of plots per
>> household ranges from 1 to 13.
>>
>> When I estimated the random intercept model using lmer, I can extract a
>> random intercept for all households, irrespective of their number of plots.
>>
>> How does lmer treat these households with just 1 plot? Also, is it
>> theoretically correct to include these observations ?
>>
>> Thank you,
>>
>> Regards,
>> Yashree
>>
>>
>>
>
> [[alternative HTML version deleted]]
>
>
>
>
> ------------------------------
>
> Message: 4
> Date: Mon, 29 Oct 2018 17:23:16 -0400
> From: Ben Bolker <bbolker using gmail.com>
> To: r-sig-mixed-models using r-project.org
> Subject: Re: [R-sig-ME] Random intercept model- unbalanced cluster
> Message-ID: <1a88c90a-8869-3b01-072d-57513f5293c4 using gmail.com>
> Content-Type: text/plain; charset="utf-8"
>
>
> In principle lme4 shouldn't have problems with a subset of groups that
> have only one observation (although clearly the model will get more
> fragile/unreliable the less information is available about within vs
> among group variation ...). I'd expect the random effects for groups
> with only one observation to be strongly shrunk toward the population
> mean ... if in doubt, it can be very useful to simulate a situation
> similar to your real data set to see what happens in cases where you
> know the real answer ...
>
>
>
> On 2018-10-29 5:13 p.m., Yashree Mehta wrote:
>> Or is there an alternative method of modeling this subset of households who
>> only own one plot?
>>
>> thank you,
>>
>> Regards,
>> Yashree
>>
>> On Thu, Oct 25, 2018 at 6:00 PM Yashree Mehta <yashree19 using gmail.com> wrote:
>>
>>> Hi,
>>>
>>> I am working with a random intercept model on a cluster dataset (Repeated
>>> measurements of plots per household). I have the usual "X" vector
>>> of covariates and one id variable which will make up the random
>>> intercept. For example,
>>>
>>> Response variable: Production of maize
>>> Covariates: Size, input quantities, soil fertility dummies etc..
>>> ID variable: Household_ID
>>>
>>> However, about 40% of the households own one plot. The number of plots per
>>> household ranges from 1 to 13.
>>>
>>> When I estimated the random intercept model using lmer, I can extract a
>>> random intercept for all households, irrespective of their number of plots.
>>>
>>> How does lmer treat these households with just 1 plot? Also, is it
>>> theoretically correct to include these observations ?
>>>
>>> Thank you,
>>>
>>> Regards,
>>> Yashree
>>>
>>>
>>>
>>
>> [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-mixed-models using r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>
>
>
>
> ------------------------------
>
> Message: 5
> Date: Mon, 29 Oct 2018 20:58:09 -0300
> From: John Wilson <jhwilson.nb using gmail.com>
> To: r-sig-mixed-models using r-project.org
> Subject: [R-sig-ME] zero-inflated glmmTMB with poly() - confidence
> band
> Message-ID:
> <CABdA5Q0o2nJTRFo+mj6Hcz_UbZ+6DM0Rxfzmyk2aOwzovjc92g using mail.gmail.com>
> Content-Type: text/plain; charset="utf-8"
>
> Hello,
>
> I've been using the newly documented predict() with group = NA to predict
> population-level values, as per the thread here (
> https://stat.ethz.ch/pipermail/r-sig-mixed-models/2018q4/027305.html). A
> follow-up question: in the case of a zero-inflated model, how would I go
> about to get the 95% CIs for response predictions? Obviously, I can get
> them for the count and the zero-component, and without poly(), I would just
> follow the example in Brooks et al (2017).
>
> However, I have a poly() predictor; how do I get the 95% CIs if I can't use
> model.matrix naively when I have a poly() in the model? The models take a
> while to converge, so I don't want to run a full bootstrapping either, if
> at all possible.
>
> Thank you!
> John
>
> Here's a toy dataset:
>
> library(ggplot2)
> library(glmmTMB)
>
> set.seed(0)
> x <- 1:20
> z <- sample(c("a", "b"), length(x), replace = TRUE)
> y <- round(5 * 2*x + 3 * x^2 + 0.1 * x^3 + rnbinom(length(x), 10, 0.03))
> y[c(2, 3, 5, 11, 13, 19)] <- 0
> group <- sample(c("i", "ii"), length(x), replace = TRUE)
> df <- data.frame(x = x, y = y, z = z, group = group)
>
> ggplot(df) +
> geom_point(aes(x = x, y = y, colour = z))
>
> m <- glmmTMB(y ~ poly(x, 3) * z +
> (1 | group),
> zi = ~ z,
> family = nbinom1,
> data = df)
> # prediction on a new grid
> newdata <- expand.grid(x = 1:20, z = unique(df$z), group = NA)
> newdata$Pred <- predict(m, type = "response")
> ### now how to add CIs?
>
> ggplot(df) +
> geom_point(aes(x = x, y = y, colour = z)) +
> geom_line(data = newdata, aes(x = x, y = Pred, colour = z, group = z), size
> = 2) +
> facet_wrap(~ z)
>
> [[alternative HTML version deleted]]
>
>
>
>
> ------------------------------
>
> Subject: Digest Footer
>
> _______________________________________________
> R-sig-mixed-models mailing list
> R-sig-mixed-models using r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
> ------------------------------
>
> End of R-sig-mixed-models Digest, Vol 142, Issue 28
> ***************************************************
>
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
More information about the R-sig-mixed-models
mailing list