National audienceIn this study, we propose a method called sequence breaking framework to search high local maximum of the likelihood by providing starting values based on the observations for the Expectation-Maximization algorithm, for Hidden semi-Markov model parameter estimation. The method is shown to be efficient on several datasets with multiple categorical sequences

Durand, Jean-Baptiste

Guérin-Dugué, Anne

Olivier, Brice

Hal - Université Grenoble Alpes

Assessment of various initialization strategies for the Expectation-Maximization algorithm for Hidden Semi-Markov Models with multiple categorical sequences

INRIA a CCSD electronic archive server

HAL Id: hal-02129122https://hal.inria.fr/hal-02129122Submitted on 14 May 2019HAL is a multi-disciplinary open accessarchive for the deposit and dissemination of sci-entific research documents, whether they are pub-lished or not. The documents may come fromteaching and research institutions in France orabroad, or from public or private research centers.L’archive ouverte pluridisciplinaire HAL, estdestinée au dépôt et à la diffusion de documentsscientifiques de niveau recherche, publiés ou non,émanant des établissements d’enseignement et derecherche français ou étrangers, des laboratoirespublics ou privés.Assessment of various initialization strategies for theExpectation-Maximization algorithm for HiddenSemi-Markov Models with multiple categoricalsequencesBrice Olivier, Anne Guérin-Dugué, Jean-Baptiste DurandTo cite this version:Brice Olivier, Anne Guérin-Dugué, Jean-Baptiste Durand. Assessment of various initialization strate-gies for the Expectation-Maximization algorithm for Hidden Semi-Markov Models with multiple cate-gorical sequences. JdS 2019 - 51èmes Journées de Statistique, Jun 2019, Vandœuvre-lès-Nancy, France.pp.1-7. ￿hal-02129122￿Assessment of various initialization strategiesfor the Expectation-Maximization algorithm forHidden Semi-Markov Models with multiplecategorical sequencesBrice Olivier 1, Anne Guérin-Dugué 2 & Jean-Baptiste Durand 11 Université Grenoble Alpes, Laboratoire Jean Kuntzmann, Inria Mistis, 51 rue desMathématiques, B.P. 53, F-38041 Grenoble cedex 9, France; brice.olivier@inria.fr,jean-baptiste.durand@univ-grenoble-alpes.fr2 Université Grenoble Alpes, CNRS, Gipsa-lab, 11 rue des Mathématiques, B.P. 46,F-38040 Grenoble, France; anne.guerin@gipsa-lab.grenoble-inp.frRésumé. Dans cette étude, nous proposons une nouvelle stratégie d’initialisationde l’algorithme d’Espérence-Maximisation afin d’augmenter la vraisemblance en sortiede l’algorithme, appliqué aux modèles semi-Markoviens cachés en fournissant des valeursinitiales des paramètres calculées à partir des observations. Cette stratégie se révèle efficacesur certains jeux de données à séquences multiples et catégorielles.Mots-clés. Modèles semi-Markoviens cachés, choix des valeurs initiales pour l’algo-rithme d’Espérance-Maximisation, séquences multiples et à valeurs catégorielles.Abstract. In this study, we propose a method called sequence breaking frameworkto search high local maximum of the likelihood by providing starting values based onthe observations for the Expectation-Maximization algorithm, for Hidden semi-Markovmodel parameters estimation. The method is shown to be efficient on several datasetswith multiple categorical sequences.Keywords. Hidden semi-Markov models, choosing starting values for the EM algo-rithm, multiple categorical sequences.1 IntroductionSince introduced by Ferguson (1980), hidden semi-Markov models (HSMM) have receiveda lot of attention in the literature as a natural and more flexible extension of the well-known Hidden Markov model (HMM) by relaxing the geometric state sojourn time as-sumption. The HSMM framework being very generic, modeling assumptions have beenproposed, mainly focused around the dependencies between the past state and currentsojourn time, see Yu (2010) for a summary, as well as the latent process time censoringGuédon (2003), leading to a wide variety of inference algorithms. Parameter estimationis generally performed using the Expectation-Maximization algorithm, McLachlan andKrishnan (2007). Barbu and Limnios (2009) have proved the asymptotic convergence andnormality of the estimators, but did not provide any detail on the convergence speed or1the multiple sequence framework. The Expectation-Maximization algorithm finds a localmaximum of the likelihood and is known to be extremely sensitive to starting values.In this context, we propose a method that we call sequence breaking framework, whichaims at finding high local maxima of the likelihood by choosing starting values for HSMM’sEM, for which the randomness is controlled by the observed sequences in order to restrictthe search space. To our knowledge, most of the contributions of this kind have been donearound the independent Mixture Models, see Biernacki et al. (2003) for Gaussian or Juanet al. (2004) for Bernoulli. While for HMMs, what seems to work best is a simple jitterof the parameters around their centers as implemented in the python library hmmlearn1,method which is here compared to ours.2 Hidden semi-Markov modelAn HSMM is composed of two stochastic processes. The former is a finite-state homo-geneous semi-Markov chain (SMC) which is latent. It conditions the latter which is theobserved process. An SMC is like a Markov Chain except that the within-state sojourntime is not necessarily geometric and can therefore be of any form with support includedin N (categorical or parametric).In this section, we focus on the simplest assumption on the HSMM also called theExplicit Duration HMM (EDHMM) for which the within-sojourn state duration is onlyconditioned by the state. Using the formalism proposed by Murphy (2012) which isDynamic Bayesian Network-oriented (DBN), we describe a SMC by {St, Rt, Ft}t∈J1,T K,∀t ∈J1, T K St ∈ J1,MK, Rt ∈ J1, DK and Ft ∈ {0, 1}. St is the state value while Rt encodesthe residual time in the current state St at time t. At the beginning of a new state, anew duration is randomly sampled from an arbitrary distribution and then counts downdeterministically to 1, and so on. Ft simply acts as a binary switch which is turned onwhen Rt−1 = 1 and off else. Initialization is performed s.t. F0 = 1 and R0 = 1, and wealso have FT = 1 which means that the process starts at time 1 and will end at time T .Follows the parameters associated to the conditional probability distributions associatedto a SMC:P (S1 = k) = πkwith πk ∈ J1,MK, a vector of size M representing all the initial probabilities, thenP (St = k|St−1 = k′, Ft−1 = f) ={1{k = k′} if f = 0Ak′k if f = 1with 1 denoted as the indicator function, and A ∈ M ×M , a matrix representing thetransition probabilities from one state to another with zeros on the diagonal to forbid1https://hmmlearn.readthedocs.io2self-transitions. We also haveP (Rt = d|Rt−1 = d′, St = k, Ft−1 = 1) = pk(d)with pk(d) being an arbitrary, upper-bounded by D, probability distribution representingthe sojourn distributions for each state k while entering a new state at time t and thensampling a new value d ≥ 1 for Rt. Finally,P (Rt = d|Rt−1 = d′, St = k, Ft−1 = 0) ={1{d = d′ − 1} if d > 1undefined if d = 1andP (Ft = f |Rt = d) ={1{d = 1} if f = 11{d > 1} if f = 0define the countdown process, i.e. the residual time in the current state.To sum up, the process can be described as the transition from a latent state to anotherat time t triggering the following changes: the finishing node switches on Fj−1 = 1,requiring a transition to a new state, St = k, from the previous one, St−1 = k′, withk′ 6= k, along with its new duration sampled s.t. Rt ∼ pk ≥ 1.The discrete observed process {Xt}t∈J1,T K takes values in {v1, ..., vG} and is conditionedby {St}t∈J1,T K s.t.P (Ot = vg|St = k) = bk(vg)where bk(vg) represent either a categorical distribution, i.e. a matrix of size M ×G, or Mparametric distributions. Finally, we denote θ the set of all model parameters. Note thathere, we treat the discrete observation case but it is generalizable when observations arecontinuous.3 Experimental strategyMultiple sequence framework. So far, we have written down the EDHMM consider-ing only one sequence of observation {Xt}t∈J1,T K for notation convenience. We now considerthat we have multiple observed sequences X = {{X(1)t }t∈J1,T1K, ..., {X(K)t }t∈J1,TKK}.Proposition 1. For an observed semi-Markov chain {St, Rt, Ft}t∈J1,T K with a correspond-ing output process {Xt}t∈J1,T K, the MLE of each set of parameters is given by the condi-tional empirical frequencies.Proof. Given that an EDHMM can be expressed as a DBN, see section 2, and since thelikelihood of a DBN can be expressed as a global decomposition of the local-likelihood ofeach node given its parents, see Ghahramani (2001), then the parameters of the EDHMMcan be locally optimized by MLE using the conditional empirical frequencies.3Choosing starting values with the sequence breaking framework. The proposedstrategy relies on two intertwined algorithms:• Algorithm 1 HighLikelihoodSearch: describes the global framework, it randomlychooses α observed sequences X(Qα) from X, generates the corresponding state se-quences S(Qα) using SequenceBreaking, computes the parameters θinit by MLEusing Proposition 1 and injects it as a starting value for the EM algorithm whichfinds the a local maximum of the likelihood for all data X. The goal of samplingsequences randomly from X is to generate starting values related to the observationprocess while keeping only a subset to maintain the randomness of the starting val-ues. Note that Sample(.) is a function which samples uniformly on the given set.Algorithm 1: HighLikelihoodSearch: High local maximum of the likelihoodsearch by sequence breakingInput: α ∈ J1, KK the number of sequence to sample,N , the number of initialization1 θ̂ ← ∅;2 for n← 0 to N do3 Sample X(Qα) ⊂ X observed sequences s.t. Qα ⊂ J1, KK;4 S(Qα) ← SequenceBreaking(X(Qα));5 θinit ← argmaxθ L(θ; {X(Qα),S(Qα)}); # MLE provided by Proposition 16 θ̂ ← θ̂ ∪ ExpectationMaximization(θinit,X);7 end8 θ̂∗ ← argmaxθ̂ L(θ̂;X)Output: θ̂∗, a high local maximum of the log-likelihood.• Algorithm 2 SequenceBreaking: randomly generates a hidden state sequence.Given each observed sequence X(q) ∈ X(Qα) with its length, it randomly chooses anumber of transitions J as well as transition instants I which "break" the sequencesinto pieces, and then affects a state randomly to each piece of sequence with theconstraints that two consecutive states should be different due to the EDHMMassumption on the transition matrix.Choosing starting values with jittered-center parameters. We compare the pro-posed methodology with the default strategy used for HMM which consists in selectingslightly perturbed parameters around their centers, i.e. s.t. each event has equal proba-bility, similarly to Juan et al. (2004).Example 1. With M = 2, we wish to randomly sample π s.t. π1 = Sample([ε, 1 − ε])and π2 = 1− π1, with ε ∈]0, 0.5].4Algorithm 2: SequenceBreakingInput: X(Qα), an observed sequence subset of size α1 for q ∈ Q(α) do2 J ← Sample(J1, Tq − 1K); # number of transitions3 I ← ∅;4 for j ← 0 to J do5 i← Sample(J1, TqK) s.t. i 6∈ I; # transition instant6 I ← I ∪ i7 {S(q)t }t∈JIj−1,iK ← Sample(J1,MK) s.t. S(q)Ij−16= S(q)Ij−1−1; # choose state8 end9 {S(q)t }t∈Ji,TqK ← Sample(J1,MK) s.t. S(q)i 6= S(q)Ij−1−1; # choose final state10 endOutput: S(Qα), a randomly sampled state sequencesExample 1 is a specific instance with a Bernoulli distribution. In order to generalizewith M ≥ 2, i.e. to the Multinomial distribution, one solution is to sample from itsconjugate, the Dirichlet distribution. Therefore, we apply a Dirichlet sample for eachset of parameters except the sojourn distribution which we initialize with a Geometricdistribution of parameter p = 0.1.4 ResultsNot only do the experiments consist of numerical comparison of both methods in findingthe highest likelihood, but also to compare the converge speed of EM, for three differentdatasets, artificial, artificial with noise, and real.Datasets. The first dataset D(a) is artificial, composed of 100 sequences of length100 each, with K = 5 clusters and G = 5 observed parameters, and parameters π =(0.2, 0.2, 0.2, 0.2, 0.2), p(.) = (G(0.2),G(0.05),NB(8, 0.5),P(4),NB(5, 0.1)) where G standsfor the Geometric distribution, NB Negative Binomial and P Poisson,A =0 0.5 0.3 0.1 0.10.5 0 0.1 0.3 0.10.25 0.25 0 0.25 0.250.2 0.2 0.2 0 0.40.1 0.1 0.4 0.4 0 B =0.1 0.2 0.4 0.2 0.10.25 0.2 0.1 0.2 0.250.2 0.2 0.2 0.2 0.20.35 0.3 0.2 0.1 0.050.05 0.1 0.2 0.3 0.35The second dataset D(an) is generated from the first one by replacing 20% of its observa-tions at random. The third dataset D(r) consists of the bounded number of words skippedduring an ocular saccade by different subjects for a reading tasks for which G = 5 and5we assume M = 5. There are 2390 sequences of different lengths, on average 17 with astandard deviation of 8, see Durand et al. (2016) for more information about this dataset.Global results. For a global analysis, we compute the mean and standard deviationof all optimal likelihoods for each method and dataset over 100 initializations and 1000iterations of EM each. Results are presented in table 1.D(a) D(an) D(r)Sequence Breaking −15448± 2.3 −15785± 4.5 −50567± 236Jittered-centers −15452± 2.9 −15782± 2.6 −50592± 274Table 1: Means and standard deviations of maximum likelihood. Significant (<5%) meandifferences are boldfaced.Local results. For a local analysis, we split the 100 initializations into 10 blocks of 10and compute the max per block. For D(a), sequence breaking performed better than thejittered-centers 9 times out of 10. For D(an), 9/10, and for D(r), 6/10.Convergence speed. For D(a), on average, it took the sequence breaking initialization133 less iterations to convergence than the jittered-centers. For D(an), it took 305 less,and for D(r), 44 less iterations on average.5 Discussion and perspectivesWe proposed a new strategy to search for improved maximum likelihood of HSMM withmultiple sequences categorical data which is a significant improvement considering thatthe current Python implementation of HSMM under the virtual plants library sequenceanalysis2, as well as the R package mhsmm by O’Connell et al. (2011) do not provide arandom start option.We can address a few remarks on this work. First, α, the number of sequences tosample can be seen as a control of the randomness of the initial parameters since the higherit is, the more the observed process will be close to the real one. Secondly, we would like tonote the similarity of the procedure with the stick-breaking process view of the Dirichletprocess. Thirdly, as shown by Meilă and Heckerman (2001), initialization techniques aredata dependent and we still need to inspect more datasets but also more techniques suchas Stochastic EM or Classification EM for initialization, to find out the relevance of thesequence breaking framework. Another interesting perspective is to propose a heuristicfor early detection of bad candidates in order to save processing time. Finally, we planon testing the strategy on more datasets and on continuous data as well.2https://github.com/openalea/StructureAnalysis6ReferencesBarbu, V. S. and Limnios, N. (2009). Semi-Markov chains and hidden semi-Markov modelstoward applications: their use in reliability and DNA analysis, volume 191. SpringerScience & Business Media.Biernacki, C., Celeux, G., and Govaert, G. (2003). Choosing starting values for the emalgorithm for getting the highest likelihood in multivariate gaussian mixture models.Computational Statistics & Data Analysis, 41(3-4):561–575.Durand, J.-B., Guérin-Dugué, A., and Achard, S. (2016). Analyse de séquences ocu-lométriques et d’électroencéphalogrammes par modèles markoviens cachés. In 48èmesJournées de Statistique.Ferguson, J. (1980). pp. 143–179, variable duration models for speech. In Proc. of theSymposium on the applications of hidden Markov models to text and speech, JD Fergu-son, Ed. Princeton: IDA-CRD.Ghahramani, Z. (2001). An introduction to hidden markov models and bayesian networks.International journal of pattern recognition and artificial intelligence, 15(01):9–42.Guédon, Y. (2003). Estimating hidden semi-markov chains from discrete sequences. Jour-nal of Computational and Graphical Statistics, 12(3):604–639.Juan, A., García-Hernández, J., and Vidal, E. (2004). Em initialisation for bernoullimixture learning. In Joint IAPR International Workshops on Statistical Techniques inPattern Recognition (SPR) and Structural and Syntactic Pattern Recognition (SSPR),pages 635–643. Springer.McLachlan, G. and Krishnan, T. (2007). The EM algorithm and extensions, volume 382.John Wiley & Sons.Meilă, M. and Heckerman, D. (2001). An experimental comparison of model-based clus-tering methods. Machine learning, 42(1-2):9–29.Murphy, K. P. (2012). Machine learning: a probabilistic perspective. MIT press.O’Connell, J., Højsgaard, S., et al. (2011). Hidden semi markov models for multipleobservation sequences: The mhsmm package for r. Journal of Statistical Software,39(4):1–22.Yu, S.-Z. (2010). Hidden semi-markov models. Artificial intelligence, 174(2):215–243.7

https://hal.inria.fr/hal-02129122/document

Assessment of various initialization strategies for the Expectation-Maximization algorithm for Hidden Semi-Markov Models with multiple categorical sequences

Abstract

Similar works

Full text

Available Versions

Hal - Université Grenoble Alpes

INRIA a CCSD electronic archive server