Description Usage Arguments Value Author(s) References See Also Examples

Constructor for a a multivariate hidden Markov model (HMM)
where each of the `n`

variables observed at the same time has a
(potentially different) standard univariate distribution conditionally
on the underlying state. The `n`

outcomes are independent
conditionally on the hidden state.

If a particular state in a HMM has such an outcome distribution, then
a call to `hmmMV`

is supplied as the corresponding element
of the `hmodel`

argument to `msm`

. See Example 2 below.

A multivariate HMM where multiple outcomes at the same time are
generated from the *same* distribution is specified in the same
way as the corresponding univariate model, so that `hmmMV`

is not required. The outcome data are simply supplied as a matrix
instead of a vector. See Example 1 below.

The outcome data for such models are supplied as a matrix, with
number of columns equal to the maximum number of arguments supplied to
the `hmmMV`

calls for each state. If some but not all of
the variables are missing (`NA`

) at a particular time, then the
observed data at that time still contribute to the likelihood. The
missing data are assumed to be missing at random. The Viterbi
algorithm may be used to predict the missing values given the fitted
model and the observed data.

Typically the outcome model for each state will be from the same family or set of families, but with different parameters. Theoretically, different numbers of distributions may be supplied for different states. If a particular state has fewer outcomes than the maximum, then the data for that state are taken from the first columns of the response data matrix. However this is not likely to be a useful model, since the number of observations will probably give information about the underlying state, violating the missing at random assumption.

Models with outcomes that are dependent conditionally on the hidden state (e.g. correlated multivariate normal observations) are not currently supported.

1 |

`...` |
The number of arguments supplied should equal the maximum
number of observations made at one time. Each argument represents
the univariate distribution of that outcome conditionally on the
hidden state, and should be the result of calling a univariate
hidden Markov model constructor (see |

A list of objects, each of class `hmmdist`

as returned by
the univariate HMM constructors documented in `hmm-dists`

.
The whole list has class `hmmMVdist`

, which inherits from
`hmmdist`

.

C. H. Jackson [email protected]

Jackson, C. H., Su, L., Gladman, D. D. and Farewell, V. T. (2015) On modelling minimal disease activity. Arthritis Care and Research (early view).

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 | ```
## Simulate data from a Markov model
nsubj <- 30; nobspt <- 5
sim.df <- data.frame(subject = rep(1:nsubj, each=nobspt),
time = seq(0, 20, length=nobspt))
set.seed(1)
two.q <- rbind(c(-0.1, 0.1), c(0, 0))
dat <- simmulti.msm(sim.df[,1:2], qmatrix=two.q, drop.absorb=FALSE)
### EXAMPLE 1
## Generate two observations at each time from the same outcome
## distribution:
## Bin(40, 0.1) for state 1, Bin(40, 0.5) for state 2
dat$obs1[dat$state==1] <- rbinom(sum(dat$state==1), 40, 0.1)
dat$obs2[dat$state==1] <- rbinom(sum(dat$state==1), 40, 0.1)
dat$obs1[dat$state==2] <- rbinom(sum(dat$state==2), 40, 0.5)
dat$obs2[dat$state==2] <- rbinom(sum(dat$state==2), 40, 0.5)
dat$obs <- cbind(obs1 = dat$obs1, obs2 = dat$obs2)
## Fitted model should approximately recover true parameters
msm(obs ~ time, subject=subject, data=dat, qmatrix=two.q,
hmodel = list(hmmBinom(size=40, prob=0.2),
hmmBinom(size=40, prob=0.2)))
### EXAMPLE 2
## Generate two observations at each time from different
## outcome distributions:
## Bin(40, 0.1) and Bin(40, 0.2) for state 1,
dat$obs1 <- dat$obs2 <- NA
dat$obs1[dat$state==1] <- rbinom(sum(dat$state==1), 40, 0.1)
dat$obs2[dat$state==1] <- rbinom(sum(dat$state==1), 40, 0.2)
## Bin(40, 0.5) and Bin(40, 0.6) for state 2
dat$obs1[dat$state==2] <- rbinom(sum(dat$state==2), 40, 0.6)
dat$obs2[dat$state==2] <- rbinom(sum(dat$state==2), 40, 0.5)
dat$obs <- cbind(obs1 = dat$obs1, obs2 = dat$obs2)
## Fitted model should approximately recover true parameters
msm(obs ~ time, subject=subject, data=dat, qmatrix=two.q,
hmodel = list(hmmMV(hmmBinom(size=40, prob=0.3),
hmmBinom(size=40, prob=0.3)),
hmmMV(hmmBinom(size=40, prob=0.3),
hmmBinom(size=40, prob=0.3))),
control=list(maxit=10000))
``` |

```
Call:
msm(formula = obs ~ time, subject = subject, data = dat, qmatrix = two.q, hmodel = list(hmmBinom(size = 40, prob = 0.2), hmmBinom(size = 40, prob = 0.2)))
Optimisation probably not converged to the maximum likelihood.
optim() reported convergence but estimated Hessian not positive-definite.
-2 * log-likelihood: 3706.02
[Note, to obtain old print format, use "printold.msm"]
Warning message:
In msm(obs ~ time, subject = subject, data = dat, qmatrix = two.q, :
Optimisation has probably not converged to the maximum likelihood - Hessian is not positive definite.
Call:
msm(formula = obs ~ time, subject = subject, data = dat, qmatrix = two.q, hmodel = list(hmmMV(hmmBinom(size = 40, prob = 0.3), hmmBinom(size = 40, prob = 0.3)), hmmMV(hmmBinom(size = 40, prob = 0.3), hmmBinom(size = 40, prob = 0.3))), control = list(maxit = 10000))
Maximum likelihood estimates
Transition intensities
Baseline
State 1 - State 1 -0.09458 (-0.13940,-0.06416)
State 1 - State 2 0.09458 ( 0.06416, 0.13940)
Hidden Markov model, 2 states
State 1
Outcome 1 - binomial distribution
Estimate LCL UCL
size 40.0000000 NA NA
prob 0.1054793 0.0948457 0.1171507
Outcome 2 - binomial distribution
Estimate LCL UCL
size 40.0000000 NA NA
prob 0.1958908 0.1818942 0.2106871
State 2
Outcome 1 - binomial distribution
Estimate LCL UCL
size 40.0000000 NA NA
prob 0.5954564 0.5780108 0.612664
Outcome 2 - binomial distribution
Estimate LCL UCL
size 40.0000000 NA NA
prob 0.5110363 0.4933762 0.528669
-2 * log-likelihood: 1523.652
[Note, to obtain old print format, use "printold.msm"]
```

msm documentation built on Feb. 3, 2018, 1:04 a.m.

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.