Ticket #3726: sage-3726-part5.patch

File sage-3726-part5.patch, 8.8 KB (added by was, 11 years ago)
  • sage/stats/hmm/all.py

    # HG changeset patch
    # User William Stein <wstein@gmail.com>
    # Date 1217660625 25200
    # Node ID 0ec7e93289ae6ddd0ec6de054b4748ca25d32355
    # Parent  ba54a820d1368e08358889763dee4206933c1cbe
    trac #3726 -- more work on continuous HMM models
    
    diff -r ba54a820d136 -r 0ec7e93289ae sage/stats/hmm/all.py
    a b from hmm import DiscreteHiddenMarkovMod 
    11from hmm  import DiscreteHiddenMarkovModel
    2 from chmm import NormalHiddenMarkovModel
     2from chmm import GaussianHiddenMarkovModel
  • sage/stats/hmm/chmm.pxd

    diff -r ba54a820d136 -r 0ec7e93289ae sage/stats/hmm/chmm.pxd
    a b cdef extern from "ghmm/smodel.h": 
    7474        # IDs of predecessor states
    7575        int *in_id
    7676        # transition probs to successor states. It is a
    77         # matrix in case of mult. transition matrices (COS > 1)
     77        # matrix in case of mult. transition matrices (cos > 1)
    7878        double **out_a
    7979        # transition probs from predecessor states. It is a
    80         # matrix in case of mult. transition matrices (COS > 1)
     80        # matrix in case of mult. transition matrices (cos > 1)
    8181        double **in_a
    8282        # number of  successor states
    8383        int out_states
    cdef extern from "ghmm/smodel.h": 
    9797        int xPosition
    9898        # y coordinate position for graph representation plotting
    9999        int yPosition
     100
     101
     102    double **ighmm_cmatrix_alloc (int nrows, int ncols)
    100103
    101104    #############################################################
    102105    # Continous Hidden Markov Model
    cdef class ContinuousHiddenMarkovModel(H 
    196199    cdef ghmm_cmodel* m
    197200    cdef bint initialized
    198201
    199 cdef class NormalHiddenMarkovModel(ContinuousHiddenMarkovModel):
     202cdef class GaussianHiddenMarkovModel(ContinuousHiddenMarkovModel):
    200203    pass
    201204
    202205
  • sage/stats/hmm/chmm.pyx

    diff -r ba54a820d136 -r 0ec7e93289ae sage/stats/hmm/chmm.pyx
    a b cdef class ContinuousHiddenMarkovModel(H 
    1212        # Set number of states
    1313        self.m.N = self.A.nrows()
    1414
    15 cdef class NormalHiddenMarkovModel(ContinuousHiddenMarkovModel):
     15cdef class GaussianHiddenMarkovModel(ContinuousHiddenMarkovModel):
    1616    """
     17    Create a Gaussian Hidden Markov Model.  The probability
     18    distribution associated with each state is a Gaussian
     19    distribution.
     20
     21    GaussianHiddenMarkovModel(A, B, pi, name)
     22
     23    INPUT:
     24        A  -- matrix; the transition matrix (n x n)
     25        B  -- list of n pairs (mu, sigma) that define the
     26              Gaussian distributions associated to each state
     27        pi -- list of floats that sums to 1.0; these are
     28              the initial probabilities of each hidden state
     29        name -- (default: None);
     30
    1731    EXAMPLES:
    18     The transition matrix:
     32    Define the transition matrix:
    1933        sage: A = [[0,1,0],[0.5,0,0.5],[0.3,0.3,0.4]]
    2034
    2135    Parameters of the normal emission distributions in pairs of (mu, sigma):
    cdef class NormalHiddenMarkovModel(Conti 
    2438    The initial probabilities per state:
    2539        sage: pi = [1,0,0]
    2640
    27     Create the continuous HMM:
    28         sage: m = hmm.NormalHiddenMarkovModel(A, B, pi); m
     41    Create the continuous Gaussian hidden Markov model:
     42        sage: m = hmm.GaussianHiddenMarkovModel(A, B, pi); m
    2943    """
    3044    def __init__(self, A, B, pi, name=None):
    3145        ContinuousHiddenMarkovModel.__init__(self, A, B, pi)
    3246
    33         # Set number of outputs
     47        # Set number of outputs.  This is 1 here because each
     48        # output is a single Gaussian distribution.
    3449        self.m.M = 1
    3550
    3651        # Set the model type to continuous
    cdef class NormalHiddenMarkovModel(Conti 
    5772
    5873        for i in range(self.m.N):
    5974            # Parameters of normal distribution
    60             mu, sigma = self.B[i]
     75            mu, sigma   = self.B[i]
    6176            # Get a reference to the i-th state for convenience of the notation below.
    6277            state = &(states[i])
    63             state.desc = NULL
    64             state.M = 1
     78            state.M     = 1
     79            state.pi    = pi[i]
     80            state.desc  = NULL
     81            state.out_states = 0
     82            state.in_states = 0
    6583            e = <ghmm_c_emission*> safe_malloc(sizeof(ghmm_c_emission))
    66             e.type = 0  # normal
     84            e.type      = 0  # normal
    6785            e.dimension = 1
    68             e.mean.val = mu
     86            e.mean.val  = mu
    6987            e.variance.val = sigma
    7088            # fixing of emissions is deactivated by default           
    71             e.fixed = 0
    72             e.sigmacd = NULL
    73             e.sigmainv = NULL
    74             state.e = e
    75             state.c = to_double_array([1.0])
     89            e.fixed     = 0
     90            e.sigmacd   = NULL
     91            e.sigmainv  = NULL
     92            state.e     = e
     93            state.c     = to_double_array([1.0])
     94            state.in_a  = ighmm_cmatrix_alloc(1, self.m.N)
     95            state.out_a = ighmm_cmatrix_alloc(1, self.m.N)
    7696
    7797        # Set states
    7898        self.m.s = states
    7999
     100        self.m.class_change = NULL
     101
    80102        self.initialized = True
    81103
    82104    def __dealloc__(self):
    83         return
    84105        if self.initialized:
    85106            ghmm_cmodel_free(&self.m)
    86107
    cdef class NormalHiddenMarkovModel(Conti 
    92113            string
    93114
    94115        EXAMPLES:
    95             sage: m = hmm.NormalHiddenMarkovModel([[0.0,1.0,0],[0.5,0.0,0.5],[0.3,0.3,0.4]], [(0.0,1.0), (-1.0,0.5), (1.0,0.2)], [1,0,0])
     116            sage: m = hmm.GaussianHiddenMarkovModel([[0.0,1.0,0],[0.5,0.0,0.5],[0.3,0.3,0.4]], [(0.0,1.0), (-1.0,0.5), (1.0,0.2)], [1,0,0])
    96117            sage: a.__repr__()
    97118            "Discrete Hidden Markov Model (2 states, 2 outputs)\nInitial probabilities: [0.5, 0.5]\nTransition matrix:\n[0.1 0.9]\n[0.1 0.9]\nEmission matrix:\n[0.9 0.1]\n[0.1 0.9]\nEmission symbols: [3/4, 'abc']"
    98119        """
    99         s = "Normal Hidden Markov Model%s (%s states, %s outputs)"%(
     120        s = "Gaussian Hidden Markov Model%s (%s states, %s outputs)"%(
    100121            ' ' + self.m.name if self.m.name else '',
    101122            self.m.N, self.m.M)
    102         #s += '\nInitial probabilities: %s'%self.initial_probabilities()
    103         #s += '\nTransition matrix:\n%s'%self.transition_matrix()
    104         #s += '\nEmission matrix:\n%s'%self.emission_matrix()
     123        s += '\nInitial probabilities: %s'%self.initial_probabilities()
     124        s += '\nTransition matrix:\n%s'%self.transition_matrix()
     125        s += '\nEmission parameters:\n%s'%self.emission_parameters()
    105126        return s
    106127
     128    def initial_probabilities(self):
     129        """
     130        Return the list of initial state probabilities.
     131
     132        OUTPUT:
     133            list of floats
     134       
     135        EXAMPLES:
     136            sage: m = hmm.GaussianHiddenMarkovModel([[0.0,1.0,0],[0.5,0.0,0.5],[0.3,0.3,0.4]], [(0.0,1.0), (-1.0,0.5), (1.0,0.2)], [0.4,0.3,0.3])
     137            sage: m.initial_probabilities()
     138            [0.4, 0.3, 0.3]
     139        """
     140        cdef Py_ssize_t i
     141        return [self.m.s[i].pi for i in range(self.m.N)]
     142   
     143    def transition_matrix(self, list_only=True):
     144        """
     145        Return the hidden state transition matrix.
     146       
     147        EXAMPLES:
     148            sage: m = hmm.GaussianHiddenMarkovModel([[0.0,1.0,0],[0.5,0.0,0.5],[0.3,0.3,0.4]], [(0.0,1.0), (-1.0,0.5), (1.0,0.2)], [1,0,0])
     149            sage: m.transition_matrix()
     150            [0.9 0.1]
     151            [0.9 0.1]
     152        """
     153        cdef Py_ssize_t i, j
     154        for i from 0 <= i < self.m.N:
     155            for j from 0 <= j < self.m.s[i].out_states:
     156                self.A.set_unsafe_double(i,j,self.m.s[i].out_a[0][j])
     157        return self.A
     158
     159    def emission_parameters(self):
     160        """
     161        Return the emission probability matrix.
     162       
     163        EXAMPLES:
     164            sage: m = hmm.GaussianHiddenMarkovModel([[0.0,1.0,0],[0.5,0.0,0.5],[0.3,0.3,0.4]], [(0.0,1.0), (-1.0,0.5), (1.0,0.2)], [0.1,0.4,0.5])
     165            sage: m.emission_parameters()
     166            [(0.0, 1.0), (-1.0, 0.5), (1.0, 0.20000...)]
     167        """
     168        cdef Py_ssize_t i, j
     169        v = []
     170        for i from 0 <= i < self.m.N:
     171            v.append((self.m.s[i].e.mean.val, self.m.s[i].e.variance.val))
     172        return v
  • sage/stats/hmm/hmm.pyx

    diff -r ba54a820d136 -r 0ec7e93289ae sage/stats/hmm/hmm.pyx
    a b cdef class DiscreteHiddenMarkovModel(Hid 
    200200        """
    201201        Return the list of initial state probabilities.
    202202       
     203        OUTPUT:
     204            list of floats
     205
    203206        EXAMPLES:
    204207            sage: a = hmm.DiscreteHiddenMarkovModel([[0.9,0.1],[0.9,0.1]], [[0.5,0.5,0,0],[0,0,.5,.5]], [0.5,0.5], [1,-1,3,-3])
    205208            sage: a.initial_probabilities()