[Depmix-commits] r459 - in tags: . release-1.0-3 release-1.0-3/R release-1.0-3/data release-1.0-3/inst release-1.0-3/inst/doc release-1.0-3/man release-1.0-3/tests

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Wed Jun 15 12:00:09 CEST 2011


Author: ingmarvisser
Date: 2011-06-15 12:00:05 +0200 (Wed, 15 Jun 2011)
New Revision: 459

Added:
   tags/release-1.0-3/
   tags/release-1.0-3/DESCRIPTION
   tags/release-1.0-3/NAMESPACE
   tags/release-1.0-3/NEWS
   tags/release-1.0-3/R/
   tags/release-1.0-3/R/EM.R
   tags/release-1.0-3/R/allGenerics.R
   tags/release-1.0-3/R/depmix-class.R
   tags/release-1.0-3/R/depmix.R
   tags/release-1.0-3/R/depmixfit-class.R
   tags/release-1.0-3/R/depmixfit.R
   tags/release-1.0-3/R/depmixsim-class.R
   tags/release-1.0-3/R/em.control.R
   tags/release-1.0-3/R/fb.R
   tags/release-1.0-3/R/forwardbackward.R
   tags/release-1.0-3/R/freepars.R
   tags/release-1.0-3/R/getpars.R
   tags/release-1.0-3/R/llratio.R
   tags/release-1.0-3/R/logLik.R
   tags/release-1.0-3/R/lystig.R
   tags/release-1.0-3/R/makeDepmix.R
   tags/release-1.0-3/R/makePriorModel.R
   tags/release-1.0-3/R/makeResponseModels.R
   tags/release-1.0-3/R/makeTransModels.R
   tags/release-1.0-3/R/mlogit.R
   tags/release-1.0-3/R/multinomial.R
   tags/release-1.0-3/R/nobs.R
   tags/release-1.0-3/R/pa2conr.R
   tags/release-1.0-3/R/response-class.R
   tags/release-1.0-3/R/responseGLM.R
   tags/release-1.0-3/R/responseGLMBINOM.R
   tags/release-1.0-3/R/responseGLMGAMMA.R
   tags/release-1.0-3/R/responseGLMMULTINOM.R
   tags/release-1.0-3/R/responseGLMPOISSON.R
   tags/release-1.0-3/R/responseMVN.R
   tags/release-1.0-3/R/responseNORM.R
   tags/release-1.0-3/R/setpars.R
   tags/release-1.0-3/R/stationary.R
   tags/release-1.0-3/R/transInit.R
   tags/release-1.0-3/R/viterbi.R
   tags/release-1.0-3/README
   tags/release-1.0-3/data/
   tags/release-1.0-3/data/balance.rda
   tags/release-1.0-3/data/speed.rda
   tags/release-1.0-3/inst/
   tags/release-1.0-3/inst/CITATION
   tags/release-1.0-3/inst/doc/
   tags/release-1.0-3/inst/doc/baldist.pdf
   tags/release-1.0-3/inst/doc/depmixS4.Rnw
   tags/release-1.0-3/inst/doc/depmixS4.bib
   tags/release-1.0-3/inst/doc/depmixS4.pdf
   tags/release-1.0-3/man/
   tags/release-1.0-3/man/GLMresponse.Rd
   tags/release-1.0-3/man/balance.Rd
   tags/release-1.0-3/man/depmix-class.Rd
   tags/release-1.0-3/man/depmix-internal.Rd
   tags/release-1.0-3/man/depmix-methods.Rd
   tags/release-1.0-3/man/depmix.Rd
   tags/release-1.0-3/man/depmix.fit.Rd
   tags/release-1.0-3/man/depmix.fitted-class.Rd
   tags/release-1.0-3/man/depmix.sim-class.Rd
   tags/release-1.0-3/man/depmixS4-package.Rd
   tags/release-1.0-3/man/em.control.Rd
   tags/release-1.0-3/man/forwardbackward.Rd
   tags/release-1.0-3/man/llratio.Rd
   tags/release-1.0-3/man/makeDepmix.Rd
   tags/release-1.0-3/man/mix-class.Rd
   tags/release-1.0-3/man/mix.Rd
   tags/release-1.0-3/man/mix.fitted-class.Rd
   tags/release-1.0-3/man/mix.sim-class.Rd
   tags/release-1.0-3/man/posterior.Rd
   tags/release-1.0-3/man/response-class.Rd
   tags/release-1.0-3/man/response-classes.Rd
   tags/release-1.0-3/man/responses.Rd
   tags/release-1.0-3/man/simulate.Rd
   tags/release-1.0-3/man/speed.Rd
   tags/release-1.0-3/man/transInit.Rd
   tags/release-1.0-3/tests/
   tags/release-1.0-3/tests/test1speed.R
   tags/release-1.0-3/tests/test1speed.Rout.save
   tags/release-1.0-3/tests/test2getsetpars.R
   tags/release-1.0-3/tests/test2getsetpars.Rout.save
   tags/release-1.0-3/tests/test3responses.R
Log:
Added tag for release 1.0-3 (rev 458 on CRAN).

Added: tags/release-1.0-3/DESCRIPTION
===================================================================
--- tags/release-1.0-3/DESCRIPTION	                        (rev 0)
+++ tags/release-1.0-3/DESCRIPTION	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,12 @@
+Package: depmixS4
+Version: 1.0-3
+Date: 2011-03-22
+Title: Dependent Mixture Models - Hidden Markov Models of GLMs and Other Distributions in S4
+Author: Ingmar Visser <i.visser at uva.nl>, Maarten Speekenbrink <m.speekenbrink at ucl.ac.uk>
+Maintainer: Ingmar Visser <i.visser at uva.nl>
+Depends: R (>= 2.12.2), stats, nnet, methods, MASS, Rsolnp, stats4
+Suggests: Rdonlp2, gamlss, gamlss.dist
+Description: Fit latent (hidden) Markov models on mixed categorical and continuous (timeseries)
+   data, otherwise known as dependent mixture models
+License: GPL (>=2)
+URL: http://depmix.r-forge.r-project.org/


Property changes on: tags/release-1.0-3/DESCRIPTION
___________________________________________________________________
Added: svn:executable
   + *
Added: svn:eol-style
   + native

Added: tags/release-1.0-3/NAMESPACE
===================================================================
--- tags/release-1.0-3/NAMESPACE	                        (rev 0)
+++ tags/release-1.0-3/NAMESPACE	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,61 @@
+import(methods)
+
+importFrom(stats, predict, simulate)
+
+importFrom(stats4, AIC, BIC, logLik)
+
+export(	
+	makeDepmix,
+	makeMix,
+	lystig,
+	fb,
+	forwardbackward,
+	MVNresponse,
+	llratio,
+	multinomial,
+	em,
+	em.control,
+	viterbi,
+	mlogit,
+	logLik
+)
+
+exportClasses(
+	depmix,
+	depmix.sim,
+	mix,
+	mix.sim,
+	depmix.fitted,
+	mix.fitted,
+	response,
+	GLMresponse,
+	MVNresponse,
+	transInit
+)
+
+exportMethods(
+	fit,
+	getConstraints,
+	npar,
+	freepars,
+	nlin,
+	getdf,
+	nobs,
+	nresp,
+	ntimes,
+	nstates,
+	depmix,
+	mix,
+	posterior,
+	GLMresponse,
+	MVNresponse,
+	transInit,
+	setpars,
+	getpars,
+	predict,
+	dens,
+	show,
+	simulate,
+	summary,
+	logLik
+)

Added: tags/release-1.0-3/NEWS
===================================================================
--- tags/release-1.0-3/NEWS	                        (rev 0)
+++ tags/release-1.0-3/NEWS	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,175 @@
+
+Changes in depmixS4 version 1.0-4
+
+  o Fixed a bug in simulation of gaussian response model, which was 
+    returning NaNs due to an error in assignment of the sd parameter
+   (introduced in version x). Thanks to Jeffrey Arnold for reporting
+   this (bug #1365). 
+
+Changes in depmixS4 version 1.0-3
+
+  o Using AIC/BIC/logLik/nobs generics from stats/stats4 rather than 
+    defining them anew (which gave clashes with other packages that did 
+    the same).
+
+Changes in depmixS4 version 1.0-2
+
+  o fixed a bug in simulation of binomial response model data (the response
+    consists of the number of successes, and the number of failures; in 
+    simulation, the number of failures was an exact copy of the number of
+    successes). 
+
+  o added a meaningful error message in the EM algorithm for lca/mixture 
+    models in case the initial log likelihood is NA (thanks to Matthias
+    Ihrke for pointing this out). 
+
+Changes in depmixS4 version 1.0-1
+
+  o minor changes in documentation to conform to R 2.12.0 standards. 
+
+  o fixed a bug concerning random start values (the argument to specify
+    this was not passed to the EM algorithm and hence was completely 
+    ineffective ...). 
+
+  o changed the emcontrol argument to the fit function; it now calls 
+    a function em.control which returns the list of control parameters, which
+    now also includes maxit, the max number of iterations of the EM algorithm. 
+    This makes future additions to EM control parameters easier. 
+
+  o random parameter initialization is now the default when using EM 
+    to fit models. 
+  
+  o fixed a bug in multinomial models with n>1; the parameters are now 
+    normalized such that they sum to unity (this bug was introduced in 
+    version 0.9-0 in multinomial models with identity link). 
+
+  o added an error message for multinomial response models with n>1 and 
+    link='mlogit' as this case is not handled; n>1 multinomial can use the
+    'identity' link function. 
+
+Changes in depmixS4 version 1.0-0
+
+  o added a vignette to the package and upped the version number 1.0-0 to 
+    celebrate publication in the Journal of Statistical Software. 
+
+Changes in depmixS4 version 0.9-1
+
+  o fixed a bug in setting the lower and upper bounds for GLMresponse
+    models (the number of bounds was wrong for models with covariates/
+    predictors; these bounds are only used in constrained optimization in
+    which case they produced an error immediately; in EM optimization these
+    bounds are not used).
+
+Changes in depmixS4 version 0.9-0
+
+  o added optimization using Rsolnp, which can be invoked by using 
+    method="rsolnp" in calling fit on (dep-)mix objects. Note that 
+    this is meant for fitting models with additional constraints. 
+    Method="rsolnp" is now the default when fitting constrained 
+    models, method="donlp" is still supported. 
+
+  o added documentation for control arguments that can be passed to 
+    em algorithm, particularly for controlling the tolerance in 
+    optimization.
+
+  o added multinomial models with identity link for transition and prior 
+    probabilities. These are now the default when no covariates are 
+    present. 
+
+  o added bounds and constraints for multinomial identity models such
+    that these constraints are satisfied when fitting models with
+    method="rsolnp" or "donlp".  Also, variance and sd parameters in
+    gaussian and multivariate normal models are given bounds to prevent
+    warnings and errors in optimization of such models using rsolnp or 
+    donlp.
+
+  o added option to generate starting values as part of the EM 
+    algorithm. 
+
+  o fixed a bug in multinomial response models with n>1; the response for
+    these models can be specified as a k-column matrix with the number of
+    observed responses for each category; the loglikelihood for these
+    models in which there was more than 1 observation per row was
+    incorrect; note that such models may lead to some numerical
+    instabilities when n is large.
+
+Changes in depmixS4 version 0.3-0
+  
+  o added multinomial response function with identity link (no covariates
+    allowed in such a model); useful when (many) boundary values occur; 
+    currently no constraints are used for such models, and hence only EM
+    can be used for optimization, or alternatively, if and when Rdonlp2
+    is used, sum constraints need to be added when fitting the model.
+    See ?GLMresponse for details. 
+
+  o added an example of how to specify a model with multivariate normal
+    responses (and fixed a bug in MVNresponse that prevented such models
+    from being specified in the first place). See ?makeDepmix for an 
+    example. 
+
+Changes in depmixS4 version 0.2-2
+
+  o fixed a warning produced when specifying conrows.upper and .lower in
+    the fit function
+
+  o added error message in case the initial log likelihood is infeasible
+
+  o fixed a bug in the fit function for multinomial response models with 
+    covariates (thanks to Gilles Dutilh for spotting this)
+
+Changes in depmixS4 version 0.2-1
+
+  o fixed a bug in the Viterbi algorithm used to compute posterior states
+    (this bug was introduced in version 0.2-0)
+  
+  o restructured test files somewhat
+
+  o fixed a bug in the use of the conrows argument in the fit function (a 
+    missing drop=FALSE statement)
+
+  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 added newx argument to predict function of transInit objects, to be used
+    for predicting probabilities depending on covariates (useful in eg plotting
+    transition probabilities as function of a covariate)
+
+  o added example of the use of conrows argument in fitting functions and other 
+    minor updates in documentation
+  
+Changes in depmixS4 version 0.2-0
+
+  o restructured R and Rd (help) files; added depmixS4 help with a short
+    overview of the package and links to appropriate help files
+  
+  o added function 'simulate' to generate new data from a (fitted) model
+  
+  o added function 'forwardbackward' to access the forward and backward 
+    variables as well as the smoothed transition and state variables
+  
+  o added new glm distributions: gamma, poisson
+  
+  o added multivariate normal distribution
+  
+  o freepars now works correctly on both depmix and depmix.fitted objects
+  
+  o added function 'nlin' to compute the number of linear constraints in 
+    a fitted model object
+
+  o added mix class for mixture and latent class models; the depmix class 
+    extends this mix class and adds a transition model to it
+  
+  o added help file for makeDepmix to provide full control in specifying 
+    models
+  
+  o minor changes to make depmixS4 compatible with R 2.7.1
+  
+
+Changes in depmixS4 version 0.1-1
+
+  o adjusted for R 2.7.0
+
+First version released on CRAN: 0.1-0


Property changes on: tags/release-1.0-3/NEWS
___________________________________________________________________
Added: svn:eol-style
   + native

Added: tags/release-1.0-3/R/EM.R
===================================================================
--- tags/release-1.0-3/R/EM.R	                        (rev 0)
+++ tags/release-1.0-3/R/EM.R	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,249 @@
+# 
+# Maarten Speekenbrink 23-3-2008
+# 
+
+em <- function(object,...) {
+	if(!is(object,"mix")) stop("object is not of class '(dep)mix'")
+	call <- match.call()
+	if(is(object,"depmix")) {
+		call[[1]] <- as.name("em.depmix")
+	} else {
+		call[[1]] <- as.name("em.mix")
+	}
+	object <- eval(call, parent.frame())
+	object
+}
+
+# em for lca and mixture models
+em.mix <- function(object,maxit=100,tol=1e-8,crit="relative",random.start=TRUE,verbose=FALSE,...) {
+	
+	if(!is(object,"mix")) stop("object is not of class 'mix'")
+		
+	ns <- nstates(object)
+	ntimes <- ntimes(object)
+	lt <- length(ntimes)
+	et <- cumsum(ntimes)
+	bt <- c(1,et[-lt]+1)
+	
+	converge <- FALSE
+	j <- 0
+	
+	if(random.start) {
+				
+		nr <- sum(ntimes(object))
+		gamma <- matrix(runif(nr*ns,min=.0001,max=.9999),nr=nr,nc=ns)
+		gamma <- gamma/rowSums(gamma)
+		LL <- -1e10
+		
+		for(i in 1:ns) {
+			for(k in 1:nresp(object)) {
+				object at response[[i]][[k]] <- fit(object at response[[i]][[k]],w=gamma[,i])
+				# update dens slot of the model
+				object at dens[,k,i] <- dens(object at response[[i]][[k]])
+			}
+		}
+		
+		# initial expectation
+		fbo <- fb(init=object at init,matrix(0,1,1),B=object at dens,ntimes=ntimes(object),stationary=object at stationary)
+		LL <- fbo$logLike
+		
+		if(is.nan(LL)) stop("Cannot find suitable starting values; please provide them.")
+		
+	} else {
+		# initial expectation
+		B <- apply(object at dens,c(1,3),prod)
+		gamma <- object at init*B
+		LL <- sum(log(rowSums(gamma)))
+		if(is.nan(LL)) stop("Starting values not feasible; please provide them.")
+		gamma <- gamma/rowSums(gamma)
+	}
+	
+	LL.old <- LL + 1
+	
+	while(j <= maxit & !converge) {
+		
+		# maximization
+		
+		# should become object at prior <- fit(object at prior)
+		object at prior@y <- gamma[bt,,drop=FALSE]
+		object at prior <- fit(object at prior, w=NULL,ntimes=NULL)
+		object at init <- dens(object at prior)
+		
+		for(i in 1:ns) {
+			for(k in 1:nresp(object)) {
+				object at response[[i]][[k]] <- fit(object at response[[i]][[k]],w=gamma[,i])
+				# update dens slot of the model
+				object at dens[,k,i] <- dens(object at response[[i]][[k]])
+			}
+		}
+		
+		# expectation
+		B <- apply(object at dens,c(1,3),prod)
+		gamma <- object at init*B
+		LL <- sum(log(rowSums(gamma)))
+
+		# normalize
+		gamma <- gamma/rowSums(gamma)
+		
+		# print stuff
+		if(verbose&((j%%5)==0)) {
+			cat("iteration",j,"logLik:",LL,"\n")
+		}
+		
+		if(LL >= LL.old) {
+		  if((crit == "absolute" &&  LL - LL.old < tol) || (crit == "relative" && (LL.old - LL)/LL.old  < tol)) {
+			  cat("iteration",j,"logLik:",LL,"\n")
+			  converge <- TRUE
+			}
+		} else {
+		  # this should not really happen...
+		  if(j > 0) warning("likelihood decreased on iteration",j)
+		}
+
+		LL.old <- LL
+		j <- j+1
+
+	}
+
+	class(object) <- "mix.fitted"
+
+	if(converge) {
+		object at message <- switch(crit,
+			relative = "Log likelihood converged to within tol. (relative change)",
+			absolute = "Log likelihood converged to within tol. (absolute change)"
+		)
+	} else object at message <- "'maxit' iterations reached in EM without convergence."
+
+	# no constraints in EM, except for the standard constraints ...
+	# which are produced by the following (only necessary for getting df right in logLik and such)
+	constraints <- getConstraints(object)
+	object at conMat <- constraints$lincon
+	object at lin.lower <- constraints$lin.l
+	object at lin.upper <- constraints$lin.u
+	
+	object
+	
+}
+
+# em for hidden markov models
+em.depmix <- function(object,maxit=100,tol=1e-8,crit="relative",random.start=TRUE,verbose=FALSE,...) {
+	
+	if(!is(object,"depmix")) stop("object is not of class 'depmix'")
+	
+	ns <- nstates(object)
+	
+	ntimes <- ntimes(object)
+	lt <- length(ntimes)
+	et <- cumsum(ntimes)
+	bt <- c(1,et[-lt]+1)
+	
+	converge <- FALSE
+	j <- 0
+	
+	if(random.start) {
+				
+		nr <- sum(ntimes(object))
+		gamma <- matrix(runif(nr*ns,min=.0001,max=.9999),nr=nr,nc=ns)
+		gamma <- gamma/rowSums(gamma)
+		LL <- -1e10
+		
+		for(i in 1:ns) {
+			for(k in 1:nresp(object)) {
+				object at response[[i]][[k]] <- fit(object at response[[i]][[k]],w=gamma[,i])
+				# update dens slot of the model
+				object at dens[,k,i] <- dens(object at response[[i]][[k]])
+			}
+		}
+		
+		# initial expectation
+		fbo <- fb(init=object at init,A=object at trDens,B=object at dens,ntimes=ntimes(object),stationary=object at stationary)
+		LL <- fbo$logLike
+		
+		if(is.nan(LL)) stop("Cannot find suitable starting values; please provide them.")
+		
+	} else {
+		# initial expectation
+		fbo <- fb(init=object at init,A=object at trDens,B=object at dens,ntimes=ntimes(object),stationary=object at stationary)
+		LL <- fbo$logLike
+		if(is.nan(LL)) stop("Starting values not feasible; please provide them.")
+	}
+	
+	LL.old <- LL + 1
+	
+	while(j <= maxit & !converge) {
+		
+		# maximization
+				
+		# should become object at prior <- fit(object at prior, gamma)
+		object at prior@y <- fbo$gamma[bt,,drop=FALSE]
+		object at prior <- fit(object at prior, w=NULL, ntimes=NULL)
+		object at init <- dens(object at prior)
+				
+		trm <- matrix(0,ns,ns)
+		for(i in 1:ns) {
+			if(!object at stationary) {
+				object at transition[[i]]@y <- fbo$xi[,,i]/fbo$gamma[,i]
+				object at transition[[i]] <- fit(object at 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])
+				}
+				# FIX THIS; it will only work with specific trinModels
+				# should become object at transition = fit(object at transition, xi, gamma)
+				object at transition[[i]]@parameters$coefficients <- switch(object at transition[[i]]@family$link,
+					identity = object at transition[[i]]@family$linkfun(trm[i,]),
+					mlogit = object at transition[[i]]@family$linkfun(trm[i,],base=object at transition[[i]]@family$base),
+					object at transition[[i]]@family$linkfun(trm[i,])
+				)
+			}
+			# update trDens slot of the model
+			object at trDens[,,i] <- dens(object at transition[[i]])
+		}
+		
+		for(i in 1:ns) {
+			for(k in 1:nresp(object)) {
+				object at response[[i]][[k]] <- fit(object at response[[i]][[k]],w=fbo$gamma[,i])
+				# update dens slot of the model
+				object at dens[,k,i] <- dens(object at response[[i]][[k]])
+			}
+		}
+		
+		# expectation
+		fbo <- fb(init=object at init,A=object at trDens,B=object at dens,ntimes=ntimes(object),stationary=object at stationary)
+		LL <- fbo$logLike
+				
+		if(verbose&((j%%5)==0)) cat("iteration",j,"logLik:",LL,"\n")
+		
+		if( (LL >= LL.old)) {
+		  if((crit == "absolute" &&  LL - LL.old < tol) || (crit == "relative" && (LL.old - LL)/LL.old  < tol)) {
+			  cat("iteration",j,"logLik:",LL,"\n")
+			  converge <- TRUE
+			}
+		} else {
+		  # this should not really happen...
+		  if(j > 0) warning("likelihood decreased on iteration",j)
+		}
+		
+		LL.old <- LL
+		j <- j+1
+		
+	}
+		
+	class(object) <- "depmix.fitted"
+	
+	if(converge) {
+		object at message <- switch(crit,
+			relative = "Log likelihood converged to within tol. (relative change)",
+			absolute = "Log likelihood converged to within tol. (absolute change)"
+		)
+	} else object at message <- "'maxit' iterations reached in EM without convergence."
+	
+	# no constraints in EM, except for the standard constraints ...
+	# which are produced by the following (only necessary for getting df right in logLik and such)
+	constraints <- getConstraints(object)
+	object at conMat <- constraints$lincon
+	object at lin.lower <- constraints$lin.l
+	object at lin.upper <- constraints$lin.u
+	
+	object
+}


Property changes on: tags/release-1.0-3/R/EM.R
___________________________________________________________________
Added: svn:eol-style
   + native

Added: tags/release-1.0-3/R/allGenerics.R
===================================================================
--- tags/release-1.0-3/R/allGenerics.R	                        (rev 0)
+++ tags/release-1.0-3/R/allGenerics.R	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,80 @@
+
+# 
+# Ingmar Visser, 23-3-2008
+# 
+
+.First.lib <- function(lib, pkg) { 
+	require(stats)
+	require(methods)
+	require(MASS)
+ 	require(nnet)
+	require(Rsolnp)
+	require(stats4)	
+}
+
+.Last.lib <- function(libpath) {}
+
+# Guess what: all generics
+
+setGeneric("depmix", function(response,data=NULL,nstates,transition=~1,family=gaussian(),prior=~1,initdata=NULL,
+		respstart=NULL,trstart=NULL,instart=NULL,ntimes=NULL, ...) standardGeneric("depmix"))
+
+setGeneric("GLMresponse", function(formula, data = NULL, family = gaussian(), pstart =
+                 NULL, fixed = NULL, prob=TRUE, ...) standardGeneric("GLMresponse"))
+                 
+setGeneric("MVNresponse", function(formula, data = NULL,pstart=NULL,fixed=NULL,...) standardGeneric("MVNresponse"))
+
+setGeneric("transInit", function(formula, nstates, data = NULL, family = multinomial(),
+                 pstart = NULL, fixed = NULL, prob=TRUE, ...) standardGeneric("transInit"))
+
+# extractors/set functions
+
+setGeneric("npar", function(object, ...) standardGeneric("npar"))
+
+setGeneric("ntimes", function(object, ...) standardGeneric("ntimes"))
+
+setGeneric("nstates", function(object, ...) standardGeneric("nstates"))
+
+setGeneric("nresp", function(object, ...) standardGeneric("nresp"))
+
+setGeneric("freepars", function(object, ...) standardGeneric("freepars"))
+
+setGeneric("getdf",function(object) standardGeneric("getdf"))
+
+setGeneric("nlin", function(object, ...) standardGeneric("nlin"))
+
+setGeneric("nobs", function(object, ...) standardGeneric("nobs"))
+
+setGeneric("getConstraints", function(object, ...) standardGeneric("getConstraints"))
+
+setGeneric("is.stationary", function(object,...) standardGeneric("is.stationary"))
+
+setGeneric("setpars", function(object,values,which="pars",...) standardGeneric("setpars"))
+
+setGeneric("getpars", function(object,which="pars",...) standardGeneric("getpars"))
+
+
+# functions 
+setGeneric("fit", function(object, ...) standardGeneric("fit"))
+
+setGeneric("posterior", function(object, ...) standardGeneric("posterior"))
+
+setGeneric("forwardbackward", function(object, ...) standardGeneric("forwardbackward"))
+
+setGeneric("simulate", function(object,nsim=1,seed=NULL, ...) standardGeneric("simulate"))
+
+setGeneric("logDens",function(object,...) standardGeneric("logDens"))
+
+setGeneric("dens",function(object,...) standardGeneric("dens"))
+
+setGeneric("predict", function(object, ...) standardGeneric("predict"))
+
+
+# redundant??
+
+# setGeneric("getModel", function(object, ...) standardGeneric("getModel"))
+
+# these are imported from stats4
+# setGeneric("logLik", function(object, ...) standardGeneric("logLik"))
+# setGeneric("AIC", function(object, ..., k=2) standardGeneric("AIC"))
+# setGeneric("BIC", function(object, ...) standardGeneric("BIC"))


Property changes on: tags/release-1.0-3/R/allGenerics.R
___________________________________________________________________
Added: svn:eol-style
   + native

Added: tags/release-1.0-3/R/depmix-class.R
===================================================================
--- tags/release-1.0-3/R/depmix-class.R	                        (rev 0)
+++ tags/release-1.0-3/R/depmix-class.R	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,306 @@
+
+# 
+# Ingmar Visser, 11-6-2008
+# 
+
+# 
+# DEPMIX CLASS BELOW THE MIX CLASS
+# 
+
+# 
+# Class definition, accessor functions, print and summary methods
+# 
+
+# 
+# MIX CLASS
+# 
+
+setClass("mix",
+	representation(response="list", # response models
+		prior="ANY", # the prior model (multinomial)
+		dens="array", # response densities (B)
+		init="array", # usually called pi 
+		nstates="numeric",
+		nresp="numeric",
+		ntimes="numeric",
+		npars="numeric" # number of parameters
+	)
+)
+
+# accessor functions
+setMethod("npar","mix",
+	function(object) return(object at npars)
+)
+
+setMethod("ntimes","mix",
+	function(object) return(object at ntimes)
+)
+
+setMethod("nstates","mix",
+	function(object) return(object at nstates)
+)
+
+setMethod("nresp","mix",
+	function(object) return(object at nresp)
+)
+
+setMethod("is.stationary",signature(object="mix"),
+  function(object) {
+		return(TRUE)
+	}
+)
+
+setMethod("simulate",signature(object="mix"),
+	function(object,nsim=1,seed=NULL,...) {
+		
+		if(!is.null(seed)) set.seed(seed)
+		
+		ntim <- ntimes(object)
+		nt <- sum(ntim)
+		bt <- 1:nt
+		
+		nr <- nresp(object)
+		ns <- nstates(object)
+		
+		# simulate state sequences first, then observations
+		
+		# random generation is slow when done separately for each t, so first draw
+		# variates for all t, and then determine state sequences iteratively
+		states <- array(,dim=c(nt,nsim))
+		states[bt,] <- simulate(object at prior,n=nsim,is.prior=T)
+		sims <- array(,dim=c(nt,ns,nsim))
+				
+		states <- as.vector(states)
+		responses <- list(length=nr)
+		#responses <- array(,dim=c(nt,nr,nsim))
+		for(i in 1:nr) {
+			tmp <- matrix(,nrow=nt*nsim,ncol=NCOL(object at response[[1]][[i]]@y))
+			for(j in 1:ns) {
+				tmp[states==j,] <- simulate(object at response[[j]][[i]],nsim=nsim)[states==j,]
+			}
+			responses[[i]] <- tmp
+		}
+		
+		# generate new mix.sim object
+		class(object) <- c("mix.sim")
+		object at states <- as.matrix(states)
+		
+		object at prior@x <- as.matrix(apply(object at prior@x,2,rep,nsim))
+		for(j in 1:ns) {
+			for(i in 1:nr) {
+				object at response[[j]][[i]]@y <- as.matrix(responses[[i]])
+				object at response[[j]][[i]]@x <- as.matrix(apply(object at response[[j]][[i]]@x,2,rep,nsim))
+			}
+		}
+		object at ntimes <- rep(object at ntimes,nsim)
+		
+		# make appropriate array for transition densities
+		nt <- sum(object at ntimes)
+		
+		# make appropriate array for response densities
+		dns <- array(,c(nt,nr,ns))
+		
+		# compute observation and transition densities
+		for(i in 1:ns) {
+			for(j in 1:nr) {
+				dns[,j,i] <- dens(object at response[[i]][[j]]) # remove this response as an argument from the call to setpars
+			}
+		}
+		
+		# compute initial state probabilties
+		object at init <- dens(object at prior)
+		object at dens <- dns
+		
+		return(object)
+	}
+)
+
+# setMethod("getModel",signature(object="mix"),
+# 	function(object,which="response",...) {
+# 		res <- switch(which,
+# 			"prior"=object at prior,
+# 			"response"=object at response)
+# 		res
+# 	}
+# )
+
+# 
+# PRINT method
+# 
+
+setMethod("show","mix",
+	function(object) {
+		cat("Initial state probabilties model \n")
+		print(object at prior)
+		cat("\n")
+		for(i in 1:object at nstates) {
+			cat("Response model(s) for state", i,"\n\n")
+			for(j in 1:object at nresp) {
+				cat("Response model for response",j,"\n")
+				print(object at response[[i]][[j]])
+				cat("\n")
+			}
+			cat("\n")
+		}
+	}
+)
+
+# 
+# SUMMARY method: to do
+# 
+
+
+# 
+# Ingmar Visser, 23-3-2008
+# 
+
+# 
+# Class definition, accessor functions, print and summary methods
+# 
+
+# 
+# DEPMIX CLASS
+# 
+
+setClass("depmix",
+	representation(transition="list", # transition models (multinomial logistic)
+		trDens="array", # transition densities (A)
+		stationary="logical"
+	),
+	contains="mix"
+)
+
+# 
+# PRINT method
+# 
+
+setMethod("show","depmix",
+	function(object) {
+		cat("Initial state probabilties model \n")
+		print(object at prior)
+		cat("\n")
+		for(i in 1:object at nstates) {
+			cat("Transition model for state (component)", i,"\n")
+			print(object at transition[[i]])
+			cat("\n")
+		}
+		cat("\n")
+		for(i in 1:object at nstates) {
+			cat("Response model(s) for state", i,"\n\n")
+			for(j in 1:object at nresp) {
+				cat("Response model for response",j,"\n")
+				print(object at response[[i]][[j]])
+				cat("\n")
+			}
+			cat("\n")
+		}
+	}
+)
+
+setMethod("is.stationary",signature(object="depmix"),
+  function(object) {
+		return(object at stationary)
+	}
+)
+
+setMethod("simulate",signature(object="depmix"),
+	function(object,nsim=1,seed=NULL,...) {
+		
+		if(!is.null(seed)) set.seed(seed)
+		
+		ntim <- ntimes(object)
+		nt <- sum(ntim)
+		lt <- length(ntim)
+		et <- cumsum(ntim)
+		bt <- c(1,et[-lt]+1)
+		
+		nr <- nresp(object)
+		ns <- nstates(object)
+		
+		# simulate state sequences first, then observations
+		
+		# random generation is slow when done separately for each t, so first draw
+		#   variates for all t, and then determine state sequences iteratively
+		states <- array(,dim=c(nt,nsim))
+		states[bt,] <- simulate(object at prior,n=nsim,is.prior=T)
+		sims <- array(,dim=c(nt,ns,nsim))
+		for(i in 1:ns) {
+			if(is.stationary(object)) {
+				# TODO: this is a temporary fix!!! 
+				sims[,i,] <- simulate(object at transition[[i]],nsim=nsim,times=rep(1,nt))
+			} else {
+				sims[,i,] <- simulate(object at transition[[i]],nsim=nsim)
+			}
+		}
+		# track states
+		for(case in 1:lt) {
+			for(i in (bt[case]+1):et[case]) {
+				states[i,] <- sims[cbind(i,states[i-1,],1:nsim)]
+			}
+		}
+		
+		states <- as.vector(states)
+		responses <- list(length=nr)
+		#responses <- array(,dim=c(nt,nr,nsim))
+		for(i in 1:nr) {
+			tmp <- matrix(,nrow=nt*nsim,ncol=NCOL(object at response[[1]][[i]]@y))
+			for(j in 1:ns) {
+				tmp[states==j,] <- simulate(object at response[[j]][[i]],nsim=nsim)[states==j,]
+			}
+			responses[[i]] <- tmp
+		}
+		
+		# generate new depmix.sim object
+		class(object) <- c("depmix.sim")
+		object at states <- as.matrix(states)
+		
+		object at prior@x <- as.matrix(apply(object at prior@x,2,rep,nsim))
+		for(j in 1:ns) {
+			if(!is.stationary(object)) object at transition[[j]]@x <- as.matrix(apply(object at transition[[j]]@x,2,rep,nsim))
+			for(i in 1:nr) {
+				object at response[[j]][[i]]@y <- as.matrix(responses[[i]])
+				object at response[[j]][[i]]@x <- as.matrix(apply(object at response[[j]][[i]]@x,2,rep,nsim))
+			}
+		}
+		object at ntimes <- rep(object at ntimes,nsim)
+		
+		# make appropriate array for transition densities
+		nt <- sum(object at ntimes)
+		if(is.stationary(object)) trDens <- array(0,c(1,ns,ns)) else trDens <- array(0,c(nt,ns,ns))
+		
+		# make appropriate array for response densities
+		dns <- array(,c(nt,nr,ns))
+		
+		# compute observation and transition densities
+		for(i in 1:ns) {
+			for(j in 1:nr) {
+				dns[,j,i] <- dens(object at response[[i]][[j]]) # remove this response as an argument from the call to setpars
+			}
+			trDens[,,i] <- dens(object at transition[[i]])
+		}
+		
+		# compute initial state probabilties
+		object at init <- dens(object at prior)
+		object at trDens <- trDens
+		object at dens <- dns
+		
+		return(object)
+	}
+)
+
+# setMethod("getModel",signature(object="depmix"),
+# 	function(object,which="response",...) {
+# 		res <- switch(which,
+# 			"prior"=object at prior,
+# 			"response"=object at response,
+# 			"transition"=object at transition)
+# 		res
+# 	}
+# )
+
+# 
+# SUMMARY method: to do
+# 
+
+
+


Property changes on: tags/release-1.0-3/R/depmix-class.R
___________________________________________________________________
Added: svn:eol-style
   + native

Added: tags/release-1.0-3/R/depmix.R
===================================================================
--- tags/release-1.0-3/R/depmix.R	                        (rev 0)
+++ tags/release-1.0-3/R/depmix.R	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,94 @@
+#
+# Ingmar Visser, 11-6-2008
+#
+
+#
+# Main function to construct mix models
+#
+
+#
+# UNIVARIATE AND MULTIVARIATE MIXTURE OF GLM'S
+#
+
+
+setGeneric("mix", function(response, data = NULL, 
+    nstates, family = gaussian(), prior = ~1, initdata = NULL, 
+    respstart = NULL, instart = NULL, ...) standardGeneric("mix"))
+
+
+setMethod("mix", signature(response = "ANY"), function(response, 
+    data = NULL, nstates, family = gaussian(), prior = ~1, initdata = NULL, 
+    respstart = NULL, instart = NULL, ...) {
+    
+    # make response models
+    response <- makeResponseModels(response = response, data = data, 
+        nstates = nstates, family = family, values = respstart)
+    
+    # FIX ME: this only works if data are actually provided ... 
+	# (maybe make this obligatory ...)
+    ntimes <- rep(1, nrow(data))
+    
+    # make prior model
+    prior <- makePriorModel(nstates = nstates, ncases = length(ntimes), 
+        formula = prior, data = initdata, values = instart)
+    
+    # call main depmix with all these models, ntimes and stationary
+    model <- makeMix(response = response, prior = prior)
+        
+    return(model)
+})
+
+#
+# Ingmar Visser, 23-3-2008
+#
+
+#
+# Main function to construct depmix models
+#
+
+#
+# UNIVARIATE AND MULTIVARIATE MARKOV MIXTURE OF GLM'S
+#
+
+setMethod("depmix", signature(response = "ANY"), function(response, 
+    data = NULL, nstates, transition = ~1, family = gaussian(), 
+    prior = ~1, initdata = NULL, respstart = NULL, trstart = NULL, 
+    instart = NULL, ntimes = NULL, ...) {
+    
+    if (is.null(data)) {
+        if (is.null(ntimes)) 
+            stop("'ntimes' must be provided if not in the data")
+    } else {
+        if (is.null(attr(data, "ntimes"))) {
+            if (is.null(ntimes)) 
+                ntimes <- nrow(data)
+        } else {
+            ntimes <- attr(data, "ntimes")
+        }
+        if (sum(ntimes) != nrow(data)) 
+            stop("'ntimes' and data do not match")
+    }
+    
+    # make response models
+    response <- makeResponseModels(response = response, data = data, 
+        nstates = nstates, family = family, values = respstart)
+    
+    # make transition models
+    stationary = FALSE
+    if (transition == ~1) 
+        stationary = TRUE
+    transition <- makeTransModels(nstates = nstates, formula = transition, 
+        data = data, stationary = stationary, values = trstart)
+    
+    # make prior model
+    prior <- makePriorModel(nstates = nstates, ncases = length(ntimes), 
+        formula = prior, data = initdata, values = instart)
+    
+    # call main depmix with all these models, ntimes and stationary
+    model <- makeDepmix(response = response, transition = transition, 
+        prior = prior, ntimes = ntimes, stationary = stationary)
+    
+    # deal with starting values here!!!!!!
+    
+    return(model)
+})


Property changes on: tags/release-1.0-3/R/depmix.R
___________________________________________________________________
Added: svn:eol-style
   + native

Added: tags/release-1.0-3/R/depmixfit-class.R
===================================================================
--- tags/release-1.0-3/R/depmixfit-class.R	                        (rev 0)
+++ tags/release-1.0-3/R/depmixfit-class.R	2011-06-15 10:00:05 UTC (rev 459)
@@ -0,0 +1,138 @@
+
+# 
+# Ingmar Visser, 11-6-2008
+# 
+
+# Changes
+# - added lin.upper and lin.lower slots to these objects
+
+# 
+# MIX.FITTED CLASS
+# 
+
+setClass("mix.fitted",
+	representation(message="character", # convergence information
+		conMat="matrix", # constraint matrix on the parameters for general linear constraints
+		lin.upper="numeric", # upper bounds for linear constraint
+		lin.lower="numeric", # lower bounds for linear constraints
+		posterior="data.frame" # posterior probabilities for the states
+	),
+	contains="mix"
+)
+
+# accessor functions
+
+setMethod("posterior","mix.fitted",
+	function(object) {
+		return(object at posterior)
+	}
+)
+
+setMethod("show","mix.fitted",
+	function(object) {
+		cat("Convergence info:",object at message,"\n")
+		print(logLik(object))
+		cat("AIC: ", AIC(object),"\n")
+		cat("BIC: ", BIC(object),"\n")
+	}
+)
+
+setMethod("summary","mix.fitted",
+	function(object,which="all") {
+		ans=switch(which,
+			"all" = 1,
+			"response" = 2,
+			"prior" = 3,
+			stop("Invalid 'which' argument in summary of fitted mix model")
+		)
+		if(ans==1|ans==3) {
+			cat("Mixture probabilities model \n")
+			print(object at prior)
+			cat("\n")
+		}
+		if(ans==1|ans==2) {
+			for(i in 1:object at nstates) {
+				cat("Response model(s) for state", i,"\n\n")
+				for(j in 1:object at nresp) {
+					cat("Response model for response",j,"\n")
+					print(object at response[[i]][[j]])
+					cat("\n")
+				}
+				cat("\n")
+			}
+		}
+	}	
+)
+
+# 
+# Ingmar Visser, 23-3-2008
+# 
+
+# 
+# DEPMIX.FITTED CLASS
+# 
+
+setClass("depmix.fitted",
+	representation(message="character", # convergence information
+		conMat="matrix", # constraint matrix on the parameters for general linear constraints
+		lin.upper="numeric", # upper bounds for linear constraints
+		lin.lower="numeric", # lower bounds for linear constraints
+		posterior="data.frame" # posterior probabilities for the states
+	),
+	contains="depmix"
+)
+
+# accessor functions
+
+setMethod("posterior","depmix.fitted",
+	function(object) {
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/depmix -r 459


More information about the depmix-commits mailing list