草庐IT

最大熵原理及代码

晓柒NLP与药物设计 2023-09-28 原文

一.最大熵原理

最大熵的思想很朴素,即将已知事实以外的未知部分看做“等可能”的,而熵是描述“等可能”大小很合适的量化指标,熵的公式如下:

这里分布的取值有种情况,每种情况的概率为,下图绘制了二值随机变量的熵:

p=np.linspace(0.1,0.9,90)
def entropy(p):
    return -np.log(p)*p-np.log(1-p)*(1-p)
plt.plot(p,entropy(p))
[<matplotlib.lines.Line2D at 0x245a3d6d278>]

当两者概率均为0.5时,熵取得最大值,通过最大化熵,可以使得分布更“等可能”;另外,熵还有优秀的性质,它是一个凹函数,所以最大化熵其实是一个凸问题。

对于“已知事实”,可以用约束条件来描述,比如4个值的随机变量分布,其中已知,它的求解可以表述如下:


显然,最优解为:

二.最大熵模型

最大熵模型是最大熵原理在分类问题上的应用,它假设分类模型是一个条件概率分布,即对于给定的输入,以概率输出,这时最大熵模型的目标函数定义为条件熵:

这里,表示边缘分布的经验分布,表示训练样本中输入出现的次数,表示训练样本的总数。

而最大熵模型的“已知事实”可以通过如下等式来约束:

为了方便,左边式子记着,右边式子记着,等式描述的是某函数关于模型与经验分布的期望与函数关于经验分布的期望相同。(这里)
所以重要的约束信息将由来表示,它的定义如下:

故最大熵模型可以理解为,模型在某些事实发生的期望和训练集相同的条件下,使得条件熵最大化。所以,对于有个约束条件的最大熵模型可以表示为:

按照优化问题的习惯,可以改写为如下:

由于目标函数为凸函数,约束条件为仿射,所以我们可以通过求解对偶问题,得到原始问题的最优解,首先引入拉格朗日乘子,定义拉格朗日函数

所以原问题等价于:

它的对偶问题:

首先,解里面的 ,其实对于都是关于的凸函数,因为是关于的凸函数,而后面的是关于的仿射函数,所以求的偏导数,并令其等于0,即可解得最优的,记为,即:

在训练集中对任意样本,都有,显然(本来就是训练集中的一个样本,自然概率大于0),所以,所以:

这就是最大熵模型的表达式(最后一步变换用到了),这里即是模型的参数,聪明的童鞋其实已经发现,最大熵模型其实就是一个线性函数外面套了一个softmax函数,它大概就是如下图所示的这么回事:
[图片上传失败...(image-fb29de-1658666093359)]

接下来,将带入外层的函数,即可求解最优的参数

推导一下模型的梯度更新公式:

这里,倒数第三步到倒数第二步用到了,最后一步中,所以:

所以,自然的更新公式:

这里,是学习率

三.对特征函数的进一步理解

上面推导出了最大熵模型的梯度更新公式,想必大家对还是有点疑惑,“满足某一事实”这句话该如何理解?这其实与我们的学习目的相关,学习目的决定了我们的“事实”,比如有这样一个任务,判断“打”这个词是量词还是动词,我们收集了如下的语料:

句子/ 目标/
一打火柴 量词
三打啤酒 量词
打电话 动词
打篮球 动词

通过观察,我们可以设计如下的两个特征函数来分别识别"量词"和"动词"任务:

当然,你也可以设计这样的特征函数来做识别“量词”的任务:


只是,这样的特征函数设计会使得模型学习能力变弱,比如遇到“三打火柴”,采用后面的特征函数设计就识别不出“打”是量词,而采用第一种特征函数设计就能很好的识别出来,所以要使模型具有更好的泛化能力,就需要设计更好的特征函数,而这往往依赖于人工经验,对于自然语言处理这类任务(比如上面的例子),我们可以较容易的归纳总结出一些有用的经验知识,但是对于其他情况,人工往往难以总结出一般性的规律,所以对于这些问题,我们需要设计更“一般”的特征函数。

一种简单的特征函数设计

我们可以简单考虑的某个特征取某个值和取某个类的组合做特征函数(对于连续型特征,可以采用分箱操作),所以我们可以设计这样两类特征函数:

(1)离散型:

(2)连续型:

四.代码实现

为了方便演示,首先构建训练数据和测试数据

# 测试
from sklearn import datasets
from sklearn import model_selection
from sklearn.metrics import f1_score

iris = datasets.load_iris()
data = iris['data']
target = iris['target']
X_train, X_test, y_train, y_test = model_selection.train_test_split(data, target, test_size=0.2,random_state=0)

为了方便对数据进行分箱操作,封装一个DataBinWrapper类,并对X_train和X_test进行转换(该类放到ml_models.wrapper_models中)

class DataBinWrapper(object):
    def __init__(self, max_bins=10):
        # 分段数
        self.max_bins = max_bins
        # 记录x各个特征的分段区间
        self.XrangeMap = None

    def fit(self, x):
        n_sample, n_feature = x.shape
        # 构建分段数据
        self.XrangeMap = [[] for _ in range(0, n_feature)]
        for index in range(0, n_feature):
            tmp = x[:, index]
            for percent in range(1, self.max_bins):
                percent_value = np.percentile(tmp, (1.0 * percent / self.max_bins) * 100.0 // 1)
                self.XrangeMap[index].append(percent_value)

    def transform(self, x):
        """
        抽取x_bin_index
        :param x:
        :return:
        """
        if x.ndim == 1:
            return np.asarray([np.digitize(x[i], self.XrangeMap[i]) for i in range(0, x.size)])
        else:
            return np.asarray([np.digitize(x[:, i], self.XrangeMap[i]) for i in range(0, x.shape[1])]).T
data_bin_wrapper=DataBinWrapper(max_bins=10)
data_bin_wrapper.fit(X_train)
X_train=data_bin_wrapper.transform(X_train)
X_test=data_bin_wrapper.transform(X_test)
X_train[:5,:]
array([[7, 6, 8, 7],
       [3, 5, 5, 6],
       [2, 8, 2, 2],
       [6, 5, 6, 7],
       [7, 2, 8, 8]], dtype=int64)
X_test[:5,:]
array([[5, 2, 7, 9],
       [5, 0, 4, 3],
       [3, 9, 1, 2],
       [9, 3, 9, 7],
       [1, 8, 2, 2]], dtype=int64)

由于特征函数可以有不同的形式,这里我们将特征函数解耦出来,构造一个SimpleFeatureFunction类(后续构造其他复杂的特征函数,需要定义和该类相同的函数名,该类放置到ml_models.linear_model中)

class SimpleFeatureFunction(object):
    def __init__(self):
        """
        记录特征函数
        {
            (x_index,x_value,y_index)
        }
        """
        self.feature_funcs = set()

    # 构建特征函数
    def build_feature_funcs(self, X, y):
        n_sample, _ = X.shape
        for index in range(0, n_sample):
            x = X[index, :].tolist()
            for feature_index in range(0, len(x)):
                self.feature_funcs.add(tuple([feature_index, x[feature_index], y[index]]))

    # 获取特征函数总数
    def get_feature_funcs_num(self):
        return len(self.feature_funcs)

    # 分别命中了那几个特征函数
    def match_feature_funcs_indices(self, x, y):
        match_indices = []
        index = 0
        for feature_index, feature_value, feature_y in self.feature_funcs:
            if feature_y == y and x[feature_index] == feature_value:
                match_indices.append(index)
            index += 1
        return match_indices

接下来对MaxEnt类进行实现,首先实现一个softmax函数的功能(ml_models.utils)

def softmax(x):
    if x.ndim == 1:
        return np.exp(x) / np.exp(x).sum()
    else:
        return np.exp(x) / np.exp(x).sum(axis=1, keepdims=True)

进行MaxEnt类的具体实现(ml_models.linear_model)

from ml_models import utils
class MaxEnt(object):
    def __init__(self, feature_func, epochs=5, eta=0.01):
        self.feature_func = feature_func
        self.epochs = epochs
        self.eta = eta

        self.class_num = None
        """
        记录联合概率分布:
        {
            (x_0,x_1,...,x_p,y_index):p
        }
        """
        self.Pxy = {}
        """
        记录边缘概率分布:
        {
            (x_0,x_1,...,x_p):p
        }
        """
        self.Px = {}

        """
        w[i]-->feature_func[i]
        """
        self.w = None

    def init_params(self, X, y):
        """
        初始化相应的数据
        :return:
        """
        n_sample, n_feature = X.shape
        self.class_num = np.max(y) + 1

        # 初始化联合概率分布、边缘概率分布、特征函数
        for index in range(0, n_sample):
            range_indices = X[index, :].tolist()

            if self.Px.get(tuple(range_indices)) is None:
                self.Px[tuple(range_indices)] = 1
            else:
                self.Px[tuple(range_indices)] += 1

            if self.Pxy.get(tuple(range_indices + [y[index]])) is None:
                self.Pxy[tuple(range_indices + [y[index]])] = 1
            else:
                self.Pxy[tuple(range_indices + [y[index]])] = 1

        for key, value in self.Pxy.items():
            self.Pxy[key] = 1.0 * self.Pxy[key] / n_sample
        for key, value in self.Px.items():
            self.Px[key] = 1.0 * self.Px[key] / n_sample

        # 初始化参数权重
        self.w = np.zeros(self.feature_func.get_feature_funcs_num())

    def _sum_exp_w_on_all_y(self, x):
        """
        sum_y exp(self._sum_w_on_feature_funcs(x))
        :param x:
        :return:
        """
        sum_w = 0
        for y in range(0, self.class_num):
            tmp_w = self._sum_exp_w_on_y(x, y)
            sum_w += np.exp(tmp_w)
        return sum_w

    def _sum_exp_w_on_y(self, x, y):
        tmp_w = 0
        match_feature_func_indices = self.feature_func.match_feature_funcs_indices(x, y)
        for match_feature_func_index in match_feature_func_indices:
            tmp_w += self.w[match_feature_func_index]
        return tmp_w

    def fit(self, X, y):
        self.eta = max(1.0 / np.sqrt(X.shape[0]), self.eta)
        self.init_params(X, y)
        x_y = np.c_[X, y]
        for epoch in range(self.epochs):
            count = 0
            np.random.shuffle(x_y)
            for index in range(x_y.shape[0]):
                count += 1
                x_point = x_y[index, :-1]
                y_point = x_y[index, -1:][0]
                # 获取联合概率分布
                p_xy = self.Pxy.get(tuple(x_point.tolist() + [y_point]))
                # 获取边缘概率分布
                p_x = self.Px.get(tuple(x_point))
                # 更新w
                dw = np.zeros(shape=self.w.shape)
                match_feature_func_indices = self.feature_func.match_feature_funcs_indices(x_point, y_point)
                if len(match_feature_func_indices) == 0:
                    continue
                if p_xy is not None:
                    for match_feature_func_index in match_feature_func_indices:
                        dw[match_feature_func_index] = p_xy
                if p_x is not None:
                    sum_w = self._sum_exp_w_on_all_y(x_point)
                    for match_feature_func_index in match_feature_func_indices:
                        dw[match_feature_func_index] -= p_x * np.exp(self._sum_exp_w_on_y(x_point, y_point)) / (
                            1e-7 + sum_w)
                # 更新
                self.w += self.eta * dw
                # 打印训练进度
                if count % (X.shape[0] // 4) == 0:
                    print("processing:\tepoch:" + str(epoch + 1) + "/" + str(self.epochs) + ",percent:" + str(
                        count) + "/" + str(X.shape[0]))

    def predict_proba(self, x):
        """
        预测为y的概率分布
        :param x:
        :return:
        """
        y = []
        for x_point in x:
            y_tmp = []
            for y_index in range(0, self.class_num):
                match_feature_func_indices = self.feature_func.match_feature_funcs_indices(x_point, y_index)
                tmp = 0
                for match_feature_func_index in match_feature_func_indices:
                    tmp += self.w[match_feature_func_index]
                y_tmp.append(tmp)
            y.append(y_tmp)
        return utils.softmax(np.asarray(y))

    def predict(self, x):
        return np.argmax(self.predict_proba(x), axis=1)
# 构建特征函数类
feature_func=SimpleFeatureFunction()
feature_func.build_feature_funcs(X_train,y_train)

maxEnt = MaxEnt(feature_func=feature_func)
maxEnt.fit(X_train, y_train)
y = maxEnt.predict(X_test)

print('f1:', f1_score(y_test, y, average='macro'))
processing: epoch:1/5,percent:30/120
processing: epoch:1/5,percent:60/120
processing: epoch:1/5,percent:90/120
processing: epoch:1/5,percent:120/120
processing: epoch:2/5,percent:30/120
processing: epoch:2/5,percent:60/120
processing: epoch:2/5,percent:90/120
processing: epoch:2/5,percent:120/120
processing: epoch:3/5,percent:30/120
processing: epoch:3/5,percent:60/120
processing: epoch:3/5,percent:90/120
processing: epoch:3/5,percent:120/120
processing: epoch:4/5,percent:30/120
processing: epoch:4/5,percent:60/120
processing: epoch:4/5,percent:90/120
processing: epoch:4/5,percent:120/120
processing: epoch:5/5,percent:30/120
processing: epoch:5/5,percent:60/120
processing: epoch:5/5,percent:90/120
processing: epoch:5/5,percent:120/120
f1: 0.9295631904327557

通过前面的分析,我们知道特征函数的复杂程度决定了模型的复杂度,下面我们添加更复杂的特征函数来增强MaxEnt的效果,上面的特征函数仅考虑了单个特征与目标的关系,我们进一步考虑二个特征与目标的关系,即:

如此,我们可以定义一个新的UserDefineFeatureFunction类(注意:类中的方法名称要和SimpleFeatureFunction一样

class UserDefineFeatureFunction(object):
    def __init__(self):
        """
        记录特征函数
        {
            (x_index1,x_value1,x_index2,x_value2,y_index)
        }
        """
        self.feature_funcs = set()

    # 构建特征函数
    def build_feature_funcs(self, X, y):
        n_sample, _ = X.shape
        for index in range(0, n_sample):
            x = X[index, :].tolist()
            for feature_index in range(0, len(x)):
                self.feature_funcs.add(tuple([feature_index, x[feature_index], y[index]]))
                for new_feature_index in range(0,len(x)):
                    if feature_index!=new_feature_index:
                        self.feature_funcs.add(tuple([feature_index, x[feature_index],new_feature_index,x[new_feature_index],y[index]]))

    # 获取特征函数总数
    def get_feature_funcs_num(self):
        return len(self.feature_funcs)

    # 分别命中了那几个特征函数
    def match_feature_funcs_indices(self, x, y):
        match_indices = []
        index = 0
        for item in self.feature_funcs:
            if len(item)==5:
                feature_index1, feature_value1,feature_index2,feature_value2, feature_y=item
                if feature_y == y and x[feature_index1] == feature_value1 and x[feature_index2]==feature_value2:
                    match_indices.append(index)
            else:
                feature_index1, feature_value1, feature_y=item
                if feature_y == y and x[feature_index1] == feature_value1:
                    match_indices.append(index)
            index += 1
        return match_indices
# 检验
feature_func=UserDefineFeatureFunction()
feature_func.build_feature_funcs(X_train,y_train)

maxEnt = MaxEnt(feature_func=feature_func)
maxEnt.fit(X_train, y_train)
y = maxEnt.predict(X_test)

print('f1:', f1_score(y_test, y, average='macro'))
processing: epoch:1/5,percent:30/120
processing: epoch:1/5,percent:60/120
processing: epoch:1/5,percent:90/120
processing: epoch:1/5,percent:120/120
processing: epoch:2/5,percent:30/120
processing: epoch:2/5,percent:60/120
processing: epoch:2/5,percent:90/120
processing: epoch:2/5,percent:120/120
processing: epoch:3/5,percent:30/120
processing: epoch:3/5,percent:60/120
processing: epoch:3/5,percent:90/120
processing: epoch:3/5,percent:120/120
processing: epoch:4/5,percent:30/120
processing: epoch:4/5,percent:60/120
processing: epoch:4/5,percent:90/120
processing: epoch:4/5,percent:120/120
processing: epoch:5/5,percent:30/120
processing: epoch:5/5,percent:60/120
processing: epoch:5/5,percent:90/120
processing: epoch:5/5,percent:120/120
f1: 0.957351290684624

我们可以根据自己对数据的认识,不断为模型添加一些新特征函数去增强模型的效果,只需要修改build_feature_funcsmatch_feature_funcs_indices这两个函数即可(但注意控制函数的数量规模
简单总结一下MaxEnt的优缺点,优点很明显:我们可以diy任意复杂的特征函数进去,缺点也很明显:训练很耗时,而且特征函数的设计好坏需要先验知识,对于某些任务很难直观获取


有关最大熵原理及代码的更多相关文章

  1. ruby - 如何在 buildr 项目中使用 Ruby 代码? - 2

    如何在buildr项目中使用Ruby?我在很多不同的项目中使用过Ruby、JRuby、Java和Clojure。我目前正在使用我的标准Ruby开发一个模拟应用程序,我想尝试使用Clojure后端(我确实喜欢功能代码)以及JRubygui和测试套件。我还可以看到在未来的不同项目中使用Scala作为后端。我想我要为我的项目尝试一下buildr(http://buildr.apache.org/),但我注意到buildr似乎没有设置为在项目中使用JRuby代码本身!这看起来有点傻,因为该工具旨在统一通用的JVM语言并且是在ruby中构建的。除了将输出的jar包含在一个独特的、仅限ruby​​

  2. ruby-on-rails - Rails 源代码 : initialize hash in a weird way? - 2

    在rails源中:https://github.com/rails/rails/blob/master/activesupport/lib/active_support/lazy_load_hooks.rb可以看到以下内容@load_hooks=Hash.new{|h,k|h[k]=[]}在IRB中,它只是初始化一个空哈希。和做有什么区别@load_hooks=Hash.new 最佳答案 查看rubydocumentationforHashnew→new_hashclicktotogglesourcenew(obj)→new_has

  3. ruby-on-rails - 浏览 Ruby 源代码 - 2

    我的主要目标是能够完全理解我正在使用的库/gem。我尝试在Github上从头到尾阅读源代码,但这真的很难。我认为更有趣、更温和的踏脚石就是在使用时阅读每个库/gem方法的源代码。例如,我想知道RubyonRails中的redirect_to方法是如何工作的:如何查找redirect_to方法的源代码?我知道在pry中我可以执行类似show-methodmethod的操作,但我如何才能对Rails框架中的方法执行此操作?您对我如何更好地理解Gem及其API有什么建议吗?仅仅阅读源代码似乎真的很难,尤其是对于框架。谢谢! 最佳答案 Ru

  4. ruby - 模块嵌套代码风格偏好 - 2

    我的假设是moduleAmoduleBendend和moduleA::Bend是一样的。我能够从thisblog找到解决方案,thisSOthread和andthisSOthread.为什么以及什么时候应该更喜欢紧凑语法A::B而不是另一个,因为它显然有一个缺点?我有一种直觉,它可能与性能有关,因为在更多命名空间中查找常量需要更多计算。但是我无法通过对普通类进行基准测试来验证这一点。 最佳答案 这两种写作方法经常被混淆。首先要说的是,据我所知,没有可衡量的性能差异。(在下面的书面示例中不断查找)最明显的区别,可能也是最著名的,是你的

  5. ruby - 寻找通过阅读代码确定编程语言的ruby gem? - 2

    几个月前,我读了一篇关于ruby​​gem的博客文章,它可以通过阅读代码本身来确定编程语言。对于我的生活,我不记得博客或gem的名称。谷歌搜索“ruby编程语言猜测”及其变体也无济于事。有人碰巧知道相关gem的名称吗? 最佳答案 是这个吗:http://github.com/chrislo/sourceclassifier/tree/master 关于ruby-寻找通过阅读代码确定编程语言的rubygem?,我们在StackOverflow上找到一个类似的问题:

  6. ruby - Net::HTTP 获取源代码和状态 - 2

    我目前正在使用以下方法获取页面的源代码:Net::HTTP.get(URI.parse(page.url))我还想获取HTTP状态,而无需发出第二个请求。有没有办法用另一种方法做到这一点?我一直在查看文档,但似乎找不到我要找的东西。 最佳答案 在我看来,除非您需要一些真正的低级访问或控制,否则最好使用Ruby的内置Open::URI模块:require'open-uri'io=open('http://www.example.org/')#=>#body=io.read[0,50]#=>"["200","OK"]io.base_ur

  7. 程序员如何提高代码能力? - 2

    前言作为一名程序员,自己的本质工作就是做程序开发,那么程序开发的时候最直接的体现就是代码,检验一个程序员技术水平的一个核心环节就是开发时候的代码能力。众所周知,程序开发的水平提升是一个循序渐进的过程,每一位程序员都是从“菜鸟”变成“大神”的,所以程序员在程序开发过程中的代码能力也是根据平时开发中的业务实践来积累和提升的。提高代码能力核心要素程序员要想提高自身代码能力,尤其是新晋程序员的代码能力有很大的提升空间的时候,需要针对性的去提高自己的代码能力。提高代码能力其实有几个比较关键的点,只要把握住这些方面,就能很好的、快速的提高自己的一部分代码能力。1、多去阅读开源项目,如有机会可以亲自参与开源

  8. 7个大一C语言必学的程序 / C语言经典代码大全 - 2

    嗨~大家好,这里是可莉!今天给大家带来的是7个C语言的经典基础代码~那一起往下看下去把【程序一】打印100到200之间的素数#includeintmain(){ inti; for(i=100;i 【程序二】输出乘法口诀表#includeintmain(){inti;for(i=1;i 【程序三】判断1000年---2000年之间的闰年#includeintmain(){intyear;for(year=1000;year 【程序四】给定两个整形变量的值,将两个值的内容进行交换。这里提供两种方法来进行交换,第一种为创建临时变量来进行交换,第二种是不创建临时变量而直接进行交换。1.创建临时变量来

  9. git使用常见问题(提交代码,合并冲突) - 2

    文章目录git常用命令(简介,详细参数往下看)Git提交代码步骤gitpullgitstatusgitaddgitcommitgitpushgit代码冲突合并问题方法一:放弃本地代码方法二:合并代码常用命令以及详细参数gitadd将文件添加到仓库:gitdiff比较文件异同gitlog查看历史记录gitreset代码回滚版本库相关操作远程仓库相关操作分支相关操作创建分支查看分支:gitbranch合并分支:gitmerge删除分支:gitbranch-ddev查看分支合并图:gitlog–graph–pretty=oneline–abbrev-commit撤消某次提交git用户名密码相关配置g

  10. ruby - 这两段代码有什么区别? - 2

    打印1:defsum(i)i=i+[2]end$x=[1]sum($x)print$x打印12:defsum(i)i.push(2)end$x=[1]sum($x)print$x后者是修改全局变量$x。为什么它在第二个例子中被修改而不是在第一个例子中?类Array的任何方法(不仅是push)都会发生这种情况吗? 最佳答案 变量范围在这里无关紧要。在第一段代码中,您仅使用赋值运算符=为变量i赋值,而在第二段代码中,您正在修改$x(也称为i)使用破坏性方法push。赋值从不修改任何对象。它只是提供一个名称来引用一个对象。方法要么是破坏性

随机推荐