Hi all.
I seem to have managed (not quite sure how! :-) ) with considerable
assistance from the maintainer of the JuliaCall package, to get that
package running.
I am however having a struggle to figure out how to fit the sort of
model that I am interested in. The vignette that Doug Bates provided a
while back gives an example of fitting a linear mixed model, and I have
managed to reproduce that example. (My results differed from those that
Doug's vignette showed, ever so slightly in the values of some of the
variance components, but I'm not going to worry my pretty little head
about that.)
Now I want to proceed to my "real problem" which is to fit a
*generalised* linear mixed model, of the binomial persuasion, and I
cannot find my way to documentation on how to do this.
The syntax for the sort of fit that I want would be, in the glmer()
context, of the form:
fit <- glmer(cbind(Dead,Alive) ~ (Trt+0)/prd + (prd | Rep),data=Dat,
family=binomial(link=logit))
where "Trt" is a factor (fixed effect), "prd" is a numeric predictor,
and "Rep" is a factor (random effect).
Can anyone instruct me (in very simple terms please, I'm slow!) as to
what the analogous syntax would be using the MixedModels package in
Julia? Or point me at an elementary tutorial on doing this?
I have the impression, from one item that I saw on StackExchange, that
MixedModels doesn't really handle general binomial models, only
Bernoulli models. Consequently one has to expand out one's data set
creating one row for each success ("Dead" in my case) and one row for
each failure ("Alive"). Is this correct, or is there an easier way to
go about it?
Thanks for any words of wisdom.
cheers,
Rolf Turner
Technical Editor ANZJS
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276
Binomial generalised linear mixed model using JuliaCall.
3 messages · Rolf Turner, Ben Bolker
Don't have time to dig into this right now (I'm not an expert), but the docs at https://github.com/dmbates/MixedModels.jl/blob/master/docs/jmd/constructors.jmd#L104 say:
Note that, in keeping with convention in the [`GLM` package](https://github.com/JuliaStats/GLM.jl), the distribution family for a binary (i.e. 0/1) response is the `Bernoulli` distribution.
The `Binomial` distribution is only used when the response is the fraction of trials returning a positive, in which case the number of trials must be specified as the case weights. ... which implies that you can do what you want. https://github.com/dmbates/MixedModels.jl/blob/fd22daf226938a4f98ce3d4b228e5a99437e3218/src/pls.jl#L61 seem to give 'weights' as the third argument to the LinearMixedModel function ... ?
On Tue, Apr 17, 2018 at 8:16 PM, Rolf Turner <r.turner at auckland.ac.nz> wrote:
Hi all.
I seem to have managed (not quite sure how! :-) ) with considerable
assistance from the maintainer of the JuliaCall package, to get that package
running.
I am however having a struggle to figure out how to fit the sort of model
that I am interested in. The vignette that Doug Bates provided a while back
gives an example of fitting a linear mixed model, and I have managed to
reproduce that example. (My results differed from those that Doug's
vignette showed, ever so slightly in the values of some of the variance
components, but I'm not going to worry my pretty little head about that.)
Now I want to proceed to my "real problem" which is to fit a *generalised*
linear mixed model, of the binomial persuasion, and I cannot find my way to
documentation on how to do this.
The syntax for the sort of fit that I want would be, in the glmer() context,
of the form:
fit <- glmer(cbind(Dead,Alive) ~ (Trt+0)/prd + (prd | Rep),data=Dat,
family=binomial(link=logit))
where "Trt" is a factor (fixed effect), "prd" is a numeric predictor,
and "Rep" is a factor (random effect).
Can anyone instruct me (in very simple terms please, I'm slow!) as to what
the analogous syntax would be using the MixedModels package in Julia? Or
point me at an elementary tutorial on doing this?
I have the impression, from one item that I saw on StackExchange, that
MixedModels doesn't really handle general binomial models, only Bernoulli
models. Consequently one has to expand out one's data set
creating one row for each success ("Dead" in my case) and one row for each
failure ("Alive"). Is this correct, or is there an easier way to go about
it?
Thanks for any words of wisdom.
cheers,
Rolf Turner
Technical Editor ANZJS
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
On 18/04/18 15:58, Ben Bolker wrote:
Don't have time to dig into this right now (I'm not an expert), but the docs at https://github.com/dmbates/MixedModels.jl/blob/master/docs/jmd/constructors.jmd#L104 say:
Note that, in keeping with convention in the [`GLM` package](https://github.com/JuliaStats/GLM.jl), the distribution family for a binary (i.e. 0/1) response is the `Bernoulli` distribution.
The `Binomial` distribution is only used when the response is the
fraction of trials returning a positive, in which case the number of
trials must be specified as the case weights.
... which implies that you can do what you want.
https://github.com/dmbates/MixedModels.jl/blob/fd22daf226938a4f98ce3d4b228e5a99437e3218/src/pls.jl#L61
seem to give 'weights' as the third argument to the LinearMixedModel
function ... ?
Thanks Ben. That looks like it holds out some hope. I'll try to pursue it. I still cannot discern aspects of the syntax however. The docs say that the "canonical link, which is `GLM.LogitLink` for the `Bernoulli` distribution, is used if no explicit link is specified." But I can't find how to specify other links (e.g. cloglog). I shall keep searching. Thanks again. cheers, Rolf
Technical Editor ANZJS Department of Statistics University of Auckland Phone: +64-9-373-7599 ext. 88276