[机器学习] 使用极大似然估计与贝叶斯估计拟合线性基函数

本文链接：https://blog.csdn.net/weixin_43720225/article/details/130793485

文章介绍了概率模型的概念，包括二项分布、正态分布等，并重点讲述了极大似然估计（MLE）和贝叶斯估计。通过数学公式展示了如何计算模型参数，并以线性基函数模型为例，解释了如何拟合非线性数据。最后，通过实验部分，使用Python生成非线性数据并应用MLE和贝叶斯估计进行模型训练，可视化结果展示了不同数据集和参数K对模型的影响。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

一、概率模型

概率模型可以写作：
$P(·|\theta)$
指模型的分布类型已知（二项分布、正态/高斯分布，Bate分布等），但是模型的参数未知。基于独立同分布的样本数据D对模型参数进行估计。常用的方法有极大似然估计（MLE）、极大后验估计（MAP）、贝叶斯估计。

二、极大似然估计Maximum Likelihood Estimation

MLE的实质为计算在参数下D出现的概率，需要计算D中每个样本出现概率的乘积，即似然函数。表达式如下：
$\{x_i\}_{i=1}^n \\ Likelihood(\theta) = P(D|\theta) = \Pi_{i = 1}^nP(x_i|\theta) \quad 似然函数 \\ lnL(\theta) = \Sigma_{i=1}^nln P(x_i|\theta)\quad 对数似然函数$
对似然函数取对数可以将概率乘法变为加法，将参数估计问题变为函数优化问题，从而可以得到能够使观察数据D出现概率最大的概率模型P参数来。

三、贝叶斯估计Bayesian

贝叶斯估计表达式如下。其中p(A|B)为B条件下A的概率，w为参数，D为独立同分布数据样本：
$\frac{p(D|w)p(w)}{p(D)}$
贝叶斯模型生成参数模型的概率分布（后验概率），其实质为基于样本数据调整先验概率，得到更合理的后验概率。

四、线性基函数模型

用于拟合非线性数据，普通模型的计算为：
$y = x^t·w$
线性基函数即对输入x进行线性变换后，再与参数进行计算。一般是为了对输入特征的维度进行变换，参数维度K决定了模型的描述能力：K=0，模型为一条水平线；K=1，模型为一条直线…然而K的值不是越高越好，过大的K会使得模型容易过拟合。

$\phi^t(x)·w+ \epsilon \quad(x\in R^D,w \in R^K,\epsilon \sim N(0,\sigma^2)) \\\phi(x) = \begin{vmatrix} \phi_0(x)\\ \phi_1(x)\\ \vdots\\ \phi_{K-1}(x)\\ \end{vmatrix} = \begin{vmatrix} 1\\ x\\ x^2\\ \vdots\\ x^{K-1}\\ \end{vmatrix}\in R^K \\\Longrightarrow \quad y = w_01 + w_1x + \cdots + w_{k-1}x^{k-1} +\epsilon$

五、实验

1、训练集生成

使用sin函数和方差为0.2高斯噪声生成S、M、L三个数据集，大小分别为10、30、60。

# 使用sin函数和方差为2高斯噪声生成非线性数据样本
import numpy as np
import random

# n 生成样本数目
def data_gen(num):
    # 在[-5,5]内生成num个样本
    X = np.linspace(-5,5,num,endpoint=True)
    Y = np.sin(X)

    # 加入高斯噪声
    for i in range(X.size):
        X[i] += random.gauss(0,0.2)
        Y[i] += random.gauss(0,0.2)
    return X,Y

2、模型训练

使用线性基函数作为回归模型，针对不同的K，使用MLE在不同数据集上进行拟合。

对于观测数据X和Y，有：
$w_{MLE} = (\Phi^T\Phi)^{-1}\Phi^Ty \\ \sigma_{MLE}^2 = \frac1N \sum_{i=1}^N(y_i - \Phi^T(x_i)w)^2$
其中：
$\Phi =\begin{vmatrix} \phi_0(x_1) & \cdots & \phi_{K-1}(x_1)\\ \phi_0(x_2) & \cdots & \phi_{K-1}(x_2)\\ \vdots & \ddots & \vdots\\ \phi_0(x_N) & \cdots & \phi_{K-1}(x_N)\\ \end{vmatrix} \in R^{N*K}$

import math

# 线性基函数模型，使用MLE拟合
class model_MLE:
    # K为w的维数，w为一个K维向量，sigma为噪声方差，XY为训练数据
    def __init__(self,K,X,Y):
        self.K = K
        
        self.X = X
        self.Y = Y
        
        self.Phi = np.array(self.Phi_gen())
        self.w = self.w_gen()
        self.sigma = self.sigma_gen()
    
    # 极大似然估计计算w
    def w_gen(self):
        w = np.dot(np.dot(np.linalg.inv(np.dot(self.Phi.T,self.Phi)),self.Phi.T),self.Y)
        return w
    
    # 极大似然估计计算sigma
    def sigma_gen(self):
        sigma = 0
        for i in range(self.X.size):
            sigma += math.pow((self.Y[i] - np.dot(np.array(self.Phi_func(self.X[i])).T, self.w)),2)
        sigma /= self.X.size
        return sigma
    
    # 依据输入X和K生成Phi矩阵
    def Phi_gen(self):
        Phi = []
        for i in range(self.X.size):
            Phi_x = self.Phi_func(self.X[i])
            Phi.append(Phi_x)
        return Phi
    
    # Phi函数，对单个x进行处理
    def Phi_func(self,x):
        # 对x进行维度变换
        Phi_x =[]
        item = 1
        for i in range(self.K):
            Phi_x.append(item)
            item *= x
#         Phi_x = np.array(Phi_x)
        return Phi_x

    def forcast(self,x):
        Phi_x = self.Phi_func(x)
#         print(Phi_x,self.w)
        y = np.dot(Phi_x,self.w)
        return y

使用线性基函数作为回归模型，针对不同的K，使用贝叶斯估计在不同数据集上进行拟合。

贝叶斯估计将模型参数w和噪声视为一个参数未知的概率分布，通过求该概率分布的均值和方差就可以得到模型w的分布概率。设模型的参数的先验分布概率服从高斯分布：（其中I为单位矩阵）
$\sim N(\mu_0,\varSigma_0) = N(0,\frac{1}{\lambda}I) \\ \epsilon \sim N(0,\sigma^2)$
有：
$P(w|y,X,\sigma^2) \sim N(\mu_{new}|\Sigma_{new}) \\ \varSigma_{new} = (\frac{1}{\sigma^2}\Phi^T\Phi + \varSigma_0^{-1})^{-1} = (\frac{1}{\sigma^2}\Phi^T\Phi + \lambda I)^{-1} \\ \mu_{new} = \varSigma_{new}(\frac{1}{\sigma^2}\Phi^Ty+\varSigma_0^{-1}\mu_0) = \varSigma_{new}(\frac{1}{\sigma^2}\Phi^Ty)$

# 线性基函数模型参数概率分布，使用贝叶斯估计拟合，得到模型参数的后验概率分布
class distribution_BYS:
    # K为w的维数，w为一个K维向量，sigma为噪声的方差
    def __init__(self,K,X,Y):
        self.K = K
        
        self.X = X
        self.Y = Y
        
        # 模型参数先验分布
        self.mu = 0
#       I为单位矩阵
        self.I = np.identity(self.K, dtype = int)
        self.sigma = 0.25*self.I
        
        # 模型参数后验分布
        self.Phi = self.Phi_gen()
        self.sigma_new = self.sigma_new_gen()
        self.mu_new = self.mu_new_gen()
        
    def sigma_new_gen(self):
        sigma_new = np.linalg.inv(np.dot(self.Phi.T,self.Phi)/0.2+4*self.I)
        return sigma_new
    
    def mu_new_gen(self):
        mu_new = np.dot(self.sigma_new,(np.dot(self.Phi.T,self.Y)/0.2))
        return mu_new
    
     # 依据输入X和K生成Phi矩阵
    def Phi_gen(self):
        Phi = []
        for i in range(self.X.size):
            Phi_x = self.Phi_func(self.X[i])
            Phi.append(Phi_x)
        Phi = np.array(Phi)
        return Phi
    
    # Phi函数，对单个x进行处理
    def Phi_func(self,x):
        # 对x进行维度变换
        Phi_x =[]
        item = 1
        for i in range(self.K):
            Phi_x.append(item)
            item *= x
#         Phi_x = np.array(Phi_x)
        return Phi_x
    
    # 基于模型参数后验概率分布随机生成w，实现模型的随机取样
    def model_gen(self):
        w = []
        for i in range(self.K):
            w.append(random.gauss(self.mu_new[i],self.sigma_new[i][i]))
        w = np.array(w)
        return w

类实例化，使用不同的K在不同大小的数据集上进行拟合，代码略。

3、可视化

样本数据可视化：

# 样本数据可视化
import matplotlib.pyplot as plt

plt.figure(figsize = (12,8))

plt.subplot(2,2,1)
plt.plot(S_X,S_Y,linestyle='',marker='.',label = 'S')
plt.plot(np.linspace(-5,5),np.sin(np.linspace(-5,5)),label = 'Sin(x)')
plt.legend(loc = 1)

plt.subplot(2,2,2)
plt.plot(M_X,M_Y,linestyle='',marker='.',label = 'M')
plt.plot(np.linspace(-5,5),np.sin(np.linspace(-5,5)),label = 'Sin(x)')
plt.legend(loc = 1)

plt.subplot(2,2,3)
plt.plot(L_X,L_Y,linestyle='',marker='.',label = 'L')
plt.plot(np.linspace(-5,5),np.sin(np.linspace(-5,5)),label = 'Sin(x)')
plt.legend(loc = 1)

plt.show()

MLE拟合结果可视化：

# print(model_MLE_1.forcast(X[0]))
# MLE拟合结果可视化
def model_show(model,lab):
    x = np.linspace(-5,5,100)
    y = []
    for i in range(x.size):
        y.append(model.forcast(x[i]))
    plt.plot(x,y,label = lab)

plt.figure(figsize = (12,8))

plt.subplot(2,2,1)
plt.plot(S_X,S_Y,linestyle='',marker='.',label = 'S')
model_show(model_MLE_1,"K=2")
model_show(model_MLE_2,"K=4")
model_show(model_MLE_3,"K=8")
plt.legend(loc = 1)

plt.subplot(2,2,2)
plt.plot(M_X,M_Y,linestyle='',marker='.',label = 'M')
model_show(model_MLE_4,"K=2")
model_show(model_MLE_5,"K=4")
model_show(model_MLE_6,"K=8")
plt.legend(loc = 1)

plt.subplot(2,2,3)
plt.plot(L_X,L_Y,linestyle='',marker='.',label = 'L')
model_show(model_MLE_7,"K=2")
model_show(model_MLE_8,"K=4")
model_show(model_MLE_9,"K=8")
plt.legend(loc = 1)

plt.show()

贝叶斯估计拟合结果可视化：

# 贝叶斯拟合结果可视化
def dis_show(dis):
    x = np.linspace(-5,5,100)
    for i in range(20):
        w = dis.model_gen()
#         print(w)
        y = []
        for j in range(x.size):
            y.append(np.dot(dis.Phi_func(x[j]),w.T))
        plt.plot(x,y)


def draw(loc,dis,data_lab,tit):
    plt.subplot(3,3,loc)
    if data_lab == 'S':
        plt.plot(S_X,S_Y,linestyle='',marker='.',label = data_lab)
    elif data_lab == 'M':
        plt.plot(M_X,M_Y,linestyle='',marker='.',label = data_lab)
    else:
        plt.plot(L_X,L_Y,linestyle='',marker='.',label = data_lab)
    dis_show(dis)
    plt.title(tit)
    plt.legend(loc = 1)

plt.figure(figsize = (12,9))
draw(1,dis_BYS_1,'S','K=2')
draw(2,dis_BYS_2,'S','K=4')
draw(3,dis_BYS_3,'S','K=8')

draw(4,dis_BYS_4,'M','K=2')
draw(5,dis_BYS_5,'M','K=4')
draw(6,dis_BYS_6,'M','K=8')

draw(7,dis_BYS_7,'L','K=2')
draw(8,dis_BYS_8,'L','K=4')
draw(9,dis_BYS_9,'L','K=8')

plt.show()