SCM

SCM Repository

[depmix] View of /trunk/depmixNew-test2.R
ViewVC logotype

View of /trunk/depmixNew-test2.R

Parent Directory Parent Directory | Revision Log Revision Log


Revision 68 - (download) (annotate)
Fri Mar 7 14:14:55 2008 UTC (11 years, 10 months ago) by ingmarvisser
File size: 2614 byte(s)
Deleted all redundant files
# 
# Started by Ingmar Visser & Maarten Speekenbrink, march 2008
# 
# Usage: go to trunk directory and source this file in R, if the program
# still works it should return TRUE at every test (or make immediate sense
# otherwise)


setwd("/Users/ivisser/Documents/projects/depmixProject/depmixNew/rforge/depmix/trunk/")

source("R/responses.R")
source("R/depmix.R")
source("R/depmix.fitted.R")
source("R/llratio.R")
source("R/EM.R")
source("R/lystig.R")
source("R/fb.R")

load("data/speed.Rda")

set.seed(1)
mod <- depmix(rt~1,data=speed,nstates=2,trstart=runif(4))
logLik(mod)
mod1 <- fit(mod)

ll <- logLik(mod1)




# 
# test model with EM optimization, no covariates
# 

trstart=c(0.899,0.101,0.084,0.916)
instart=c(0.5,0.5)
resp <- c(5.52,0.202,0.472,0.528,6.39,0.24,0.098,0.902)

mod <- depmix(list(rt~1,corr~1),data=speed,nstates=2,family=list(gaussian(),multinomial()),trstart=trstart)
# 	respstart=resp,trstart=trstart,instart=instart)

logLik(mod)

mod1 <- fit(mod)
ll <- logLik(mod1)


# 
# Test optimization using Rdonlp2
# 

trstart=c(0.899,0.101,0,0.01,0.084,0.916,0,0)
instart=c(0.5,0.5)
resp <- c(5.52,0.202,0.472,0.528,6.39,0.24,0.098,0.902)

mod <- depmix(list(rt~1,corr~1),data=speed,transition=~Pacc,nstates=2,family=list(gaussian(),multinomial()),
	respstart=resp,trstart=trstart,instart=instart)

logLik(mod)

mod1 <- fit(mod)
ll <- logLik(mod1)

cat("Test 1: ", ll, "(loglike of speed data with covariate, hopefully better than 300.5701) \n")

# 
# now constrain the initial state probs to be 0 and 1
# also constrain the guessing probs to be 0.5 and 0.5 
# 

# set the starting values of this model to the optimized 
# values of the previously fitted model to speed optimization
pars <- c(unlist(getpars(mod1)))
# set some to other values, ie the ones that we want to constrain
pars[2]=+Inf
pars[14]=0

mod <- setpars(mod,pars)

# fix some parameters
free <- c(0,0,rep(c(0,1),4),1,1,0,0,1,1,0,1)

mod2 <- fit(mod,fixed=!free)

llratio(mod1,mod2)

cat("Test 2: likelihood ratio not significant, hence mod2 better than mod1.\n")

# 
# 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(mod2)))
mod <- setpars(mod,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

mod3 <- fit(mod,equal=conpat)

llratio(mod2,mod3)

cat("Test3: log Likelihood ratio significant, constraints not tenable.")



root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge