jEM <-
function(n,Y,idx,nj,p,Covars,Alphaopt,Wopt,Sigopt,Iq,repV,eps,Np,Nvp,C2p,Ip,picon,V){
Covarsj=Covars[,-idx[n]]
J=Y[-idx[n],]
muhat=colMeans(J)
tJc=t(J)-muhat
Jc=t(tJc)
sumdiSj=mdiag(fM(tJc,Jc)/nj,p)
Wj=Wopt
Sigj=Sigopt
iM=chol.inv(crossprod(Wj)+Sigj*Iq)
iMtW=tcrossprod(iM,Wj)
Alphaj=Alphaopt
delj=Alphaj%*%Covarsj
Xproj=crossprod(Covarsj,chol.inv(tcrossprod(Covarsj)))
ll=lla=repV
tol=eps+1; v=0
while(tol>eps){
v=v+1
u=fM(iMtW,tJc)+iM%*%(Sigj*delj)
tu=t(u)
sumEuu=nj*Sigj*iM+u%*%tu
Alphaj=u%*%Xproj
delj=Alphaj%*%Covarsj
Wj=(tJc%*%tu)%*%chol.inv(sumEuu)
ww=crossprod(Wj)
MLESigj=(nj*sumdiSj+sum(ww*sumEuu)-2*sum(fM(Jc,Wj)*tu))/(nj*p)
Sigj=(Np*MLESigj+C2p)/Nvp
SIGj=Sigj*Iq
M=ww+SIGj
iM=chol.inv(M)
iMtW=tcrossprod(iM,Wj)
InvS=(Ip-Wj%*%iMtW)/Sigj
tJcc=tJc-Wj%*%delj
MhalDist=colSums(tJcc*fM(InvS,tJcc))
logdetD=p*log(Sigj)+log(det(M/Sigj))
ll[v]=sum(-picon-.5*logdetD-.5*MhalDist)
converge=Aitkenf(ll,lla,v,eps)
tol=converge[[1]]
lla[v]=converge[[2]]
if(v==V){tol=eps-1}}
Alphaj}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.