[R-sig-ME] R: Re: R: Re: pMCMC and HPD in MCMCglmm
m.fenati at libero.it
m.fenati at libero.it
Thu Aug 25 18:52:43 CEST 2011
Hi,
I was sure to have fixed indipendent effect between intercept and sex!!
Yes, I have mistyped "()".......
thanks a lot!
ps MCMCglmm is a great package!
Massimo
------------------
Massimo Fenati
DVM
Padova - Italy
>----Messaggio originale----
>Da: j.hadfield at ed.ac.uk
>Data: 25/08/2011 18.17
>A: "m.fenati at libero.it"<m.fenati at libero.it>
>Cc: <r-sig-mixed-models at r-project.org>
>Ogg: Re: R: Re: [R-sig-ME] pMCMC and HPD in MCMCglmm
>
>Hi,
>
>You have specified a strong prior correlation (~0.5) between the
>intercept and sex effect:
>
>> prior$B
>$mu
>[1] 0 0
>
>$V
> [,1] [,2]
>[1,] 6.289868 3.289868
>[2,] 3.289868 6.289868
>
>Since the data have 50 0's and one 1 the data give a lot of weight to
>the intercept being negative. As you believe, the data are not very
>informative about sex-differences in this example but with your prior
>specification we expect a priori that the intercept and sexM effect
>are positively correlated. Hence, the "signifcant" negative sex effect.
>
>Perhaps you had mistyped the prior specification, and had intended
>
>V=diag(2)*(1+pi^2/3) rather than V=diag(2)+pi^2/3 ?
>
>as in Section 2.6 of the CourseNotes. Treating the effects as
>independent in the prior gives results closer to what you would hope
>for.
>
>Jarrod
>
>
>
>Quoting "m.fenati at libero.it" <m.fenati at libero.it> on Thu, 25 Aug 2011
>17:54:21 +0200 (CEST):
>
>> Hi Jarrod,
>> In the past example, where HPD and pMCMC were slightly different, I tested
an
>> extreme dataset: 1 positive event on 51 sample of 34 animals. In this
>> circumstance, even assuming not repeated data, the posterior distribution
of
>> the sex beta coefficient MCMC estimates seems to suggest a possible effect
of
>> “sex” on the response “dis”. But if I perform the analysis on the
>> same dataset
>> under frequentist approach this fails (using glm for perfect separation) or
>> returns with high p-value (using aalysis of frequency table via fisher
exact
>> test). See the following example:
>>
>> sex<-c(rep("F",21),rep("M",30))
>> dis<-c(1,rep(0,50))
>> dat<-data.frame(sex,dis)
>> prior<-list(R=list(V=1,fix=1),G=list(G1=list(V=1,nu=0.002)),B=list(mu=c(rep
>> (0,2)),V=diag(2)*3+pi^2/3))
>> m.1<-MCMCglmm(dis~sex,slice=T,prior=priorS,data=dat,nitt=800000,thin=100,
>> burnin=250000,family="categorical",verbose=FALSE)
>> summary(m.1)
>>
>> fisher.test(dat$dis,dat$sex)
>>
>> summary(glm(dis~sex,data=dat,family=binomial))
>>
>>
>> How can I interpret the differences between Bayesian (MCMCglmm) and
>> Frequentist approaches in these circumstances?
>>
>> Sorry for the basic question, but I am new in Bayesian world!
>>
>> Thanks
>>
>> Massimo
>>
>>
>>
>>
>> -----------------------
>> Massimo Fenati
>> DVM
>> Padova - Italy
>>
>>
>>
>>
>>> ----Messaggio originale----
>>> Da: j.hadfield at ed.ac.uk
>>> Data: 24/08/2011 18.11
>>> A: "m.fenati at libero.it"<m.fenati at libero.it>
>>> Cc: <r-sig-mixed-models at r-project.org>
>>> Ogg: Re: [R-sig-ME] pMCMC and HPD in MCMCglmm
>>>
>>> Hi,
>>>
>>> pMCMC is the two times the smaller of the two quantities: MCMC
>>> estimates of i) the probability that a<0 or ii) the probability that
>>> a>0, where a is the parameter value. Its not a p-value as such, and
>>> better ways of obtaining Bayesian p-values exist.
>>>
>>> HPDinterval finds the closest points (c and d) for which Fa(d)-Fa(c) =
>>> 0.95 (If prob=0.95 in HPDinterval) and Fa is the empirical cumulative
>>> distribution of a.
>>>
>>> Cheers,
>>>
>>> Jarrod
>>>
>>>
>>>
>>>
>>>
>>> On 24 Aug 2011, at 16:14, m.fenati at libero.it wrote:
>>>
>>>> Hi Jarrod,
>>>> thanks for your answer, but I have again a lot of confusion. If
>>>> possible,
>>>> could you explain to me the definition of pMCMC?
>>>> Maybe, knowing the right definition of pMCMC I will be able to
>>>> understand
>>>> completely your answer.
>>>>
>>>> Thank a lot!
>>>>
>>>> Massimo
>>>>
>>>> -----------------------
>>>> Massimo Fenati
>>>> DVM
>>>> Padova - Italy
>>>>
>>>>
>>>>
>>>>> ----Messaggio originale----
>>>>> Da: j.hadfield at ed.ac.uk
>>>>> Data: 24/08/2011 13.24
>>>>> A: "m.fenati at libero.it"<m.fenati at libero.it>
>>>>> Cc: <ndjido at gmail.com>, <r-sig-mixed-models at r-project.org>
>>>>> Ogg: Re: [R-sig-ME] pMCMC and HPD in MCMCglmm
>>>>>
>>>>> Hi Massimo,
>>>>>
>>>>> They only need to be slightly skewed (even up to Monte Carlo error
>>>>> probably) - conclusions drawn from HPDinterval and pMCMC are in
>>>>> reality almost identical in your example, it is the consequences of
>>>>> the (arbitrary) distinction between <0.05 and >0.05 that makes them
>>>>> "feel" different. Lets say we used the cutoff <0.06 and >0.06. Does
>>>>> HPDinterval(m1$Sol[,3], prob=0.94) overlap zero? If not then
>>>>> HPDinterval and pMCMC "agree" with respect to which side of the
>>>>> cutoff
>>>>> the probability lies ? It may make us happier, but it shouldn't.
>>>>>
>>>>> Jarrod
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On 24 Aug 2011, at 11:45, m.fenati at libero.it wrote:
>>>>>
>>>>>> The posterior distribution seem to be only slightly skewed.
>>>>>> However the question remains: what is the sense of the discrepancy
>>>>>> between HPD
>>>>>> and pMCMC?
>>>>>>
>>>>>> Thanks
>>>>>>
>>>>>> Massimo
>>>>>>
>>>>>>
>>>>>>
>>>>>> ----Messaggio originale----
>>>>>> Da: ndjido at gmail.com
>>>>>> Data: 24/08/2011 11.43
>>>>>> A: "m.fenati at libero.it"<m.fenati at libero.it>
>>>>>> Cc: <r-sig-mixed-models at r-project.org>
>>>>>> Ogg: Re: [R-sig-ME] pMCMC and HPD in MCMCglmm
>>>>>>
>>>>>> Check your posterior distributions, the one corresponding to GENDER
>>>>>> seems to
>>>>>> be skewed.
>>>>>> Ardo.
>>>>>> On Wed, Aug 24, 2011 at 11:33 AM, m.fenati at libero.it <m.fenati at libero.
it
>>>>>>>
>>>>>> wrote:
>>>>>> As suggested by Ben Bolker, I re-post the following question in this
>>>>>> list.
>>>>>> Thanks
>>>>>>
>>>>>>> Dear R users,
>>>>>>> I’d like to pose aquestion about pMCMC and HDP.
>>>>>>> I have performed a mixed logistic regression by MCMCglmm (a very
>>>>>>> good
>>>>>> package)
>>>>>>> obtaining the following results:
>>>>>>>
>>>>>>> Iterations = 250001:799901
>>>>>>> Thinning interval = 100
>>>>>>> Sample size = 5500
>>>>>>>
>>>>>>> DIC: 10.17416
>>>>>>>
>>>>>>> G-structure: ~ID_an
>>>>>>>
>>>>>>> post.mean l-95% CI u-95% CIeff.samp
>>>>>>> ID_an 0.7023 0.0001367 3.678 2126
>>>>>>>
>>>>>>> R-structure: ~units
>>>>>>>
>>>>>>> post.mean l-95% CIu-95% CI eff.samp
>>>>>>> units 1 1 1 0
>>>>>>>
>>>>>>> Location effects: febbreq~ as.factor(sex)
>>>>>>>
>>>>>>> post.mean l-95% CIu-95% CI eff.samp pMCMC
>>>>>>> (Intercept) -3.6332 -5.6136 -1.7719 3045 <2e-04 ***
>>>>>>> as.factor(sex)M -2.9959 -6.0690 0.1969 2628 0.0455 *
>>>>>>> ---
>>>>>>> Signif. codes: 0 ‘***’0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>>>>>>>
>>>>>>>
>>>>>>> As you can see, pMCMC for gender is just less than 5%, but the
>>>>>>> credible
>>>>>>> interval (HPD) is wide and includes the 0 value.
>>>>>>> How can I interpret these different results?
>>>>>>>
>>>>>>> Thank you in advance
>>>>>>>
>>>>>>> Massimo
>>>>>>>
>>>>>>> -----------------------
>>>>>>> Massimo Fenati
>>>>>>> DVM
>>>>>>> Padova - Italy
>>>>>>
>>>>>> _______________________________________________
>>>>>> 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.
>>>>>
>>>>>
>>>>
>>>>
>>>>
>>>
>>>
>>> --
>>> 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