Statistical-Learning-Method.../HMM/HMM.py

277 lines
18 KiB
Python
Raw Normal View History

2019-12-08 03:11:58 +08:00
#coding=utf-8
#Author:Dodo
#Date:2018-12-10
#Email:lvtengchao@pku.edu.cn
#Blog:www.pkudodo.com
'''
数据集人民日报1998年中文标注语料库
------------------------------
运行结果
-------------------原文----------------------
深圳有个打工者阅览室
去年12月我在广东深圳市出差听说南山区工商分局为打工者建了个免费图书阅览室这件新鲜事引起了我的兴趣
月18日下午我来到了这个阅览室阅览室位于桂庙临南油大道是一间轻体房面积约有40平方米内部装修得整洁干净四周的书架上摆满了书并按政治哲学法律法规文化教育经济科技艺术中国文学外国文学等分类屋中央有两排书架上面也摆满了图书和杂志一些打工青年或站或蹲认真地阅读不时有人到借阅台前办理借书或还书手续南山区在深圳市西边地处城乡结合部外来打工者较多去年2月南山区工商分局局长王安全发现分局对面的公园里常有不少打工者业余时间闲逛有时还滋扰生事为了给这些打工者提供一个充实自己的场所他提议由全分局工作人员捐款兴建一个免费阅览室领导带头群众响应大家捐款1万元购买了近千册图书月6日建在南头繁华的南新路和金鸡路交叉口的阅览室开放了从此这里每天都吸引了众多借书看书的人们其中不仅有打工者还有机关干部公司职员和个体户到了夏天由于阅览室所在地被工程征用南山区工商分局便把阅览室迁到了桂庙阅览室的管理人员是两名青年男的叫张攀女的叫赵阳张攀自己就是湖北来的打工者听说南山区工商分局办免费阅览室便主动应聘来服务阅览室每天从早9时开到晚10时夜里张攀就住在这里他谈起阅览室里的图书翻着一本本的借阅名册如数家珍对图书和工作的挚爱之情溢于言表我在这里碰到南山区华英大厦一位叫聂煜的女青年她说她也是个打工者由于春节探家回来后就要去市内工作很留恋这里的这个免费阅览室想抓紧时间多看些书她还把自己买的几本杂志捐给了阅览室在阅览室的捐书登记簿上记录着这样的数字工商系统内部捐书3550册社会各界捐书250册我在阅览室读到了这样几封感谢信深圳瑞兴光学厂的王志明写道我们这些年轻人远离了家乡来到繁华紧张的都市打工辛劳之余能有机会看书读报感到特别充实深圳文光灯泡厂的江虹说南山区工商分局的干部职工捐款捐书给我们打工者提供良好的学习环境鼓励我们求知上进真是办了一件大好事他们是我们打工者的知音本报记者罗华
-------------------分词后----------------------
深圳|有个||工者|阅览室
去年|||||广东|深圳|市出||||||山区|工商|分局|||工者|建了||免费|图书|阅览室||这件||鲜事|引起||||兴趣|
||下午||我来|||这个|阅览室||阅览|室位||桂庙||临南油|大道|||一间||体房||面积||||平方||||部装|修得|整洁|干净|||||书架||摆满||||||政治||哲学||法律|法规||文化|教育||经济||科技||艺术||中国|文学||外国|文学||分类|||中央||两排|书架||上面||摆满||图书||杂志||一些|打工|青年|或站||||认真|地阅|||不时||人到|借阅|台前|办理|借书||还书|手续|||山区||深圳||西边||地处|城乡|结合部||外来||工者||||去年||||山区|工商|分局|局长||安全|发现|分局|对面||公园||常有|不少||工者|业余|时间|闲逛||有时||滋扰|生事|||||这些||工者|提供|一个|充实|自己||场所|||提议|由全|分局|工作|人员|捐款||兴建|一个|免费|阅览室||领导|带头||群众|响应||大家|捐款||||购买|||千册|图书|||||建在|南头|繁华|||新路||金鸡|路交|叉口||阅览室|开放|||从此||这里|每天||吸引||众多|借书||看书||人们||其中|不仅|||工者||||机关|干部||公司|职员|||体户||||夏天||由于|阅览室|||||工程|征用|||山区|工商|分局|便||阅览|室迁|||桂庙||阅览室||管理|人员||||青年|||||||||||赵阳||||自己|||湖北||||工者|||||山区|工商|分局||免费|阅览室||便|主动||聘来|服务||阅览室|每天||早9时|开到|||||夜里||||||这里||他谈||阅览|室里||图书||翻着||本本||借阅|名册||如数|家珍|||图书||工作||挚爱||情溢||言表||||这里|碰到||山区|华英|大厦|一位||||||青年|||||||||工者||由于|春节|探家|回来||||去市||工作|||留恋|这里||这个|免费|阅览室|||抓紧|时间||看些||||||自己|||几本|杂志|捐给||阅览室|||阅览室||捐书||记簿|||记录||这样||数字||工商|系统|内部|捐书|||社会|各界|捐书|||||阅览室|读到||这样||封感|谢信||深圳|瑞兴|光学||||志明|写道|||我们|这些||轻人|远离||家乡||来到|繁华|紧张||||打工||辛劳|之余||||机会|看书|读报||感到|特别|充实|||深圳|文光||泡厂||江虹|||||山区|工商|分局||干部|职工|捐款||捐书|||我们||工者|提供|良好||学习|环境||鼓励|我们|求知|上进||真是|||一件|大好|||他们||我们||工者||知音||||本报|记者|罗华|
运行时长3.6s
'''
import numpy as np
import time
def trainParameter(fileName):
'''
依据训练文本统计PIAB
:param fileName: 训练文本
:return: 三个参数
'''
#定义一个查询字典,用于映射四种标记在数组中对应的位置,方便查询
# B词语的开头
# M一个词语的中间词
# E一个词语的结果
# S非词语单个词
statuDict = {'B':0, 'M':1, 'E':2, 'S':3}
#每个字只有四种状态所以下方的各类初始化中大小的参数均为4
#初始化PI的一维数组因为对应四种状态大小为4
PI = np.zeros(4)
#初始化状态转移矩阵A涉及到四种状态各自到四种状态的转移因为大小为4x4
A = np.zeros((4, 4))
#初始化观测概率矩阵,分别为四种状态到每个字的发射概率
#因为是中文分词使用ord(汉字)即可找到其对应编码这里用一个65536的空间来保证对于所有的汉字都能
#找到对应的位置来存储
B = np.zeros((4, 65536))
#去读训练文本
fr = open(fileName, encoding='utf-8')
#文本中的每一行认为是一个训练样本
#在统计上三个参数依据“10.3.2” Baum-Welch算法内描述的统计
#PI依据式10.35
#A依据10.37
#B依据10.38
#注并没有使用Baum-Welch算法只是借助了其内部的三个参数生成公式其实
#公式并不是Baum-Welch特有的只是在那一节正好有描述
for line in fr.readlines():
#---------------------训练集单行样例--------------------
#深圳 有 个 打工者 阅览室
#------------------------------------------------------
#可以看到训练样本已经分词完毕,词语之间空格隔开,因此我们在生成统计时主要借助以下思路:
# 1.先将句子按照空格隔开例如例句中5个词语隔开后变成一个长度为5的列表每个元素为一个词语
# 2.对每个词语长度进行判断:
# 如果为1认为该词语是S即单个字
# 如果为2则第一个是B表开头第二个为E表结束
# 如果大于2则第一个为B最后一个为E中间全部标为M表中间词
# 3.统计PI该句第一个字的词性对应的PI中位置加1
# 例如PI = [0 0 0 0]当本行第一个字是B即表示开头时PI中B对应位置为0
# 则PI = [1 0 0 0],全部统计结束后,按照计数值再除以总数得到概率
# 统计A对状态链中位置t和t-1的状态进行统计在矩阵中相应位置加1全部结束后生成概率
# 统计B对于每个字的状态以及字内容生成状态到字的发射计数全部结束后生成概率
# 注可以看一下“10.1.1 隐马尔可夫模型的定义”一节中三个参数的定义,会有更清晰一点的认识
#-------------------------------------------------------
#对单行句子按空格进行切割
curLine = line.strip().split()
#对词性的标记放在该列表中
wordLabel = []
#对每一个单词进行遍历
for i in range(len(curLine)):
#如果长度为1则直接将该字标记为S即单个词
if len(curLine[i]) == 1:
label = 'S'
else:
#如果长度不为1开头为B最后为E中间添加长度-2个M
#如果长度刚好为2长度-2=0也就不添加了反之添加对应个数的M
label = 'B' + 'M' * (len(curLine[i]) - 2) + 'E'
#如果是单行开头第一个字PI中对应位置加1,
if i == 0: PI[statuDict[label[0]]] += 1
#对于该单词中的每一个字在生成的状态链中统计B
for j in range(len(label)):
#遍历状态链中每一个状态并找到对应的中文汉字在B中
#对应位置加1
B[statuDict[label[j]]][ord(curLine[i][j])] += 1
#在整行的状态链中添加该单词的状态链
#注意extend表直接在原先元素的后方添加
#可以百度一下extend和append的区别
wordLabel.extend(label)
#单行所有单词都结束后统计A信息
#因为A涉及到前一个状态因此需要等整条状态链都生成了才能开始统计
for i in range(1, len(wordLabel)):
#统计t时刻状态和t-1时刻状态的所有状态组合的出现次数
A[statuDict[wordLabel[i - 1]]][statuDict[wordLabel[i]]] += 1
#上面代码在统计上全部是统计的次数,实际运算需要使用概率,
#下方代码是将三个参数的次数转换为概率
#----------------------------------------
#对PI求和概率生成中的分母
sum = np.sum(PI)
#遍历PI中每一个元素元素出现的次数/总次数即为概率
for i in range(len(PI)):
#如果某元素没有出现过该位置为0在后续的计算中这是不被允许的
#比如说某个汉字在训练集中没有出现过,那在后续不同概率相乘中只要有
#一项为0其他都是0了此外整条链很长的情况下太多0-1的概率相乘
#不管怎样最后的结果都会很小,很容易下溢出
#所以在概率上我们习惯将其转换为log对数形式这在书上是没有讲的
#x大的时候log也大x小的时候log也相应小我们最后比较的是不同
#概率的大小所以使用log没有问题
#那么当单向概率为0的时候log没有定义因此需要单独判断
#如果该项为0则手动赋予一个极小值
if PI[i] == 0: PI[i] = -3.14e+100
#如果不为0则计算概率再对概率求log
else: PI[i] = np.log(PI[i] / sum)
#与上方PI思路一样求得A的概率对数
for i in range(len(A)):
sum = np.sum(A[i])
for j in range(len(A[i])):
if A[i][j] == 0: A[i][j] = -3.14e+100
else: A[i][j] = np.log(A[i][j] / sum)
#与上方PI思路一样求得B的概率对数
for i in range(len(B)):
sum = np.sum(B[i])
for j in range(len(B[i])):
if B[i][j] == 0: B[i][j] = -3.14e+100
else:B[i][j] = np.log(B[i][j] / sum)
#返回统计得到的三个参数
return PI, A, B
def loadArticle(fileName):
'''
加载文章
:param fileName:文件路径
:return: 文章内容
'''
#初始化文章列表
artical = []
#打开文件
fr = open(fileName, encoding='utf-8')
#按行读取文件
for line in fr.readlines():
#读到的每行最后都有一个\n使用strip将最后的回车符去掉
line = line.strip()
#将该行放入文章列表中
artical.append(line)
#将文章返回
return artical
def participle(artical, PI, A, B):
'''
分词
算法依据10.4.2 维特比算法
:param artical:要分词的文章
:param PI: 初始状态概率向量PI
:param A: 状态转移矩阵
:param B: 观测概率矩阵
:return: 分词后的文章
'''
#初始化分词后的文章列表
retArtical = []
#对文章按行读取
for line in artical:
#初始化δ,δ存放四种状态的概率值,因为状态链中每个状态都有
#四种概率值,因此长度时该行的长度
delta = [[0 for i in range(4)] for i in range(len(line))]
#依据算法10.5 第一步:初始化
for i in range(4):
#初始化δ状态链中第一个状态的四种状态概率
delta[0][i] = PI[i] + B[i][ord(line[0])]
#初始化ψ初始时为0
psi = [[0 for i in range(4)] for i in range(len(line))]
#算法10.5中的第二步:递推
#for循环的符号与书中公式一致可以对比着看来理解
#依次处理整条链
for t in range(1, len(line)):
#对于链中的米格状态,求四种状态概率
for i in range(4):
#初始化一个临时列表,用于存放四种概率
tmpDelta = [0] * 4
for j in range(4):
# 计算第二步中的δ该部分只计算max内部不涉及后面的bi(o)
# 计算得到四个结果以后再去求那个max即可
# 注bi(Ot)并不在max的式子中是求出max以后再乘b的
# 此外读者可能注意到书中的乘法在这里变成了加法,这是由于原先是概率
# 直接相乘但我们在求得概率时同时取了log取完log以后概率的乘法
# 也就转换为加法了,同时也简化了运算
# 所以log优点还是很多的对不
tmpDelta[j] = delta[t - 1][j] + A[j][i]
#找到最大的那个δ * a
maxDelta = max(tmpDelta)
#记录最大值对应的状态
maxDeltaIndex = tmpDelta.index(maxDelta)
#将找到的最大值乘以b放入
#注意这里同样因为log变成了加法
delta[t][i] = maxDelta + B[i][ord(line[t])]
#在ψ中记录对应的最大状态索引
psi[t][i] = maxDeltaIndex
#建立一个状态链列表,开始生成状态链
sequence = []
#算法10.5 第三步:终止
#在上面for循环全部结束后很明显就到了第三步了
#获取最后一个状态的最大状态概率对应的索引
i_opt = delta[len(line) - 1].index(max(delta[len(line) - 1]))
#在状态链中添加索引
#注状态链应该是B、M、E、S这里图方便用了0、1、2、3其实一样的
sequence.append(i_opt)
#算法10.5 第四步:最优路径回溯
#从后往前遍历整条链
for t in range(len(line) - 1, 0, -1):
#不断地从当前时刻t的ψ列表中读取到t-1的最优状态
i_opt = psi[t][i_opt]
#将状态放入列表中
sequence.append(i_opt)
#因为是从后往前将状态放入的列表,所以这里需要翻转一下,变成了从前往后
sequence.reverse()
#开始对该行分词
curLine = ''
#遍历该行每一个字
for i in range(len(line)):
#在列表中放入该字
curLine += line[i]
#如果该字是3S->单个词 或 2:E->结尾词 ,则在该字后面加上分隔符 |
#此外如果改行的最后一个字了,也就不需要加 |
if (sequence[i] == 3 or sequence[i] == 2) and i != (len(line) - 1):
curLine += '|'
#在返回列表中添加分词后的该行
retArtical.append(curLine)
#返回分词后的文章
return retArtical
if __name__ == '__main__':
# 开始时间
start = time.time()
#依据现有训练集统计PI、A、B
PI, A, B = trainParameter('HMMTrainSet.txt')
#读取测试文章
artical = loadArticle('testArtical.txt')
#打印原文
print('-------------------原文----------------------')
for line in artical:
print(line)
#进行分词
partiArtical = participle(artical, PI, A, B)
#打印分词结果
print('-------------------分词后----------------------')
for line in partiArtical:
print(line)
#结束时间
print('time span:', time.time() - start)