[R-sig-ME] Error message in ordinal() package running cumulative link mixed models function clmm()
Sverre Stausland
johnsen at fas.harvard.edu
Wed Jan 12 21:49:21 CET 2011
Dear Rune,
1) I appended control=clmm.control(innerCtrl="noWarn") to my clmm
call, as you suggested, and indeed, it says:
Warning message:
clmm may not have converged:
optimizer 'ucminf' terminated with max|gradient|: 0.000952031664607097
2) The codings 1.5 and 6.5 for the Response thresholds were intended.
One subject interpreted the explanations "Not at all" and "Perfectly"
as the endpoints of the scale rather than as comments for the intended
endpoints 1 and 7. Instead of expanding the endpoints of the scale for
this subject, I encoded "Not at all" as 1 and "1" as 1.5, and
"Perfectly" as 7 and "7" as 6.5. Since this is true only for one
subject, I collapsed these into 1 and 7, as you suggested. It reduces
the Condition number of Hessian in the clm() call from 342793.76 to
242545.48, but clmm() still does not converge, now with optimizer
'ucminf' terminated with max|gradient|: 0.00111011989404809
2.2) It sounds like a good idea to fit the clmm() function without
fixed effects to see if it still doesn't converge, but I haven't been
able to construct the right call for clmm() to do that. How would I do
it?
I will send you an e-mail with my data set.
Thanks for you help
Sverre
On Wed, Jan 12, 2011 at 7:44 AM, Rune Haubo <rhbc at imm.dtu.dk> wrote:
> 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