| 22 | |
| 23 | |
| 24 | class HMM: |
| 25 | def __init__(self, M): |
| 26 | self.M = M # number of hidden states |
| 27 | |
| 28 | def fit(self, X, learning_rate=0.001, max_iter=10, V=None, print_period=1): |
| 29 | # train the HMM model using stochastic gradient descent |
| 30 | # print "X to train:", X |
| 31 | |
| 32 | # determine V, the vocabulary size |
| 33 | # assume observables are already integers from 0..V-1 |
| 34 | # X is a jagged array of observed sequences |
| 35 | if V is None: |
| 36 | V = max(max(x) for x in X) + 1 |
| 37 | N = len(X) |
| 38 | print("number of train samples:", N) |
| 39 | |
| 40 | preSoftmaxPi0 = np.zeros(self.M) # initial state distribution |
| 41 | preSoftmaxA0 = np.random.randn(self.M, self.M) # state transition matrix |
| 42 | preSoftmaxB0 = np.random.randn(self.M, V) # output distribution |
| 43 | |
| 44 | thx, cost = self.set(preSoftmaxPi0, preSoftmaxA0, preSoftmaxB0) |
| 45 | |
| 46 | pi_update = self.preSoftmaxPi - learning_rate*T.grad(cost, self.preSoftmaxPi) |
| 47 | A_update = self.preSoftmaxA - learning_rate*T.grad(cost, self.preSoftmaxA) |
| 48 | B_update = self.preSoftmaxB - learning_rate*T.grad(cost, self.preSoftmaxB) |
| 49 | |
| 50 | updates = [ |
| 51 | (self.preSoftmaxPi, pi_update), |
| 52 | (self.preSoftmaxA, A_update), |
| 53 | (self.preSoftmaxB, B_update), |
| 54 | ] |
| 55 | |
| 56 | train_op = theano.function( |
| 57 | inputs=[thx], |
| 58 | updates=updates, |
| 59 | allow_input_downcast=True, |
| 60 | ) |
| 61 | |
| 62 | costs = [] |
| 63 | for it in range(max_iter): |
| 64 | if it % print_period == 0: |
| 65 | print("it:", it) |
| 66 | |
| 67 | for n in range(N): |
| 68 | # this would of course be much faster if we didn't do this on |
| 69 | # every iteration of the loop |
| 70 | c = self.get_cost_multi(X).sum() |
| 71 | costs.append(c) |
| 72 | train_op(X[n]) |
| 73 | |
| 74 | # print "A:", self.A.get_value() |
| 75 | # print "B:", self.B.get_value() |
| 76 | # print "pi:", self.pi.get_value() |
| 77 | plt.plot(costs) |
| 78 | plt.show() |
| 79 | |
| 80 | def get_cost(self, x): |
| 81 | # returns log P(x | model) |