version 0.2-2

ChangeLog

@@ -285,3 +285,37 @@ procedure).
 Various other small infelicities fixed.
 
 Uploaded to CRAN 31 January 2014.
+
+Version 0.2-1 --> 0.2-2
+
+Fixed typo (an omitted parenthesis) in fitted.hmm.discnp.Rd.
+(02/02/2014)
+
+Augmented the details in fitted.hmm.discnp.Rd to explain the fitted
+value concept.  (02/02/2014)
+
+Revised code of fitted.hmm.discnp() to (a) make sure that "y" is
+a list and to check that all components of the list are numeric.
+(02/02/2014)
+
+Fixed up the description of "ispd" in logLikeHmm.Rd and in mps.Rd
+to cover the case where ispd is a matrix.  (02/02/2014)
+
+Fixed bug which occurred when the observations contained missing
+values causing Rho[x,1:ncol(Rho)] to throw a "subscript out of
+bounds error. (18/05/2014)
+
+Revised code to make sure that Rho *always* has row names and that
+the observations *always* get treated as *character*. (18/05/2014)
+
+Added some example data sets:  The Sydney coliform count data and
+the Albert et al. multiple sclerosis lesion count data. (18/05/2014)
+
+Changed mat2list() to charList() which does the conversion from
+matrix to list (and the whinge-producing) and coerces the y-values
+to character. (19/05/2014)
+
+Changed check.yval() to return an emission probability matrix Rho
+if the given Rho doesn't have row names.  The returned Rho is
+equipped with appropriate row names if possible.  Otherwise an
+error is thrown. (19/05/2014)

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: hmm.discnp
-Version: 0.2-1
-Date: 2014-01-31
+Version: 0.2-2
+Date: 2014-05-19
 Title: Hidden Markov models with discrete non-parametric observation
         distributions.
 Author: Rolf Turner <r.turner@auckland.ac.nz> and Limin Liu.
@@ -12,9 +12,10 @@ Description: Fits hidden Markov models with discrete non-parametric
 	states, the most probable sequences of such states, and the
 	log likelihood of a collection of observations given the
 	parameters of the model.
+LazyData: true
 License: GPL (>= 2)
 URL: http://www.stat.auckland.ac.nz/~rolf/
-Packaged: 2014-01-31 03:12:30 UTC; rolf
+Packaged: 2014-05-19 07:00:28 UTC; rolf
 NeedsCompilation: yes
 Repository: CRAN
-Date/Publication: 2014-01-31 08:24:18
+Date/Publication: 2014-05-19 09:21:01

MD5

@@ -1,5 +1,5 @@
-a628dc9137e7a9200a19d13ed016fa56 *ChangeLog
-962e3601155784f197faf1d8146a981f *DESCRIPTION
+b9e485939776a83e288b2a8c9b59173c *ChangeLog
+0e8d3e7715647a85572d5dc34b7ee9cb *DESCRIPTION
 e97c7b6e0d7f0eb3b8e1303451431001 *NAMESPACE
 bb6f6344934b965a455a21347756e847 *R/First.R
 33b7623305e837aa24a6407e2484ab9f *R/RCS/First.R,v
@@ -20,22 +20,24 @@ cb8b173f762a5c70002fa1fd88d8b8b6 *R/RCS/revise.rho.R,v
 7bec57092db176335935eac4c4e81726 *R/RCS/sp.R,v
 7138d81efb770e4b4d3c27538f2476e6 *R/RCS/tidyup.R,v
 830058e16c6db708863ab4ebc518a58d *R/RCS/viterbi.R,v
-90355144d471e120fc01fd1685cf445f *R/check.yval.R
-ab1a44feb53f6bb9d961ef3b2d731638 *R/ffun.R
-d9b85d50148cd101bed871410ac73a0d *R/fitted.hmm.discnp.R
-5b21b23625b351b6ec9361db4e161eba *R/hmm.R
+c9e5308e5060f990e650c5cedd4589f4 *R/charList.R
+6c8678115742b85e6ecfef77c5a7b7d8 *R/check.yval.R
+ae212cc9cb39095ad627a889b3786014 *R/ffun.R
+135777d19d976b4ba921d83c8b092b27 *R/fitted.hmm.discnp.R
+0092335aa6f73da636741dd8acfe49fd *R/hmm.R
 eb22d9be687321b5243b12879f272cee *R/init.all.R
-da3f83ad5ee8b659b6630ca9fcbf380b *R/logLikHmm.R
-de0b8ed3091fa489ada6c6e5355b1856 *R/mat2list.R
-b92e3051eafb9511f4cae92645158ba8 *R/mps.R
-0dfefa748fd129b0dbe03289e8c45f35 *R/pr.R
+927cbdf60ee88f3bab5fbc46866aead5 *R/logLikHmm.R
+d9d18d44d2bea87b6e9bb428bece75bb *R/mps.R
+8bccab6ebff4d659e00deab28535f46b *R/pr.R
 dc97c65ddc21eacd17cacdf3b7c93993 *R/recurse.R
 fced2f094ff5a74046a258e8463b590d *R/revise.ispd.R
-453a3658f3432a9b1b1feda4400134d7 *R/revise.rho.R
+6003184b9be779bbdd42cecc7412b034 *R/revise.rho.R
 f180f9d05c14a25f746e32f579875ff1 *R/revise.tpm.R
 1343dfe5864d155602c1457e62891944 *R/sim.hmm.R
-4c2d789f7ee7157a3f4c66fadbe74a7a *R/sp.R
-1195febe2765cea0cc092eff470fa11b *R/viterbi.R
+76c005d4a086f5a735d0a4ffb92b4eca *R/sp.R
+10e4606077501a01eb153b505f2fe808 *R/viterbi.R
+edd536ed6ec8824fe8abfee411b640a3 *data/colifCount.rda
+ddd1dbb89d36f705ee69dae474e9c072 *data/lesionCount.rda
 2656ac33803f83f5d131053efc7c4b6b *inst/READ_ME
 dd41961396b612dbb4afc7fdcf73e55e *inst/Ratfor/RCS/afun.r,v
 75869a842561e24d282d8ed2fc73eab8 *inst/Ratfor/RCS/bfun.r,v
@@ -49,11 +51,13 @@ f1d36db9f8e9fda32f12221e55c6dbdd *inst/Ratfor/RCS/gfun.r,v
 40201e8652101f856140d17279cf3f4c *inst/Ratfor/makefor
 41395d57b572e9fac85b6616aa3853e0 *inst/Ratfor/recurse.r
 b97264d24cd7880ce7fbad579da6953b *inst/Ratfor/xfun.r
-e7adcaa1233f0f2fb31efc55874a7edb *man/fitted.hmm.discnp.Rd
-bd10aec7610def882c11a6bd29eb6e88 *man/hmm.Rd
-1803bd4405ea01fd41683f62c18e447c *man/hmm.discnp-internal.Rd
-a7e7164497cf88e5f2338028b7949b74 *man/logLikHmm.Rd
-d9194a8fc5f40590705782296fd9dafc *man/mps.Rd
+1bdd1b4d40f5938a949e43adf2985c3b *man/colifCount.Rd
+5df37caa85124c656e2ed3ba747602b3 *man/fitted.hmm.discnp.Rd
+bfcea043173e4ba815cc40a7ba867247 *man/hmm.Rd
+d78a62e221d2d12092c9421137c40727 *man/hmm.discnp-internal.Rd
+380b1a58283d1143ada76b6e97f0af93 *man/lesionCount.Rd
+34b9623d3bfada4508602ef8089926ee *man/logLikHmm.Rd
+9bc1f1e14d8aca22b2abda07aff194da *man/mps.Rd
 56d004e4901d847eaf69c2359e302c5e *man/pr.Rd
 8b52a1ca1d73021d3f0f7f228ab8ed9e *man/sim.hmm.Rd
 92f7052aee5abd5d239040a33db46cde *man/sp.Rd

R/mat2list.R → R/charList.R

@@ -1,4 +1,4 @@
-mat2list <- function(y) {
+charList <- function(y) {
 if(is.matrix(y)) {
 	warning(paste("Presenting \"y\" as a matrix is deprecated.\n",
                       "Change to presenting \"y\" either as a vector\n",
@@ -12,5 +12,8 @@ if(is.matrix(y)) {
 		      })
 }
 if(!is.list(y)) y <- list(y)
+uval <- as.character(sort(unique(unlist(y))))
+y    <- lapply(y,as.character)
+attr(y,"uval") <- uval
 y
 }

R/check.yval.R

@@ -1,17 +1,20 @@
-check.yval <- function(y,Rho) {
+check.yval <- function(y,Rho,warn=TRUE) {
 yval  <- unique(unlist(y))
+rn    <- rownames(Rho)
 fname <- as.character(sys.call(-1))[1]
 if(is.na(fname)) fname <- "call from the command line"
-if(is.null(row.names(Rho))) {
-	yval <- as.numeric(yval)
-	OK   <- all(yval%in%(1:nrow(Rho)))
-	if(!OK) stop(paste("In ",fname,".  The values of \"y\" must be in ",
-                     "\"1:nrow(Rho)\".\n",sep=""),call.=FALSE)
-} else {
-	yval <- as.character(yval)
-	OK   <- all(yval%in%row.names(Rho))
-	if(!OK) stop(paste("In ",fname,".  The values of \"y\" must be in ",
-                     "\"row.names(Rho)\".\n",sep=""),call.=FALSE)
+if(!is.null(rn)) {
+  if(all(yval %in% rn)) return(Rho)
+  stop(paste("In ",fname," some y values do not match the row names",
+             "of \"Rho\".\n",sep=""),call.=FALSE)
 }
-return(invisible())
+if(length(yval) != nrow(Rho))
+   stop(paste("In ",fname," wrong number of rows in \"Rho\".\n",sep=""),call.=FALSE)
+whinge <- paste("Matrix \"Rho\" has no row names.  Assuming that the\n",
+                "rows of Rho correspond to the sorted unique values of \"y\".\n")
+if(warn) warning(whinge)
+nval <- as.numeric(yval)
+yval <- if(!any(is.na(nval))) yval[order(nval)] else sort(yval)
+rownames(Rho) <- yval
+return(Rho)
 }

R/ffun.R

@@ -1,4 +1,16 @@
-ffun <- function(y,Rho)
+ffun <- local({
+
+fixNA <- function(x,Rho) {
+  if(any(is.na(x))) {
+    rn <- c(rownames(Rho),"missing")
+    Rho  <- rbind(Rho,rep(1,ncol(Rho)))
+    rownames(Rho) <- rn
+    x[is.na(x)] <- "missing"
+  }
+  Rho[x,1:ncol(Rho)]
+}
+
+function(y,Rho)
 {
 #
 # Function ffun to calculate f(x) = Pr(Y=x | the model parameters)
@@ -8,13 +20,8 @@ ffun <- function(y,Rho)
 # correspond to the observations y.
 #
 
-if(is.null(rownames(Rho))) {
-	y <- lapply(y,as.numeric)
-} else {
-	y <- lapply(y,as.character)
-}
-fy <- lapply(y,function(x,Rho){Rho[x,1:ncol(Rho)]},Rho=Rho)
+fy <- lapply(y,fixNA,Rho=Rho)
 fy <- do.call(rbind,fy)
-fy[is.na(fy)] <- 1
 t(fy)
 }
+})

R/fitted.hmm.discnp.R

@@ -1,7 +1,10 @@
 fitted.hmm.discnp <- function(object,...) {
 y <- object$y
 if(is.null(y)) stop("Observations \"y\" were not kept.\n")
+if(!is.list(y)) y <- list(y)
 if(!is.numeric(y[[1]]))
-	stop("Observations are not numeric; fitted values make no sense.\n")
+if(!all(sapply(y,is.numeric)))
+	stop(paste("Some observations are not numeric;\n",
+                   "fitted values make no sense.\n"))
 sp(y,object,means=TRUE)$means
 }

R/hmm.R

@@ -27,13 +27,13 @@ if(stationary & !cis)
 # Put together a data name tag for the output.
 if(is.null(data.name)) data.name <- deparse(substitute(y))
 
-# If y is a matrix, change it to a list, and put out a
-# snarky message to the user.
-y <- mat2list(y)
+# Change y into *character* data.  If y is a matrix, change it
+# to a list, and put out a snarky message to the user.
+y <- charList(y)
 
 # Check that the observation values are compatible
 # with yval if it is specified.
-uval <- sort(unique(unlist(y)))
+uval <- attr(y,"uval")
 if(is.null(yval)) yval <- uval
 if(!all(uval%in%yval))
         stop("Specified y values do not include all observed y values.\n")
@@ -54,7 +54,15 @@ else {
 		stop(paste("Row dimension of \"Rho\" not equal to\n",
                            "the number of distinct y-values.\n"))
 }
-row.names(par0$Rho) <- yval
+if(is.null(row.names(par0$Rho))) {
+    row.names(par0$Rho) <- yval
+} else {
+    if(!all.equal(yval,row.names(par0$Rho))) {
+        whinge <- paste("The row names of the initial value of \"Rho\" are not\n",
+                        "equal to the list of possible y-values, \"yval\".\n")
+        stop(whinge)
+    }
+}
 
 # If K=1 do the triv thing:
 if(K==1) {

R/logLikHmm.R

@@ -9,7 +9,7 @@ logLikHmm <- function(y,par) {
 
 # If y is a matrix, change it to a list, and put out a
 # snarky message to the user.
-y <- mat2list(y)
+y <- charList(y)
 
 # Get the parameters.
 Rho  <- par$Rho
@@ -21,7 +21,7 @@ if(is.null(ispd)) {
 
 # Make sure that the entries of the vectors in y correspond
 # to the row names of Rho.
-check.yval(y,Rho)
+Rho <- check.yval(y,Rho)
 
 # If K=1 do the triv thing:
 K <- length(ispd)

R/mps.R

@@ -13,8 +13,8 @@ if(missing(y)) {
 	y <- if(!is.null(object)) object$y else NULL
 	if(is.null(y)) stop("No observation sequence supplied.\n")
 }
-y <- mat2list(y)
-check.yval(y,Rho)
+y    <- charList(y)
+Rho  <- check.yval(y,Rho)
 lns  <- sapply(y,length)
 nseq <- length(y)
 fy   <- ffun(y,Rho)

R/pr.R

@@ -10,8 +10,8 @@ if(missing(y)) {
 	if(is.null(y)) stop("No observation sequence supplied.\n")
 }
 if(is.null(ispd)) ispd <- revise.ispd(tpm)
-y <- mat2list(y)
-check.yval(y,Rho)
+y <- charList(y)
+Rho <- check.yval(y,Rho)
 
 if(!is.list(s)) s <- list(s)
 nseq <- length(s)

R/revise.rho.R

@@ -4,5 +4,7 @@ revise.rho <- function(y,gamma,yval) {
 	t1 <- apply(gamma[,there],1,
 		function(x,index){tapply(x,index,sum)},y[there])
 	t1[is.na(t1)] <- 0
-	t(t(t1)/apply(t1,2,sum))
+	Rho <- t(t(t1)/apply(t1,2,sum))
+        rownames(Rho) <- yval
+        Rho
 }

R/sp.R

@@ -10,20 +10,18 @@ sp <- function (y, object = NULL, tpm, Rho, ispd=NULL, means=FALSE)
 	y <- if(!is.null(object)) object$y else NULL
 	if(is.null(y)) stop("No observation sequence supplied.\n")
     }
-    y <- mat2list(y)
-    check.yval(y,Rho)
-    if(is.null(row.names(Rho))) row.names(Rho) <- 1:nrow(Rho)
+    y    <- charList(y)
+    Rho  <- check.yval(y,Rho)
     lns  <- sapply(y,length)
     fy   <- ffun(y, Rho)
     rp   <- recurse(fy, tpm, ispd, lns)
     prbs <- rp$gamma
     if(means) {
-	yval <- ( if(!is.null(row.names(Rho)))
-                      as.numeric(row.names(Rho)) else 1:nrow(Rho) )
+	yval <- as.numeric(row.names(Rho))
 	if(any(is.na(yval)))
 		stop("Non-numeric y-values; means make no sense.\n")
         cmns <- apply(yval*Rho,2,sum)
-	mns <- apply(cmns*prbs,2,sum)
+	mns  <- apply(cmns*prbs,2,sum)
     }
     nseq <- length(lns)
     if (nseq == 1) {

R/viterbi.R

@@ -16,13 +16,13 @@ if(missing(y)) {
 	y <- if(!is.null(object)) object$y else NULL
 	if(is.null(y)) stop("No observation sequence supplied.\n")
 }
-y <- mat2list(y)
+y <- charList(y)
 
 # Build ispd if it was given as NULL
 if(is.null(ispd)) ispd <- revise.ispd(tpm)
 
 # Make sure that the y-values are compatible with Rho.
-check.yval(y,Rho)
+Rho <- check.yval(y,Rho)
 
 # Make sure y is a list, and get the number of sequences and
 # lengths of these sequences.