Skip to content
Prev 8624 / 20628 Next

MCMCglmm: priors for ordinal regression

Dear Massimo,

The function below (prior.scale) may be useful. It returns a  
covariance matrix that can be passed to prior$B$V.  If the predictors  
had been scaled according to the procedure outlined in Gelman et al.  
(2006) then this would induce an identity prior covariance matrix on  
the new regression coefficients (i.e. they're iid a priori). Gelman  
recommends a scaled Cauchy or scaled t prior on these new regression  
ceofficients. I don't think this is possible in MCMCglmm* but a normal  
with variance equal to the sum of the variance components + 1 (probit)  
  pi^2/3 (logit) I think is not too unreasonable - but it does put  
less weight on very extreme probabilities than the Cauchy and t. To  
achieve this just multiply the output of prior.scale by the variance  
you want to use.

Bear in mind that the inputs are scaled, not the columns of the design  
matrix - this means that interactions are penalised more than main  
effects. In my field just-so stories involving 'significant'  
interactions  are often used to resuscitate a `failed' experiment, so  
penalising them by default may be no bad thing.

Cheers,

Jarrod

*setting a weaker prior on the non-identified residual variance rather  
than fixing it may achieve a t-prior - but really not sure

Gelman et al. A WEAKLY INFORMATIVE DEFAULT PRIOR DISTRIBUTION FOR
LOGISTIC AND OTHER REGRESSION MODELS The Annals of Applied Statistics
2008, Vol. 2, No. 4, 1360?1383

# formula = the fixed formula in the model
# data = the data set

prior.scale<-function(formula, data){
   X1<-model.matrix(formula, data)
   X2<-get_all_vars(formula, data)
   X2<-as.data.frame(lapply(X2, function(x){if(is.numeric(x)){scale(x,  
scale=sd(x)*2)}else{x}}))
   X2<-model.matrix(formula, data=X2)
   X2[,-1]<-apply(X2[,-1], 2,  
function(x){if(any(!x%in%c(0,1))){x}else{scale(x,  
center=sum(x)/length(x), scale=1)}})
   crossprod(t(solve(t(X1)%*%X1, t(X1)%*%X2)))
}


Quoting "m.fenati at libero.it" <m.fenati at libero.it> on Mon, 16 Jul 2012  
11:04:45 +0200 (CEST):