PRML:二元变量分布

版权声明:本文为博主原创文章,遵循 CC 4.0 by-sa 版权协议,转载请附上原文出处链接和本声明。
本文链接:https://blog.csdn.net/qilixuening/article/details/76472597

伯努利分布

考虑二元随机变量 x{0,1}(抛硬币,正面为 1,反面为 0),其概率分布由参数 μ 决定:

p(x=1)=μ

其中 (0μ1),并且有 p(x=0)=1μ。这就是伯努利分布(Bernoulli distribution),其概率分布可以写成:

Bern(x|μ)=μx(1μ)1x

均值和方差为:

E[x]var[x]=μ=μ(1μ)

伯努利分布的最大似然估计

考虑一组 x 的观测数据 D={x1,,xN},在独立同分布的假设下,其似然函数为

p(D|μ)=n=1Np(xn|μ)=n=1Nμxn(1μ)1xn

对数似然为

lnp(D|μ)=n=1Nlnp(xn|μ)=n=1N{xnlnμ+(1xn)ln(1μ)}

对数似然值只依赖于 Nn=1xn 的取值,而事实上 Nn=1xi 就是伯努利分布的一个充分统计量,它可以提供参数 μ 的全部信息。

μ 最大化对数似然,我们很容易得到

μML=1Nn=1Nxn

即最大似然估计值为样本均值(sample mean),若样本中 x=1 的数目为 m 则:

μML=mN

考虑抛三次硬币出现了三次正面的情况,此时 N=m=3,μML=1。在这种情况下,最大似然估计会得到每次都是正面的结果,这显然违背了我们的正常认知。事实上,这是一种过拟合的典型表现。

为了解决这个问题,我们可以考虑引入先验知识。

二项分布

给定数据总数 Nx=1 的总次数 m 满足一定的分布,这个分布叫做二项分布(binomial distribution)。

从伯努利分布的似然函数中可以看出它应该正比于 μm(1μ)Nm,事实上它可以写成:

Bin(m | N,μ)=(Nm)μm(1μ)Nm

其中

(Nm)N!(Nm)!m!

是组合数。

验证它是一个概率分布,二项式定理给出:

m=0N(Nm)μm(1μ)Nm=(μ+1μ)N=1

例如下面 N=20,μ=0.6 的一个分布直方图。

import matplotlib.pyplot as plt
import numpy as np
import scipy as sp
%matplotlib inline
from scipy.stats import binom

n = 20
mu = 0.6

X = binom.rvs(n, mu, size=10000)
fig, ax = plt.subplots()
ax.hist(X, bins=range(21), rwidth=0.7)
ax.set_xlabel("$m$", fontsize='x-large')
# ax.set_xlim(0, 21)
# ax.set_yticks(np.arange(0, 0.31, 0.1))
# ax.set_xticks(np.arange(0.5, 10.6, 1))
# ax.set_xticklabels(range(21))
ax.set_title(r'$N = 20, \mu=0.6$', fontsize='x-large')
plt.show()

output_7_0.png
计算均值时考虑下式对 μ 的导数,方差考虑对 μ 的二阶导数,其均值和方差分别为:

E[m]var[m]=Nμ=Nμ(1μ)

beta 分布

之前看到,当数据量较少时,最大似然的结果很可能会过拟合。为了减少这样的情况,从 Bayes 概率的观点出发,我们引入一个关于 μ 的先验分布 p(μ)

我们观察到似然函数是一系列 μx(1μ)1x 形式的乘积,如果我们选择一个正比于 μ 的某个幂次和 1μ 的某个幂次的先验分布,那么对 μ 来说,后验分布应当满足同样的形式。这样的性质叫做共轭性(conjugacy)。

在这里,我们引入的是 01 间的 beta 分布:

Beta(μ | a,b)=Γ(a+b)Γ(a)Γ(b)μa1(1μ)b1

其中:

Γ(x)0ux1eudu

满足如下性质:

Γ(x+1)Γ(1)=0uxeudu=[euux]0+x0ux1eudu=0+xΓ(x)=xΓ(x)=0eudu=[eu]0=1

验证它是一个概率分布:

由定义

Γ(a)Γ(b)=0xa1exdx0yb1eydy

t=y+x,dt=dy,则有:

Γ(a)Γ(b)=0xa1{x(tx)b1etdt}dx

交换积分次序,原来 t 是从 x 积分到 ,现在 x 是从 0 积分到 t

Γ(a)Γ(b)=0t0xa1(tx)b1etdxdt

x=tμ,dx=tdμ,则有

Γ(a)Γ(b)=0etta1tb1tdt10μa1(1μ)b1dμ=Γ(a+b)10μa1(1μ)b1dμ

于是我们有:

10Beta(μ | a,b)dμ=1

其均值和方差为

E[μ]var[μ]=aa+b=ab(a+b)2(a+b+1)

求均值:

从归一化我们知道:

10μa1(1μ)b1dμ=Γ(a)Γ(b)Γ(a+b)

从而,利用 Γ(x+1)=xΓ(x)

E[μ]=Γ(a+b)Γ(a)Γ(b)μa+11(1μ)b1dμ=Γ(a+b)Γ(a)Γ(b)Γ(a+1)Γ(b)Γ(a+b+1)=aa+b

类似可以得到:

E[μ2]=Γ(a+b)Γ(a)Γ(b)μa+21(1μ)b1dμ=Γ(a+b)Γ(a)Γ(b)Γ(a+2)Γ(b)Γ(a+b+2)=a(a+1)(a+b)(a+b+1)

从而可以计算出方差。

ab 叫做超参数,因为它们控制分布的参数 μ

from scipy.stats import beta

fig, axes = plt.subplots(2, 2,figsize=(10, 7))

axes = axes.flatten()

A = (0.1, 1, 2, 8)
B = (0.1, 1, 3, 4)

xx = np.linspace(0, 1, 100)

for a, b, ax in zip(A, B, axes):
    yy = beta.pdf(xx, a, b)
    ax.plot(xx, yy, 'r')
    ax.set_ylim(0, 3)

    ax.set_xticks([0, 0.5, 1])
    ax.set_xticklabels(["$0$", "$0.5$", "$1$"], fontsize="large")
    ax.set_yticks([0, 1, 2, 3])
    ax.set_yticklabels(["$0$", "$1$", "$2$", "$3$"], fontsize="large")
    ax.set_xlabel("$\mu$", fontsize="x-large")

    ax.text(0.1, 2.5, r"$a={}$".format(a), fontsize="x-large")
    ax.text(0.1, 2.2, r"$b={}$".format(b), fontsize="x-large")

output_11_0.png
有了先验分布,我们的后验分布为

p(μ | m,l,a,b)μm+a+1(1μ)l+b1

其中 l=Nm

由共轭性,我们知道后验分布还是一个 beta 分布,从而有

p(μ | m,l,a,b)Beta(μ | a+m,b+l)

至此,我们看到,如果我们观测到了一组 mx=1lx=0 的数据,那么超参由 a,b 变成 a+m,b+l。因此,超参 a,b 可以看成是 x=1x=0 的有效观测次数(注意:a,b 可以不是整数)。

更进一步,当有新的数据到来时,后验概率可以看成新的先验概率,因此这个过程可以序列化进行。下图表示观测到一个新的 x=1 数据时,先验到后验的变化。

xx = np.linspace(0, 1, 100)

fig, axes = plt.subplots(1, 3, figsize=(10, 2))

axes = axes.flatten()

axes[0].plot(xx, beta.pdf(xx, 2, 2), 'r')
axes[0].set_ylim(0, 2)
axes[0].text(0.1, 1.6, "prior", fontsize="x-large")
axes[0].set_xlabel("$\mu$", fontsize="x-large")
axes[0].set_xticks([0, 0.5, 1])
axes[0].set_yticks([0, 1, 2])

axes[1].plot(xx, xx)
axes[1].set_ylim(0, 2)
axes[1].text(0.1, 1.6, "likelihood", fontsize="x-large")
axes[1].set_xlabel("$\mu$", fontsize="x-large")
axes[1].set_xticks([0, 0.5, 1])
axes[1].set_yticks([0, 1, 2])

axes[2].plot(xx, beta.pdf(xx, 3, 2), 'r')
axes[2].set_ylim(0, 2)
axes[2].text(0.1, 1.6, "posterior", fontsize="x-large")
axes[2].set_xlabel("$\mu$", fontsize="x-large")
axes[2].set_xticks([0, 0.5, 1])
axes[2].set_yticks([0, 1, 2])

plt.show()

output_13_0.png
如果我们的目的是预测下一次实验的结果,那么我们有

p(x=1|D)=10p(x=1|μ)p(μ|D)dμ10μp(μ|D)=E[μ|D]

而我们知道后验概率的分布,所以可以计算出其均值:

p(x=1|D)=m+am+a+l+b=m+aN+a+b

m,l 时,有

p(x=1|D)=m+am+a+l+b=m+aN+a+bmN

即趋向于最大似然的解。当 N 有限时,我们的解总在先验均值和最大似然解之间。

另外我们观察到,随着 (a,b) 的增加,函数分布变得越来越尖,即方差越来越小。事实上,随着观测数据的增加,我们对于后验分布的不确定性也在不断减小。

另外考虑条件期望和方差的性质,我们有:

Eθ[θ]=ED[Eθ[θ|D]]

varθ[θ]=ED[varθ[θ|D]]+varD[Eθ[θ|D]]

从而我们知道,从平均意义上来说,后验分布的方差要比先验分布的方差要小。后验分布在数据分布上的均值等于先验分布的均值

展开阅读全文

c 内存 变量分布

12-25

rn学了这么久的c/c++。连他们的内存分布都还不清楚。想请教一下!rnrn我写了这样的一段代码来测试变量在虚拟内存的分布情况:rn[code=C/C++]rn int ia(0);rn int ib(0);rn rn static int sa(0);rn static int sb(0);rnrn int *paa = new int(0);rn int *pbb = new int(0);rnrn const int ca(0);rn const int cb(0);rn rn //输出变量的内存分布情况(虚拟内存)rn cout << "local variable:\n" rn << "ia's address = " << &ia << '\n'rn << "ib's address = " << &ib << "\n\nstatic variable:\n"rnrn << "sa's address = " << &sa << '\n'rn << "sb's address = " << &sb << "\n\ndynamic variable:\n"rnrn << "paa point address = " << paa << '\n'rn << "pbb point address = " << pbb << "\n\nconstant:\n"rnrn << "ca's address = " << &ca << '\n'rn << "cb's address = " << &cb << endl;rnrn delete paa;rn paa = 0;rn delete pbb;rn pbb = 0;rn rn[/code]rnrn我是在vs2003,window xp 测试的,结果如下:rnlocal variable:rnia's address = 0012FED4rnib's address = 0012FEC8rnrnstatic variable:rnsa's address = 004585E0rnsb's address = 004585E4rnrndynamic variable:rnpaa point address = 003711C0rnpbb point address = 003711F0rnrnconstant:rnca's address = 0012FEA4rncb's address = 0012FE98rn==============rnrn局部变量: 高地址 --> 低地址 rn静态变量: 低地址 --> 高地址rn动态变量: 不明确rn常量: 不明确rnrn请问我的这个断言正确吗?rn那每两个局部变量之间为什么会差距8个字节?这之间是什么内容?rnrnrn谢谢!!!!rn 论坛

没有更多推荐了,返回首页