-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathdmr.py
103 lines (89 loc) · 3.43 KB
/
dmr.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
#!/usr/bin/env python
# -*- coding: utf-8 -*-
import numpy as np
import scipy.special as special
import scipy.optimize as optimize
import scipy.misc as misc
from lda import LDA
class DMR(LDA):
'''
Topic Model with Dirichlet Multinomial Regression
'''
def __init__(self, K, sigma, beta, docs, vecs, V, trained=None):
super(DMR, self).__init__(K, 0.0, beta, docs, V, trained)
self.L = vecs.shape[1]
self.vecs = vecs
self.sigma = sigma
self.Lambda = np.random.multivariate_normal(np.zeros(self.L),
(self.sigma ** 2) * np.identity(self.L), size=self.K)
if self.trained is None:
self.alpha = self.get_alpha()
else:
self.alpha = self.get_alpha(self.trained.Lambda)
def hyperparameter_learning(self):
'''
update alpha (overwrite)
'''
if self.trained is None:
self.bfgs()
self.alpha = self.get_alpha()
def get_alpha_n_m_z(self, idx=None):
if idx is None:
return self.n_m_z + self.alpha
else:
return self.n_m_z[idx] + self.alpha[idx]
def get_alpha(self, Lambda=None):
'''
alpha = exp(Lambda^T x)
'''
if Lambda is None:
if self.trained is None:
Lambda = self.Lambda
else:
Lambda = self.trained.Lambda
return np.exp(np.dot(self.vecs, Lambda.T))
def bfgs(self):
def ll(x):
x = x.reshape((self.K, self.L))
return self._ll(x)
def dll(x):
x = x.reshape((self.K, self.L))
result = self._dll(x)
result = result.reshape(self.K * self.L)
return result
Lambda = np.random.multivariate_normal(np.zeros(self.L),
(self.sigma ** 2) * np.identity(self.L), size=self.K)
Lambda = Lambda.reshape(self.K * self.L)
newLambda, fmin, res = optimize.fmin_l_bfgs_b(ll, Lambda, dll)
self.Lambda = newLambda.reshape((self.K, self.L))
def _ll(self, x):
result = 0.0
# P(w|z)
result += self.K * special.gammaln(self.beta * self.K)
result += -np.sum(special.gammaln(np.sum(self.n_z_w, axis=1)))
result += np.sum(special.gammaln(self.n_z_w))
result += -self.V * special.gammaln(self.beta)
# P(z|Lambda)
alpha = self.get_alpha(x)
result += np.sum(special.gammaln(np.sum(alpha, axis=1)))
result += -np.sum(special.gammaln(
np.sum(self.n_m_z+alpha, axis=1)))
result += np.sum(special.gammaln(self.n_m_z+alpha))
result += -np.sum(special.gammaln(alpha))
# P(Lambda)
result += -self.K / 2.0 * np.log(2.0 * np.pi * (self.sigma ** 2))
result += -1.0 / (2.0 * (self.sigma ** 2)) * np.sum(x ** 2)
result = -result
return result
def _dll(self, x):
alpha = self.get_alpha(x)
result = np.sum(self.vecs[:,np.newaxis,:] * alpha[:,:,np.newaxis]\
* (special.digamma(np.sum(alpha, axis=1))[:,np.newaxis,np.newaxis]\
- special.digamma(np.sum(self.n_m_z+alpha, axis=1))[:,np.newaxis,np.newaxis]\
+ special.digamma(self.n_m_z+alpha)[:,:,np.newaxis]\
- special.digamma(alpha)[:,:,np.newaxis]), axis=0)\
- x / (self.sigma ** 2)
result = -result
return result
def params(self):
return '''K=%d, sigma=%s, beta=%s''' % (self.K, self.sigma, self.beta)