# R/GMI.R In hmmm: Hierarchical Multinomial Marginal Models

```GMI<-function(freq,marg,lev,names,mflag="M"){
marg<-marg.list(marg, mflag="m")
interf<-hmmm.model(marg=marg,lev=lev,names=names)
interfun<-interf\$functions\$L.fct
freq<-as.matrix(freq)
gmi<-apply(freq,2,interfun)
x<-hmmm.model.summary(interf,printflag=FALSE)
#print(interf,quote=FALSE)
na<-rep(x[,2],times=x[,6])
for(i in 1:dim(x)[1]){
na[x[i,7]:x[i,8]]<-paste(x[i,2],1:x[i,6],sep="")}
rownames(gmi)<-na
colnames(gmi)<-paste("F",1:dim(gmi)[2],sep="")
gmi<-list(marginals=interf,gmi=gmi)
}

inv_GMI<-function(etpar,mod,start=rep(0,prod(mod\$modello\$livelli))){
#given a vector of complete hyerarchical generalized marginal interactions  computes the joint probabilities
#etpar=vettorog generalized marginal interactions
#modello marginale oggetto classe hmmmmod
#vector of starting values for the loglinear parameters zero vector is the default
#if(is.null(start)) start<-rep(0,prod(mod\$modello\$livelli))
interf<-mod\$functions\$L.fct
interfder<-mod\$functions\$derLt.fct

myfun<-function(x){x<-as.matrix(exp(Zlist\$DMAT%*%x))
x<-x/sum(x)
interf(x)-etpar}
myder<-function(x){
x<-as.matrix(exp(Zlist\$DMAT%*%x))
x<-x/sum(x)
(
t( t(Zlist\$DMAT)%*%
diag(c(x+1e-08)) %*%interfder(x)
)) }
r<-nleqslv(rep(0,15),myfun,myder)\$x
r<-as.matrix(exp(Zlist\$DMAT%*%r))
r<-r/sum(r)
}
```

## Try the hmmm package in your browser

Any scripts or data that you put into this service are public.

hmmm documentation built on May 2, 2019, 12:27 p.m.