[Depmix-commits] r624 - pkg/depmixS4/R

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

Author: maarten
Date: 2014-03-28 00:14:42 +0100 (Fri, 28 Mar 2014)
New Revision: 624

Modified:
   pkg/depmixS4/R/EM.R
Log:
solves  some problems in classification EM when a state/component does not have any assignments (keeps parameters the same as on the previous iteration)

Modified: pkg/depmixS4/R/EM.R
===================================================================
--- pkg/depmixS4/R/EM.R	2014-03-26 18:28:09 UTC (rev 623)
+++ pkg/depmixS4/R/EM.R	2014-03-27 23:14:42 UTC (rev 624)
@@ -145,9 +145,10 @@
 		if(clsf == "hard") {
 		  fbo <- list()
 		  vstate <- apply(gamma,1,which.max)
+      fbo$gamma <- t(apply(gamma,1,ind.max))
 		  B <- dens
 		  if(na.allow) B[is.na(B)] <- 1
-		  fbo$gamma <- t(apply(gamma,1,ind.max))
+		  #fbo$gamma <- t(apply(gamma,1,ind.max))
 		  fbo$logLike <- sum(log((apply(B,c(1,3),prod))[cbind(1:sum(ntimes),vstate)])) + sum(log(init[cbind(1:lt,vstate)]))
 		} else {
 		  fbo <- fb(init=init,matrix(0,1,1),B=dens,ntimes=ntimes(object))
@@ -182,9 +183,11 @@
 		
 		for(i in 1:ns) {
 			for(k in 1:nresp(object)) {
-				response[[i]][[k]] <- fit(response[[i]][[k]],w=fbo$gamma[,i])
-				# update dens slot of the model
-				dens[,k,i] <- dens(response[[i]][[k]])
+				if(sum(fbo$gamma[,i]) > 0) {
+          response[[i]][[k]] <- fit(response[[i]][[k]],w=fbo$gamma[,i])
+				  # update dens slot of the model
+				  dens[,k,i] <- dens(response[[i]][[k]])
+				}
 			}
 		}
 		
@@ -325,12 +328,18 @@
 		trm <- matrix(0,ns,ns)
 		for(i in 1:ns) {
 			if(!object at homogeneous) {
+        # TODO: check whether fbo$gamma > 0, otherwise set to previous value....
 				transition[[i]]@y <- fbo$xi[,,i]/fbo$gamma[,i]
 				transition[[i]] <- fit(transition[[i]],w=as.matrix(fbo$gamma[,i]),ntimes=ntimes(object)) # check this
 			} else {
-				for(k in 1:ns) {
-					trm[i,k] <- sum(fbo$xi[-c(et),k,i])/sum(fbo$gamma[-c(et),i])
-				}
+			  if(sum(fbo$gamma[-c(et),i]) == 0) {
+          # set unidentified transition probs to previous value
+          trm[i,] <- trDens[1,,i]
+			  } else {
+  				for(k in 1:ns) {
+  					trm[i,k] <- sum(fbo$xi[-c(et),k,i])/sum(fbo$gamma[-c(et),i])
+  				}
+			  }
 				# FIX THIS; it will only work with specific trinModels
 				# should become object at transition = fit(object at transition, xi, gamma)
 				transition[[i]]@parameters$coefficients <- switch(transition[[i]]@family$link,