[Depmix-commits] r250 - in trunk: . R tests

[noreply at r-forge.r-project.org]

Author: ingmarvisser
Date: 2009-02-17 15:22:50 +0100 (Tue, 17 Feb 2009)
New Revision: 250

Modified:
   trunk/NEWS
   trunk/R/responseGLM.R
   trunk/R/responseGLMGAMMA.R
   trunk/R/transInit.R
   trunk/tests/test3responses.R
Log:
Fixed a bug in setting starting values of transInit models

Modified: trunk/NEWS
===================================================================
--- trunk/NEWS	2009-02-17 10:35:00 UTC (rev 249)
+++ trunk/NEWS	2009-02-17 14:22:50 UTC (rev 250)
@@ -10,9 +10,12 @@
 
   o Updated help files for mix classes
 
+  o Fixed a bug in setting the starting values of regression coefficients in 
+    prior and transInit models with covariates (thanks to Verena Schmittmann 
+    for reporting this)
+
   o 
 
-
 Changes in depmixS4 version 0.2-0
 
   o restructured R and Rd (help) files; added depmixS4 help with a short

Modified: trunk/R/responseGLM.R
===================================================================
--- trunk/R/responseGLM.R	2009-02-17 10:35:00 UTC (rev 249)
+++ trunk/R/responseGLM.R	2009-02-17 14:22:50 UTC (rev 250)
@@ -120,7 +120,6 @@
 				else {
 					cat(object at family$linkinv(object at parameters$coefficients,base=object at family$base),"\n")
 				}
-				
 			}
 		}
 		if(object at family$family=="binomial") {
Modified: trunk/R/responseGLMGAMMA.R
===================================================================
--- trunk/R/responseGLMGAMMA.R	2009-02-17 10:35:00 UTC (rev 249)
+++ trunk/R/responseGLMGAMMA.R	2009-02-17 14:22:50 UTC (rev 250)
@@ -9,13 +9,13 @@
 # returns: matrix with log(p(y|x,parameters))
 setMethod("logDens","GAMMAresponse",
 	function(object) {
-		dpois(x=object at y,shape=predict(object),log=TRUE)
+		dgamma(x=object at y,shape=predict(object),log=TRUE)
 	}
 )
 
 setMethod("dens","GAMMAresponse",
 	function(object,log=FALSE) {
-		dpois(x=object at y,shape=predict(object),log=log)
+		dgamma(x=object at y,shape=predict(object),log=log)
 	}
 )
 

Modified: trunk/R/transInit.R
===================================================================
--- trunk/R/transInit.R	2009-02-17 10:35:00 UTC (rev 249)
+++ trunk/R/transInit.R	2009-02-17 14:22:50 UTC (rev 250)
@@ -44,7 +44,7 @@
 				} else {
 					parameters$coefficients[1,] <- pstart[1:ncol(parameters$coefficients)]
 				}
-				pstart <- matrix(pstart,,ncol(x),byrow=TRUE)
+				pstart <- matrix(pstart,ncol(x),byrow=TRUE)
 				if(ncol(x)>1) parameters$coefficients[2:ncol(x),] <- pstart[2:ncol(x),]
 			} else {
 				if(family$link=="identity") parameters$coefficients <- family$linkfun(pstart[1:length(parameters$coefficients)])
Modified: trunk/tests/test3responses.R
===================================================================
--- trunk/tests/test3responses.R	2009-02-17 10:35:00 UTC (rev 249)
+++ trunk/tests/test3responses.R	2009-02-17 14:22:50 UTC (rev 250)
@@ -4,79 +4,54 @@
 
 require(depmixS4)
 
-
-# gaussian response model
-y <- rnorm(1000)
-mod <- GLMresponse(y~1)
-fm <- fit(mod)
-cat("Test gaussian fit: ", all.equal(getpars(fm),c(mean(y),sd(y)),check=FALSE))
-
-
-# multinomial response model
-y<-sample(1:3,1000,rep=TRUE)
-mod <- GLMresponse(y~1,family=multinomial())
-fm <- fit(mod)
-cat("Test multinomial fit: ", all.equal(getpars(fm),c(table(y)/1000),check=FALSE))
-
-
-x <- sample(0:1,1000,rep=TRUE)
-
-mod <- GLMresponse(sample(1:3,1000,rep=TRUE)~x,family=multinomial(),pstart=c(0.33,0.33,0.33,0,0,1))
-
-mod at y <- simulate(mod)
-
-fit(mod)
-
-colSums(mod at y[which(x==0),])/length(which(x==0))
-
-colSums(mod at y[which(x==1),])/length(which(x==1))
-
+# binomial response model
 x <- rnorm(1000)
-
 library(boot)
-
 p <- inv.logit(x)
-
 ss <- rbinom(1000,1,p)
-
 mod <- GLMresponse(cbind(ss,1-ss)~x,family=binomial())
-
 fit(mod)
-
 glm(cbind(ss,1-ss)~x, family=binomial)
 
-x <- abs(rnorm(1000,2))
-
-res <- rpois(1000,x)
-
-mod <- GLMresponse(res~x,family=poisson())
-
-fit(mod)
-
-glm(res~x, family=poisson)
-
+# gamma response model
 x=runif(1000,1,5)
-
 res <- rgamma(1000,x)
-
 ## note that gamma needs proper starting values which are not
 ## provided by depmixS4 (even with them, this may produce warnings)
 mod <- GLMresponse(res~x,family=Gamma(),pst=c(0.8,1/0.8))
-
 fit(mod)
-
 glm(res~x,family=Gamma)
 
-mn <- c(1,2,3)
+# multinomial response model
+x <- sample(0:1,1000,rep=TRUE)
+mod <- GLMresponse(sample(1:3,1000,rep=TRUE)~x,family=multinomial(),pstart=c(0.33,0.33,0.33,0,0,1))
+mod at y <- simulate(mod)
+fit(mod)
+colSums(mod at y[which(x==0),])/length(which(x==0))
+colSums(mod at y[which(x==1),])/length(which(x==1))
 
+# multivariate normal response model
+mn <- c(1,2,3)
 sig <- matrix(c(1,.5,0,.5,1,0,0,0,2),3,3)
-
 y <- mvrnorm(1000,mn,sig)
-
 mod <- MVNresponse(y~1)
+fit(mod)
+colMeans(y)
+var(y)
 
+# normal (gaussian) response model
+y <- rnorm(1000)
+mod <- GLMresponse(y~1)
+fm <- fit(mod)
+cat("Test gaussian fit: ", all.equal(getpars(fm),c(mean(y),sd(y)),check=FALSE))
+
+
+# poisson response model
+x <- abs(rnorm(1000,2))
+res <- rpois(1000,x)
+mod <- GLMresponse(res~x,family=poisson())
 fit(mod)
+glm(res~x, family=poisson)
 
-colMeans(y)
 
-var(y)
+