[Depmix-commits] r167 - trunk/R

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

Author: maarten
Date: 2008-06-17 12:06:32 +0200 (Tue, 17 Jun 2008)
New Revision: 167

Modified:
   trunk/R/responseGLM.R
   trunk/R/responseGLMBINOM.R
Log:
changed BINOMresponse so that n>1 is now contained in y matrix as a second column.

Modified: trunk/R/responseGLM.R
===================================================================
--- trunk/R/responseGLM.R	2008-06-17 09:34:30 UTC (rev 166)
+++ trunk/R/responseGLM.R	2008-06-17 10:06:32 UTC (rev 167)
@@ -38,19 +38,19 @@
 			switch(is(y)[1],
     		factor = {
     		  y <- as.matrix(as.numeric(as.numeric(y)==1))
-    		  n <- matrix(1,nrow=nrow(y))
+    		  #n <- matrix(1,nrow=nrow(y))
     		},
     		matrix = {
     		  if(ncol(y) == 2) {
-    			n <- as.matrix(rowSums(y))
-    			y <- as.matrix(y[,1])
+    			#n <- as.matrix(rowSums(y))
+    			#y <- as.matrix(y[,1])
     		  } else {
     			stop("model response not valid for binomial model")
     		  }
     		},
     		numeric = {
     		  if(sum(y %in% c(0,1)) != length(y)) stop("model response not valid for binomial model")
-    		  n <- matrix(1,nrow=length(y))
+    		  #n <- matrix(1,nrow=length(y))
     		  y <- as.matrix(y)
     		},
     		stop("model response not valid for binomial model")
@@ -89,7 +89,7 @@
 		}
 		mod <- switch(family$family,
 			gaussian = new("NORMresponse",formula=formula,family=family,parameters=parameters,fixed=fixed,x=x,y=y,npar=npar),
-			binomial = new("BINOMresponse",formula=formula,family=family,parameters=parameters,fixed=fixed,x=x,y=y,npar=npar,n=n),
+			binomial = new("BINOMresponse",formula=formula,family=family,parameters=parameters,fixed=fixed,x=x,y=y,npar=npar),
 			multinomial = new("MULTINOMresponse",formula=formula,family=family,parameters=parameters,fixed=fixed,x=x,y=y,npar=npar),
 			poisson = new("POISSONresponse",formula=formula,family=family,parameters=parameters,fixed=fixed,x=x,y=y,npar=npar),
 			Gamma = new("GAMMAresponse",formula=formula,family=family,parameters=parameters,fixed=fixed,x=x,y=y,npar=npar),
Modified: trunk/R/responseGLMBINOM.R
===================================================================
--- trunk/R/responseGLMBINOM.R	2008-06-17 09:34:30 UTC (rev 166)
+++ trunk/R/responseGLMBINOM.R	2008-06-17 10:06:32 UTC (rev 167)
@@ -1,5 +1,5 @@
 
-setClass("BINOMresponse",representation(n="matrix"),contains="GLMresponse")
+setClass("BINOMresponse",contains="GLMresponse")
 
 # method 'fit'
 # use: in EM (M step)
@@ -10,13 +10,23 @@
 # returns: matrix with log(p(y|x,parameters))
 setMethod("logDens","BINOMresponse",
 	function(object) {
-		dbinom(x=object at y,size=object at n,prob=predict(object),log=TRUE)
+    if(NCOL(object at y) == 2) {
+      dbinom(x=object at y[,1],size=object at y[,2],prob=predict(object),log=TRUE)
+    } else {
+      if(!NCOL(object at y==1)) stop("not a valid response matrix for BINOMresponse")
+		  dbinom(x=object at y,prob=predict(object),log=TRUE)
+    }
 	}
 )
 
 setMethod("dens","BINOMresponse",
 	function(object,log=FALSE) {
-		dbinom(x=object at y,size=object at n,prob=predict(object),log=log)
+    if(NCOL(object at y) == 2) {
+      dbinom(x=object at y[,1],size=object at y[,2],prob=predict(object),log=log)
+    } else {
+      if(!NCOL(object at y==1)) stop("not a valid response matrix for BINOMresponse")
+		  dbinom(x=object at y,prob=predict(object),log=log)
+    }
 	}
 )
 
@@ -27,9 +37,13 @@
       pr <- predict(object)
     } else {
       pr <- predict(object)[time,]
-    }  
+    }
     nt <- nrow(pr)
-    response <- rbinom(nt*nsim,size=object at n,prob=pr)
+    if(NCOL(object at y) == 2) {
+      response <- rbinom(nt*nsim,size=object at y[,2],prob=pr)
+    } else {
+      response <- rbinom(nt*nsim,size=1,prob=pr)
+    }
     #if(nsim > 1) response <- matrix(response,ncol=nsim)
     response <- as.matrix(response)
     return(response)