data(speed)
# 2-state model on rt and corr from speed data set with Pacc as covariate on the transition matrix
# starting values for the transition pars (without those EM does not get off the ground)
set.seed(1)
tr=runif(6)
trst=c(tr[1:2],0,tr[3:5],0,tr[6])
mod1 <- depmix(list(rt~1,corr~1),data=speed,transition=~Pacc,nstates=2,family=list(gaussian(),multinomial()),
trstart=trst)
# fit the model
fmod1 <- fit(mod1)
fmod1 # to see the logLik and optimization information
# to see the parameters
summary(fmod1)
# NOTE: this requires Rdonlp2 package to be installed
# FIX SOME PARAMETERS
# get the starting values of this model to the optimized
# values of the previously fitted model to speed optimization
pars <- c(unlist(getpars(fmod1)))
# constrain the initial state probs to be 0 and 1
# also constrain the guessing probs to be 0.5 and 0.5
# (ie the probabilities of corr in state 1)
# change the ones that we want to constrain
pars[2]=+Inf # this means the process will always start in state 2
pars[14]=0 # the corr parameters in state 1 are now both 0, corresponding the 0.5 prob
mod2 <- setpars(mod1,pars)
# fix the parameters by setting:
free <- c(0,0,rep(c(0,1),4),1,1,0,0,1,1,0,1)
# fit the model
fmod2 <- fit(mod2,fixed=!free)
# likelihood ratio insignificant, hence fmod2 better than fmod1
llratio(fmod1,fmod2)
# NOW ADD SOME GENERAL LINEAR CONSTRAINTS
# set the starting values of this model to the optimized
# values of the previously fitted model to speed optimization
pars <- c(unlist(getpars(fmod2)))
mod3 <- setpars(mod2,pars)
# start with fixed and free parameters
conpat <- c(0,0,rep(c(0,1),4),1,1,0,0,1,1,0,1)
# constrain the beta's on the transition parameters to be equal
conpat[4] <- conpat[8] <- 2
conpat[6] <- conpat[10] <- 3
fmod3 <- fit(mod3,equal=conpat)
llratio(fmod2,fmod3)
# above constraints can also be specified using the conrows argument as follows
conr <- matrix(0,2,18)
# pars 4 and 8 have to be equal, otherwise stated, their diffence should be zero
conr[1,4] <- 1
conr[1,8] <- -1
conr[2,6] <- 1
conr[2,10] <- -1
# note here that we use the fitted model fmod2 as that has appropriate
# starting values
fmod3b <- fit(fmod2,conrows=conr,fixed=!free) # using free defined above
data(balance)
# four binary items on the balance scale task
instart=c(0.5,0.5)
set.seed(1)
respstart=runif(16)
# note that ntimes argument is used to make this a mixture model
mod4 <- mix(list(d1~1,d2~1,d3~1,d4~1), data=balance, nstates=2,
family=list(multinomial(),multinomial(),multinomial(),multinomial()),
respstart=respstart,instart=instart)
fmod4 <- fit(mod4)
# add age as covariate on class membership by using the prior argument
instart=c(0.5,0.5,0,0) # we need the initial probs and the coefficients of age
set.seed(2)
respstart=runif(16)
mod5 <- mix(list(d1~1,d2~1,d3~1,d4~1), data=balance, nstates=2,
family=list(multinomial(),multinomial(),multinomial(),multinomial()),
instart=instart, respstart=respstart, prior=~age, initdata=balance)
fmod5 <- fit(mod5)
# check the likelihood ratio; adding age significantly improves the goodness-of-fit
llratio(fmod5,fmod4)Run the code above in your browser using DataLab