inst/examples/hmm-multivariate-gaussian-priors.R

In luisdamiano/BayesHMM: Full Bayesian Inference for Hidden Markov Models

mySpec <- hmm(
  K = 3, R = 2,
  observation = # Priors vary per state
    MVGaussian(
      mu    = MVGaussian(mu = c(-10, -10), sigma = matrix(c(100, 0, 0, 100), 2, 2)),
      sigma = matrix(c(1, 0, 0, 1), 2, 2)
    ) +
    MVGaussian(
      mu    = MVGaussian(mu = c(  0,   0), sigma = matrix(c(100, 0, 0, 100), 2, 2)),
      sigma = matrix(c(1, 0, 0, 1), 2, 2)
    ) +
    MVGaussian(
      mu    = MVGaussian(mu = c( 10,  10), sigma = matrix(c(100, 0, 0, 100), 2, 2)),
      sigma = matrix(c(1, 0, 0, 1), 2, 2)
    ),
  initial     = Dirichlet(alpha = c(1, 1, 1)),
  transition  = Dirichlet(alpha = c(1, 1, 1)),
  name = "Multivariate Gaussian"
)

browse_model(mySpec)

set.seed(9000)
y = rbind(
  MASS::mvrnorm(n = 100, mu = c(-10, -10), Sigma = matrix(c(1, 0, 0, 1), 2, 2)),
  MASS::mvrnorm(n = 100, mu = c(  0,   0), Sigma = matrix(c(1, 0, 0, 1), 2, 2)),
  MASS::mvrnorm(n = 100, mu = c( 10,  10), Sigma = matrix(c(1, 0, 0, 1), 2, 2))
)

myModel <- compile(mySpec)
myFit   <- draw_samples(mySpec, stanModel = myModel, y = y, chains = 1, iter = 500, seed = 9000)
print_obs(myFit)