user12702896
user12702896

Reputation:

depmixs4 list of state estimates

I have fitted a time series of data to a Hidden Markov model using the fit(depmix()) functions in the package depmixs4. I want to obtain a time series of estimated states, where the estimate for state x is the mean value assigned to state x. Currently, I'm only getting a times series of values that say the index of the state (e.g. 1,3,5,2,5...).

This is my current code:

set.seed(9)

hmm9 <- depmix(volume ~ 1, data=data.frame(volume), nstates=9)

fitted_hmm9 <- fit(hmm9)

summary(fitted_hmm9)

state_ests9 <- posterior(fitted_hmm9)

state_ests9[,1]

The last part, state_ests9[,1], is the time series of state indices, while the actual expected values for states are stored somewhere in summary(fitted_hmm9).

Upvotes: 1

Views: 655

Answers (1)

mspeek
mspeek

Reputation: 176

By "estimated state" I assume you want the predicted response (i.e. "volume") according to each state. To get state-wise predicted values for responses, you can use the predict method, which however needs to be called on the response sub-model, which is somewhat hidden in the fitted depmix object.

In a fitted depmix object, there is a slot called "response", which is a list of lists, structured as

fitted_model@response[[state_id]][[response_id]]

In your case, I think "volume" is univariate, so that there is a single response and response_id is always 1. For a model with 9 states, state_id would vary from 1 to 9.

The following code (including randomly generated values for "volume" to make it reproducible) gives you what you want (I think):

set.seed(123)
volume <- rnorm(10000)

hmm9 <- depmix(volume ~ 1, data=data.frame(volume), nstates=9)

fitted_hmm9 <- fit(hmm9)

summary(fitted_hmm9)

state_ests9 <- posterior(fitted_hmm9)

state_ests9[,1]

# construct matrix for state-dependent predictions
pred_resp9 <- matrix(0.0,ncol=9,nrow=nrow(state_ests9))
# fill matrix column-wise by using the "predict" method of
# corresponding response model
for(i in 1:9) {
  pred_resp9[,i] <- predict(fitted_hmm9@response[[i]][[1]])
}
## use MAP state estimates to extract a single time-series
## with predictions
pred_resp9[cbind(1:nrow(pred_resp9),state_ests9[,1])]

Upvotes: 2

Related Questions