[R-sig-ME] MCMCglmm
Jarrod Hadfield
j.hadfield at ed.ac.uk
Fri Apr 9 14:38:35 CEST 2010
Hi,
In order to update the covariance matrix it is much easier if every
combination of phase and marker exist. MCMCglmm generates any missing
combinations and treats the unknown responses as missing data. It is
just a computational strategy and the warning message can be ignored
(I may suppress it in future versions).
The error message is unrelated. I presume phase is a factor (?) with n
levels. At some iteration the nXn covariance matrix associated with
us(0+phase):marker becomes singular (or close to). If n=2 the error
message implies that a variance has hit zero, or a correlation has hit
-1 and 1. If n>2 then this implies that one (or more) eigenvalues of
the covariance matrix has hit zero. Numerical problems arise when
these conditions occur so MCMCglmm terminates. In your analysis you
have used the default flat priors, but if a proper prior is specified
these conditions do not generally arise. You can choose from the
standard inverse-Wishart prior or from the parameter expanded non-
central F prior (see Gelman 2006 Bayesian Analysis 1 515-533, or the
CourseNotes). The latter is particularly useful if the variances are
close to zero because you can get dramatic improvements in mixing.
Cheers,
Jarrod
On 9 Apr 2010, at 13:03, Iasonas Lamprianou wrote:
> Dear Jarrod Hadfield
> I followed your advice but got this message:
>
> m5_06.mcmc<- MCMCglmm(score ~ 1, random=~us(0+phase):marker+candidate
> +batch, data=mg2006_sub)
>
> Warning in MCMCglmm(score ~ 1, random = ~us(0 + phase):marker +
> candidate + :
> some combinations in us(0 + phase):marker do not exist and 118
> missing records have been generated
> Error in MCMCglmm(score ~ 1, random = ~us(0 + phase):marker +
> candidate + :
> ill-conditioned G/R structure: use proper priors if you haven't or
> rescale data if you have
>
>
>
> could you please help me?
>
> Dr. Iasonas Lamprianou
>
>
> Assistant Professor (Educational Research and Evaluation)
> Department of Education Sciences
> European University-Cyprus
> P.O. Box 22006
> 1516 Nicosia
> Cyprus
> Tel.: +357-22-713178
> Fax: +357-22-590539
>
>
> Honorary Research Fellow
> Department of Education
> The University of Manchester
> Oxford Road, Manchester M13 9PL, UK
> Tel. 0044 161 275 3485
> iasonas.lamprianou at manchester.ac.uk
>
>
> --- On Thu, 8/4/10, Jarrod Hadfield <j.hadfield at ed.ac.uk> wrote:
>
>> From: Jarrod Hadfield <j.hadfield at ed.ac.uk>
>> Subject: Re: [R-sig-ME] MCMCglmm
>> To: "Iasonas Lamprianou" <lamprianou at yahoo.com>
>> Cc: r-sig-mixed-models at r-project.org
>> Date: Thursday, 8 April, 2010, 15:10
>> Dear Jason,
>>
>> These should work:
>>
>> m4.mcmc <- MCMCglmm(score ~ 1,
>> random=~marker+marker:day+candidate
>> +batch, data=mg2006)
>> m5 .mcmc<- MCMCglmm(score ~ 1,
>> random=~us(1+day):marker+candidate
>> +batch, data=mg2006)
>>
>> You may need to code day as a factor for m4, and as numeric
>> for m5
>> depending on the model you actually want to fit.
>>
>> Cheers,
>>
>> Jarrod
>>
>>
>> On 8 Apr 2010, at 14:58, Iasonas Lamprianou wrote:
>>
>>> Dear all,
>>> I recently experimented with MCMCglmm and I loved
>> (really loved) the
>>> fact that it will give me confidence intervals for the
>> variance of
>>> the random effects. It seems that MCMC is a reasonable
>> method to do
>>> so, in contrast to REML which seems to have problems
>> on this front.
>>> However, MCMCglmm is painfully slower than lmer which
>> is more
>>> familiar to me. The good news is that the point
>> estimates of lmer
>>> are near the centre of the confidence intervals by
>> MCMCglmm.
>>>
>>> I reduced my sample size a bit and managed to fit
>> those two models
>>> with lmer (the second would not fit because it needed
>> 1.5GB or RAM).
>>> Both seem to have a reasonable fit (at least at first
>> look).
>>>
>>> m4 <- lmer(score ~
>> 1+(1|marker/day)+(1|candidate)+(1|batch), mg2006)
>>>
>>> m5 <- lmer(score ~
>> 1+(1+day|marker)+(1|candidate)+(1|batch), mg2006)
>>>
>>> I would like to run these two models above with
>> MCMCglmm. Does
>>> anyone know how to do it?
>>>
>>> Thank you for the help
>>>
>>> Jason
>>>
>>> Dr. Iasonas Lamprianou
>>>
>>>
>>> Assistant Professor (Educational Research and
>> Evaluation)
>>> Department of Education Sciences
>>> European University-Cyprus
>>> P.O. Box 22006
>>> 1516 Nicosia
>>> Cyprus
>>> Tel.: +357-22-713178
>>> Fax: +357-22-590539
>>>
>>>
>>> Honorary Research Fellow
>>> Department of Education
>>> The University of Manchester
>>> Oxford Road, Manchester M13 9PL, UK
>>> Tel. 0044 161 275 3485
>>> iasonas.lamprianou at manchester.ac.uk
>>>
>>>
>>> --- On Thu, 8/4/10, r-sig-mixed-models-request at r-project.org
>> <r-sig-mixed-models-request at r-project.org
>>
>>>> wrote:
>>>
>>>> From: r-sig-mixed-models-request at r-project.org
>> <r-sig-mixed-models-request at r-project.org
>>
>>>>>
>>>> Subject: R-sig-mixed-models Digest, Vol 40, Issue
>> 15
>>>> To: r-sig-mixed-models at r-project.org
>>>> Date: Thursday, 8 April, 2010, 12:13
>>>> Send R-sig-mixed-models mailing list
>>>> submissions to
>>>> r-sig-mixed-models at 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 at r-project.org
>>>>
>>>> You can reach the person managing the list at
>>>> r-sig-mixed-models-owner at 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: Multi-level models Odds ratio
>> (E
>>>> T)
>>>> 2. Re: Multi-level models Odds ratio
>>>> (ONKELINX, Thierry)
>>>> 3. Re: Multi-level models Odds ratio
>> (E
>>>> T)
>>>> 4. Re: Multi-level models Odds ratio
>>>> (Andy Fugard (Work))
>>>> 5. Re: Multi-level models Odds ratio
>> (E
>>>> T)
>>>>
>>>>
>>>>
>> ----------------------------------------------------------------------
>>>>
>>>> Message: 1
>>>> Date: Thu, 8 Apr 2010 11:27:49 +0100
>>>> From: E T <2nuzzbot at gmail.com>
>>>> To: Daniel Ezra Johnson <danielezrajohnson at gmail.com>
>>>> Cc: "r-sig-mixed-models at r-project.org"
>>>> <r-sig-mixed-models at r-project.org>
>>>> Subject: Re: [R-sig-ME] Multi-level models Odds
>> ratio
>>>> Message-ID:
>>>> <l2w706f8d1f1004080327re8708f46mc7d334f41ae19a10 at mail.gmail.com
>>>> >
>>>> Content-Type: text/plain
>>>>
>>>> odds.ratios = exp(coefs(model))
>>>>
>>>> Thanks, however unfortunately when I try the above
>> command
>>>> I receive the
>>>> following error:
>>>>
>>>> Error: could not find function "coefs"
>>>>
>>>> Regards
>>>>
>>>> Et
>>>>
>>>>
>>>>
>>>> On Wed, Apr 7, 2010 at 5:47 PM, Daniel Ezra
>> Johnson <
>>>> danielezrajohnson at gmail.com>
>>>> wrote:
>>>>
>>>>> something like odds.ratios =
>> exp(coefs(model))
>>>>>
>>>>>
>>>>> On Apr 7, 2010, at 12:28 PM, E T <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>
>>>>> Hi all,
>>>>>>
>>>>>> Apologies for the simplicity of my
>> question....
>>>> however any advice is
>>>>>> greatly appreciated. Thanks
>>>>>>
>>>>>> Is there a specific command available to
>> obtain
>>>> the odds ratios produced
>>>>>> from a multilevel logistic model?
>>>>>>
>>>>>> I have estimated a multi-level logistic
>> model
>>>> using the lme4 package. I
>>>>>> can
>>>>>> obtain results using the 'summary'
>> command,
>>>> however I would like to obtain
>>>>>> the computed odds ratios.
>>>>>> (Similar to the output that can be
>> produced for
>>>> logistic GLM using the
>>>>>> logistic.display command from the epicalc
>>>> package).
>>>>>>
>>>>>> [[alternative HTML version
>>>> deleted]]
>>>>>>
>>>>>>
>> _______________________________________________
>>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>>
>>>>>
>>>>
>>>> [[alternative HTML version
>> deleted]]
>>>>
>>>>
>>>>
>>>> ------------------------------
>>>>
>>>> Message: 2
>>>> Date: Thu, 8 Apr 2010 12:32:01 +0200
>>>> From: "ONKELINX, Thierry" <Thierry.ONKELINX at inbo.be>
>>>> To: "E T" <2nuzzbot at gmail.com>,
>>>> "Daniel Ezra Johnson"
>>>> <danielezrajohnson at gmail.com>
>>>> Cc: r-sig-mixed-models at r-project.org
>>>> Subject: Re: [R-sig-ME] Multi-level models Odds
>> ratio
>>>> Message-ID:
>>>> <2E9C414912813E4EB981326983E0A104071B69A6 at inboexch.inbo.be>
>>>> Content-Type: text/plain;
>>>> charset="us-ascii"
>>>>
>>>> It should be
>>>>
>>>> exp(coef(model))
>>>>
>>>> Without the "s"
>>>>
>>>> HTH,
>>>>
>>>> Thierry
>>>>
>> ------------------------------------------------------------------------
>>>> ----
>>>> ir. Thierry Onkelinx
>>>> Instituut voor natuur- en bosonderzoek
>>>> team Biometrie & Kwaliteitszorg
>>>> Gaverstraat 4
>>>> 9500 Geraardsbergen
>>>> Belgium
>>>>
>>>> Research Institute for Nature and Forest
>>>> team Biometrics & Quality Assurance
>>>> Gaverstraat 4
>>>> 9500 Geraardsbergen
>>>> Belgium
>>>>
>>>> tel. + 32 54/436 185
>>>> Thierry.Onkelinx at inbo.be
>>>> www.inbo.be
>>>>
>>>> To call in the statistician after the experiment
>> is done
>>>> may be no more
>>>> than asking him to perform a post-mortem
>> examination: he
>>>> may be able to
>>>> say what the experiment died of.
>>>> ~ Sir Ronald Aylmer Fisher
>>>>
>>>> The plural of anecdote is not data.
>>>> ~ Roger Brinner
>>>>
>>>> The combination of some data and an aching desire
>> for an
>>>> answer does not
>>>> ensure that a reasonable answer can be extracted
>> from a
>>>> given body of
>>>> data.
>>>> ~ John Tukey
>>>>
>>>>
>>>>> -----Oorspronkelijk bericht-----
>>>>> Van: r-sig-mixed-models-bounces at r-project.org
>>>>
>>>>> [mailto:r-sig-mixed-models-bounces at r-project.org]
>>>> Namens E T
>>>>> Verzonden: donderdag 8 april 2010 12:28
>>>>> Aan: Daniel Ezra Johnson
>>>>> CC: r-sig-mixed-models at r-project.org
>>>>> Onderwerp: Re: [R-sig-ME] Multi-level models
>> Odds
>>>> ratio
>>>>>
>>>>> odds.ratios = exp(coefs(model))
>>>>>
>>>>> Thanks, however unfortunately when I try the
>> above
>>>> command I
>>>>> receive the following error:
>>>>>
>>>>> Error: could not find function "coefs"
>>>>>
>>>>> Regards
>>>>>
>>>>> Et
>>>>>
>>>>>
>>>>>
>>>>> On Wed, Apr 7, 2010 at 5:47 PM, Daniel Ezra
>> Johnson
>>>> <
>>>>> danielezrajohnson at gmail.com>
>>>> wrote:
>>>>>
>>>>>> something like odds.ratios =
>> exp(coefs(model))
>>>>>>
>>>>>>
>>>>>> On Apr 7, 2010, at 12:28 PM, E T <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>>
>>>>>> Hi all,
>>>>>>>
>>>>>>> Apologies for the simplicity of my
>>>> question.... however
>>>>> any advice is
>>>>>>> greatly appreciated. Thanks
>>>>>>>
>>>>>>> Is there a specific command available
>> to
>>>> obtain the odds ratios
>>>>>>> produced from a multilevel logistic
>> model?
>>>>>>>
>>>>>>> I have estimated a multi-level
>> logistic model
>>>> using the
>>>>> lme4 package.
>>>>>>> I can obtain results using the
>> 'summary'
>>>> command, however I would
>>>>>>> like to obtain the computed odds
>> ratios.
>>>>>>> (Similar to the output that can be
>> produced
>>>> for logistic GLM using
>>>>>>> the logistic.display command from the
>> epicalc
>>>> package).
>>>>>>>
>>>>>>> [[alternative HTML
>> version
>>>> deleted]]
>>>>>>>
>>>>>>>
>>>> _______________________________________________
>>>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>>>
>>>>>>
>>>>>
>>>>> [[alternative HTML
>> version
>>>> deleted]]
>>>>>
>>>>>
>> _______________________________________________
>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>
>>>>
>>>> Druk dit bericht a.u.b. niet onnodig af.
>>>> Please do not print this message unnecessarily.
>>>>
>>>> Dit bericht en eventuele bijlagen geven enkel de
>> visie van
>>>> de schrijver weer
>>>> en binden het INBO onder geen enkel beding, zolang
>> dit
>>>> bericht niet bevestigd is
>>>> door een geldig ondertekend document. The views
>> expressed
>>>> in this message
>>>> and any annex are purely those of the writer and
>> may not be
>>>> regarded as stating
>>>> an official position of INBO, as long as the
>> message is not
>>>> confirmed by a duly
>>>> signed document.
>>>>
>>>>
>>>>
>>>> ------------------------------
>>>>
>>>> Message: 3
>>>> Date: Thu, 8 Apr 2010 11:35:03 +0100
>>>> From: E T <2nuzzbot at gmail.com>
>>>> To: Daniel Ezra Johnson <danielezrajohnson at gmail.com>
>>>> Cc: "r-sig-mixed-models at r-project.org"
>>>> <r-sig-mixed-models at r-project.org>
>>>> Subject: Re: [R-sig-ME] Multi-level models Odds
>> ratio
>>>> Message-ID:
>>>> <i2y706f8d1f1004080335q84f61b78u8c7b656b67a08a8e at mail.gmail.com
>>>> >
>>>> Content-Type: text/plain
>>>>
>>>> If I use the command coef(model) this extracts
>> the
>>>> coefficients in the
>>>> model, however if I try exp(coef(model)) I receive
>> an
>>>> error:
>>>>
>>>> Error in exp(coef(model)) : Non-numeric argument
>> to
>>>> mathematical function
>>>>
>>>> I could manually get the exp of each factor in my
>>>> model..... but as I have a
>>>> large model (and also have numerous other models
>> to
>>>> produce), I was
>>>> wondering if there was an automated method
>>>>
>>>> Regards
>>>>
>>>> Et
>>>>
>>>> On Thu, Apr 8, 2010 at 11:27 AM, E T <2nuzzbot at gmail.com>
>>>> wrote:
>>>>
>>>>> odds.ratios = exp(coefs(model))
>>>>>
>>>>> Thanks, however unfortunately when I try the
>> above
>>>> command I receive the
>>>>> following error:
>>>>>
>>>>> Error: could not find function "coefs"
>>>>>
>>>>> Regards
>>>>>
>>>>> Et
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On Wed, Apr 7, 2010 at 5:47 PM, Daniel Ezra
>> Johnson
>>>> <
>>>>> danielezrajohnson at gmail.com>
>>>> wrote:
>>>>>
>>>>>> something like odds.ratios =
>> exp(coefs(model))
>>>>>>
>>>>>>
>>>>>> On Apr 7, 2010, at 12:28 PM, E T <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>>
>>>>>> Hi all,
>>>>>>>
>>>>>>> Apologies for the simplicity of my
>>>> question.... however any advice is
>>>>>>> greatly appreciated. Thanks
>>>>>>>
>>>>>>> Is there a specific command available
>> to
>>>> obtain the odds ratios produced
>>>>>>> from a multilevel logistic model?
>>>>>>>
>>>>>>> I have estimated a multi-level
>> logistic model
>>>> using the lme4 package. I
>>>>>>> can
>>>>>>> obtain results using the 'summary'
>> command,
>>>> however I would like to
>>>>>>> obtain
>>>>>>> the computed odds ratios.
>>>>>>> (Similar to the output that can be
>> produced
>>>> for logistic GLM using the
>>>>>>> logistic.display command from the
>> epicalc
>>>> package).
>>>>>>>
>>>>>>> [[alternative HTML
>> version
>>>> deleted]]
>>>>>>>
>>>>>>>
>>>> _______________________________________________
>>>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>>>
>>>>>>
>>>>>
>>>>
>>>> [[alternative HTML version
>> deleted]]
>>>>
>>>>
>>>>
>>>> ------------------------------
>>>>
>>>> Message: 4
>>>> Date: Thu, 08 Apr 2010 12:48:30 +0200
>>>> From: "Andy Fugard (Work)" <andy.fugard at sbg.ac.at>
>>>> To: E T <2nuzzbot at gmail.com>
>>>> Cc: "r-sig-mixed-models at r-project.org"
>>>> <r-sig-mixed-models at r-project.org>
>>>> Subject: Re: [R-sig-ME] Multi-level models Odds
>> ratio
>>>> Message-ID: <4BBDB47E.8030305 at sbg.ac.at>
>>>> Content-Type: text/plain; charset=ISO-8859-1
>>>>
>>>> Here's another example, borrowed from the help for
>> "lmer":
>>>>
>>>>> gm1 <- glmer(cbind(incidence, size -
>> incidence) ~
>>>> period + (1 | herd),
>>>>
>> family =
>>>> binomial, data = cbpp)
>>>>
>>>> As you say, coef works:
>>>>
>>>>> coef(gm1)
>>>> $herd
>>>> (Intercept)
>>>>
>> period2 period3 period4
>>>> 1 -0.8085096 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 2 -1.6974292 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 3 -0.9922697 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 4 -1.3592525 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 5 -1.5885461 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 6 -1.7987950 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 7 -0.5091313 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 8 -0.7991613 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 9 -1.6361848 -0.9923347
>> -1.128675
>>>> -1.580374
>>>> 10 -1.9394614 -0.9923347 -1.128675
>> -1.580374
>>>> 11 -1.4831632 -0.9923347 -1.128675
>> -1.580374
>>>> 12 -1.4633469 -0.9923347 -1.128675
>> -1.580374
>>>> 13 -2.0884474 -0.9923347 -1.128675
>> -1.580374
>>>> 14 -0.4278151 -0.9923347 -1.128675
>> -1.580374
>>>> 15 -1.9290041 -0.9923347 -1.128675
>> -1.580374
>>>>
>>>> But note the "$herd" bit. Since this model
>> has a
>>>> varying intercept by
>>>> herd, you get a column in the resulting data frame
>> called
>>>> "herd".
>>>>
>>>> So you could try, for this example:
>>>>
>>>>> exp(coef(gm1)$herd)
>>>>
>> (Intercept) period2 period3 period4
>>>> 1 0.4455216 0.3707102 0.3234614
>> 0.2058981
>>>> 2 0.1831538 0.3707102 0.3234614
>> 0.2058981
>>>> 3 0.3707343 0.3707102 0.3234614
>> 0.2058981
>>>> 4 0.2568527 0.3707102 0.3234614
>> 0.2058981
>>>> 5 0.2042223 0.3707102 0.3234614
>> 0.2058981
>>>> 6 0.1654982 0.3707102 0.3234614
>> 0.2058981
>>>> 7 0.6010174 0.3707102 0.3234614
>> 0.2058981
>>>> 8 0.4497060 0.3707102 0.3234614
>> 0.2058981
>>>> 9 0.1947215 0.3707102 0.3234614
>> 0.2058981
>>>> 10 0.1437814 0.3707102 0.3234614
>>>> 0.2058981
>>>> 11 0.2269188 0.3707102 0.3234614
>>>> 0.2058981
>>>> 12 0.2314603 0.3707102 0.3234614
>>>> 0.2058981
>>>> 13 0.1238793 0.3707102 0.3234614
>>>> 0.2058981
>>>> 14 0.6519320 0.3707102 0.3234614
>>>> 0.2058981
>>>> 15 0.1452928 0.3707102 0.3234614
>>>> 0.2058981
>>>>
>>>> Since the slopes don't vary by herd, you might
>> also want
>>>> just the fixed
>>>> effects:
>>>>
>>>>> exp(fixef(gm1))
>>>> (Intercept) period2
>>>> period3
>> period4
>>>>
>>>>
>> 0.2469585 0.3707102 0.3234614 0.2058981
>>>>
>>>> HTH,
>>>>
>>>> Andy
>>>>
>>>>
>>>> E T wrote:
>>>>> If I use the command coef(model) this extracts
>> the
>>>> coefficients in the
>>>>> model, however if I try exp(coef(model)) I
>> receive an
>>>> error:
>>>>>
>>>>> Error in exp(coef(model)) : Non-numeric
>> argument to
>>>> mathematical function
>>>>>
>>>>> I could manually get the exp of each factor in
>> my
>>>> model..... but as I have a
>>>>> large model (and also have numerous other
>> models to
>>>> produce), I was
>>>>> wondering if there was an automated method
>>>>>
>>>>> Regards
>>>>>
>>>>> Et
>>>>>
>>>>> On Thu, Apr 8, 2010 at 11:27 AM, E T <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>
>>>>>> odds.ratios = exp(coefs(model))
>>>>>>
>>>>>> Thanks, however unfortunately when I try
>> the above
>>>> command I receive the
>>>>>> following error:
>>>>>>
>>>>>> Error: could not find function "coefs"
>>>>>>
>>>>>> Regards
>>>>>>
>>>>>> Et
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> On Wed, Apr 7, 2010 at 5:47 PM, Daniel
>> Ezra
>>>> Johnson <
>>>>>> danielezrajohnson at gmail.com>
>>>> wrote:
>>>>>>
>>>>>>> something like odds.ratios =
>>>> exp(coefs(model))
>>>>>>>
>>>>>>>
>>>>>>> On Apr 7, 2010, at 12:28 PM, E T
>> <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>>>
>>>>>>> Hi all,
>>>>>>>> Apologies for the simplicity of
>> my
>>>> question.... however any advice is
>>>>>>>> greatly appreciated. Thanks
>>>>>>>>
>>>>>>>> Is there a specific command
>> available to
>>>> obtain the odds ratios produced
>>>>>>>> from a multilevel logistic model?
>>>>>>>>
>>>>>>>> I have estimated a multi-level
>> logistic
>>>> model using the lme4 package. I
>>>>>>>> can
>>>>>>>> obtain results using the
>> 'summary'
>>>> command, however I would like to
>>>>>>>> obtain
>>>>>>>> the computed odds ratios.
>>>>>>>> (Similar to the output that can
>> be
>>>> produced for logistic GLM using the
>>>>>>>> logistic.display command from the
>> epicalc
>>>> package).
>>>>>>>>
>>>>>>>> [[alternative HTML
>>>> version deleted]]
>>>>>>>>
>>>>>>>>
>>>> _______________________________________________
>>>>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>>>>
>>>>>
>>>>> [[alternative HTML
>> version
>>>> deleted]]
>>>>>
>>>>>
>> _______________________________________________
>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>
>>>>
>>>> --
>>>> Andy Fugard, Postdoctoral researcher, ESF LogICCC
>> project
>>>> "Modeling human inference within the framework of
>>>> probability logic"
>>>> Department of Psychology, University of Salzburg,
>> Austria
>>>> http://www.andyfugard.info
>>>>
>>>>
>>>>
>>>> ------------------------------
>>>>
>>>> Message: 5
>>>> Date: Thu, 8 Apr 2010 12:13:23 +0100
>>>> From: E T <2nuzzbot at gmail.com>
>>>> To: "Andy Fugard (Work)" <andy.fugard at sbg.ac.at>
>>>> Cc: "r-sig-mixed-models at r-project.org"
>>>> <r-sig-mixed-models at r-project.org>
>>>> Subject: Re: [R-sig-ME] Multi-level models Odds
>> ratio
>>>> Message-ID:
>>>> <p2p706f8d1f1004080413k3014fe10q444382c927c2f90e at mail.gmail.com
>>>> >
>>>> Content-Type: text/plain
>>>>
>>>> exp(coef(model)$group)
>>>>
>>>> exp(fixef(model))
>>>>
>>>> Thanks.... yes this worked successfully :o)
>>>>
>>>> Et
>>>>
>>>> On Thu, Apr 8, 2010 at 11:48 AM, Andy Fugard
>> (Work)
>>>> <andy.fugard at sbg.ac.at>wrote:
>>>>
>>>>> Here's another example, borrowed from the help
>> for
>>>> "lmer":
>>>>>
>>>>>> gm1 <- glmer(cbind(incidence, size -
>>>> incidence) ~ period + (1 | herd),
>>>>>
>> family
>>>> = binomial, data = cbpp)
>>>>>
>>>>> As you say, coef works:
>>>>>
>>>>>> coef(gm1)
>>>>> $herd
>>>>> (Intercept)
>>>>
>> period2 period3 period4
>>>>> 1 -0.8085096 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 2 -1.6974292 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 3 -0.9922697 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 4 -1.3592525 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 5 -1.5885461 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 6 -1.7987950 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 7 -0.5091313 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 8 -0.7991613 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 9 -1.6361848 -0.9923347
>> -1.128675
>>>> -1.580374
>>>>> 10 -1.9394614 -0.9923347 -1.128675
>> -1.580374
>>>>> 11 -1.4831632 -0.9923347 -1.128675
>> -1.580374
>>>>> 12 -1.4633469 -0.9923347 -1.128675
>> -1.580374
>>>>> 13 -2.0884474 -0.9923347 -1.128675
>> -1.580374
>>>>> 14 -0.4278151 -0.9923347 -1.128675
>> -1.580374
>>>>> 15 -1.9290041 -0.9923347 -1.128675
>> -1.580374
>>>>>
>>>>> But note the "$herd" bit. Since this
>> model has a
>>>> varying intercept by
>>>>> herd, you get a column in the resulting data
>> frame
>>>> called "herd".
>>>>>
>>>>> So you could try, for this example:
>>>>>
>>>>>> exp(coef(gm1)$herd)
>>>>>
>> (Intercept) period2 period3 period4
>>>>> 1 0.4455216 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 2 0.1831538 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 3 0.3707343 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 4 0.2568527 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 5 0.2042223 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 6 0.1654982 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 7 0.6010174 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 8 0.4497060 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 9 0.1947215 0.3707102 0.3234614
>>>> 0.2058981
>>>>> 10 0.1437814 0.3707102
>> 0.3234614
>>>> 0.2058981
>>>>> 11 0.2269188 0.3707102
>> 0.3234614
>>>> 0.2058981
>>>>> 12 0.2314603 0.3707102
>> 0.3234614
>>>> 0.2058981
>>>>> 13 0.1238793 0.3707102
>> 0.3234614
>>>> 0.2058981
>>>>> 14 0.6519320 0.3707102
>> 0.3234614
>>>> 0.2058981
>>>>> 15 0.1452928 0.3707102
>> 0.3234614
>>>> 0.2058981
>>>>>
>>>>> Since the slopes don't vary by herd, you might
>> also
>>>> want just the fixed
>>>>> effects:
>>>>>
>>>>>> exp(fixef(gm1))
>>>>> (Intercept) period2
>>>> period3
>> period4
>>>>>
>>>>
>> 0.2469585 0.3707102 0.3234614 0.2058981
>>>>>
>>>>> HTH,
>>>>>
>>>>> Andy
>>>>>
>>>>>
>>>>> E T wrote:
>>>>>> If I use the command coef(model) this
>> extracts
>>>> the coefficients in the
>>>>>> model, however if I try exp(coef(model))
>> I
>>>> receive an error:
>>>>>>
>>>>>> Error in exp(coef(model)) : Non-numeric
>> argument
>>>> to mathematical function
>>>>>>
>>>>>> I could manually get the exp of each
>> factor in my
>>>> model..... but as I
>>>>> have a
>>>>>> large model (and also have numerous other
>> models
>>>> to produce), I was
>>>>>> wondering if there was an automated
>> method
>>>>>>
>>>>>> Regards
>>>>>>
>>>>>> Et
>>>>>>
>>>>>> On Thu, Apr 8, 2010 at 11:27 AM, E T
>> <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>>
>>>>>>> odds.ratios = exp(coefs(model))
>>>>>>>
>>>>>>> Thanks, however unfortunately when I
>> try the
>>>> above command I receive the
>>>>>>> following error:
>>>>>>>
>>>>>>> Error: could not find function
>> "coefs"
>>>>>>>
>>>>>>> Regards
>>>>>>>
>>>>>>> Et
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On Wed, Apr 7, 2010 at 5:47 PM, Daniel
>> Ezra
>>>> Johnson <
>>>>>>> danielezrajohnson at gmail.com>
>>>> wrote:
>>>>>>>
>>>>>>>> something like odds.ratios =
>>>> exp(coefs(model))
>>>>>>>>
>>>>>>>>
>>>>>>>> On Apr 7, 2010, at 12:28 PM, E T
>> <2nuzzbot at gmail.com>
>>>> wrote:
>>>>>>>>
>>>>>>>> Hi all,
>>>>>>>>> Apologies for the simplicity
>> of my
>>>> question.... however any advice is
>>>>>>>>> greatly appreciated. Thanks
>>>>>>>>>
>>>>>>>>> Is there a specific command
>> available
>>>> to obtain the odds ratios
>>>>> produced
>>>>>>>>> from a multilevel logistic
>> model?
>>>>>>>>>
>>>>>>>>> I have estimated a
>> multi-level
>>>> logistic model using the lme4 package.
>>>>> I
>>>>>>>>> can
>>>>>>>>> obtain results using the
>> 'summary'
>>>> command, however I would like to
>>>>>>>>> obtain
>>>>>>>>> the computed odds ratios.
>>>>>>>>> (Similar to the output that
>> can be
>>>> produced for logistic GLM using the
>>>>>>>>> logistic.display command from
>> the
>>>> epicalc package).
>>>>>>>>>
>>>>>>>>> [[alternative
>> HTML
>>>> version deleted]]
>>>>>>>>>
>>>>>>>>>
>>>> _______________________________________________
>>>>>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>>>>>
>>>>>>
>>>>>> [[alternative
>> HTML
>>>> version deleted]]
>>>>>>
>>>>>>
>> _______________________________________________
>>>>>> R-sig-mixed-models at r-project.org
>>>> mailing list
>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>
>>>>>
>>>>> --
>>>>> Andy Fugard, Postdoctoral researcher, ESF
>> LogICCC
>>>> project
>>>>> "Modeling human inference within the framework
>> of
>>>> probability logic"
>>>>> Department of Psychology, University of
>> Salzburg,
>>>> Austria
>>>>> http://www.andyfugard.info
>>>>>
>>>>
>>>> [[alternative HTML version
>> deleted]]
>>>>
>>>>
>>>>
>>>> ------------------------------
>>>>
>>>> _______________________________________________
>>>> R-sig-mixed-models mailing list
>>>> R-sig-mixed-models at r-project.org
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>
>>>>
>>>> End of R-sig-mixed-models Digest, Vol 40, Issue
>> 15
>>>>
>> **************************************************
>>>>
>>>
>>>
>>>
>>>
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org
>> mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>
>>
>> --
>> The University of Edinburgh is a charitable body,
>> registered in
>> Scotland, with registration number SC005336.
>>
>>
>
>
>
>
--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.
More information about the R-sig-mixed-models
mailing list