Estimation of mixed hidden Markov models with SAEM. Application to daily seizures data.
M. Delattre (1), R. Savic (2), R. Miller (3), M. O. Karlsson (4), M. Lavielle (5)
(1) University of Paris-Sud; (2) INSERM U738; (3) Pfizer Global R & D; (4) Uppsala University; (5) INRIA Saclay ╬le-de-France
Objectives: In some specific medical contexts, the values of biological markers at successive time points are the only informations available to assess the seriousness of a given pathology in patients. Considering a unique sequence of observations, hidden Markov models (HMM) are thus a particularly relevant modeling tool. In those models, the different presupposed disease stages are treated as a Markov process with finite state space and memory one. Such models also allow a correct handle on the dependency between consecutive observations.
When the data to be described include several individuals, specific care is needed to account correctly for the between-subjects heterogeneity. Mixed effects hidden Markov models (MHMM) have been recently developped  as an extension of hidden Markov models to population studies. In our area, mixed hidden Markov models would provide an accurate description of longitudinal data collected during certain clinical trials, especially when distinct (hidden) disease stages are supposed to condition the distribution of some biological markers. Those particular models are quite easily interpretable and could even show similarities in the biological process that governs certain pathologies.
Mixed hidden Markov models include several levels of definition. Assume we have at our disposal observations from n subjects, which respective distributions could reasonnably be supposed to be driven by an underlying Markov chain. First, a hidden Markov model is put on the observations of each of the n subjects. Each individual model is thus specified by its own transition probabilities and its own emission probabilities. Second, those individual parameters are given a common probability distribution. The parameters of this shared distribution, also called population parameters, give access to the mean tendency of the examined phenomenon and capture the potential heterogeneity in the population studied.
Our work mainly aimed at developping and evaluating a complete methodology for estimating parameters in those new models. Our algorithms were applied in the clinical context of epilepsy, to model daily seizure counts in epileptic patients and to assess the effects of a given anti-epileptic drug on the evolution of the epileptic symptoms.
Methods: Making inference on mixed hidden Markov models is a challenging issue. We need to interest in three consecutive angles. The MHMM's population parameters have to be estimated to allow next the estimation of the individual parameters and the decoding of the the most likely sequence of hidden states for each subject.
The maximum likelihood approach is often chosen in practice to estimate the population parameters. However, in addition to their highly non linear structure, mixed hidden Markov models show similarities with incomplete data models. Indeed, both the (random) individual parameters and the hidden sequences of visited states could be considered as "missing" data. As a consequence, the likelihood has a complex expression, and locating its maximum is directly intractable. In a classical HMM, where only emissions are given, the likelihood is difficult to express also, but the Baum-Welch algorithm makes us able to compute it quickly. We consequently suggest estimating the population parameters of mixed hidden Markov models by combining the SAEM algorithm with the Baum-Welch algorithm. Then, the individual parameters for each subject's HMM are estimated using the MAP (Maximum A Posteriori) approach. The estimates for the individual parameters incorporate all the prior information on the data. Therefore, each individual HMM can be considered separately, and the Viterbi algorithm can finally be computed to decode the optimal sequence of hidden states for each subject.
The evaluation of the estimation properties was based on Monte Carlo studies, especially focusing on the performances of the SAEM algorithm.
An application on a real dataset followed. The data coming from a double-blind, placebo-controlled, parallel-group and multicenter study consisted of daily seizure counts collected in epileptic patients during 12 weeks screening phase and 12 weeks treatment phase. A placebo/drug model was suggested using mixed hidden Markov models. For that purpose, a two state Poisson MHMM was built, assuming the epileptic patients go through periods of low and high seizure susceptibility . The treatment dose was included as a covariate at both transition and emission levels in the model to identify clearly the treatment effects on epileptic symptoms.
Our analysis were performed using Monolix and Matlab programs.
Results: First, the good behavior of the SAEM algorithm was a very encouraging result. The convergence was clear and fast. Then, based on the Monte Carlo studies, the population estimates were close to the true values. Indeed, the relative estimation errors (REE) were computed and showed small ranges for the estimates and very little bias. This suggested our algorithm would estimate parameters with a certain accuracy in large databases. Then, the estimated standard errors for each parameter were low.
A first application of mixed hidden Markov models on real data gave good results also. Based on the 788 individual sets of daily seizures in screening phase, a two state Poisson MHMM provided a good description of daily seizures' evolution over time. According to the BIC criteria, Poisson mixed hidden Markov models appeared to be better candidates than Poisson models and mixtures of Poisson for describing epilepsy data. In particular, MHMMs pretty well described the characteristic overdispersion of the data. Moreover, our models mainly showed the drug had a non negligible effect on the Poisson parameters describing the daily seizure counts in each state. To be more precise, the estimations suggested the drug reduces the number of daily seizures in both states of epileptic activity. The estimations also revealed a large interpatient variability at both transition and emission levels.
Conclusion: The algorithms developped for estimating parameters in mixed hidden Markov models appeared to be performant and fast. Based on Monte Carlo studies, the Baum-Welch-SAEM algorithm was shown to provide accurate estimates. The consistency of the maximum likelihood estimates is thus expected, but this point keeps to be studied rigorously by the following.
More generally, mixed hidden Markov models offer very promising statistical applications. In some cases, their structure could even help better understand some disease mechanisms and provide a new way to analyze some drugs' pharmacodynamics. Those new models should thus offer improvement in the analysis of some clinical trials, by envisaging a given treatment could influence not only the mean disease symptoms but the time spent in each disease stage too.
 Altman, Mixed hidden markov models : an extension of the hidden markov model to the longitudinal data setting, Journal of the American Statistical Association (2007).
 Albert, A two state markov mixture model for a time series of epileptic seizure counts, Biometrics 37 (1991).
 Kuhn, Lavielle, Maximum likelihood estimation in nonlinear mixed effects models, Computational Statistics and Data Analysis 49 (2005).
 Savic, Lavielle, Performance in population models for count data : a new SAEM algorithm, Journal Pharmacokinetics Pharmacodynamics 36 (2009).