贝叶斯线性回归模型 R语言贝叶斯线性回归的预测

转载

mob6454cc65110a 2024-01-15 22:47:37

文章标签 贝叶斯线性回归模型 R语言数据集线性回归数据 文章分类 R语言后端开发

贝叶斯推断

贝叶斯定理：通过观察到的数据D D ，把先验概率p(θ)p(θ)转化为后验概率p(θ∣D) p ( θ ∣ D )

p(θ∣D)=p(D∣θ)p(θ)∫p(D∣θ)p(θ)dθ=p(D∣θ)p(θ)p(D) p ( θ ∣ D ) = p ( D ∣ θ ) p ( θ ) ∫ p ( D ∣ θ ) p ( θ ) d θ = p ( D ∣ θ ) p ( θ ) p ( D )

显然，分母是一个归一化常数。故有p(θ∣D)∝p(D∣θ)p(θ) p ( θ ∣ D ) ∝ p ( D ∣ θ ) p ( θ ) 即后验∝似然×先验后验 ∝ 似然 × 先验。

贝叶斯线性回归

问题是这样的，不能够一次性接收到整个数据集，而是不断接收到小的数据集Di,i=1,2,...,n D i , i = 1 , 2 , . . . , n ，同时由于存储的限制不能存储已经接收到的所有数据集,每次可以处理的数据仅为Di D i 。这就导致不能对所有数据做线性回归，但是可以通过贝叶斯线性回归达到同样的效果。

第i i 个数据集Di D i 中有m m 个训练样本，构成 (X(i),y(i))(X(i),y(i))

p(y(i)∣X(i),θ)=N(y(i);X(i)θ,I)∝exp(−12(y(i)−X(i)θ)T(y(i)−X(i)θ)) p ( y ( i ) ∣ X ( i ) , θ ) = N ( y ( i ) ; X ( i ) θ , I ) ∝ e x p ( − 1 2 ( y ( i ) − X ( i ) θ ) T ( y ( i ) − X ( i ) θ ) )

为了确定模型参数向量 θ θ 的后验分布
假设其先验分布

p(θ)=N(θ;μ0,Λ0)∝exp(−12(θ−μ0)TΛ−10(θ−μ0)) p ( θ ) = N ( θ ; μ 0 , Λ 0 ) ∝ e x p ( − 1 2 ( θ − μ 0 ) T Λ 0 − 1 ( θ − μ 0 ) )

其中μ0 μ 0 ,Λ0 Λ 0 分别是先验分布的均值向量和协方差矩阵。通过贝叶斯回归得到的目标为θ θ 的期望。
模型参数的后验分布：

p(θ∣X(i),y(i)) ∝p(y(i)∣X(i),θ)p(θ) ∝exp(−12(y(i)−X(i)θ)T(y(i)−X(i)θ))exp(−12(θ−μ0)TΛ−10(θ−μ0)) ∝exp(−12(−2y(i)TX(i)θ+θTX(i)TX(i)θ+θTΛ−10θ−2μT0Λ−10θ)) p ( θ ∣ X ( i ) , y ( i ) ) ∝ p ( y ( i ) ∣ X ( i ) , θ ) p ( θ ) ∝ e x p ( − 1 2 ( y ( i ) − X ( i ) θ ) T ( y ( i ) − X ( i ) θ ) ) e x p ( − 1 2 ( θ − μ 0 ) T Λ 0 − 1 ( θ − μ 0 ) ) ∝ e x p ( − 1 2 ( − 2 y ( i ) T X ( i ) θ + θ T X ( i ) T X ( i ) θ + θ T Λ 0 − 1 θ − 2 μ 0 T Λ 0 − 1 θ ) )

Λi=(X(i)TX(i)+Λ−10)−1,μi=Λi(X(i)Ty(i)+Λ−10μ0) Λ i = ( X ( i ) T X ( i ) + Λ 0 − 1 ) − 1 , μ i = Λ i ( X ( i ) T y ( i ) + Λ 0 − 1 μ 0 )

p(θ∣X(i),y(i))∝exp(−12(θ−μi)TΛ−1i(θ−μi)) p ( θ ∣ X ( i ) , y ( i ) ) ∝ e x p ( − 1 2 ( θ − μ i ) T Λ i − 1 ( θ − μ i ) )

缺点：

1 参数先验分布的不同假设形式，可能会带来计算上的不便。
2 参数先验分布的假设有偏，对于小数据会有较大的影响。

解决方法：

1 参数的先验分布假设为数据分布假设的共轭先验
共轭先验：对于一个给定的概率分布p(x∣∣w) p ( x ∣∣ w ) ,能够寻找一个先验 p(w) p ( w )

Bern(x∣μ)=μx(1−μ)1−xBeta(μ∣a,b)=Γ(a+b)Γ(a)Γ(b)μa−1(1−μ)b−1 B e r n ( x ∣ μ ) = μ x ( 1 − μ ) 1 − x B e t a ( μ ∣ a , b ) = Γ ( a + b ) Γ ( a ) Γ ( b ) μ a − 1 ( 1 − μ ) b − 1

N(x∣m,Λ)=12π|Λ|exp(−12(x−m)TΛ−1(x−m))N(m∣μ,Λ′)=12π|Λ′|exp(−12(m−μ)T(Λ′)−1(m−μ))W(Λ∣W,v)=B|Λ|v−D−12exp(−12Tr(W−1Λ)) N ( x ∣ m , Λ ) = 1 2 π | Λ | e x p ( − 1 2 ( x − m ) T Λ − 1 ( x − m ) ) N ( m ∣ μ , Λ ′ ) = 1 2 π | Λ ′ | e x p ( − 1 2 ( m − μ ) T ( Λ ′ ) − 1 ( m − μ ) ) W ( Λ ∣ W , v ) = B | Λ | v − D − 1 2 e x p ( − 1 2 T r ( W − 1 Λ ) )

2 合理初始化，迭代求解
对于接收到的第1个数据集有：

Λ1=(X(1)TX(1)+Λ−10)−1,μ1=Λ1(X(1)Ty(1)+Λ−10μ0) Λ 1 = ( X ( 1 ) T X ( 1 ) + Λ 0 − 1 ) − 1 , μ 1 = Λ 1 ( X ( 1 ) T y ( 1 ) + Λ 0 − 1 μ 0 )

p(θ∣X(1),y(1))∝exp(−12(θ−μ1)TΛ−11(θ−μ1)) p ( θ ∣ X ( 1 ) , y ( 1 ) ) ∝ e x p ( − 1 2 ( θ − μ 1 ) T Λ 1 − 1 ( θ − μ 1 ) )

这里根据极大似然估计得到的解θ=(X(1)TX(1))−1X(1)Ty(1) θ = ( X ( 1 ) T X ( 1 ) ) − 1 X ( 1 ) T y ( 1 ) ，所以假设Λ−10=O Λ 0 − 1 = O , 此时极大似然的解和贝叶斯回归的参数期望一致。

对于接收到的第i i 个数据集Di(i>1)Di(i>1)，将第i−1 i − 1

Λi=(X(i)TX(i)+Λ−1i−1)−1,μi=Λi(X(i)Ty(i)+Λ−1i−1μi−1) Λ i = ( X ( i ) T X ( i ) + Λ i − 1 − 1 ) − 1 , μ i = Λ i ( X ( i ) T y ( i ) + Λ i − 1 − 1 μ i − 1 )

p(θ∣X(i),y(i))∝exp(−12(θ−μi)TΛ−1i(θ−μi)) p ( θ ∣ X ( i ) , y ( i ) ) ∝ e x p ( − 1 2 ( θ − μ i ) T Λ i − 1 ( θ − μ i ) )

具体算法
输入：D1,D2,D3,...,Dn D 1 , D 2 , D 3 , . . . , D n 其中 Di=(X(i),y(i)) D i = ( X ( i ) , y ( i ) )
输出：μn μ n
初始化

Λ1=(X(1)TX(1))−1μ1=Λ1(X(1)Ty(1))i+=1 Λ 1 = ( X ( 1 ) T X ( 1 ) ) − 1 μ 1 = Λ 1 ( X ( 1 ) T y ( 1 ) ) i + = 1

while i<=n i <= n

Λi=(X(i)TX(i)+Λ−1i−1)−1μi=Λi(X(i)Ty(i)+Λ−1i−1μi−1)i+=1 Λ i = ( X ( i ) T X ( i ) + Λ i − 1 − 1 ) − 1 μ i = Λ i ( X ( i ) T y ( i ) + Λ i − 1 − 1 μ i − 1 ) i + = 1

代码：

def BayesLR(path):
    la=10
    mu=np.mat(np.zeros(3)).T
    gama=np.mat(np.eye(3)*la)
    for i in range(n):
        fileName = path + "%d.csv" % i
        x0,y0 = loadDataFromFile(fileName)#从文件中加载数据
        X, y = data2Mat(x0,y0)#将数据转换成np.mat的格式
        mu0 = mu
        gama0 = gama
        if i==1:
            gama = (X.T*X).I
            mu = gama*(X.T*y)

        else:
            gama = (X.T*X+gama0.I).I
            mu = gama*(X.T*y+gama0.I*mu0)

    return np.array(mu)