[R-sig-ME] predictions for MCMCglmm

Fri Mar 22 18:54:56 CET 2013

Hi,

In the original code I sent you I marginalised the random effects. If  
you use predict on your new model with the missing data you should get  
the same answer if you have marginal=your_model$Random$formula.

Cheers,

Jarrod

Quoting "Antonio P. Ramos" <ramos.grad.student at gmail.com> on Fri, 22  
Mar 2013 09:00:56 -0700:

> Thanks Jarrod, but isn't it true that the random effects I get for
> predictions will be different in each scenario?
>
>
> On Fri, Mar 22, 2013 at 2:54 AM, Jarrod Hadfield <j.hadfield at ed.ac.uk>wrote:
>
>> cc-ing back to the list ...
>>
>> Hi,
>>
>> Yes you can do it that way, and it should give the same answer.  It will
>> slow down the MCMCing in terms of time per iteration and mixing though.
>>
>> Cheers,
>>
>> Jarrod
>>
>>
>>
>> Quoting "Antonio P. Ramos" <ramos.grad.student at gmail.com> on Thu, 21 Mar
>> 2013 18:18:41 -0700:
>>
>>   Hi all,
>>>
>>> I have a follow up question: what is different between getting predictions
>>> used this code - which works fine - or, instead, by adding rows in the
>>> data
>>> matrix with NA values in the response and the desired covariates in the X
>>> matrix, such as often done by people using JAGS/BUGS? I am not sure at
>>> all,
>>> but it seems to be that by running this code, I am assuming that the
>>> predictions I want have some kind of random structure that might not be
>>> the
>>> case - by marginalizing over them all. Thanks for any further
>>> clarification. Antonio Pedro.
>>>
>>>
>>> On Tue, Mar 19, 2013 at 11:35 AM, Antonio P. Ramos <
>>> ramos.grad.student at gmail.com> wrote:
>>>
>>>  sorry: the full code here
>>>> > pred.data.pred <-
>>>> data.frame(maternal_age_c=rep(**20,25),wealth=rep("Lowest quintile",25),
>>>> +
>>>>  birth_year=1971:1995,birth_**order=rep(1,25),
>>>> +
>>>>  sex=rep("Female",25),**residence=rep("Rural",25),
>>>> +                                    maternal_educ=rep("Primary",**25))
>>>> > # creating predictions "by hand"
>>>> > X <- model.matrix(~ maternal_age_c + I(maternal_age_c^2)  +
>>>> as.factor(birth_year) +
>>>> +                   as.factor(birth_order) + residence + sex + wealth,
>>>> data=pred.data.pred)
>>>> Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
>>>>   contrasts can be applied only to factors with 2 or more levels
>>>>
>>>>
>>>>
>>>>
>>>> On Tue, Mar 19, 2013 at 11:35 AM, Antonio P. Ramos <
>>>> ramos.grad.student at gmail.com> wrote:
>>>>
>>>>  Thanks for you reply Jarrod. The problem is that model.matrix doesn't
>>>>> allow me to created the predication I need as factors cannot have just
>>>>> one
>>>>> level. Any ways around data? thanks a bunch
>>>>>
>>>>>
>>>>> pred.data.pred <- data.frame(maternal_age_c=rep(**
>>>>> 20,25),wealth=rep("Lowest
>>>>> quintile",25),
>>>>>
>>>>>  birth_year=1971:1995,birth_**order=rep(1,25),
>>>>>
>>>>>  sex=rep("Female",25),**residence=rep("Rural",25),
>>>>>                                    maternal_educ=rep("Primary",**25))
>>>>> # creating predictions "by hand"
>>>>> X <- model.matrix(~ maternal_age_c + I(maternal_age_c^2)  +
>>>>> as.factor(birth_year) +
>>>>>                   as.factor(birth_order) + residence + sex + wealth,
>>>>> data=pred.data.pred)
>>>>>
>>>>>
>>>>>
>>>>> On Tue, Mar 19, 2013 at 2:27 AM, Jarrod Hadfield <j.hadfield at ed.ac.uk
>>>>> >wrote:
>>>>>
>>>>>  Hi Antonio,
>>>>>>
>>>>>> With (simple) random effects marginalised:
>>>>>>
>>>>>> X<-model.matrix(~ maternal_age_c + I(maternal_age_c^2)  +
>>>>>> as.factor(birth_year) + residence + sex + wealth, data=newdata)
>>>>>>
>>>>>> V<-rowSums(glm.MC.2$VCV)
>>>>>>
>>>>>> beta<-glm.MC.2$Sol
>>>>>>
>>>>>> c2 <- (16 * sqrt(3)/(15 * pi))^2
>>>>>>
>>>>>> pred<-t(plogis(t(beta%*%t(X)/****sqrt(1+c2*V))))
>>>>>>
>>>>>>
>>>>>> pred[i,j] is the prediction for the jth new data point for the ith MCMC
>>>>>> sample. colSums(pred)  should be equivalent to the output from
>>>>>> predict.MCMCglmm.
>>>>>>
>>>>>> Cheers,
>>>>>>
>>>>>> Jarrod
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> Quoting "Antonio P. Ramos" <ramos.grad.student at gmail.com> on Mon, 18
>>>>>> Mar 2013 20:04:07 -0700:
>>>>>>
>>>>>>  Hi all.
>>>>>>
>>>>>>>
>>>>>>> As far as I can tell newdata is still not implemented for this nice
>>>>>>> package. Thus I wonder what would be the best way to get predictions
>>>>>>> "by
>>>>>>> hand". My model is actually very simple. Still I need to marginalize
>>>>>>> the
>>>>>>> random effects. Any hints? Thanks in advance, Antonio Pedro.
>>>>>>>
>>>>>>>
>>>>>>> glm.MC.2 <- MCMCglmm(mortality.under.2 ~ maternal_age_c +
>>>>>>> I(maternal_age_c^2)  +
>>>>>>>                        as.factor(birth_year) + residence +
>>>>>>>                        sex + wealth,
>>>>>>>                      nitt=20000, thin=10, burnin=1000,
>>>>>>>                      random= ~CASEID, prior=prior.2,data=egypt2,
>>>>>>> family='categorical')
>>>>>>>
>>>>>>>         [[alternative HTML version deleted]]
>>>>>>>
>>>>>>> ______________________________****_________________
>>>>>>> R-sig-mixed-models at r-project.****org <R-sig-mixed-models at r-project.**
>>>>>>> org <R-sig-mixed-models at r-project.org>>mailing list
>>>>>>> https://stat.ethz.ch/mailman/****listinfo/r-sig-mixed-models<https://stat.ethz.ch/mailman/**listinfo/r-sig-mixed-models>
>>>>>>> <h**ttps://stat.ethz.ch/mailman/**listinfo/r-sig-mixed-models<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.