Initializing HMM parameters

[maforinw]

Hi Shaun,

I recently came across HMM to improve outputs from ensemble learners in order to infer animal behaviours. I am particularly interested in your 'aphid' package because I would like to use a list of different sequences (e.g. obtained from direct observations) for parameter tuning (using the Baum-Welch algorithm) before using Viterbi algorithm in a routine process.

The first step is to initialize the parameters of the HMM :

• emission probabilities using classifier performance metrics,
• initial state and state transition probabilities from direct observations (e.g. discrete data).

By using the R package HMM, each parameter initialization is quite obvious. The final object returns a list of 5 elements containing all the parameters.

library(HMM) symbols = paste0(levelBehav[2:7],"E") startProbs = c(0.040, 0.407, 0.077, 0.063, 0.360, 0.053)

transProbs = matrix(c(0.976,    0.009,    0.008,    0.006,    0.000,    0.000,
                      0.001,    0.976,    0.003,    0.003,    0.014,    0.003,
                      0.016,    0.050,    0.895,    0.008,    0.030,    0.002,
                      0.004,    0.040,    0.008,    0.931,    0.015,    0.002,
                      0.000,    0.040,    0.021,    0.007,    0.929,    0.003,
                      0.000,    0.039,    0.001,    0.002,    0.012,    0.945),6)

emissionProbs = matrix(c(0.846,    0.006,    0.068,    0.042,    0.011,    0.027,
                         0.012,    0.807,    0.111,    0.047,    0.008,    0.016,
                         0.103,    0.087,    0.625,    0.083,    0.035,    0.067,
                         0.068,    0.063,    0.104,    0.613,    0.070,    0.082,
                         0.015,    0.002,    0.038,    0.082,    0.843,    0.020,
                         0.095,    0.015,    0.085,    0.105,    0.029,    0.672),6)

hmm<-initHMM(levelBehav[2:7], symbols, startProbs, transProbs, emissionProbs) hmm

$States [1] "event1" "event2" "event3" "event4" "event5" "event6"

$Symbols [1] "event1E" "event2E" "event3E" "event4E" "event5E" "event6E"

$startProbs event1 event2 event3 event4 event5 event6 0.040 0.407 0.077 0.063 0.360 0.053

$transProbs to from event1 event2 event3 event4 event5 event6 event1 0.976 0.001 0.016 0.004 0.000 0.000 event2 0.009 0.976 0.050 0.040 0.040 0.039 event3 0.008 0.003 0.895 0.008 0.021 0.001 event4 0.006 0.003 0.008 0.931 0.007 0.002 event5 0.000 0.014 0.030 0.015 0.929 0.012 event6 0.000 0.003 0.002 0.002 0.003 0.945

$emissionProbs symbols states event1E event2E event3E event4E event5E event6E event1 0.846 0.012 0.103 0.068 0.015 0.095 event2 0.006 0.807 0.087 0.063 0.002 0.015 event3 0.068 0.111 0.625 0.104 0.038 0.085 event4 0.042 0.047 0.083 0.613 0.082 0.105 event5 0.011 0.008 0.035 0.070 0.843 0.029 event6 0.027 0.016 0.067 0.082 0.020 0.672

So my question is how porting objects between the HMM and aphid packages ? Or, how initializing all the parameters, including initial states probabilities in aphid (given that 'residues' means 'symbols') ?

Thank you.

[shaunpwilkinson]

Hi @maforinw Thanks for your question. There are a couple of small differences between the HMM objects, the main one being that the transitions from the begin state and to the end state are included in the first row and column of the transition matrix, respectively. In your case the end state is not modeled, so the first column will just be zeros. Also, the states and residues (symbols) are included as the row names and column names of the emission matrix, respectively. This means there are only two objects in the list, the transition matrix (A) and the emission matrix (E). So you can initialize and plot your model like this:

states <- c("Begin", "feeding", "inactive", "exploring", "others", "moving", "grooming")
residues <- paste0(states[2:7],"E")
transProbs <- rbind(startProbs, transProbs)
transProbs <- cbind(0, transProbs)
dimnames(transProbs) <- list(from = states, to = states)
residues <- c("F", "I", "E", "O", "M", "G") # abbreviated for plotting
dimnames(emissionProbs) <- list(states = states[-1], residues = residues)
model <- list(A = transProbs, E = emissionProbs)
class(model) <- "HMM"
plot(model, begin = TRUE) 

[maforinw]

Thank you for taking the time to tackle this question.