[R-sig-ME] bivariate logistic/probit with hierarchical structure

[Jarrod Hadfield]

Hi Arthur,

Example code for fitting bivariate logistic (m1) / probit (m2) models  
with hierarchical structures are given below. The function will not  
give you the likelihood but it will sample from the posterior.  
rcov=~cor(trait) fixes the variances of E1 and E2 at 1 but estimates  
the correlation.

I believe MLwin is capable of doing this also and is free to academic users.

Cheers,

Jarrod



m1<-MCMCglm(cbind(Y1,Y2)~trait-1, random=~us(trait):rfactor,  
rcov=~cor(trait), family=c("categorical", "categorical"),  
data=your.data, prior=your.prior)

m2<-MCMCglm(cbind(Y1,Y2)~trait-1, random=~us(trait):rfactor,  
rcov=~cor(trait), family=c("ordinal", "ordinal"), data=your.data,  
prior=your.prior)



Quoting Arthur Charpentier <arthur.charpentier at gmail.com> on Fri, 6  
May 2011 14:10:36 -0400:

> I would like to fit a bivariate probit model, i.e. I observe Y=(Y1,Y2),
> where Y1 and Y2 take values 0 or 1,
> and I assume that there is a latent bivariate Gaussian vector Z=(Z1,Z2), and
> * if (Z1<s1,Z2<s2) then Y=(0,0)
> * if (Z1>s1,Z2<s2) then Y=(1,0)
> * if (Z1<s1,Z2>s2) then Y=(0,1)
> * if (Z1>s1,Z2>s2) then Y=(1,1)
> The point here is that I assume that the noise E=(E1,E2) can be correlated
> between the two components...
> I have some covariates X1,...Xk, and I assume that there is a hierarchical
> structure (and random effects), i.e. Y_i is actually a Y_{i,j} where j is a
> subgroup index.
> So this yield two correlations levels,
> - within subgroup correlation between Y1_i's (if they belong to the same
> subgroup)
> - between equations, i.e. Y1_i and Y2_i
> I have found packages for logistic hierarchical regressions (eg *lmer*), and
> others for bivaraite logistic (eg *zelig*)
> without the hierarchical structure, I can compute the log likelihood and
> optimize it (and it works)
> but with the hierarchical structure, I can still get the likelihood of the
> latent vector (written in dimension 2n), but I have trouble when I want to
> go from the latent model to the observed one since Poincarré equation
> involves 2^n summations (e.g. if I observe ((0,1),(1,0),(1,1)) I have to get
> the probability that ((Y1<s.Y2>s),(Y1>s.Y2<s),(Y1>s.Y2>s)) which is related
> to cdf's with a lot a summations) which cannot be computed here
> So I was wondering if someone has already computed such a thing ?
> thanks
>
> 	[[alternative HTML version deleted]]
>
>



-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.