# R/NIgroupscontrib.R In prefmod: Utilities to Fit Paired Comparison Models for Preferences

```# contribution of single cov group to likelihood
NIgroupscontrib<-function(obj,lambda,X,nobj,ENV)
{
# set up design matrix by expanding basic matrix X with covariates
if (!is.null(names(obj\$cov)))       # only in case of covariates
X<-do.call("cbind", lapply(1:length(obj\$cov), function(i) X *obj\$cov[i])) # design matrix for covariates

## MCAR
if (! ENV\$NI) {

if(ENV\$undec) X<-cbind(X,ENV\$U)

if(ENV\$ia) X<-cbind(X,ENV\$XI)

######patprob was a function once
patt <- exp(X %*% lambda)
sumpatt<-sum(patt)
p.patt<-patt/sumpatt
####
# get likelihood contribution for single cov group
ll.bl<-lapply(obj[1:(length(obj)-1)],blcontrib,p.patt)   # last element of obj is covariate structure

## MNAR
} else if(ENV\$NItest) {   # for testing ignorable missing

## first complete cases
NIcovs<-c(1,0)
XM<-do.call("cbind", lapply(1:2, function(i) X * NIcovs[i])) # design matrix for ign missings

if(ENV\$undec) XM<-cbind(XM,ENV\$U)

if(ENV\$ia) XM<-cbind(XM,ENV\$XI)

patt <- exp(XM %*% lambda)
sumpatt<-sum(patt)
p.patt<-patt/sumpatt

ll.bl<-lapply(obj[1],blcontrib,p.patt)

## now incompletes
NIcovs<-c(1,1)
XM<-do.call("cbind", lapply(1:2, function(i) X * NIcovs[i])) # design matrix for ign missings

if(ENV\$undec) XM<-cbind(XM,ENV\$U)

if(ENV\$ia) XM<-cbind(XM,ENV\$XI)

patt <- exp(XM %*% lambda)
sumpatt<-sum(patt)
p.patt<-patt/sumpatt

ll.mbl<-lapply(obj[2:(length(obj)-1)],blcontrib,p.patt)   # last element of obj is covariate structure
ll.bl<-c(ll.bl,ll.mbl)

} else { # nonresponse models

likparts<-lapply(obj[1:(length(obj)-1)],NIblcontrib,lambda,X,nobj,ENV)
ll<-sum(unlist(lapply(likparts,function(x)x\$ll)))          # now normalised part
fl<-sum(unlist(lapply(likparts,function(x)x\$fl))) # loglik for saturated model

ll.bl<-list(list(ll=ll,fl=fl))

# nur mehr zwei elemente in liste
ll.bl

############################################################################################
# this part is now in NIblcontrib
# ll.bl<-list()             # some initialisation                                          #
# ncomp<-ENV\$ncomp                                                                         #
# patt.all<-list()                                                                         #
# summ.patt.all<-0                                                                         #
# summ.cnts<-0                                                                             #
# ll1<-0                                                                                   #
# fl<-0                                                                                    #
#                                                                                          #
# for (block in 1:(length(obj)-1)) {                                                       #
#                                                                                          #
#     naidx<-1-as.numeric(obj[[block]]\$notnaidx)  # R=1 missing, R=0 observed              #
#                                                                                          #
#     R<-matrix(rep(naidx,2^ncomp),nrow=2^ncomp,byrow=T)                                     #
#     RBstar<-R %*%abs(pcdesign(nobj))  # alpha_i + alpha_j                                #
#     #RBstar<-R %*%(pcdesign(4))    # alpha_i - alpha_j                                   #
#                                                                                          #
#     YRBstar<-do.call(cbind,lapply(1:(nobj),function(i) RBstar[,i]*ENV\$Y[,i])) # betas    #
#                                                                                          #
#     ###################################################################                  #
#     #####    nonresponse model alpha_i+alpha_j/alpha_i-alpha_j ########                  #
#     XX<-X                              # only lambdas                                    #
#     if (ENV\$MISalpha) XX<-cbind(X,RBstar[,1:nobj])          # lambdas, alpha_i           #
#     if (ENV\$MISbeta)  XX<-cbind(X,RBstar[,1:nobj],YRBstar)  # lambdas, alpha_i. beta_i   #
#     ###################################################################                  #
#                                                                                          #
#                                                                                          #
#     ####################################################################                 #
#     #################    nonresponse model alpha_ij ????    ############                 #
#     #XX<-cbind(X,R)                                                                      #
#     #XX<-cbind(X,R[,1:4],RX[,1:2])                                                       #
#     #XX<-cbind(X,RX[,1:6])                                                               #
#     ####################################################################                 #
#                                                                                          #
#     # add dependencies if ia=TRUE                                                        #
#     if(ENV\$ia) XX<-cbind(XX,ENV\$XI)                                                      #
#                                                                                          #
#     patt.all[[block]] <- exp(XX %*% lambda)                                              #
#     patt<-tapply(patt.all[[block]],obj[[block]]\$s,sum)                                   #
#     ll1<-ll1+sum(obj[[block]]\$counts*log(patt))                                          #
#                                                                                          #
#     summ.patt.all<-summ.patt.all+sum(patt.all[[block]])                                  #
#     summ.cnts<-summ.cnts+sum(obj[[block]]\$counts)                                        #
#                                                                                          #
#     p.cnts<-obj[[block]]\$counts/sum(obj[[block]]\$counts)                                 #
#     fl<-fl+sum(log(p.cnts[p.cnts>0])*obj[[block]]\$counts[obj[[block]]\$counts>0])         #
#                                                                                          #
# }                                                                                        #
#                                                                                          #
# ll<-ll1-summ.cnts*log(summ.patt.all)                                                     #
# ll.bl[[1]]<-list(ll=ll,fl=fl)                                                            #
#                                                                                          #
# # nur mehr zwei elemente in liste: ll und fl mit 1 ebene dazwischen (das waeren Covs)    #
# ll.bl                                                                                    #
############################################################################################

}
}
```

## Try the prefmod package in your browser

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

prefmod documentation built on May 2, 2019, 4:59 p.m.