[R-sig-ME] Error message in ordinal() package running cumulative link mixed models function clmm()
Rune Haubo
rhbc at imm.dtu.dk
Wed Jan 12 13:44:57 CET 2011
Dear Sverre,
You seem to be using clmm as it was intended to be used and my best
guess is that the problems you observe are related to data structure
and model identifiability. I have a few ideas that might help, but for
a full diagnose I would need to see the data.
1) My guess is that the clmm did not converge and that the last
warning is exactly a warning about this. You probably cannot see that
unless you suppress the warnings in the random effects update in the
inner loop. You can do this by appending control =
clmm.control(innerCtrl = "noWarn") to your clmm call.
2) The condition number of the Hessian in the fixed effect clm
indicates that the likelihood function is ill-defined and possibly
there is a (near-) ridge. Looking at the threshold estimates, the
first two and the last two are very close and with rather peculiar
codings - is that intended? I would expect that the likelihood is more
well-defined if you collapse the first two and the last two response
categories, and this might be enough to make the clmm converge to a
well-defined optimum. Collapsing categories does not change the nature
of the model interpretations, so it is safe to do so - especially in a
case like this where the thresholds are indistinguishably for all
practical purposes.
2.2) Starting with the "full" model does not always work in mixed
effects models: You could try to fit the clmm without fixed effects at
all. If the convergence problem pertains, this means that the
thresholds are most likely the problem. Given the significance of your
fixed effects in the clm, I don't think the fixed effects are the
problem (though removing Trial could help). You could take a look at
the correlation between parameter estimates from
summary(FullModel.clm, correlation = TRUE) to get a better feeling for
which parameters are hard to distinguish.
3) Sometimes the likelihood function as defined by the Laplace
approximation is not sufficiently well-behaved on the way to the
optimum that the default optimizer used in clmm (ucminf) will converge
smoothly. This can often be rectified by using the more accurate
adaptive Gauss-Hermite quadrature (AGQ) approximation. If you append
nAGQ = 10 to your clmm call you will ask for 10 quadrature nodes which
usually gives reasonable precision.
4) In hard-to-optimize situations it should help to allow more
iterations in the inner loop (random effect update) and you can do
that by setting, e.g. control = clmm.control(maxIter = 200,
maxLineIter = 200). Decreasing gradTol to 1e-4 or 1e-5 might also
help. That being said, I have never seen a case where this was
necessary for a well-defined clmm, so i suggest that you consider the
points above first.
In CLMs model identifiability is a much more frequent issue than in
GLMs, LMs etc. and this is magnified when mixed effects versions are
in play. To provide more ideas, I would need more information, in
particular the results of str(dataset) and sessionInfo() would be
valuable. Ultimately I would need the data - if you feel
uncomfortable sharing them with the list, you can send them to me
privately.
Cheers,
Rune
Ps: You will not be able to have two random effect terms (Subject and
Pair) in clmm in its current implementation.
On 12 January 2011 03:49, Sverre Stausland <johnsen at fas.harvard.edu> wrote:
> Dear all,
>
> I have a data set where subjects gave a 1-7 rating to pairs of items,
> and there's a set of numeric independent variables. The data set has
> the following set-up:
>
> Subject Pair Trial Response A B E
> AM25Y x,y 54 1 2 4 3.13
> AR82M x,y 65 3 2 4 3.13
> BA89W b,n 27 1 5 6 5.17
> CK44O b,n 19 1 5 6 5.17
> etc.
>
> Given the ordinal dependent variable, I need to run a ordinal
> regression model. Since I want to model Subject and Pair as random
> effects, I am trying Rune Haubo's clmm() function from the ordinal()
> package.
>
> First, running a cumulative link model with clm() works fine:
>
>> clm(Response~Trial+A+B+C,data=dataset)->FullModel.clm
>> summary(FullModel.clm)
> Call:
> clm(location = Response ~ Trial + A + B + C, data = dataset)
>
> Location coefficients:
> Estimate Std. Error z value Pr(>|z|)
> Trial -0.0010 0.0009 -1.0638 0.28743
> A -0.1953 0.0247 -7.9144 2.4841e-15
> B -0.1446 0.0240 -6.0166 1.7809e-09
> C -0.2151 0.0204 -10.5287 < 2.22e-16
>
> No scale coefficients
>
> Threshold coefficients:
> Estimate Std. Error z value
> 1|1.5 -2.7595 0.0915 -30.1696
> 1.5|2 -2.7513 0.0914 -30.0925
> 2|3 -1.8590 0.0874 -21.2622
> 3|4 -1.1453 0.0853 -13.4291
> 4|5 -0.4202 0.0849 -4.9511
> 5|6 0.6346 0.0897 7.0737
> 6|6.5 2.3828 0.1289 18.4908
> 6.5|7 2.3942 0.1293 18.5119
>
> log-likelihood: -9104.93
> AIC: 18233.86
> Condition number of Hessian: 342793.76
>
>
> But when I try to run a cumulative link mixed model with Subject as a
> random effect, I get errors:
>
>> clmm(Response~Trial+A+B+C,random=Subject,data=dataset)->FullModel.clmm
> There were 50 or more warnings (use warnings() to see the first 50)
>> warnings()
> Warning messages:
> 1: In rho$updateU(rho) : Non finite negative log-likelihood
> at iteration 112
> 2: In rho$updateU(rho) : Non finite negative log-likelihood
> at iteration 112
> 3: In rho$updateU(rho) : Non finite negative log-likelihood
> at iteration 112
>
>> summary(FullModel.clmm)
> Cumulative Link Mixed Model fitted with the Laplace approximation
>
> Call:
> clmm(location = Response ~ Trial + A + B + C, random = Subject, data = dataset)
>
> Random effects:
> Var Std.Dev
> Subject 2.241108 1.497033
>
> Location coefficients:
> Estimate Std. Error z value Pr(>|z|)
> Trial -0.0015 NaN NaN NA
> A -0.2581 NaN NaN NA
> B -0.1833 NaN NaN NA
> C -0.2833 NaN NaN NA
>
> No scale coefficients
>
> Threshold coefficients:
> Estimate Std. Error z value
> 1|1.5 -3.6036 NaN NaN
> 1.5|2 -3.5921 NaN NaN
> 2|3 -2.3596 NaN NaN
> 3|4 -1.3958 NaN NaN
> 4|5 -0.4317 NaN NaN
> 5|6 0.9001 NaN NaN
> 6|6.5 2.8441 NaN NaN
> 6.5|7 2.8557 NaN NaN
>
> log-likelihood: -8022.125
> AIC: 16070.25
> Condition number of Hessian: NaN
> Warning message:
> In summary.clmm(FullModel.clmm) :
> Variance-covariance matrix of the parameters is not defined
>
> I would appreciate if anyone could point out to me what I am doing
> wrong in my use of the clmm() function.
>
> Thanks
> Sverre
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
--
Rune Haubo Bojesen Christensen
PhD Student, M.Sc. Eng.
Phone: (+45) 45 25 33 63
Mobile: (+45) 30 26 45 54
DTU Informatics, Section for Statistics
Technical University of Denmark, Build. 305, Room 122,
DK-2800 Kgs. Lyngby, Denmark
More information about the R-sig-mixed-models
mailing list