These distribution classes are designed specifically for HMM’s and not for general use in statistics. For example, they have fixed or non-fixed status, which only make sense relative to being used in a hidden Markov model.
AUTHOR:
- William Stein, 2010-03
Bases: object
A distribution.
Return a plot of the probability density function.
INPUT:
- args and kwds, passed to the Sage plot function
OUTPUT:
- a Graphics object
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.plot(-10,30)
The probability density function evaluated at x.
INPUT:
- x – object
OUTPUT:
- float
EXAMPLES:
This method must be defined in a derived class:
sage: sage.stats.hmm.distributions.Distribution().prob(0)
Traceback (most recent call last):
...
NotImplementedError
Return either a single sample (the default) or n samples from this probability distribution.
INPUT:
- n – None or a positive integer
OUTPUT:
- a single sample if n is 1; otherwise many samples
EXAMPLES:
This method must be defined in a derived class:
sage: sage.stats.hmm.distributions.Distribution().sample()
Traceback (most recent call last):
...
NotImplementedError
Bases: sage.stats.hmm.distributions.Distribution
A probability distribution defined by taking a weighted linear combination of Gaussian distributions.
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.3,1,2),(.7,-1,1)]); P
0.3*N(1.0,2.0) + 0.7*N(-1.0,1.0)
sage: P[0]
(0.29999999999999999, 1.0, 2.0)
sage: P.is_fixed()
False
sage: P.fix(1)
sage: P.is_fixed(0)
False
sage: P.is_fixed(1)
True
sage: P.unfix(1)
sage: P.is_fixed(1)
False
Set that this GaussianMixtureDistribution (or its ith component) is fixed when using Baum-Welch to update the corresponding HMM.
INPUT:
- i - None (default) or integer; if given, only fix the i-th component
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.fix(1); P.is_fixed()
False
sage: P.is_fixed(1)
True
sage: P.fix(); P.is_fixed()
True
Return whether or not this GaussianMixtureDistribution is fixed when using Baum-Welch to update the corresponding HMM.
INPUT:
- i - None (default) or integer; if given, only return whether the i-th component is fixed
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.is_fixed()
False
sage: P.is_fixed(0)
False
sage: P.fix(0); P.is_fixed()
False
sage: P.is_fixed(0)
True
sage: P.fix(); P.is_fixed()
True
Return the probability of x.
Since this is a continuous distribution, this is defined to be the limit of the p’s such that the probability of [x,x+h] is p*h.
INPUT:
- x – float
OUTPUT:
- float
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.prob(.5)
0.21123919605857971
sage: P.prob(-100)
0.0
sage: P.prob(20)
0.1595769121605731
Return the probability of x using just the m-th summand.
INPUT:
- x – float
- m – integer
OUTPUT:
- float
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.prob_m(.5, 0)
2.760811768050888...e-97
sage: P.prob_m(.5, 1)
0.21123919605857971
sage: P.prob_m(.5, 2)
0.0
Return a single sample from this distribution (by default), or if n>1, return a TimeSeries of samples.
INPUT:
- n – integer or None (default: None)
OUTPUT:
- float if n is None (default); otherwise a TimeSeries
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.sample()
19.658243610875129
sage: P.sample(1)
[-10.4683]
sage: P.sample(5)
[-0.1688, -10.3479, 1.6812, 20.1083, -9.9801]
sage: P.sample(0)
[]
sage: P.sample(-3)
Traceback (most recent call last):
...
ValueError: n must be nonnegative
Set that this GaussianMixtureDistribution (or its ith component) is not fixed when using Baum-Welch to update the corresponding HMM.
INPUT:
- i - None (default) or integer; if given, only fix the i-th component
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: P.fix(1); P.is_fixed(1)
True
sage: P.unfix(1); P.is_fixed(1)
False
sage: P.fix(); P.is_fixed()
True
sage: P.unfix(); P.is_fixed()
False
Used in unpickling GaussianMixtureDistribution’s.
EXAMPLES:
sage: P = hmm.GaussianMixtureDistribution([(.2,-10,.5),(.6,1,1),(.2,20,.5)])
sage: loads(dumps(P)) == P # indirect doctest
True
Continuous Emission Hidden Markov Models
Hidden Markov Models – Utility functions
Enter search terms or a module, class or function name.