[R] multilevel binary and ordered regression models

Rune Haubo rune.haubo at gmail.com
Thu Jun 13 21:11:07 CEST 2013


Could you share the results of sessionInfo() and str(alllev)?

Also please share the exact in- and output with relevant error
messages; for example 'cntnew:male' does not make much sense without
context.

Unfortunately I don't understand your model specification and is lost
in the interpretation of gammas, w's, x's, mu's beta's and their
indexes, so it is hard for me to say whether you are specifying the
models you seek correctly i clmm. Perhaps the package vignettes can be
an inspiration here? Ahh, I see you have discovered those.

And yes, the current implementation only allows random effects terms
on the form (1| <grouping_factor>), so unless x1 is a factor
(1|group:x1) wont work. I currently have a working implementation that
allows all kinds of random effects terms on the form
(<general_formula> | <grouping_factor>) like glmer allows for, but I
lack the time to wrap it up and test it before I can release it on
CRAN.

Others have successfully fitted models to large data sets (nobs > 1e6)
with several random effects terms, so I don't expect that is a
problem. If you there are a lot of parameters in your models, then,
yes, they can take a while to fit. Equally important is the
identifiability of the estimation problem; the optimum of a ill
defined model can be very hard to locate.

Cheers,
Rune

On 9 June 2013 16:44, Xu Jun <junxu.r at gmail.com> wrote:
> Rune,
>
> Thanks a lot for pointing me to your ordinal package. It is wonderful,
> and I tried a random intercept model and it worked well except that
> probably there is something wrong with my data (size is big), I got
> some warning messages indicating that "In sqrt(diag(vc)[1:npar]) :
> NaNs produced." Therefore, some of the coefficients do not have
> standard errors. I also have a follow-up question: if I want to
> estimate a slope as outcome model, how am I supposed to specify the
> model. I tried following the few tutorials you posted on CRAN, but was
> not able to figure it out:
>
>
> Here is my model:
>
> level 1: y*_ij = beta_0j + beta_1j*x1_ij + beta_2j*x2_ij +
> beta_3j*x3_ij + epsilon_ij
>
> where y* is a latent continuous variable and y is an observed binary
> dependent variable.
> y  is an observed ordinal variable, say ranging from 1-4, and has been
> coded as a factor variable.
>
> x's are predictors at level 1
> beta's are regression coefficients at level 1
> epsilon's are error terms in level 1 equations
>
> level 2 Eq1: beta_0j  = gamma_00 + gamma_01*w1_j + gamma_02*w2_j + mu_0j
> Level 2 Eq2: beta_1j  = gamma_10 + gamma_11*w1_j + gamma_02*w2_j + mu_1j
>
> My guess would be:
> # try 1
> clmm(y~ x1 + x2 + x3 + w1 + w2 + w1:x1 + w2:x2 + (1 + x1 | group), #
> like the syntax in lmer or glmer
>                      data=alllev, link="logit",
>                      na.action=na.omit,
>                      Hess=T)
>
> # or try 2
>
> clmm(y~ x1 + x2 + x3 + w1 + w2 + w1:x1 + w2:x2 + (1  | group) + (1 | group:x1),
>                      data=alllev, link="logit",
>                      na.action=na.omit,
>                      Hess=T)
>
> but none worked. After I issued the first try, I got the following message:
>
> Error: Matrices must have same number of columns in rbind2(..1, r)
> In addition: Warning messages:
> 1: In cntnew:male :
>   numerical expression has 177770 elements: only the first used
>
> and after the second try, simply it says that I got to following the
> (1|factor) format. I would appreciate that if you could point me to
> the right direction. Also, I know I am dealing with a relatively large
> data set, but is there any way to speed up the estimation a bit.
> Thanks a lot!
>
> Jun
>
> On Fri, Jun 7, 2013 at 1:04 AM, Rune Haubo <rune.haubo at gmail.com> wrote:
>> On 6 June 2013 00:13, Xu Jun <junxu.r at gmail.com> wrote:
>>> Dear r-helpers,
>>>
>>> I have two questions on multilevel binary and ordered regression models,
>>> respectively:
>>>
>>> 1. Is there any r function (like lmer or glmer) to run multilevel ordered
>>> regression models?
>>
>> Yes, package ordinal will fit such models.
>>
>> Cheers,
>> Rune



More information about the R-help mailing list