Continuous Emission Hidden Markov Models¶
AUTHOR:
- William Stein, 2010-03 
- class sage.stats.hmm.chmm.GaussianHiddenMarkovModel[source]¶
- Bases: - HiddenMarkovModel- Gaussian emissions Hidden Markov Model. - INPUT: - A– matrix; the \(N \times N\) transition matrix
- B– list of pairs- (mu, sigma)that define the distributions
- pi– initial state probabilities
- normalize– boolean (default:- True)
 - EXAMPLES: - We illustrate the primary functions with an example 2-state Gaussian HMM: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,1), (-1,1)], ....: [.5,.5]); m Gaussian Hidden Markov Model with 2 States Transition matrix: [0.1 0.9] [0.5 0.5] Emission parameters: [(1.0, 1.0), (-1.0, 1.0)] Initial probabilities: [0.5000, 0.5000] - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),Integer(1)), (-Integer(1),Integer(1))], ... [RealNumber('.5'),RealNumber('.5')]); m Gaussian Hidden Markov Model with 2 States Transition matrix: [0.1 0.9] [0.5 0.5] Emission parameters: [(1.0, 1.0), (-1.0, 1.0)] Initial probabilities: [0.5000, 0.5000] - We query the defining transition matrix, emission parameters, and initial state probabilities: - sage: m.transition_matrix() [0.1 0.9] [0.5 0.5] sage: m.emission_parameters() [(1.0, 1.0), (-1.0, 1.0)] sage: m.initial_probabilities() [0.5000, 0.5000] - >>> from sage.all import * >>> m.transition_matrix() [0.1 0.9] [0.5 0.5] >>> m.emission_parameters() [(1.0, 1.0), (-1.0, 1.0)] >>> m.initial_probabilities() [0.5000, 0.5000] - We obtain a sample sequence with 10 entries in it, and compute the logarithm of the probability of obtaining this sequence, given the model: - sage: obs = m.sample(5); obs # random [-1.6835, 0.0635, -2.1688, 0.3043, -0.3188] sage: log_likelihood = m.log_likelihood(obs) sage: counter = 0 sage: n = 0 sage: def add_samples(i): ....: global counter, n ....: for _ in range(i): ....: n += 1 ....: obs2 = m.sample(5) ....: if all(abs(obs2[i] - obs[i]) < 0.25 for i in range(5)): ....: counter += 1 sage: add_samples(10000) sage: while abs(log_likelihood - log(counter*1.0/n/0.5^5)) < 0.1: ....: add_samples(10000) - >>> from sage.all import * >>> obs = m.sample(Integer(5)); obs # random [-1.6835, 0.0635, -2.1688, 0.3043, -0.3188] >>> log_likelihood = m.log_likelihood(obs) >>> counter = Integer(0) >>> n = Integer(0) >>> def add_samples(i): ... global counter, n ... for _ in range(i): ... n += Integer(1) ... obs2 = m.sample(Integer(5)) ... if all(abs(obs2[i] - obs[i]) < RealNumber('0.25') for i in range(Integer(5))): ... counter += Integer(1) >>> add_samples(Integer(10000)) >>> while abs(log_likelihood - log(counter*RealNumber('1.0')/n/RealNumber('0.5')**Integer(5))) < RealNumber('0.1'): ... add_samples(Integer(10000)) - We compute the Viterbi path, and probability that the given path of states produced obs: - sage: m.viterbi(obs) # random ([1, 0, 1, 0, 1], -8.714092684611794) - >>> from sage.all import * >>> m.viterbi(obs) # random ([1, 0, 1, 0, 1], -8.714092684611794) - We use the Baum-Welch iterative algorithm to find another model for which our observation sequence is more likely: - sage: try: ....: p, s = m.baum_welch(obs) ....: assert p > log_likelihood ....: assert (1 <= s <= 500) ....: except RuntimeError: ....: pass - >>> from sage.all import * >>> try: ... p, s = m.baum_welch(obs) ... assert p > log_likelihood ... assert (Integer(1) <= s <= Integer(500)) ... except RuntimeError: ... pass - Notice that running Baum-Welch changed our model: - sage: m # random Gaussian Hidden Markov Model with 2 States Transition matrix: [ 0.4154981366185841 0.584501863381416] [ 0.9999993174253741 6.825746258991804e-07] Emission parameters: [(0.4178882427119503, 0.5173109664360919), (-1.5025208631331122, 0.5085512836055119)] Initial probabilities: [0.0000, 1.0000] - >>> from sage.all import * >>> m # random Gaussian Hidden Markov Model with 2 States Transition matrix: [ 0.4154981366185841 0.584501863381416] [ 0.9999993174253741 6.825746258991804e-07] Emission parameters: [(0.4178882427119503, 0.5173109664360919), (-1.5025208631331122, 0.5085512836055119)] Initial probabilities: [0.0000, 1.0000] - baum_welch(obs, max_iter=500, log_likelihood_cutoff=0.0001, min_sd=0.01, fix_emissions=False, v=False)[source]¶
- Given an observation sequence - obs, improve this HMM using the Baum-Welch algorithm to increase the probability of observing- obs.- INPUT: - obs– a time series of emissions
- max_iter– integer (default: 500); maximum number of Baum-Welch steps to take
- log_likelihood_cutoff– positive float (default: 1e-4); the minimal improvement in likelihood with respect to the last iteration required to continue. Relative value to log likelihood.
- min_sd– positive float (default: 0.01); when reestimating, the standard deviation of emissions is not allowed to be less than- min_sd.
- fix_emissions– boolean (default:- False); if- True, do not change emissions when updating
 - OUTPUT: - changes the model in place, and returns the log likelihood and number of iterations. - EXAMPLES: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: m.log_likelihood([-2,-1,.1,0.1]) -8.858282215986275 sage: m.baum_welch([-2,-1,.1,0.1]) (4.534646052182..., 7) sage: m.log_likelihood([-2,-1,.1,0.1]) 4.534646052182... sage: m # rel tol 3e-14 Gaussian Hidden Markov Model with 2 States Transition matrix: [ 0.9999999992430161 7.569839394440382e-10] [ 0.49998462791192644 0.5000153720880736] Emission parameters: [(0.09999999999999999, 0.01), (-1.4999508147591902, 0.5000710504895474)] Initial probabilities: [0.0000, 1.0000] - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> m.log_likelihood([-Integer(2),-Integer(1),RealNumber('.1'),RealNumber('0.1')]) -8.858282215986275 >>> m.baum_welch([-Integer(2),-Integer(1),RealNumber('.1'),RealNumber('0.1')]) (4.534646052182..., 7) >>> m.log_likelihood([-Integer(2),-Integer(1),RealNumber('.1'),RealNumber('0.1')]) 4.534646052182... >>> m # rel tol 3e-14 Gaussian Hidden Markov Model with 2 States Transition matrix: [ 0.9999999992430161 7.569839394440382e-10] [ 0.49998462791192644 0.5000153720880736] Emission parameters: [(0.09999999999999999, 0.01), (-1.4999508147591902, 0.5000710504895474)] Initial probabilities: [0.0000, 1.0000] - We illustrate bounding the standard deviation below. Note that above we had different emission parameters when the - min_sdwas the default of 0.01:- sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: m.baum_welch([-2,-1,.1,0.1], min_sd=1) (-4.07939572755..., 32) sage: m.emission_parameters() [(-0.2663018798..., 1.0), (-1.99850979..., 1.0)] - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> m.baum_welch([-Integer(2),-Integer(1),RealNumber('.1'),RealNumber('0.1')], min_sd=Integer(1)) (-4.07939572755..., 32) >>> m.emission_parameters() [(-0.2663018798..., 1.0), (-1.99850979..., 1.0)] - We watch the log likelihoods of the model converge, step by step: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: v = m.sample(10) sage: l = stats.TimeSeries([m.baum_welch(v, max_iter=1)[0] ....: for _ in range(len(v))]) sage: all(l[i] <= l[i+1] + 0.0001 for i in range(9)) True sage: l # random [-20.1167, -17.7611, -16.9814, -16.9364, -16.9314, -16.9309, -16.9309, -16.9309, -16.9309, -16.9309] - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> v = m.sample(Integer(10)) >>> l = stats.TimeSeries([m.baum_welch(v, max_iter=Integer(1))[Integer(0)] ... for _ in range(len(v))]) >>> all(l[i] <= l[i+Integer(1)] + RealNumber('0.0001') for i in range(Integer(9))) True >>> l # random [-20.1167, -17.7611, -16.9814, -16.9364, -16.9314, -16.9309, -16.9309, -16.9309, -16.9309, -16.9309] - We illustrate fixing emissions: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.9,.1]], ....: [(1,2),(-1,.5)], ....: [.3,.7]) sage: set_random_seed(0); v = m.sample(100) sage: m.baum_welch(v,fix_emissions=True) (-164.72944548204..., 23) sage: m.emission_parameters() [(1.0, 2.0), (-1.0, 0.5)] sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.9,.1]], ....: [(1,2),(-1,.5)], ....: [.3,.7]) sage: m.baum_welch(v) (-162.854370397998..., 49) sage: m.emission_parameters() # rel tol 3e-14 [(1.2722419172602375, 2.371368751761901), (-0.9486174675179113, 0.5762360385123765)] - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.9'),RealNumber('.1')]], ... [(Integer(1),Integer(2)),(-Integer(1),RealNumber('.5'))], ... [RealNumber('.3'),RealNumber('.7')]) >>> set_random_seed(Integer(0)); v = m.sample(Integer(100)) >>> m.baum_welch(v,fix_emissions=True) (-164.72944548204..., 23) >>> m.emission_parameters() [(1.0, 2.0), (-1.0, 0.5)] >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.9'),RealNumber('.1')]], ... [(Integer(1),Integer(2)),(-Integer(1),RealNumber('.5'))], ... [RealNumber('.3'),RealNumber('.7')]) >>> m.baum_welch(v) (-162.854370397998..., 49) >>> m.emission_parameters() # rel tol 3e-14 [(1.2722419172602375, 2.371368751761901), (-0.9486174675179113, 0.5762360385123765)] 
 - emission_parameters()[source]¶
- Return the parameters that define the normal distributions associated to all of the states. - OUTPUT: - a list - Bof pairs- B[i] = (mu, std), such that the distribution associated to state \(i\) is normal with mean- muand standard deviation- std.- EXAMPLES: - sage: M = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: M.emission_parameters() [(1.0, 0.5), (-1.0, 3.0)] - >>> from sage.all import * >>> M = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> M.emission_parameters() [(1.0, 0.5), (-1.0, 3.0)] 
 - generate_sequence(length, starting_state=None)[source]¶
- Return a sample of the given length from this HMM. - INPUT: - length– positive integer
- starting_state– integer (or- None); if specified then generate a sequence using this model starting with the given state instead of the initial probabilities to determine the starting state.
 - OUTPUT: - an - IntListor list of emission symbols
- TimeSeriesof emissions
 - EXAMPLES: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: m.generate_sequence(5) # random ([-3.0505, 0.5317, -4.5065, 0.6521, 1.0435], [1, 0, 1, 0, 1]) sage: m.generate_sequence(0) ([], []) sage: m.generate_sequence(-1) Traceback (most recent call last): ... ValueError: length must be nonnegative - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> m.generate_sequence(Integer(5)) # random ([-3.0505, 0.5317, -4.5065, 0.6521, 1.0435], [1, 0, 1, 0, 1]) >>> m.generate_sequence(Integer(0)) ([], []) >>> m.generate_sequence(-Integer(1)) Traceback (most recent call last): ... ValueError: length must be nonnegative - Verify numerically that the starting state is 0 with probability about 0.1: - sage: counter = 0 sage: n = 0 sage: def add_samples(i): ....: global counter, n ....: for i in range(i): ....: n += 1 ....: if m.generate_sequence(1)[1][0] == 0: ....: counter += 1 sage: add_samples(10^5) sage: while abs(counter*1.0 / n - 0.1) > 0.01: add_samples(10^5) - >>> from sage.all import * >>> counter = Integer(0) >>> n = Integer(0) >>> def add_samples(i): ... global counter, n ... for i in range(i): ... n += Integer(1) ... if m.generate_sequence(Integer(1))[Integer(1)][Integer(0)] == Integer(0): ... counter += Integer(1) >>> add_samples(Integer(10)**Integer(5)) >>> while abs(counter*RealNumber('1.0') / n - RealNumber('0.1')) > RealNumber('0.01'): add_samples(Integer(10)**Integer(5)) - Example in which the starting state is 0 (see Issue #11452): - sage: set_random_seed(23); m.generate_sequence(2) ([0.6501, -2.0151], [0, 1]) - >>> from sage.all import * >>> set_random_seed(Integer(23)); m.generate_sequence(Integer(2)) ([0.6501, -2.0151], [0, 1]) - Force a starting state of 1 even though as we saw above it would be 0: - sage: set_random_seed(23); m.generate_sequence(2, starting_state=1) ([-3.1491, -1.0244], [1, 1]) - >>> from sage.all import * >>> set_random_seed(Integer(23)); m.generate_sequence(Integer(2), starting_state=Integer(1)) ([-3.1491, -1.0244], [1, 1]) 
 - log_likelihood(obs)[source]¶
- Return the logarithm of a continuous analogue of the probability that this model produced the given observation sequence. - Note that the “continuous analogue of the probability” above can be bigger than 1, hence the logarithm can be positive. - INPUT: - obs– sequence of observations
 - OUTPUT: float - EXAMPLES: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: m.log_likelihood([1,1,1]) -4.297880766072486 sage: s = m.sample(20) sage: -80 < m.log_likelihood(s) < -20 True - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> m.log_likelihood([Integer(1),Integer(1),Integer(1)]) -4.297880766072486 >>> s = m.sample(Integer(20)) >>> -Integer(80) < m.log_likelihood(s) < -Integer(20) True 
 - viterbi(obs)[source]¶
- Determine “the” hidden sequence of states that is most likely to produce the given sequence - obsof observations, along with the probability that this hidden sequence actually produced the observation.- INPUT: - obs– sequence of emitted ints or symbols
 - OUTPUT: - list– “the” most probable sequence of hidden states, i.e., the Viterbi path
- float– log of probability that the observed sequence was produced by the Viterbi sequence of states
 - EXAMPLES: - We find the optimal state sequence for a given model: - sage: m = hmm.GaussianHiddenMarkovModel([[0.5,0.5],[0.5,0.5]], ....: [(0,1),(10,1)], ....: [0.5,0.5]) sage: m.viterbi([0,1,10,10,1]) ([0, 0, 1, 1, 0], -9.0604285688230...) - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('0.5'),RealNumber('0.5')],[RealNumber('0.5'),RealNumber('0.5')]], ... [(Integer(0),Integer(1)),(Integer(10),Integer(1))], ... [RealNumber('0.5'),RealNumber('0.5')]) >>> m.viterbi([Integer(0),Integer(1),Integer(10),Integer(10),Integer(1)]) ([0, 0, 1, 1, 0], -9.0604285688230...) - Another example in which the most likely states change based on the last observation: - sage: m = hmm.GaussianHiddenMarkovModel([[.1,.9],[.5,.5]], ....: [(1,.5), (-1,3)], ....: [.1,.9]) sage: m.viterbi([-2,-1,.1,0.1]) ([1, 1, 0, 1], -9.61823698847639...) sage: m.viterbi([-2,-1,.1,0.3]) ([1, 1, 1, 0], -9.566023653378513) - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[RealNumber('.1'),RealNumber('.9')],[RealNumber('.5'),RealNumber('.5')]], ... [(Integer(1),RealNumber('.5')), (-Integer(1),Integer(3))], ... [RealNumber('.1'),RealNumber('.9')]) >>> m.viterbi([-Integer(2),-Integer(1),RealNumber('.1'),RealNumber('0.1')]) ([1, 1, 0, 1], -9.61823698847639...) >>> m.viterbi([-Integer(2),-Integer(1),RealNumber('.1'),RealNumber('0.3')]) ([1, 1, 1, 0], -9.566023653378513) 
 
- class sage.stats.hmm.chmm.GaussianMixtureHiddenMarkovModel[source]¶
- Bases: - GaussianHiddenMarkovModel- Gaussian mixture Hidden Markov Model. - INPUT: - A– matrix; the \(N \times N\) transition matrix
- B– list of mixture definitions for each state. Each state may have a varying number of gaussians with selection probabilities that sum to 1 and encoded as- (p, (mu,sigma))
- pi– initial state probabilities
- normalize– boolean (default:- True); if given, input is normalized to define valid probability distributions, e.g., the entries of \(A\) are made nonnegative and the rows sum to 1, and the probabilities in- piare normalized.
 - EXAMPLES: - sage: A = [[0.5,0.5],[0.5,0.5]] sage: B = [[(0.9,(0.0,1.0)), (0.1,(1,10000))],[(1,(1,1)), (0,(0,0.1))]] sage: hmm.GaussianMixtureHiddenMarkovModel(A, B, [1,0]) Gaussian Mixture Hidden Markov Model with 2 States Transition matrix: [0.5 0.5] [0.5 0.5] Emission parameters: [0.9*N(0.0,1.0) + 0.1*N(1.0,10000.0), 1.0*N(1.0,1.0) + 0.0*N(0.0,0.1)] Initial probabilities: [1.0000, 0.0000] - >>> from sage.all import * >>> A = [[RealNumber('0.5'),RealNumber('0.5')],[RealNumber('0.5'),RealNumber('0.5')]] >>> B = [[(RealNumber('0.9'),(RealNumber('0.0'),RealNumber('1.0'))), (RealNumber('0.1'),(Integer(1),Integer(10000)))],[(Integer(1),(Integer(1),Integer(1))), (Integer(0),(Integer(0),RealNumber('0.1')))]] >>> hmm.GaussianMixtureHiddenMarkovModel(A, B, [Integer(1),Integer(0)]) Gaussian Mixture Hidden Markov Model with 2 States Transition matrix: [0.5 0.5] [0.5 0.5] Emission parameters: [0.9*N(0.0,1.0) + 0.1*N(1.0,10000.0), 1.0*N(1.0,1.0) + 0.0*N(0.0,0.1)] Initial probabilities: [1.0000, 0.0000] - baum_welch(obs, max_iter=1000, log_likelihood_cutoff=1e-12, min_sd=0.01, fix_emissions=False)[source]¶
- Given an observation sequence - obs, improve this HMM using the Baum-Welch algorithm to increase the probability of observing- obs.- INPUT: - obs– a time series of emissions
- max_iter– integer (default: 1000); maximum number of Baum-Welch steps to take
- log_likelihood_cutoff– positive float (default: 1e-12); the minimal improvement in likelihood with respect to the last iteration required to continue. Relative value to log likelihood.
- min_sd– positive float (default: 0.01); when reestimating, the standard deviation of emissions is not allowed to be less than- min_sd.
- fix_emissions– boolean (default:- False); if- True, do not change emissions when updating
 - OUTPUT: - changes the model in place, and returns the log likelihood and number of iterations. - EXAMPLES: - sage: m = hmm.GaussianMixtureHiddenMarkovModel( ....: [[.9,.1],[.4,.6]], ....: [[(.4,(0,1)), (.6,(1,0.1))], [(1,(0,1))]], ....: [.7,.3]) sage: set_random_seed(0); v = m.sample(10); v [0.3576, -0.9365, 0.9449, -0.6957, 1.0217, 0.9644, 0.9987, -0.5950, -1.0219, 0.6477] sage: m.log_likelihood(v) -8.31408655939536... sage: m.baum_welch(v) (2.18905068682..., 15) sage: m.log_likelihood(v) 2.18905068682... sage: m # rel tol 6e-12 Gaussian Mixture Hidden Markov Model with 2 States Transition matrix: [ 0.8746363339773399 0.12536366602266016] [ 1.0 1.451685202290174e-40] Emission parameters: [0.500161629343*N(-0.812298726239,0.173329026744) + 0.499838370657*N(0.982433690378,0.029719932009), 1.0*N(0.503260056832,0.145881515324)] Initial probabilities: [0.0000, 1.0000] - >>> from sage.all import * >>> m = hmm.GaussianMixtureHiddenMarkovModel( ... [[RealNumber('.9'),RealNumber('.1')],[RealNumber('.4'),RealNumber('.6')]], ... [[(RealNumber('.4'),(Integer(0),Integer(1))), (RealNumber('.6'),(Integer(1),RealNumber('0.1')))], [(Integer(1),(Integer(0),Integer(1)))]], ... [RealNumber('.7'),RealNumber('.3')]) >>> set_random_seed(Integer(0)); v = m.sample(Integer(10)); v [0.3576, -0.9365, 0.9449, -0.6957, 1.0217, 0.9644, 0.9987, -0.5950, -1.0219, 0.6477] >>> m.log_likelihood(v) -8.31408655939536... >>> m.baum_welch(v) (2.18905068682..., 15) >>> m.log_likelihood(v) 2.18905068682... >>> m # rel tol 6e-12 Gaussian Mixture Hidden Markov Model with 2 States Transition matrix: [ 0.8746363339773399 0.12536366602266016] [ 1.0 1.451685202290174e-40] Emission parameters: [0.500161629343*N(-0.812298726239,0.173329026744) + 0.499838370657*N(0.982433690378,0.029719932009), 1.0*N(0.503260056832,0.145881515324)] Initial probabilities: [0.0000, 1.0000] - We illustrate bounding the standard deviation below. Note that above we had different emission parameters when the min_sd was the default of 0.01: - sage: m = hmm.GaussianMixtureHiddenMarkovModel( ....: [[.9,.1],[.4,.6]], ....: [[(.4,(0,1)), (.6,(1,0.1))], [(1,(0,1))]], ....: [.7,.3]) sage: m.baum_welch(v, min_sd=1) (-12.617885761692..., 1000) sage: m.emission_parameters() # rel tol 6e-12 [0.503545634447*N(0.200166509595,1.0) + 0.496454365553*N(0.200166509595,1.0), 1.0*N(0.0543433426535,1.0)] - >>> from sage.all import * >>> m = hmm.GaussianMixtureHiddenMarkovModel( ... [[RealNumber('.9'),RealNumber('.1')],[RealNumber('.4'),RealNumber('.6')]], ... [[(RealNumber('.4'),(Integer(0),Integer(1))), (RealNumber('.6'),(Integer(1),RealNumber('0.1')))], [(Integer(1),(Integer(0),Integer(1)))]], ... [RealNumber('.7'),RealNumber('.3')]) >>> m.baum_welch(v, min_sd=Integer(1)) (-12.617885761692..., 1000) >>> m.emission_parameters() # rel tol 6e-12 [0.503545634447*N(0.200166509595,1.0) + 0.496454365553*N(0.200166509595,1.0), 1.0*N(0.0543433426535,1.0)] - We illustrate fixing all emissions: - sage: m = hmm.GaussianMixtureHiddenMarkovModel( ....: [[.9,.1],[.4,.6]], ....: [[(.4,(0,1)), (.6,(1,0.1))], [(1,(0,1))]], ....: [.7,.3]) sage: set_random_seed(0); v = m.sample(10) sage: m.baum_welch(v, fix_emissions=True) (-7.58656858997..., 36) sage: m.emission_parameters() [0.4*N(0.0,1.0) + 0.6*N(1.0,0.1), 1.0*N(0.0,1.0)] - >>> from sage.all import * >>> m = hmm.GaussianMixtureHiddenMarkovModel( ... [[RealNumber('.9'),RealNumber('.1')],[RealNumber('.4'),RealNumber('.6')]], ... [[(RealNumber('.4'),(Integer(0),Integer(1))), (RealNumber('.6'),(Integer(1),RealNumber('0.1')))], [(Integer(1),(Integer(0),Integer(1)))]], ... [RealNumber('.7'),RealNumber('.3')]) >>> set_random_seed(Integer(0)); v = m.sample(Integer(10)) >>> m.baum_welch(v, fix_emissions=True) (-7.58656858997..., 36) >>> m.emission_parameters() [0.4*N(0.0,1.0) + 0.6*N(1.0,0.1), 1.0*N(0.0,1.0)] 
 - emission_parameters()[source]¶
- Return a list of all the emission distributions. - OUTPUT: list of Gaussian mixtures - EXAMPLES: - sage: m = hmm.GaussianMixtureHiddenMarkovModel([[.9,.1],[.4,.6]], ....: [[(.4,(0,1)), (.6,(1,0.1))], [(1,(0,1))]], ....: [.7,.3]) sage: m.emission_parameters() [0.4*N(0.0,1.0) + 0.6*N(1.0,0.1), 1.0*N(0.0,1.0)] - >>> from sage.all import * >>> m = hmm.GaussianMixtureHiddenMarkovModel([[RealNumber('.9'),RealNumber('.1')],[RealNumber('.4'),RealNumber('.6')]], ... [[(RealNumber('.4'),(Integer(0),Integer(1))), (RealNumber('.6'),(Integer(1),RealNumber('0.1')))], [(Integer(1),(Integer(0),Integer(1)))]], ... [RealNumber('.7'),RealNumber('.3')]) >>> m.emission_parameters() [0.4*N(0.0,1.0) + 0.6*N(1.0,0.1), 1.0*N(0.0,1.0)] 
 
- sage.stats.hmm.chmm.unpickle_gaussian_hmm_v0(A, B, pi, name)[source]¶
- EXAMPLES: - sage: m = hmm.GaussianHiddenMarkovModel([[1]], [(0,1)], [1]) sage: sage.stats.hmm.chmm.unpickle_gaussian_hmm_v0(m.transition_matrix(), m.emission_parameters(), m.initial_probabilities(), 'test') Gaussian Hidden Markov Model with 1 States Transition matrix: [1.0] Emission parameters: [(0.0, 1.0)] Initial probabilities: [1.0000] - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[Integer(1)]], [(Integer(0),Integer(1))], [Integer(1)]) >>> sage.stats.hmm.chmm.unpickle_gaussian_hmm_v0(m.transition_matrix(), m.emission_parameters(), m.initial_probabilities(), 'test') Gaussian Hidden Markov Model with 1 States Transition matrix: [1.0] Emission parameters: [(0.0, 1.0)] Initial probabilities: [1.0000] 
- sage.stats.hmm.chmm.unpickle_gaussian_hmm_v1(A, B, pi, prob, n_out)[source]¶
- EXAMPLES: - sage: m = hmm.GaussianHiddenMarkovModel([[1]], [(0,1)], [1]) sage: loads(dumps(m)) == m # indirect test True - >>> from sage.all import * >>> m = hmm.GaussianHiddenMarkovModel([[Integer(1)]], [(Integer(0),Integer(1))], [Integer(1)]) >>> loads(dumps(m)) == m # indirect test True 
- sage.stats.hmm.chmm.unpickle_gaussian_mixture_hmm_v1(A, B, pi, mixture)[source]¶
- EXAMPLES: - sage: m = hmm.GaussianMixtureHiddenMarkovModel([[1]], [[(.4,(0,1)), (.6,(1,0.1))]], [1]) sage: loads(dumps(m)) == m # indirect test True - >>> from sage.all import * >>> m = hmm.GaussianMixtureHiddenMarkovModel([[Integer(1)]], [[(RealNumber('.4'),(Integer(0),Integer(1))), (RealNumber('.6'),(Integer(1),RealNumber('0.1')))]], [Integer(1)]) >>> loads(dumps(m)) == m # indirect test True