参考资料:http://www.52nlp.cn/lda-math-mcmc-%e5%92%8c-gibbs-sampling1
python脚本:
import math
import time
import numpy as np
import matplotlib.pylab as plt
import random
class Samples:
def __init__(self):
pass
def rand(self, num, seed = 1):
m = math.pow(2, 32)
a = 214013
c = 2531011
i = 1
x = np.zeros(num)
x[0] = seed
while(i < num):
x[i] = (a * x[i-1] + c) % m
i += 1
return x
def uniform(self, num, seed = 1):
m = math.pow(2, 32)
x = self.rand(num, seed)
return x / m
def normal(self, num):
t = int(time.time())
u1 = self.uniform(num, t)
t1 = int(time.time())+1
u2 = self.uniform(num, t1)
z0 = np.sqrt(-2 * np.log(u2)) * np.sin(2.0 * np.pi * u1)
#plt.hist(z0)
#plt.show()
return z0
if __name__=='__main__':
s = Samples()
x00 = s.normal(1000)
plt.hist(x00,100)
plt.show()
样本的直方图如下: