[R-sig-ME] [R-sig-eco] LRT tests in lmer
Chris Mcowen
cm744 at st-andrews.ac.uk
Wed Aug 11 18:44:57 CEST 2010
Hi
> are all 5's for example associated with a single fixed factor, or something like this?
The ordinal response are categorical - different levels of threat, however they can be successfully used as a continuous response ( Purvis, Mace etc) they are not associated with any of the fixed factors, i am trying to use the fixed factors (life history traits) to predict the ordinal response.
I shall have a play with the priors as
> G1=list(V=1, nu=1, alpha.mu=0, alpha.V=1000)
Has improved things, but not greatly
Thanks
Chris
Intercept) -0.23325 -2.89744 2.83429 793.1 0.884
STOStorage organ -0.04486 -0.28088 0.23706 1306.4 0.722
BSUnisexual flower 0.21329 -0.11396 0.52257 861.1 0.206
BSUnisexual plant 0.33547 -0.04818 0.75086 806.5 0.122
PDBiotic 0.28292 -0.13199 0.63020 599.1 0.184
PDMammalia -0.46017 -2.15330 1.44028 862.8 0.640
FRNon_fleshy_fruit -0.22784 -0.54680 0.10850 764.5 0.192
ENDNon_endospermous 0.44173 0.10830 0.74418 747.8 0.016 *
WOWoody -0.22039 -0.59506 0.11227 631.0 0.252
RGTwo+ -0.04816 -0.24944 0.15221 816.4 0.666
SEAHapaxanthic -1.53904 -4.55702 1.67797 688.6 0.330
SEAHapaxanthic 0.18037 -1.72087 2.27258 796.5 0.800
SEAPerennial -0.07601 -0.44810 0.33258 926.0 0.712
SEAPleonanthic -0.14699 -1.14695 0.81452 723.9 0.748
ALTHigh -0.13191 -0.46780 0.22911 725.0 0.452
ALTLow -0.17699 -0.51173 0.10969 772.8 0.292
ALTMid 0.06855 -0.21312 0.41342 882.1 0.684
BIOBoreal 1.74800 -1.18782 4.72759 782.0 0.242
BIOMediterranean-type 2.08074 -0.62533 5.05527 780.1 0.140
BIOSubantarctic 2.17686 -1.13669 5.24883 806.7 0.180
BIOSubarctic 2.39551 -0.91077 5.41454 839.1 0.138
BIOSubtropical/Tropical 2.31132 -0.36795 5.24304 791.5 0.110
BIOTemperate 2.29529 -0.41744 5.18185 795.5 0.104
SEFew-Several 1.86331 -0.57544 4.01647 732.1 0.106
SENumerous 0.20823 -0.14937 0.57547 851.4 0.226
SESeveral 0.66868 -0.13298 1.45685 894.6 0.102
SESingle 0.42408 0.07265 0.80295 872.5 0.022 *
FSZygomorphic 0.01505 -0.22554 0.27481 760.5 0.908
On 11 Aug 2010, at 17:34, Jarrod Hadfield wrote:
Hi Chris,
The model syntax looks reasonable but there seems to be some large posterior means (outside of the 95% credible range). I bet plot(model$VCV) looks pretty horrible too. You need to consider using proper priors in this instance because the chain is getting stuck at zero for long periods of time and generating numerical problems. I tend to use parameter expanded priors more and more as they improve mixing and seem to be only weakly informative. For example: G1=list(V=1, nu=1, alpha.mu=0, alpha.V=1000) .... There is also the possibility that you have complete separation as you have a lot of fixed effects and many levels in the ordinal response - are all 5's for example associated with a single fixed factor, or something like this?
Jarrod
On 11 Aug 2010, at 17:20, Chris Mcowen wrote:
> Sorry about the formatting,
>
> i was not going to use P values for model selection, rather the DIC value
>
> Iterations = 12991
> Thinning interval = 3001
> Sample size = 1000
>
> DIC: 3171.501
>
> G-structure: ~order
>
> post.mean l-95% CI u-95% CI eff.samp
> order 7720 4.023e-13 0.09208 1000
>
> ~fam:fam
>
> post.mean l-95% CI u-95% CI eff.samp
> fam:fam 4092456 2.376e-12 0.02938 1000
>
> R-structure: ~units
>
> post.mean l-95% CI u-95% CI eff.samp
> units 1 1 1 0
>
> Location effects: IUCN ~ STO + BS + PD + FR + END + WO + RG + SEA + ALT + BIO + SE + FS
>
> post.mean l-95% CI u-95% CI eff.samp pMCMC
> (Intercept) 39.065870 -3.510793 2.407406 1000.0 0.776
> STOStorage organ -0.004916 -0.299409 0.230731 757.2 0.946
> BSUnisexual flower 0.211852 -0.131660 0.548879 708.0 0.212
> BSUnisexual plant 0.370895 0.003567 0.817429 770.3 0.070 .
> PDBiotic 0.381261 0.054626 0.724368 774.4 0.040 *
> PDMammalia 26.364377 -2.139720 1.397539 1000 .0 0.724
> FRNon_fleshy_fruit -0.208198 -0.536699 0.083012 964.2 0.202
> ENDNon_endospermous 0.503829 0.200868 0.822120 591.7 0.004 **
> WOWoody -0.203632 -0.565069 0.139240 857.5 0.272
> RGTwo+ -0.052508 -0.250675 0.163811 831.8 0.588
> SEAHapaxanthic -1.344993 -4.504625 1.848373 890.4 0.406
> SEAHapaxanthic 0.223060 -1.590483 2.012970 785.9 0.800
> SEAPerennial -0.097971 -0.460607 0.304681 849.9 0.580
> SEAPleonanthic -0.069756 -0.813837 0.704066 969.4 0.872
> ALTHigh -0.129331 -0.483238 0.200436 1000.0 0.472
> ALTLow -0.171467 -0.514753 0.121200 842.9 0.316
> ALTMid 0.068307 -0.227978 0.379701 814.9 0.660
> BIOBoreal 1.785916 -1.222387 4.769563 860.2 0.254
> BIOMediterranean-type 2.105530 -0.888236 4.786029 817.9 0.156
> BIOSubantarctic 2.214561 -0.888921 5.239470 841.3 0.190
> BIOSubarctic 2.441894 -0.667793 5.677992 849.5 0.142
> BIOSubtropical/Tropical 2.336425 -0.660675 4.899198 928.3 0.124
> BIOTemperate 2.315834 -0.761101 4.826330 809.2 0.132
> SEFew-Several 146.220538 -0.620787 3.933475 1000.0 0.172
> SENumerous 0.206148 -0.117869 0.572987 734.9 0.236
> SESeveral 0.626675 -0.236956 1.456895 881.7 0.134
> SESingle 0.399690 0.030041 0.779923 709.8 0.032 *
> FSZygomorphic 0.032334 -0.215194 0.265597 355.7 0.814
> ---
> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Cutpoints:
> post.mean l-95% CI u-95% CI eff.samp
> cutpoint.traitIUCN.1 0.6593 0.5211 0.793 48.46
> cutpoint.traitIUCN.2 2.4694 2.2952 2.663 41.37
> cutpoint.traitIUCN.3 3.6258 3.4220 3.827 38.02
> cutpoint.traitIUCN.4 4.1156 3.9166 4.341 52.46
> On 11 Aug 2010, at 17:15, Jarrod Hadfield wrote:
>
> Hi,
>
> Could you give summary(model) with the new version (2.05) - it will be easier to see what is going on?
>
> Jarrod
> On 11 Aug 2010, at 17:08, Chris Mcowen wrote:
>
>> Hi Jarrord,
>>
>> I have tried using MCMCglmm, however the posterior distributions of the majority of the fixed factors straddle 0, which i have read is a problem, likely with the priors.
>>
>> HPDintervals - https://files.me.com/chrismcowen/wqq1lu
>>
>> prior=list(R=list(V=1, fix=1), G=list(G1=list(V=1, nu=0), G2=list(V=1, nu=0)))
>>
>> So i am unsure how to interpret the results, as to ascertain the importance of each factor.
>>
>> Unfortunately i don't know enough about baysian statistics or R to alter my model so the interpretations become clearer.
>>
>> An example
>>
>> lower upper
>> (Intercept) -3.510792767 2.40740650
>> STOStorage organ -0.299408836 0.23073133
>> BSUnisexual flower -0.131660436 0.54887912
>> BSUnisexual plant 0.003566637 0.81742862
>> PDBiotic 0.054625970 0.72436838
>> PDMammalia -2.139720264 1.39753939
>>
>>
>>
>> On 11 Aug 2010, at 16:37, Jarrod Hadfield wrote:
>>
>> Hi Chris,
>>
>> It is hard to say as it will depend on the fixed effects. In addition its not clear whether such a situation is diagnostic of a problem. Imagine you just have an intercept which is estimated to be exactly zero. The residuals on the data scale will be either 0.5 or -0.5, but this does not imply the model is wrong.
>>
>> Cheers,
>>
>> Jarrod
>>
>> On 11 Aug 2010, at 15:41, Chris Mcowen wrote:
>>
>>> Thats great thanks,
>>>
>>> But will this work where you have a binary response variable or will the residuals clump around 1 and 0?
>>>
>>> Chris
>>> On 11 Aug 2010, at 15:31, Ben Bolker wrote:
>>>
>>> On 10-08-11 10:21 AM, Chris Mcowen wrote:
>>>> Dear Ben/Rob.
>>>>
>>>>
>>>>> As far as I can tell, the standard advice is simply to look at the predictions of the model, compare them with the data, and try to spot any systematic patterns in the residuals.
>>>>>
>>>>
>>>> I have plotted the residuals of my model - https://files.me.com/chrismcowen/v586vx
>>>>
>>>> I have been made aware that that lmer uses the random effects in its prediction ( Jarrord Hadfield). And this is reflected in the residual plot with the the long lines of equal residuals all belonging to the same family - i.e 200 - 600 is the orchid family and 650-100 is the grass family.
>>>>
>>>> So is there a work around with a glmm?
>>>>
>>>>
>>>>
>>>> Thanks
>>>>
>>>> Chris
>>>>
>>>>
>>>
>>> If you want to do population-level predictions from a GLMM (i.e. setting all random effects to zero), the basic recipe is to (1) construct a model (design) matrix for the desired sets of predictor variables (if you want to the predict the observed data rather than some other set, you can just extract the model matrix from the fitted object); (2) multiply it by the vector of fixed effect coefficients; (3) transform it back to the scale of the observations with the inverse link function. There's an example on p. 6 of http://glmm.wdfiles.com/local--files/examples/Owls.pdf ...
>>>
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>> _______________________________________________
>>> 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.
>>
>> _______________________________________________
>> 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.
>
> _______________________________________________
> 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.
More information about the R-sig-mixed-models
mailing list