二、概率分布 (Probability Distributions)
参数方法 (Parametric method):
- 预先假设数据服从一个特定的分布
非参数方法 (Nonparametric method):
- 数据的分布依赖于数据集的大小,且仅有控制模型复杂度的参数
共轭先验 (Conjugate prior):
- 使后验分布具有与先验分布相同的函数形式,进而简化贝叶斯分析
2.1 二元变量 (Binary Variables)
伯努利分布:
- 概率质量函数:
p ( x ) = { p if x = 1 1 − p if x = 0 p(x)=\left\{\begin{aligned}& p && \text{if } x = 1 \\& 1-p && \text{if } x = 0\end{aligned}\right.\ \ \ \ \ \ \ p(x)={p1−pif x=1if x=0 - 期望推导:
E [ X ] = p + 0 ( 1 − p ) = p E[X]=p+0(1-p)=p E[X]=p+0(1−p)=p - 方差推导:
D [ X ] = E [ X 2 ] − E [ X ] 2 = ( 1 2 p + 0 2 ( 1 − p ) ) − p 2 = p ( 1 − p ) D[X]=E[X^2]-E[X]^2=(1^2p+0^2(1-p))-p^2=p(1-p) D[X]=E[X2]−E[X]2=(12p+02(1−p))−p2=p(1−p)
二项分布:
2.1.1 贝塔分布 (The beta distribution)
伽马函数 (The gamma function):
- 定义: Γ ( x ) ≡ ∫ 0 ∞ u x − 1 e − u d u \Gamma(x) \equiv \int_{0}^{\infty} u^{x-1} e^{-u} \mathrm{~d} u Γ(x)≡∫0∞ux−1e−u du
- 性质:
- Γ ( x + 1 ) = x Γ ( x ) \Gamma(x+1)=x\Gamma(x) Γ(x+1)=xΓ(x)
-
∫
0
1
μ
a
−
1
(
1
−
μ
)
b
−
1
d
μ
=
Γ
(
a
)
Γ
(
b
)
Γ
(
a
+
b
)
\int_{0}^{1} \mu^{a-1}(1-\mu)^{b-1} \mathrm{~d} \mu=\frac{\Gamma(a) \Gamma(b)}{\Gamma(a+b)}
∫01μa−1(1−μ)b−1 dμ=Γ(a+b)Γ(a)Γ(b)
贝塔分布 (beta distribution):
-
定义:
Beta ( μ ∣ a , b ) = Γ ( a + b ) Γ ( a ) Γ ( b ) μ a − 1 ( 1 − μ ) b − 1 \operatorname{Beta}(\mu \mid a, b)=\frac{\Gamma(a+b)}{\Gamma(a) \Gamma(b)} \mu^{a-1}(1-\mu)^{b-1} Beta(μ∣a,b)=Γ(a)Γ(b)Γ(a+b)μa−1(1−μ)b−1 -
性质:
E [ μ ] = a a + b var [ μ ] = a b ( a + b ) 2 ( a + b + 1 ) \begin{aligned} \mathbb{E}[\mu] &=\frac{a}{a+b} \\ \operatorname{var}[\mu] &=\frac{a b}{(a+b)^{2}(a+b+1)} \end{aligned} E[μ]var[μ]=a+ba=(a+b)2(a+b+1)ab
将贝塔分布作为先验分布,则后验分布也将为贝塔分布,这个性质为共轭 (conjugacy):
p
(
μ
∣
a
,
b
)
=
Γ
(
a
+
b
)
Γ
(
a
)
Γ
(
b
)
μ
a
−
1
(
1
−
μ
)
b
−
1
p
(
μ
∣
m
,
l
,
a
,
b
)
=
Γ
(
m
+
a
+
l
+
b
)
Γ
(
m
+
a
)
Γ
(
l
+
b
)
μ
m
+
a
−
1
(
1
−
μ
)
l
+
b
−
1
\begin{aligned} & p(\mu \mid a, b)=\frac{\Gamma(a+b)}{\Gamma(a) \Gamma(b)} \mu^{a-1}(1-\mu)^{b-1} \\ & p(\mu \mid m, l, a, b)=\frac{\Gamma(m+a+l+b)}{\Gamma(m+a) \Gamma(l+b)} \mu^{m+a-1}(1-\mu)^{l+b-1} \end{aligned}
p(μ∣a,b)=Γ(a)Γ(b)Γ(a+b)μa−1(1−μ)b−1p(μ∣m,l,a,b)=Γ(m+a)Γ(l+b)Γ(m+a+l+b)μm+a−1(1−μ)l+b−1
共轭先验推导方式:
- 根据似然概率,确定先验分布的大致函数形式(与什么成正比)
- 再通过归一化确定系数
共轭先验 (Conjugate prior)的意义:
- sequential Bayesian inference:得到一个 observation 后,可以算出 posterior;由于选取的是共轭先验,所以 posterior 和原来的 prior 形式一样,可以把该 posterior 当作新的 prior,用于下一个 observation,如此迭代下去。对于 stream of data 的情况,这种方式可以实现 real-time learning
2.2 多元变量 (Multinomial Variables)
若变量的类别有
K
K
K 个,则用
K
K
K 维向量来表示该变量,例如:
x
=
(
0
,
0
,
1
,
0
,
0
,
0
)
T
\mathbf{x}=(0,0,1,0,0,0)^{\mathrm{T}}
x=(0,0,1,0,0,0)T
数据定义:
- μ i \mu_i μi:变量为第 i i i 类的概率
- m i m_i mi:原始数据中为 i i i 类的数据数量
2.2.1 狄利克雷分布 (The Dirichlet distribution)
定义:
Dir
(
μ
∣
α
)
=
Γ
(
α
0
)
Γ
(
α
1
)
⋯
Γ
(
α
K
)
∏
k
=
1
K
μ
k
α
k
−
1
,
α
0
=
∑
k
=
1
K
α
k
\operatorname{Dir}(\boldsymbol{\mu} \mid \boldsymbol{\alpha})=\frac{\Gamma\left(\alpha_{0}\right)}{\Gamma\left(\alpha_{1}\right) \cdots \Gamma\left(\alpha_{K}\right)} \prod_{k=1}^{K} \mu_{k}^{\alpha_{k}-1},\alpha_{0}=\sum_{k=1}^{K} \alpha_{k}
Dir(μ∣α)=Γ(α1)⋯Γ(αK)Γ(α0)k=1∏Kμkαk−1,α0=k=1∑Kαk
共轭先验:
p
(
D
∣
μ
)
=
∏
n
=
1
N
∏
k
=
1
K
μ
k
x
n
k
=
∏
k
=
1
K
μ
k
(
∑
n
x
n
k
)
=
∏
k
=
1
K
μ
k
m
k
p
(
μ
∣
α
)
=
Dir
(
μ
∣
α
)
p
(
μ
∣
D
,
α
)
∝
p
(
D
∣
μ
)
p
(
μ
∣
α
)
∝
∏
k
=
1
K
μ
k
α
k
+
m
k
−
1
p
(
μ
∣
D
,
α
)
=
Dir
(
μ
∣
α
+
m
)
=
Γ
(
α
0
+
N
)
Γ
(
α
1
+
m
1
)
⋯
Γ
(
α
K
+
m
K
)
∏
k
=
1
K
μ
k
α
k
+
m
k
−
1
\begin{aligned} &p(\mathcal{D} \mid \boldsymbol{\mu})=\prod_{n=1}^{N} \prod_{k=1}^{K} \mu_{k}^{x_{n k}}=\prod_{k=1}^{K} \mu_{k}^{\left(\sum_{n} x_{n k}\right)}=\prod_{k=1}^{K} \mu_{k}^{m_{k}}\\ &p(\boldsymbol{\mu} \mid \boldsymbol{\alpha})=\operatorname{Dir}(\boldsymbol{\mu} \mid \boldsymbol{\alpha})\\ &p(\boldsymbol{\mu} \mid \mathcal{D}, \boldsymbol{\alpha}) \propto p(\mathcal{D} \mid \boldsymbol{\mu}) p(\boldsymbol{\mu} \mid \boldsymbol{\alpha}) \propto \prod_{k=1}^{K} \mu_{k}^{\alpha_{k}+m_{k}-1}\\ &\begin{aligned} p(\boldsymbol{\mu} \mid \mathcal{D}, \boldsymbol{\alpha}) &=\operatorname{Dir}(\boldsymbol{\mu} \mid \boldsymbol{\alpha}+\mathbf{m}) \\ &=\frac{\Gamma\left(\alpha_{0}+N\right)}{\Gamma\left(\alpha_{1}+m_{1}\right) \cdots \Gamma\left(\alpha_{K}+m_{K}\right)} \prod_{k=1}^{K} \mu_{k}^{\alpha_{k}+m_{k}-1} \end{aligned} \end{aligned}
p(D∣μ)=n=1∏Nk=1∏Kμkxnk=k=1∏Kμk(∑nxnk)=k=1∏Kμkmkp(μ∣α)=Dir(μ∣α)p(μ∣D,α)∝p(D∣μ)p(μ∣α)∝k=1∏Kμkαk+mk−1p(μ∣D,α)=Dir(μ∣α+m)=Γ(α1+m1)⋯Γ(αK+mK)Γ(α0+N)k=1∏Kμkαk+mk−1
2.3 高斯分布 (The Gaussian Distribution)
一元高斯:
N
(
x
∣
μ
,
σ
2
)
=
1
(
2
π
σ
2
)
1
/
2
exp
{
−
1
2
σ
2
(
x
−
μ
)
2
}
\mathcal{N}\left(x \mid \mu, \sigma^{2}\right)=\frac{1}{\left(2 \pi \sigma^{2}\right)^{1 / 2}} \exp \left\{-\frac{1}{2 \sigma^{2}}(x-\mu)^{2}\right\}
N(x∣μ,σ2)=(2πσ2)1/21exp{−2σ21(x−μ)2}
多元高斯:
N
(
x
∣
μ
,
Σ
)
=
1
(
2
π
)
D
/
2
1
∣
Σ
∣
1
/
2
exp
{
−
1
2
(
x
−
μ
)
T
Σ
−
1
(
x
−
μ
)
}
\mathcal{N}(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Sigma})=\frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\mathbf{\Sigma}|^{1 / 2}} \exp \left\{-\frac{1}{2}(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \mathbf{\Sigma}^{-1}(\mathbf{x}-\boldsymbol{\mu})\right\}
N(x∣μ,Σ)=(2π)D/21∣Σ∣1/21exp{−21(x−μ)TΣ−1(x−μ)}
接下来推导多元高斯的一阶矩、二阶矩以及协方差,推导前需要做如下设定:
-
Σ
\boldsymbol{\Sigma}
Σ 为实对称矩阵,因此可以正交对角化,其中
λ
i
,
u
i
\lambda_{i},\mathbf{u}_{i}
λi,ui 分别为特征值和特征向量
Σ = U Λ U T = ∑ i = 1 D λ i u i u i T Σ − 1 = U Λ − 1 U T = ∑ i = 1 D 1 λ i u i u i T ∣ Σ ∣ 1 / 2 = ∏ j = 1 D λ j 1 / 2 \begin{aligned} &\boldsymbol{\Sigma}=\mathbf{U}\Lambda\mathbf{U}^{\mathrm{T}}=\sum_{i=1}^{D} \lambda_{i} \mathbf{u}_{i} \mathbf{u}_{i}^{\mathrm{T}}\\ &\boldsymbol{\Sigma}^{-1}=\mathbf{U}\Lambda^{-1}\mathbf{U}^{\mathrm{T}}=\sum_{i=1}^{D} \frac{1}{\lambda_{i}} \mathbf{u}_{i} \mathbf{u}_{i}^{\mathrm{T}}\\ &|\boldsymbol{\Sigma}|^{1 / 2}=\prod_{j=1}^{D} \lambda_{j}^{1 / 2} \end{aligned} Σ=UΛUT=i=1∑DλiuiuiTΣ−1=UΛ−1UT=i=1∑Dλi1uiuiT∣Σ∣1/2=j=1∏Dλj1/2 - 为了接下来推导方便,定义
y
i
=
u
i
T
(
x
−
μ
)
y_{i}=\mathbf{u}_{i}^{\mathrm{T}}(\mathbf{x}-\boldsymbol{\mu})
yi=uiT(x−μ),由此可以得到:
Δ 2 = ( x − μ ) T Σ − 1 ( x − μ ) = ∑ i = 1 D y i 2 λ i \Delta^{2}=(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Sigma}^{-1}(\mathbf{x}-\boldsymbol{\mu})=\sum_{i=1}^{D} \frac{y_{i}^{2}}{\lambda_{i}} Δ2=(x−μ)TΣ−1(x−μ)=i=1∑Dλiyi2 - 此处深入一步,推导
p
(
y
)
p(\mathbf{y})
p(y) 表达式
y = U ( x − μ ) J i j = ∂ x i ∂ y j = U j i ∣ J ∣ 2 = ∣ U T ∣ 2 = ∣ U T ∣ ∣ U ∣ = ∣ U T U ∣ = ∣ I ∣ = 1 p ( y ) = p ( x ) ∣ J ∣ = ∏ j = 1 D 1 ( 2 π λ j ) 1 / 2 exp { − y j 2 2 λ j } \begin{aligned} &\mathbf{y}=\mathbf{U}(\mathbf{x}-\boldsymbol{\mu})\\ &J_{i j}=\frac{\partial x_{i}}{\partial y_{j}}=U_{j i}\\ &|\mathbf{J}|^{2}=\left|\mathbf{U}^{\mathrm{T}}\right|^{2}=\left|\mathbf{U}^{\mathrm{T}}\right||\mathbf{U}|=\left|\mathbf{U}^{\mathrm{T}} \mathbf{U}\right|=|\mathbf{I}|=1\\ &p(\mathbf{y})=p(\mathbf{x})|\mathbf{J}|=\prod_{j=1}^{D} \frac{1}{\left(2 \pi \lambda_{j}\right)^{1 / 2}} \exp \left\{-\frac{y_{j}^{2}}{2 \lambda_{j}}\right\} \end{aligned} y=U(x−μ)Jij=∂yj∂xi=Uji∣J∣2=∣∣UT∣∣2=∣∣UT∣∣∣U∣=∣∣UTU∣∣=∣I∣=1p(y)=p(x)∣J∣=j=1∏D(2πλj)1/21exp{−2λjyj2}
首先推导一阶矩 E [ x ] = μ \mathbb{E}[x]=\boldsymbol{\mu} E[x]=μ:
- 令 z = x − μ \mathbf{z}=\mathbf{x}-\boldsymbol{\mu} z=x−μ
- 积分中
(
z
+
μ
)
(\mathbf{z}+\boldsymbol{\mu})
(z+μ) 的部分,第一部分根据奇函数消除,第二部分积分为
1
1
1
E [ x ] = 1 ( 2 π ) D / 2 1 ∣ Σ ∣ 1 / 2 ∫ exp { − 1 2 ( x − μ ) T Σ − 1 ( x − μ ) } x d x = 1 ( 2 π ) D / 2 1 ∣ Σ ∣ 1 / 2 ∫ exp { − 1 2 z T Σ − 1 z } ( z + μ ) d z = μ \begin{aligned} \mathbb{E}[\mathbf{x}] &=\frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\mathbf{\Sigma}|^{1 / 2}} \int \exp \left\{-\frac{1}{2}(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Sigma}^{-1}(\mathbf{x}-\boldsymbol{\mu})\right\} \mathbf{x} \mathrm{d} \mathbf{x} \\ &=\frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\mathbf{\Sigma}|^{1 / 2}} \int \exp \left\{-\frac{1}{2} \mathbf{z}^{\mathrm{T}} \boldsymbol{\Sigma}^{-1} \mathbf{z}\right\}(\mathbf{z}+\boldsymbol{\mu}) \mathrm{d} \mathbf{z}\\ &=\boldsymbol{\mu} \end{aligned} E[x]=(2π)D/21∣Σ∣1/21∫exp{−21(x−μ)TΣ−1(x−μ)}xdx=(2π)D/21∣Σ∣1/21∫exp{−21zTΣ−1z}(z+μ)dz=μ
接下来推导二阶矩 E [ x x T ] = μ μ T + Σ \mathbb{E}\left[\mathbf{x} \mathbf{x}^{\mathrm{T}}\right]=\boldsymbol{\mu} \boldsymbol{\mu}^{\mathrm{T}}+\mathbf{\Sigma} E[xxT]=μμT+Σ:
- 积分式中
(
z
+
μ
)
(
z
+
μ
)
T
(\mathbf{z}+\boldsymbol{\mu})(\mathbf{z}+\boldsymbol{\mu})^{\mathrm{T}}
(z+μ)(z+μ)T 只需考虑
z
z
T
\mathbf{z}\mathbf{z}^{\mathrm{T}}
zzT
E [ x x T ] = 1 ( 2 π ) D / 2 1 ∣ Σ ∣ 1 / 2 ∫ exp { − 1 2 ( x − μ ) T Σ − 1 ( x − μ ) } x x T d x = 1 ( 2 π ) D / 2 1 ∣ Σ ∣ 1 / 2 ∫ exp { − 1 2 z T Σ − 1 z } ( z + μ ) ( z + μ ) T d z \begin{array}{c} \mathbb{E}\left[\mathbf{x x}^{\mathrm{T}}\right]&=\frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\boldsymbol{\Sigma}|^{1 / 2}} \int \exp \left\{-\frac{1}{2}(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Sigma}^{-1}(\mathbf{x}-\boldsymbol{\mu})\right\} \mathbf{x} \mathbf{x}^{\mathrm{T}} \mathrm{d} \mathbf{x} \\ &=\frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\mathbf{\Sigma}|^{1 / 2}} \int \exp \left\{-\frac{1}{2} \mathbf{z}^{\mathrm{T}} \mathbf{\Sigma}^{-1} \mathbf{z}\right\}(\mathbf{z}+\boldsymbol{\mu})(\mathbf{z}+\boldsymbol{\mu})^{\mathrm{T}} \mathrm{d} \mathbf{z} \end{array} E[xxT]=(2π)D/21∣Σ∣1/21∫exp{−21(x−μ)TΣ−1(x−μ)}xxTdx=(2π)D/21∣Σ∣1/21∫exp{−21zTΣ−1z}(z+μ)(z+μ)Tdz - 求解
z
=
∑
j
=
1
D
y
j
u
j
\mathbf{z}=\sum_{j=1}^{D} y_{j} \mathbf{u}_{j}
z=∑j=1Dyjuj,进而推出
E
[
x
x
T
]
=
μ
μ
T
+
Σ
\mathbb{E}\left[\mathbf{x} \mathbf{x}^{\mathrm{T}}\right]=\boldsymbol{\mu} \boldsymbol{\mu}^{\mathrm{T}}+\mathbf{\Sigma}
E[xxT]=μμT+Σ,下式积分坐标转换采用雅可比行列式的方法,
U
U
U 为正交矩阵,行列式绝对值为
1
1
1
1 ( 2 π ) D / 2 1 ∣ Σ ∣ 1 / 2 ∫ exp { − 1 2 z T Σ − 1 z } z z T d z = 1 ( 2 π ) D / 2 1 ∣ Σ ∣ 1 / 2 ∑ i = 1 D ∑ j = 1 D u i u j T ∫ exp { − ∑ k = 1 D y k 2 2 λ k } y i y j d y = ∑ i = 1 D u i u i T λ i = Σ \begin{array}{l} \frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\boldsymbol{\Sigma}|^{1 / 2}} \int \exp \left\{-\frac{1}{2} \mathbf{z}^{\mathrm{T}} \boldsymbol{\Sigma}^{-1} \mathbf{z}\right\} \mathbf{z z}^{\mathrm{T}} \mathrm{d} \mathbf{z} \\ =\frac{1}{(2 \pi)^{D / 2}} \frac{1}{|\boldsymbol{\Sigma}|^{1 / 2}} \sum_{i=1}^{D} \sum_{j=1}^{D} \mathbf{u}_{i} \mathbf{u}_{j}^{\mathrm{T}} \int \exp \left\{-\sum_{k=1}^{D} \frac{y_{k}^{2}}{2 \lambda_{k}}\right\} y_{i} y_{j} \mathrm{~d} \mathbf{y} \\ =\sum_{i=1}^{D} \mathbf{u}_{i} \mathbf{u}_{i}^{\mathrm{T}} \lambda_{i}=\mathbf{\Sigma} \end{array} (2π)D/21∣Σ∣1/21∫exp{−21zTΣ−1z}zzTdz=(2π)D/21∣Σ∣1/21∑i=1D∑j=1DuiujT∫exp{−∑k=1D2λkyk2}yiyj dy=∑i=1DuiuiTλi=Σ
最后可以直接求出:
cov
[
x
]
=
E
[
(
x
−
E
[
x
]
)
(
x
−
E
[
x
]
)
T
]
=
Σ
\operatorname{cov}[\mathbf{x}]=\mathbb{E}\left[(\mathbf{x}-\mathbb{E}[\mathbf{x}])(\mathbf{x}-\mathbb{E}[\mathbf{x}])^{\mathrm{T}}\right]=\mathbf{\Sigma}
cov[x]=E[(x−E[x])(x−E[x])T]=Σ
不同协方差矩阵对应的图像:
多元高斯特点:
- 参数太多,计算复杂
- 单峰函数,难以建模多峰函数,但可以进行扩展来支持多峰
- introducing discrete latent variables (Gaussian Mixtures model)
- introducing continuous latent variables (Linear Dynamic System)
2.3.1 条件高斯分布 (Conditional Gaussian distributions)
本章重点:证明 p ( x a , x b ) p(\mathbf{x_a,x_b}) p(xa,xb) 为高斯分布时, p ( x a ∣ x b ) p(\mathbf{x_a|x_b}) p(xa∣xb) 也是高斯分布
前提准备:
-
x
∼
N
(
x
∣
μ
,
Σ
)
\mathbf{x}\sim\mathcal{N}(\mathbf{x}|\boldsymbol{\mu},\boldsymbol{\Sigma})
x∼N(x∣μ,Σ),将
x
\mathbf{x}
x 如下划分:
x = ( x a x b ) \mathbf{x}=\left(\begin{array}{l} \mathbf{x}_{a} \\ \mathbf{x}_{b} \end{array}\right) x=(xaxb) - 将
μ
,
Σ
\boldsymbol{\mu},\boldsymbol{\Sigma}
μ,Σ 进行相同分解:
μ = ( μ a μ b ) Σ = ( Σ a a Σ a b Σ b a Σ b b ) \begin{aligned} &\boldsymbol{\mu}=\left(\begin{array}{l} \boldsymbol{\mu}_{a} \\ \boldsymbol{\mu}_{b} \end{array}\right)\\ &\boldsymbol{\Sigma}=\left(\begin{array}{ll} \boldsymbol{\Sigma}_{a a} & \boldsymbol{\Sigma}_{a b} \\ \boldsymbol{\Sigma}_{b a} & \boldsymbol{\Sigma}_{b b} \end{array}\right) \end{aligned} μ=(μaμb)Σ=(ΣaaΣbaΣabΣbb)
- 定义精度矩阵 (precision matrix):
Λ ≡ Σ − 1 Λ = ( Λ a a Λ a b Λ b a Λ b b ) \begin{aligned} &\mathbf{\Lambda} \equiv \mathbf{\Sigma}^{-1}\\ &\mathbf{\Lambda}=\left(\begin{array}{ll} \boldsymbol{\Lambda}_{a a} & \boldsymbol{\Lambda}_{a b} \\ \boldsymbol{\Lambda}_{b a} & \boldsymbol{\Lambda}_{b b} \end{array}\right) \end{aligned} Λ≡Σ−1Λ=(ΛaaΛbaΛabΛbb)
证明 p ( x a ∣ x b ) p(\mathbf{x_a|x_b}) p(xa∣xb) 也是高斯分布:
- 下式可以表示为
x
a
\mathbf{x_a}
xa 的二次型,由此得证
− 1 2 ( x − μ ) T Σ − 1 ( x − μ ) = − 1 2 ( x a − μ a ) T Λ a a ( x a − μ a ) − 1 2 ( x a − μ a ) T Λ a b ( x b − μ b ) = − 1 2 ( x b − μ b ) T Λ b a ( x a − μ a ) − 1 2 ( x b − μ b ) T Λ b b ( x b − μ b ) \begin{aligned} -\frac{1}{2}(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Sigma}^{-1}(\mathbf{x}-\boldsymbol{\mu})&= -\frac{1}{2}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{a a}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right)-\frac{1}{2}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{a b}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right) \\ &=-\frac{1}{2}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{b a}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right)-\frac{1}{2}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{b b}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right) \end{aligned} −21(x−μ)TΣ−1(x−μ)=−21(xa−μa)TΛaa(xa−μa)−21(xa−μa)TΛab(xb−μb)=−21(xb−μb)TΛba(xa−μa)−21(xb−μb)TΛbb(xb−μb)
接下来需要求解
p
(
x
a
∣
x
b
)
p(\mathbf{x_a|x_b})
p(xa∣xb) 的均值与协方差矩阵,这里可以使用如下技巧 (completing the square),即求解二次项与一次项的系数:
−
1
2
(
x
−
μ
)
T
Σ
−
1
(
x
−
μ
)
=
−
1
2
x
T
Σ
−
1
x
+
x
T
Σ
−
1
μ
+
const
-\frac{1}{2}(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Sigma}^{-1}(\mathbf{x}-\boldsymbol{\mu})=-\frac{1}{2} \mathbf{x}^{\mathrm{T}} \boldsymbol{\Sigma}^{-1} \mathbf{x}+\mathbf{x}^{\mathrm{T}} \boldsymbol{\Sigma}^{-1} \boldsymbol{\mu}+\text { const }
−21(x−μ)TΣ−1(x−μ)=−21xTΣ−1x+xTΣ−1μ+ const
由此求出均值与协方差矩阵如下:
Σ
a
∣
b
=
Λ
a
a
−
1
μ
a
∣
b
=
Σ
a
∣
b
{
Λ
a
a
μ
a
−
Λ
a
b
(
x
b
−
μ
b
)
}
=
μ
a
−
Λ
a
a
−
1
Λ
a
b
(
x
b
−
μ
b
)
\begin{aligned} &\boldsymbol{\Sigma}_{a \mid b}=\mathbf{\Lambda}_{a a}^{-1}\\ &\begin{aligned} \boldsymbol{\mu}_{a \mid b} &=\boldsymbol{\Sigma}_{a \mid b}\left\{\boldsymbol{\Lambda}_{a a} \boldsymbol{\mu}_{a}-\boldsymbol{\Lambda}_{a b}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right)\right\} \\ &=\boldsymbol{\mu}_{a}-\boldsymbol{\Lambda}_{a a}^{-1} \boldsymbol{\Lambda}_{a b}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right) \end{aligned} \end{aligned}
Σa∣b=Λaa−1μa∣b=Σa∣b{Λaaμa−Λab(xb−μb)}=μa−Λaa−1Λab(xb−μb)
将精度矩阵用原始协方差矩阵代替,则需要应用下述分解方式,其中
M
−
1
\mathbf{M}^{-1}
M−1 为矩阵的 Schur complement:
(
A
B
C
D
)
−
1
=
(
M
−
M
B
D
−
1
−
D
−
1
C
M
D
−
1
+
D
−
1
C
M
B
D
−
1
)
M
=
(
A
−
B
D
−
1
C
)
−
1
\begin{aligned} \left(\begin{array}{cc} \mathbf{A} & \mathbf{B} \\ \mathbf{C} & \mathbf{D} \end{array}\right)^{-1}&=\left(\begin{array}{cc} \mathbf{M} & -\mathbf{M B D}^{-1} \\ -\mathbf{D}^{-1} \mathbf{C M} & \mathbf{D}^{-1}+\mathbf{D}^{-1} \mathbf{C M B D}^{-1} \end{array}\right)\\ \mathbf{M}&=\left(\mathbf{A}-\mathbf{B D}^{-1} \mathbf{C}\right)^{-1} \end{aligned}
(ACBD)−1M=(M−D−1CM−MBD−1D−1+D−1CMBD−1)=(A−BD−1C)−1
由此最终得到:
μ
a
∣
b
=
μ
a
+
Σ
a
b
Σ
b
b
−
1
(
x
b
−
μ
b
)
Σ
a
∣
b
=
Σ
a
a
−
Σ
a
b
Σ
b
b
−
1
Σ
b
a
\begin{aligned} \boldsymbol{\mu}_{a \mid b} &=\boldsymbol{\mu}_{a}+\boldsymbol{\Sigma}_{a b} \boldsymbol{\Sigma}_{b b}^{-1}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right) \\ \boldsymbol{\Sigma}_{a \mid b} &=\boldsymbol{\Sigma}_{a a}-\boldsymbol{\Sigma}_{a b} \boldsymbol{\Sigma}_{b b}^{-1} \boldsymbol{\Sigma}_{b a} \end{aligned}
μa∣bΣa∣b=μa+ΣabΣbb−1(xb−μb)=Σaa−ΣabΣbb−1Σba
观察上述均值与方差,可得到如下结论:
- p ( x a ∣ x b ) p(\mathbf{x_a}|\mathbf{x_b}) p(xa∣xb) 的均值是 x b \mathbf{x_b} xb 的线性函数
- p ( x a ∣ x b ) p(\mathbf{x_a}|\mathbf{x_b}) p(xa∣xb) 的方差与 x b \mathbf{x_b} xb 无关
- 由此 p ( x a ∣ x b ) p(\mathbf{x_a}|\mathbf{x_b}) p(xa∣xb) 是线性高斯 (linear-Guassian) 模型
2.3.2 边缘高斯分布 (Marginal Gaussian distributions)
求解 p ( x a ) = ∫ p ( x a , x b ) d x b p(\mathbf{x_a})=\int p(\mathbf{x_a},\mathbf{x_b})d\mathbf{x_b} p(xa)=∫p(xa,xb)dxb 的概率分布。
与条件高斯分布的过程一致,先只关注
p
(
x
a
,
x
b
)
p(\mathbf{x_a},\mathbf{x_b})
p(xa,xb) 指数项中包含
x
b
\mathbf{x_b}
xb 的部分:
−
1
2
x
b
T
Λ
b
b
x
b
+
x
b
T
m
=
−
1
2
(
x
b
−
Λ
b
b
−
1
m
)
T
Λ
b
b
(
x
b
−
Λ
b
b
−
1
m
)
+
1
2
m
T
Λ
b
b
−
1
m
m
=
Λ
b
b
μ
b
−
Λ
b
a
(
x
a
−
μ
a
)
\begin{aligned} -\frac{1}{2} \mathbf{x}_{b}^{\mathrm{T}} \boldsymbol{\Lambda}_{b b} \mathbf{x}_{b}+\mathbf{x}_{b}^{T} \mathbf{m}&=-\frac{1}{2}\left(\mathbf{x}_{b}-\mathbf{\Lambda}_{b b}^{-1} \mathbf{m}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{b b}\left(\mathbf{x}_{b}-\mathbf{\Lambda}_{b b}^{-1} \mathbf{m}\right)+\frac{1}{2} \mathbf{m}^{\mathrm{T}} \mathbf{\Lambda}_{b b}^{-1} \mathbf{m}\\ \mathbf{m}&=\mathbf{\Lambda}_{b b} \boldsymbol{\mu}_{b}-\boldsymbol{\Lambda}_{b a}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right) \end{aligned}
−21xbTΛbbxb+xbTmm=−21(xb−Λbb−1m)TΛbb(xb−Λbb−1m)+21mTΛbb−1m=Λbbμb−Λba(xa−μa)
由此对
x
b
x_b
xb 积分如下 (忽略常数):
∫
exp
{
−
1
2
(
x
b
−
Λ
b
b
−
1
m
)
T
Λ
b
b
(
x
b
−
Λ
b
b
−
1
m
)
}
d
x
b
=
const
\int \exp \left\{-\frac{1}{2}\left(\mathbf{x}_{b}-\mathbf{\Lambda}_{b b}^{-1} \mathbf{m}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{b b}\left(\mathbf{x}_{b}-\mathbf{\Lambda}_{b b}^{-1} \mathbf{m}\right)\right\} \mathrm{d} \mathbf{x}_{b}=\text{const}
∫exp{−21(xb−Λbb−1m)TΛbb(xb−Λbb−1m)}dxb=const
对
x
b
x_b
xb 积分结束后,再次关注
p
(
x
a
,
x
b
)
p(\mathbf{x_a},\mathbf{x_b})
p(xa,xb) 指数项中包含
x
a
\mathbf{x_a}
xa 的部分:
1
2
[
Λ
b
b
μ
b
−
Λ
b
a
(
x
a
−
μ
a
)
]
T
Λ
b
b
−
1
[
Λ
b
b
μ
b
−
Λ
b
a
(
x
a
−
μ
a
)
]
−
1
2
x
a
T
Λ
a
a
x
a
+
x
a
T
(
Λ
a
a
μ
a
+
Λ
a
b
μ
b
)
+
c
o
n
s
t
=
−
1
2
x
a
T
(
Λ
a
a
−
Λ
a
b
Λ
b
b
−
1
Λ
b
a
)
x
a
+
x
a
T
(
Λ
a
a
−
Λ
a
b
Λ
b
b
−
1
Λ
b
a
)
−
1
μ
a
+
c
o
n
s
t
\begin{aligned} \frac{1}{2}\left[\boldsymbol{\Lambda}_{b b}\right.&\left.\boldsymbol{\mu}_{b}-\boldsymbol{\Lambda}_{b a}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right)\right]^{\mathrm{T}} \mathbf{\Lambda}_{b b}^{-1}\left[\boldsymbol{\Lambda}_{b b} \boldsymbol{\mu}_{b}-\boldsymbol{\Lambda}_{b a}\left(\mathbf{x}_{a}-\boldsymbol{\mu}_{a}\right)\right] \\ &-\frac{1}{2} \mathbf{x}_{a}^{\mathrm{T}} \boldsymbol{\Lambda}_{a a} \mathbf{x}_{a}+\mathbf{x}_{a}^{\mathrm{T}}\left(\boldsymbol{\Lambda}_{a a} \boldsymbol{\mu}_{a}+\mathbf{\Lambda}_{a b} \boldsymbol{\mu}_{b}\right)+\mathrm{const} \\ =&-\frac{1}{2} \mathbf{x}_{a}^{\mathrm{T}}\left(\boldsymbol{\Lambda}_{a a}-\boldsymbol{\Lambda}_{a b} \mathbf{\Lambda}_{b b}^{-1} \boldsymbol{\Lambda}_{b a}\right) \mathbf{x}_{a} \\ &+\mathbf{x}_{a}^{\mathrm{T}}\left(\boldsymbol{\Lambda}_{a a}-\boldsymbol{\Lambda}_{a b} \mathbf{\Lambda}_{b b}^{-1} \boldsymbol{\Lambda}_{b a}\right)^{-1} \boldsymbol{\mu}_{a}+\mathrm{const} \end{aligned}
21[Λbb=μb−Λba(xa−μa)]TΛbb−1[Λbbμb−Λba(xa−μa)]−21xaTΛaaxa+xaT(Λaaμa+Λabμb)+const−21xaT(Λaa−ΛabΛbb−1Λba)xa+xaT(Λaa−ΛabΛbb−1Λba)−1μa+const
采用与条件高斯推导一样的方式,求出均值与协方差矩阵如下:
Σ
a
=
(
Λ
a
a
−
Λ
a
b
Λ
b
b
−
1
Λ
b
a
)
−
1
μ
a
=
Σ
a
(
Λ
a
a
−
Λ
a
b
Λ
b
b
−
1
Λ
b
a
)
μ
a
=
μ
a
\begin{aligned} &\boldsymbol{\Sigma}_{a}=\left(\boldsymbol{\Lambda}_{a a}-\boldsymbol{\Lambda}_{a b} \mathbf{\Lambda}_{b b}^{-1} \boldsymbol{\Lambda}_{b a}\right)^{-1} \\ &\boldsymbol{\mu}_{a}=\boldsymbol{\Sigma}_{a}\left(\boldsymbol{\Lambda}_{a a}-\boldsymbol{\Lambda}_{a b} \mathbf{\Lambda}_{b b}^{-1} \boldsymbol{\Lambda}_{b a}\right) \boldsymbol{\mu}_{a}=\boldsymbol{\mu}_{a} \end{aligned}
Σa=(Λaa−ΛabΛbb−1Λba)−1μa=Σa(Λaa−ΛabΛbb−1Λba)μa=μa
利用矩阵的 Schur complement,将精度矩阵用原始协方差矩阵代替:
E
[
x
a
]
=
μ
a
cov
[
x
a
]
=
Σ
a
a
\begin{aligned} \mathbb{E}\left[\mathbf{x}_{a}\right] &=\boldsymbol{\mu}_{a} \\ \operatorname{cov}\left[\mathbf{x}_{a}\right] &=\mathbf{\Sigma}_{a a} \end{aligned}
E[xa]cov[xa]=μa=Σaa
总结一下:
- 条件高斯分布 (Conditional Gaussian distribution):
p ( x a ∣ x b ) = N ( x ∣ μ a ∣ b , Λ a a − 1 ) μ a ∣ b = μ a − Λ a a − 1 Λ a b ( x b − μ b ) \begin{aligned} p\left(\mathbf{x}_{a} \mid \mathbf{x}_{b}\right) &=\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}_{a \mid b}, \boldsymbol{\Lambda}_{a a}^{-1}\right) \\ \boldsymbol{\mu}_{a \mid b} &=\boldsymbol{\mu}_{a}-\boldsymbol{\Lambda}_{a a}^{-1} \boldsymbol{\Lambda}_{a b}\left(\mathbf{x}_{b}-\boldsymbol{\mu}_{b}\right) \end{aligned} p(xa∣xb)μa∣b=N(x∣μa∣b,Λaa−1)=μa−Λaa−1Λab(xb−μb) - 边缘高斯分布 (Marginal Gaussian distribution):
p ( x a ) = N ( x a ∣ μ a , Σ a a ) p\left(\mathbf{x}_{a}\right)=\mathcal{N}\left(\mathbf{x}_{a} \mid \boldsymbol{\mu}_{a}, \boldsymbol{\Sigma}_{a a}\right) p(xa)=N(xa∣μa,Σaa) - 举例图示:
2.3.3 高斯变量的贝叶斯定理 (Bayes’ theorem for Gaussian variables)
本节主要目的是根据给定的高斯分布 p ( x ) , p ( y ∣ x ) p(\mathbf{x}),p(\mathbf{y|x}) p(x),p(y∣x),求出 p ( x , y ) , p ( y ) , p ( x ∣ y ) p(\mathbf{x},\mathbf{y}),p(\mathbf{y}),p(\mathbf{x|y}) p(x,y),p(y),p(x∣y)。
首先给出已知先验与似然:
p
(
x
)
=
N
(
x
∣
μ
,
Λ
−
1
)
p
(
y
∣
x
)
=
N
(
y
∣
A
x
+
b
,
L
−
1
)
\begin{aligned} p(\mathbf{x}) &=\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Lambda}^{-1}\right) \\ p(\mathbf{y} \mid \mathbf{x}) &=\mathcal{N}\left(\mathbf{y} \mid \mathbf{A} \mathbf{x}+\mathbf{b}, \mathbf{L}^{-1}\right) \end{aligned}
p(x)p(y∣x)=N(x∣μ,Λ−1)=N(y∣Ax+b,L−1)
定义
x
,
y
\mathbf{x},\mathbf{y}
x,y 的联合
z
=
(
x
y
)
T
\mathbf{z}=(\mathbf{x}\ \mathbf{y})^{\mathrm{T}}
z=(x y)T,进行
p
(
x
,
y
)
p(\mathbf{x},\mathbf{y})
p(x,y) 的求取,依旧只关注指数内的一次项与二次项:
ln
p
(
z
)
=
ln
p
(
x
)
+
ln
p
(
y
∣
x
)
=
−
1
2
(
x
−
μ
)
T
Λ
(
x
−
μ
)
−
1
2
(
y
−
A
x
−
b
)
T
L
(
y
−
A
x
−
b
)
+
const
\begin{aligned} \ln p(\mathbf{z})=& \ln p(\mathbf{x})+\ln p(\mathbf{y} \mid \mathbf{x}) \\ =&-\frac{1}{2}(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Lambda}(\mathbf{x}-\boldsymbol{\mu}) \\ &-\frac{1}{2}(\mathbf{y}-\mathbf{A} \mathbf{x}-\mathbf{b})^{\mathrm{T}} \mathbf{L}(\mathbf{y}-\mathbf{A} \mathbf{x}-\mathbf{b})+\text { const } \end{aligned}
lnp(z)==lnp(x)+lnp(y∣x)−21(x−μ)TΛ(x−μ)−21(y−Ax−b)TL(y−Ax−b)+ const
提取指数中的二次项,得到:
−
1
2
x
T
(
Λ
+
A
T
L
A
)
x
−
1
2
y
T
L
y
+
1
2
y
T
L
A
x
+
1
2
x
T
A
T
L
y
=
−
1
2
(
x
y
)
T
(
Λ
+
A
T
L
A
−
A
T
L
−
L
A
L
)
(
x
y
)
=
−
1
2
z
T
R
z
\begin{array}{c} -\frac{1}{2} \mathbf{x}^{\mathrm{T}}\left(\boldsymbol{\Lambda}+\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{A}\right) \mathbf{x}-\frac{1}{2} \mathbf{y}^{\mathrm{T}} \mathbf{L} \mathbf{y}+\frac{1}{2} \mathbf{y}^{\mathrm{T}} \mathbf{L} \mathbf{A} \mathbf{x}+\frac{1}{2} \mathbf{x}^{\mathrm{T}} \mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{y} \\ =-\frac{1}{2}\left(\begin{array}{l} \mathbf{x} \\ \mathbf{y} \end{array}\right)^{\mathrm{T}}\left(\begin{array}{cc} \mathbf{\Lambda}+\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{A} & -\mathbf{A}^{\mathrm{T}} \mathbf{L} \\ -\mathbf{L} \mathbf{A} & \mathbf{L} \end{array}\right)\left(\begin{array}{l} \mathbf{x} \\ \mathbf{y} \end{array}\right)=-\frac{1}{2} \mathbf{z}^{\mathrm{T}} \mathbf{R} \mathbf{z} \end{array}
−21xT(Λ+ATLA)x−21yTLy+21yTLAx+21xTATLy=−21(xy)T(Λ+ATLA−LA−ATLL)(xy)=−21zTRz
提取指数中的一次项,得到:
x
T
Λ
μ
−
x
T
A
T
L
b
+
y
T
L
b
=
(
x
y
)
T
(
Λ
μ
−
A
T
L
b
L
b
)
\mathbf{x}^{\mathrm{T}} \boldsymbol{\Lambda} \boldsymbol{\mu}-\mathbf{x}^{\mathrm{T}} \mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{b}+\mathbf{y}^{\mathrm{T}} \mathbf{L} \mathbf{b}=\left(\begin{array}{l} \mathbf{x} \\ \mathbf{y} \end{array}\right)^{\mathrm{T}}\left(\begin{array}{c} \mathbf{\Lambda} \boldsymbol{\mu}-\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{b} \\ \mathbf{L b} \end{array}\right)
xTΛμ−xTATLb+yTLb=(xy)T(Λμ−ATLbLb)
因此可以得到
p
(
z
)
p(\mathbf{z})
p(z) 的分布:
R
=
(
Λ
+
A
T
L
A
−
A
T
L
−
L
A
L
)
cov
[
z
]
=
R
−
1
=
(
Λ
−
1
Λ
−
1
A
T
A
Λ
−
1
L
−
1
+
A
Λ
−
1
A
T
)
E
[
z
]
=
R
−
1
(
Λ
μ
−
A
T
L
b
L
b
)
=
(
μ
A
μ
+
b
)
\begin{aligned} &\mathbf{R}=\left(\begin{array}{cc} \mathbf{\Lambda}+\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{A} & -\mathbf{A}^{\mathrm{T}} \mathbf{L} \\ -\mathbf{L} \mathbf{A} & \mathbf{L} \end{array}\right)\\ &\operatorname{cov}[\mathbf{z}]=\mathbf{R}^{-1}=\left(\begin{array}{cc} \mathbf{\Lambda}^{-1} & \mathbf{\Lambda}^{-1} \mathbf{A}^{\mathrm{T}} \\ \mathbf{A} \mathbf{\Lambda}^{-1} & \mathbf{L}^{-1}+\mathbf{A} \mathbf{\Lambda}^{-1} \mathbf{A}^{\mathrm{T}} \end{array}\right)\\ &\mathbb{E}[\mathbf{z}]=\mathbf{R}^{-1}\left(\begin{array}{c} \mathbf{\Lambda} \boldsymbol{\mu}-\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{b} \\ \mathbf{L b} \end{array}\right)=\left(\begin{array}{c} \boldsymbol{\mu} \\ \mathbf{A} \boldsymbol{\mu}+\mathbf{b} \end{array}\right) \end{aligned}
R=(Λ+ATLA−LA−ATLL)cov[z]=R−1=(Λ−1AΛ−1Λ−1ATL−1+AΛ−1AT)E[z]=R−1(Λμ−ATLbLb)=(μAμ+b)
得到
p
(
z
)
p(\mathbf{z})
p(z) 后,可以根据之前对于条件高斯与边缘高斯的探究,得到后验与证据因子分布,总结如下:
已知
=
{
[
先验
]
p
(
x
)
=
N
(
x
∣
μ
,
Λ
−
1
)
[
似然
]
p
(
y
∣
x
)
=
N
(
y
∣
A
x
+
b
,
L
−
1
)
求解
=
{
[
证据
]
p
(
y
)
=
N
(
y
∣
A
μ
+
b
,
L
−
1
+
A
Λ
−
1
A
T
)
[
后验
]
p
(
x
∣
y
)
=
N
(
x
∣
Σ
{
A
T
L
(
y
−
b
)
+
Λ
μ
}
,
Σ
)
,
Σ
=
(
Λ
+
A
T
L
A
)
−
1
\begin{aligned} &\text{已知}=\left\{ \begin{aligned} & [\text{先验}]\ \ p(\mathbf{x}) =\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Lambda}^{-1}\right) \\ & [\text{似然}]\ \ p(\mathbf{y} \mid \mathbf{x}) =\mathcal{N}\left(\mathbf{y} \mid \mathbf{A} \mathbf{x}+\mathbf{b}, \mathbf{L}^{-1}\right) \\ \end{aligned} \right.\\ &\text{求解}=\left\{ \begin{aligned} &[\text{证据}] \ \ p(\mathbf{y})= \mathcal{N}\left(\mathbf{y} \mid \mathbf{A} \boldsymbol{\mu}+\mathbf{b}, \mathbf{L}^{-1}+\mathbf{A} \mathbf{\Lambda}^{-1} \mathbf{A}^{\mathrm{T}}\right) \\ & [\text{后验}] \ \ p(\mathbf{x} \mid \mathbf{y})= \mathcal{N}\left(\mathbf{x} \mid \mathbf{\Sigma}\left\{\mathbf{A}^{\mathrm{T}} \mathbf{L}(\mathbf{y}-\mathbf{b})+\mathbf{\Lambda} \boldsymbol{\mu}\right\}, \mathbf{\Sigma}\right), \boldsymbol{\Sigma} =\left(\mathbf{\Lambda}+\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{A}\right)^{-1} \\ \end{aligned} \right. \end{aligned}
已知={[先验] p(x)=N(x∣μ,Λ−1)[似然] p(y∣x)=N(y∣Ax+b,L−1)求解={[证据] p(y)=N(y∣Aμ+b,L−1+AΛ−1AT)[后验] p(x∣y)=N(x∣Σ{ATL(y−b)+Λμ},Σ),Σ=(Λ+ATLA)−1
2.3.4 高斯极大似然估计 (Maximum likelihood for the Gaussian)
似然函数:
ln
p
(
X
∣
μ
,
Σ
)
=
−
N
D
2
ln
(
2
π
)
−
N
2
ln
∣
Σ
∣
−
1
2
∑
n
=
1
N
(
x
n
−
μ
)
T
Σ
−
1
(
x
n
−
μ
)
\ln p(\mathbf{X} \mid \boldsymbol{\mu}, \mathbf{\Sigma})=-\frac{N D}{2} \ln (2 \pi)-\frac{N}{2} \ln |\boldsymbol{\Sigma}|-\frac{1}{2} \sum_{n=1}^{N}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)^{\mathrm{T}} \boldsymbol{\Sigma}^{-1}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)
lnp(X∣μ,Σ)=−2NDln(2π)−2Nln∣Σ∣−21n=1∑N(xn−μ)TΣ−1(xn−μ)
极大似然估计:
∂
∂
μ
ln
p
(
X
∣
μ
,
Σ
)
=
∑
n
=
1
N
Σ
−
1
(
x
n
−
μ
)
=
0
⇒
μ
M
L
=
1
N
∑
n
=
1
N
x
n
∂
∂
Σ
ln
p
(
X
∣
μ
,
Σ
)
=
0
⇒
Σ
M
L
=
1
N
∑
n
=
1
N
(
x
n
−
μ
M
L
)
(
x
n
−
μ
M
L
)
T
\begin{aligned} &\frac{\partial}{\partial \boldsymbol{\mu}} \ln p(\mathbf{X} \mid \boldsymbol{\mu}, \boldsymbol{\Sigma})=\sum_{n=1}^{N} \boldsymbol{\Sigma}^{-1}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)=0\Rightarrow \boldsymbol{\mu}_{\mathrm{ML}}=\frac{1}{N} \sum_{n=1}^{N} \mathbf{x}_{n}\\ &\frac{\partial}{\partial \boldsymbol{\Sigma}} \ln p(\mathbf{X} \mid \boldsymbol{\mu}, \boldsymbol{\Sigma})=0\Rightarrow \boldsymbol{\Sigma}_{\mathrm{ML}}=\frac{1}{N} \sum_{n=1}^{N}\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{\mathrm{ML}}\right)\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{\mathrm{ML}}\right)^{\mathrm{T}} \end{aligned}
∂μ∂lnp(X∣μ,Σ)=n=1∑NΣ−1(xn−μ)=0⇒μML=N1n=1∑Nxn∂Σ∂lnp(X∣μ,Σ)=0⇒ΣML=N1n=1∑N(xn−μML)(xn−μML)T
检验极大似然估计的期望:
E
[
μ
M
L
]
=
μ
E
[
Σ
M
L
]
=
N
−
1
N
Σ
\begin{aligned} \mathbb{E}\left[\boldsymbol{\mu}_{\mathrm{ML}}\right] &=\boldsymbol{\mu} \\ \mathbb{E}\left[\boldsymbol{\Sigma}_{\mathrm{ML}}\right] &=\frac{N-1}{N} \boldsymbol{\Sigma} \end{aligned}
E[μML]E[ΣML]=μ=NN−1Σ
极大似然估计的协方差矩阵与真实情况有偏差 (bias),因此进行修正:
Σ
~
=
1
N
−
1
∑
n
=
1
N
(
x
n
−
μ
M
L
)
(
x
n
−
μ
M
L
)
T
\widetilde{\mathbf{\Sigma}}=\frac{1}{N-1} \sum_{n=1}^{N}\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{\mathrm{ML}}\right)\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{\mathrm{ML}}\right)^{\mathrm{T}}
Σ
=N−11n=1∑N(xn−μML)(xn−μML)T
2.3.5 序列估计 (Sequential estimation)
序列估计:
- 基于在线学习,实时更新模型
μ M L ( N ) = 1 N ∑ n = 1 N x n = 1 N x N + 1 N ∑ n = 1 N − 1 x n = 1 N x N + N − 1 N μ M L ( N − 1 ) = μ M L ( N − 1 ) + 1 N ( x N − μ M L ( N − 1 ) ) \begin{aligned} \boldsymbol{\mu}_{\mathrm{ML}}^{(N)} &=\frac{1}{N} \sum_{n=1}^{N} \mathbf{x}_{n} \\ &=\frac{1}{N} \mathbf{x}_{N}+\frac{1}{N} \sum_{n=1}^{N-1} \mathbf{x}_{n} \\ &=\frac{1}{N} \mathbf{x}_{N}+\frac{N-1}{N} \boldsymbol{\mu}_{\mathrm{ML}}^{(N-1)} \\ &=\boldsymbol{\mu}_{\mathrm{ML}}^{(N-1)}+\frac{1}{N}\left(\mathbf{x}_{N}-\boldsymbol{\mu}_{\mathrm{ML}}^{(N-1)}\right) \end{aligned} μML(N)=N1n=1∑Nxn=N1xN+N1n=1∑N−1xn=N1xN+NN−1μML(N−1)=μML(N−1)+N1(xN−μML(N−1))
Robbins-Monro algorithm (更通用的方法):
-
对于下述式子,寻找到根 (root) 使 f ( θ ∗ ) = 0 f(\theta^*)=0 f(θ∗)=0
f ( θ ) ≡ E [ z ∣ θ ] = ∫ z p ( z ∣ θ ) d z f(\theta) \equiv \mathbb{E}[z \mid \theta]=\int z p(z \mid \theta) \mathrm{d} z f(θ)≡E[z∣θ]=∫zp(z∣θ)dz -
假设 E [ ( z − f ) 2 ∣ θ ] < ∞ \mathbb{E}[(z-f)^2|\theta]<\infty E[(z−f)2∣θ]<∞,且 f ( θ ) > 0 f(\theta)>0 f(θ)>0 当 θ > θ ∗ \theta>\theta^* θ>θ∗; f ( θ ) < 0 f(\theta)<0 f(θ)<0 当 θ < θ ∗ \theta<\theta^* θ<θ∗
-
因此可以使用如下序列迭代求解 θ ∗ \theta^* θ∗:
θ ( N ) = θ ( N − 1 ) + α N − 1 z ( θ N − 1 ) \theta^{(N)}=\theta^{(N-1)}+\alpha_{N-1}z(\theta^{N-1}) θ(N)=θ(N−1)+αN−1z(θN−1) -
其中 { a N } \{a_N\} {aN} 需要满足如下条件:
[保证能收敛] lim N → ∞ a N = 0 [能收敛到根] ∑ N = 1 ∞ a N = ∞ [累计噪声方差有限,不会影响收敛] ∑ N = 1 ∞ a N 2 < ∞ \begin{aligned} \text{[保证能收敛]} \ \ & \lim_{N\rightarrow \infty} a_N=0 \\ \text{[能收敛到根]} \ \ & \sum_{N=1}^{\infty}a_N=\infty \\ \text{[累计噪声方差有限,不会影响收敛]} \ \ & \sum_{N=1}^{\infty}a_N^2<\infty \\ \end{aligned} [保证能收敛] [能收敛到根] [累计噪声方差有限,不会影响收敛] N→∞limaN=0N=1∑∞aN=∞N=1∑∞aN2<∞
2.3.6 高斯的贝叶斯推断 (Bayesian inference for the Gaussian)
设随机变量 x ∼ N ( μ , σ 2 ) x\sim\mathcal{N}(\mu,\sigma^2) x∼N(μ,σ2),进行高斯的贝叶斯推断。
第一种情况,已知 σ 2 \sigma^2 σ2,求 μ \mu μ 的概率分布 p ( μ ∣ X ) p(\mu \mid \mathbf{X}) p(μ∣X)。
使用共轭先验的常规做法:
- 写出似然函数
p ( X ∣ μ ) = ∏ n = 1 N p ( x n ∣ μ ) = 1 ( 2 π σ 2 ) N / 2 exp { − 1 2 σ 2 ∑ n = 1 N ( x n − μ ) 2 } p(\mathbf{X} \mid \mu)=\prod_{n=1}^{N} p\left(x_{n} \mid \mu\right)=\frac{1}{\left(2 \pi \sigma^{2}\right)^{N / 2}} \exp \left\{-\frac{1}{2 \sigma^{2}} \sum_{n=1}^{N}\left(x_{n}-\mu\right)^{2}\right\} p(X∣μ)=n=1∏Np(xn∣μ)=(2πσ2)N/21exp{−2σ21n=1∑N(xn−μ)2} - 仿照似然函数的形式,给出共轭先验的形式,即高斯分布
p ( μ ) = N ( μ ∣ μ 0 , σ 0 2 ) p(\mu)=\mathcal{N}\left(\mu \mid \mu_{0}, \sigma_{0}^{2}\right) p(μ)=N(μ∣μ0,σ02) - 给出后验的形式
p ( μ ∣ X ) ∝ p ( X ∣ μ ) p ( μ ) p ( μ ∣ X ) = N ( μ ∣ μ N , σ N 2 ) μ N = σ 2 N σ 0 2 + σ 2 μ 0 + N σ 0 2 N σ 0 2 + σ 2 μ M L , μ M L = 1 N ∑ n = 1 N x n 1 σ N 2 = 1 σ 0 2 + N σ 2 \begin{aligned} & p(\mu \mid \mathbf{X}) \propto p(\mathbf{X} \mid \mu) p(\mu)\\ & p(\mu \mid \mathbf{X})=\mathcal{N}\left(\mu \mid \mu_{N}, \sigma_{N}^{2}\right)\\ & \mu_{N}=\frac{\sigma^{2}}{N \sigma_{0}^{2}+\sigma^{2}} \mu_{0}+\frac{N \sigma_{0}^{2}}{N \sigma_{0}^{2}+\sigma^{2}} \mu_{\mathrm{ML}},\mu_{\mathrm{ML}}=\frac{1}{N} \sum_{n=1}^{N} x_{n}\\ & \frac{1}{\sigma_{N}^{2}}=\frac{1}{\sigma_{0}^{2}}+\frac{N}{\sigma^{2}} \end{aligned} p(μ∣X)∝p(X∣μ)p(μ)p(μ∣X)=N(μ∣μN,σN2)μN=Nσ02+σ2σ2μ0+Nσ02+σ2Nσ02μML,μML=N1n=1∑NxnσN21=σ021+σ2N
此处有一些性质值得我们进行观察,当观测数据
N
=
0
N=0
N=0 时,后验的
μ
N
\mu_N
μN 等于先验的
μ
0
\mu_0
μ0,后验的
σ
N
\sigma_N
σN 也等于先验的
σ
0
\sigma_0
σ0;而当观察数据
N
→
∞
N \rightarrow \infty
N→∞ 时,后验的
μ
N
\mu_N
μN 将等于极大似然估计得到的
μ
ML
\mu_{\text{ML}}
μML,后验的
σ
N
\sigma_N
σN 将等于 0,即后验分布变为
p
(
μ
=
μ
ML
∣
X
)
=
1
p(\mu=\mu_{\text{ML}} \mid \mathbf{X})=1
p(μ=μML∣X)=1。具体情况如下图所示,其中
μ
0
=
0
,
μ
ML
=
0.8
\mu_0=0,\mu_{\text{ML}}=0.8
μ0=0,μML=0.8。
另外,由于似然函数假设数据点之间独立,因此上述贝叶斯推断可以很容易地扩展成在线学习:
This sequential view of Bayesian inference is very general and applies to any problem in which the observed data are assumed to be independent and identically distributed.
接下来是第二种情况,已知 μ \mu μ,求 σ 2 \sigma^2 σ2 的分布,此处令 λ ≡ 1 / σ 2 \lambda \equiv 1 / \sigma^{2} λ≡1/σ2,其中 λ \lambda λ 为精度 (precision)。
依然是相同的套路:
-
写出似然函数
p ( X ∣ λ ) = ∏ n = 1 N N ( x n ∣ μ , λ − 1 ) ∝ λ N / 2 exp { − λ 2 ∑ n = 1 N ( x n − μ ) 2 } p(\mathbf{X} \mid \lambda)=\prod_{n=1}^{N} \mathcal{N}\left(x_{n} \mid \mu, \lambda^{-1}\right) \propto \lambda^{N / 2} \exp \left\{-\frac{\lambda}{2} \sum_{n=1}^{N}\left(x_{n}-\mu\right)^{2}\right\} p(X∣λ)=n=1∏NN(xn∣μ,λ−1)∝λN/2exp{−2λn=1∑N(xn−μ)2} -
因为似然函数中 exp 里只有 λ \lambda λ 的一次项,因此用 gamma 分布来表示先验,其中 gamma 分布的期望和方差可以通过积分时将 b λ b\lambda bλ 换元为 t t t 得到;另外若 a > 0 a>0 a>0,则函数积分值有限,若 a ≥ 1 a\geq 1 a≥1,则分布有限,且概率存在极大值
Gam ( λ ∣ a , b ) = 1 Γ ( a ) b a λ a − 1 exp ( − b λ ) E [ λ ] = a b var [ λ ] = a b 2 \begin{aligned} & \operatorname{Gam}(\lambda \mid a, b)=\frac{1}{\Gamma(a)} b^{a} \lambda^{a-1} \exp (-b \lambda) \\ &\mathbb{E}[\lambda] =\frac{a}{b} \\ &\operatorname{var}[\lambda] =\frac{a}{b^{2}} \end{aligned} Gam(λ∣a,b)=Γ(a)1baλa−1exp(−bλ)E[λ]=bavar[λ]=b2a -
给出后验形式 Gam ( λ ∣ a N , b N ) \operatorname{Gam}\left(\lambda \mid a_{N}, b_{N}\right) Gam(λ∣aN,bN),其中 σ M L 2 \sigma_{\mathrm{ML}}^{2} σML2 可以通过令 p ( X ∣ λ ) p(\mathbf{X} \mid \lambda) p(X∣λ) 导数为 0 得到
p ( λ ∣ X ) ∝ λ a 0 − 1 λ N / 2 exp { − b 0 λ − λ 2 ∑ n = 1 N ( x n − μ ) 2 } a N = a 0 + N 2 b N = b 0 + 1 2 ∑ n = 1 N ( x n − μ ) 2 = b 0 + N 2 σ M L 2 \begin{aligned} & p(\lambda \mid \mathbf{X}) \propto \lambda^{a_{0}-1} \lambda^{N / 2} \exp \left\{-b_{0} \lambda-\frac{\lambda}{2} \sum_{n=1}^{N}\left(x_{n}-\mu\right)^{2}\right\} \\ & a_{N}=a_{0}+\frac{N}{2}\\ & b_{N}=b_{0}+\frac{1}{2} \sum_{n=1}^{N}\left(x_{n}-\mu\right)^{2}=b_{0}+\frac{N}{2} \sigma_{\mathrm{ML}}^{2} \end{aligned} p(λ∣X)∝λa0−1λN/2exp{−b0λ−2λn=1∑N(xn−μ)2}aN=a0+2NbN=b0+21n=1∑N(xn−μ)2=b0+2NσML2
然后是第三种情况, μ \mu μ 与 λ \lambda λ 均未知时, p ( μ , λ ) p(\mu,\lambda) p(μ,λ) 的分布。此时略为复杂,因此只推导共轭先验的形式,后验的形式也可根据之前的方法继续推出。
- 写出似然函数
p ( X ∣ μ , λ ) = ∏ n = 1 N ( λ 2 π ) 1 / 2 exp { − λ 2 ( x n − μ ) 2 } ∝ [ λ 1 / 2 exp ( − λ μ 2 2 ) ] N exp { λ μ ∑ n = 1 N x n − λ 2 ∑ n = 1 N x n 2 } \begin{array}{c} p(\mathbf{X} \mid \mu, \lambda)=\prod_{n=1}^{N}\left(\frac{\lambda}{2 \pi}\right)^{1 / 2} \exp \left\{-\frac{\lambda}{2}\left(x_{n}-\mu\right)^{2}\right\} \\ \propto\left[\lambda^{1 / 2} \exp \left(-\frac{\lambda \mu^{2}}{2}\right)\right]^{N} \exp \left\{\lambda \mu \sum_{n=1}^{N} x_{n}-\frac{\lambda}{2} \sum_{n=1}^{N} x_{n}^{2}\right\} \end{array} p(X∣μ,λ)=∏n=1N(2πλ)1/2exp{−2λ(xn−μ)2}∝[λ1/2exp(−2λμ2)]Nexp{λμ∑n=1Nxn−2λ∑n=1Nxn2} - 将与
μ
\mu
μ 和
λ
\lambda
λ 无关的变量均变为常数,得到共轭先验的形式
p ( μ , λ ) ∝ [ λ 1 / 2 exp ( − λ μ 2 2 ) ] β exp { c λ μ − d λ } = exp { − β λ 2 ( μ − c / β ) 2 } λ β / 2 exp { − ( d − c 2 2 β ) λ } \begin{array}{l} p(\mu, \lambda) \propto\left[\lambda^{1 / 2} \exp \left(-\frac{\lambda \mu^{2}}{2}\right)\right]^{\beta} \exp \{c \lambda \mu-d \lambda\} \\ \quad=\exp \left\{-\frac{\beta \lambda}{2}(\mu-c / \beta)^{2}\right\} \lambda^{\beta / 2} \exp \left\{-\left(d-\frac{c^{2}}{2 \beta}\right) \lambda\right\} \end{array} p(μ,λ)∝[λ1/2exp(−2λμ2)]βexp{cλμ−dλ}=exp{−2βλ(μ−c/β)2}λβ/2exp{−(d−2βc2)λ}
由于
p
(
μ
,
λ
)
=
p
(
μ
∣
λ
)
p
(
λ
)
p(\mu, \lambda)=p(\mu \mid \lambda) p(\lambda)
p(μ,λ)=p(μ∣λ)p(λ),刚好对应表达式中的一、二部分,因此整理后得到:
p
(
μ
,
λ
)
=
N
(
μ
∣
μ
0
,
(
β
λ
)
−
1
)
Gam
(
λ
∣
a
,
b
)
μ
0
=
c
/
β
,
a
=
1
+
β
/
2
,
b
=
d
−
c
2
/
2
β
\begin{aligned} & p(\mu, \lambda)=\mathcal{N}\left(\mu \mid \mu_{0},(\beta \lambda)^{-1}\right) \operatorname{Gam}(\lambda \mid a, b) \\ &\mu_{0}=c / \beta, a=1+\beta / 2, b=d-c^{2} / 2 \beta \end{aligned}
p(μ,λ)=N(μ∣μ0,(βλ)−1)Gam(λ∣a,b)μ0=c/β,a=1+β/2,b=d−c2/2β
此时
p
(
μ
,
λ
)
p(\mu,\lambda)
p(μ,λ) 服从 normal-gamma / Gaussian-gamma 分布:
当问题扩展成多变量时 N ( x ∣ μ , Λ − 1 ) \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Lambda}^{-1}\right) N(x∣μ,Λ−1),已知 Λ − 1 \mathbf{\Lambda}^{-1} Λ−1, μ \boldsymbol{\mu} μ 分布依然为高斯;但其它两种情况发生了改变。
已知
μ
\boldsymbol{\mu}
μ,
Λ
−
1
\mathbf{\Lambda}^{-1}
Λ−1 的共轭先验为 Wishart distribution,其中
ν
\nu
ν 为分布*度:
W
(
Λ
∣
W
,
ν
)
=
B
∣
Λ
∣
(
ν
−
D
−
1
)
/
2
exp
(
−
1
2
Tr
(
W
−
1
Λ
)
)
,
B
is a constant
B
(
W
,
ν
)
=
∣
W
∣
−
ν
/
2
(
2
ν
D
/
2
π
D
(
D
−
1
)
/
4
∏
i
=
1
D
Γ
(
ν
+
1
−
i
2
)
)
−
1
\begin{aligned} & \mathcal{W}(\boldsymbol{\Lambda} \mid \mathbf{W}, \nu)=B|\boldsymbol{\Lambda}|^{(\nu-D-1) / 2} \exp \left(-\frac{1}{2} \operatorname{Tr}\left(\mathbf{W}^{-1} \boldsymbol{\Lambda}\right)\right),B \ \text{is a constant} \\ & B(\mathbf{W}, \nu)=|\mathbf{W}|^{-\nu / 2}\left(2^{\nu D / 2} \pi^{D(D-1) / 4} \prod_{i=1}^{D} \Gamma\left(\frac{\nu+1-i}{2}\right)\right)^{-1} \end{aligned}
W(Λ∣W,ν)=B∣Λ∣(ν−D−1)/2exp(−21Tr(W−1Λ)),B is a constantB(W,ν)=∣W∣−ν/2(2νD/2πD(D−1)/4i=1∏DΓ(2ν+1−i))−1
另外,
μ
\boldsymbol{\mu}
μ,
Λ
−
1
\mathbf{\Lambda}^{-1}
Λ−1 均未知时,共轭先验为 normal-Wishart or Gaussian-Wishart distribution:
p
(
μ
,
Λ
∣
μ
0
,
β
,
W
,
ν
)
=
N
(
μ
∣
μ
0
,
(
β
Λ
)
−
1
)
W
(
Λ
∣
W
,
ν
)
p\left(\boldsymbol{\mu}, \boldsymbol{\Lambda} \mid \boldsymbol{\mu}_{0}, \beta, \mathbf{W}, \nu\right)=\mathcal{N}\left(\boldsymbol{\mu} \mid \boldsymbol{\mu}_{0},(\beta \mathbf{\Lambda})^{-1}\right) \mathcal{W}(\boldsymbol{\Lambda} \mid \mathbf{W}, \nu)
p(μ,Λ∣μ0,β,W,ν)=N(μ∣μ0,(βΛ)−1)W(Λ∣W,ν)
最后总结一下:
2.3.7 t-分布 (Student’s t-distribution)
假设
x
∼
N
(
x
∣
μ
,
τ
−
1
)
x\sim \mathcal{N}\left(x \mid \mu, \tau^{-1}\right)
x∼N(x∣μ,τ−1),其中
μ
\mu
μ 已知,而
τ
−
1
\tau^{-1}
τ−1 服从
Gam
(
τ
∣
a
,
b
)
\operatorname{Gam}(\tau \mid a, b)
Gam(τ∣a,b) 的分布,则
x
x
x 的边缘分布如下:
p
(
x
∣
μ
,
a
,
b
)
=
∫
0
∞
N
(
x
∣
μ
,
τ
−
1
)
Gam
(
τ
∣
a
,
b
)
d
τ
=
∫
0
∞
b
a
e
(
−
b
τ
)
τ
a
−
1
Γ
(
a
)
(
τ
2
π
)
1
/
2
exp
{
−
τ
2
(
x
−
μ
)
2
}
d
τ
=
b
a
Γ
(
a
)
(
1
2
π
)
1
/
2
[
b
+
(
x
−
μ
)
2
2
]
−
a
−
1
/
2
Γ
(
a
+
1
/
2
)
\begin{aligned} p(x \mid \mu, a, b) &=\int_{0}^{\infty} \mathcal{N}\left(x \mid \mu, \tau^{-1}\right) \operatorname{Gam}(\tau \mid a, b) \mathrm{d} \tau \\ &=\int_{0}^{\infty} \frac{b^{a} e^{(-b \tau)} \tau^{a-1}}{\Gamma(a)}\left(\frac{\tau}{2 \pi}\right)^{1 / 2} \exp \left\{-\frac{\tau}{2}(x-\mu)^{2}\right\} \mathrm{d} \tau \\ &=\frac{b^{a}}{\Gamma(a)}\left(\frac{1}{2 \pi}\right)^{1 / 2}\left[b+\frac{(x-\mu)^{2}}{2}\right]^{-a-1 / 2} \Gamma(a+1 / 2) \end{aligned}
p(x∣μ,a,b)=∫0∞N(x∣μ,τ−1)Gam(τ∣a,b)dτ=∫0∞Γ(a)bae(−bτ)τa−1(2πτ)1/2exp{−2τ(x−μ)2}dτ=Γ(a)ba(2π1)1/2[b+2(x−μ)2]−a−1/2Γ(a+1/2)
其中第二步的推导如下:
令
ν
=
2
a
,
λ
=
a
/
b
\nu=2 a, \lambda=a / b
ν=2a,λ=a/b,替换上式中的
a
,
b
a,b
a,b 后,即可推导出学生分布 (Student’s t-distribution):
St
(
x
∣
μ
,
λ
,
ν
)
=
Γ
(
ν
/
2
+
1
/
2
)
Γ
(
ν
/
2
)
(
λ
π
ν
)
1
/
2
[
1
+
λ
(
x
−
μ
)
2
ν
]
−
ν
/
2
−
1
/
2
\operatorname{St}(x \mid \mu, \lambda, \nu)=\frac{\Gamma(\nu / 2+1 / 2)}{\Gamma(\nu / 2)}\left(\frac{\lambda}{\pi \nu}\right)^{1 / 2}\left[1+\frac{\lambda(x-\mu)^{2}}{\nu}\right]^{-\nu / 2-1 / 2}
St(x∣μ,λ,ν)=Γ(ν/2)Γ(ν/2+1/2)(πνλ)1/2[1+νλ(x−μ)2]−ν/2−1/2
其中
λ
\lambda
λ 被称为 t-分布的精度 (precision),而
ν
\nu
ν 被称为*度,其对分布的影响如下:
当
ν
=
1
\nu=1
ν=1 时,t-分布变为柯西分布 (Cauchy);当
ν
→
∞
\nu\rightarrow\infty
ν→∞ 时,t-分布变为
N
(
x
∣
μ
,
λ
−
1
)
\mathcal{N}\left(x \mid \mu, \lambda^{-1}\right)
N(x∣μ,λ−1)。
t-分布可以看作是无限个高斯分布的组合,这使得 t-分布比起高斯分布在分布图像上 tail 更长,这也使得 t-分布更加鲁棒。由此,基于高斯分布下极大似然的最小二乘法在面对离群点时可能效果较差,可以考虑换成 t-分布。
接下来进行多元变量的 t-分布推导。
回到
p
(
x
∣
μ
,
a
,
b
)
=
∫
0
∞
N
(
x
∣
μ
,
τ
−
1
)
Gam
(
τ
∣
a
,
b
)
d
τ
p(x \mid \mu, a, b)=\int_{0}^{\infty} \mathcal{N}\left(x \mid \mu, \tau^{-1}\right) \operatorname{Gam}(\tau \mid a, b) \mathrm{d} \tau
p(x∣μ,a,b)=∫0∞N(x∣μ,τ−1)Gam(τ∣a,b)dτ 式子中,令
ν
=
2
a
,
λ
=
a
/
b
,
η
=
τ
b
/
a
\nu=2a,\lambda=a/b,\eta=\tau b / a
ν=2a,λ=a/b,η=τb/a,目的是在
p
(
x
∣
μ
,
a
,
b
)
p(x \mid \mu, a, b)
p(x∣μ,a,b) 中加入高斯分布的
λ
\lambda
λ,转换后得到:
St
(
x
∣
μ
,
λ
,
ν
)
=
∫
0
∞
N
(
x
∣
μ
,
(
η
λ
)
−
1
)
Gam
(
η
∣
ν
/
2
,
ν
/
2
)
d
η
\operatorname{St}(x \mid \mu, \lambda, \nu)=\int_{0}^{\infty} \mathcal{N}\left(x \mid \mu,(\eta \lambda)^{-1}\right) \operatorname{Gam}(\eta \mid \nu / 2, \nu / 2) \mathrm{d} \eta
St(x∣μ,λ,ν)=∫0∞N(x∣μ,(ηλ)−1)Gam(η∣ν/2,ν/2)dη
将
x
x
x 扩展成
x
\mathbf{x}
x,得到:
St
(
x
∣
μ
,
Λ
,
ν
)
=
∫
0
∞
N
(
x
∣
μ
,
(
η
Λ
)
−
1
)
Gam
(
η
∣
ν
/
2
,
ν
/
2
)
d
η
\operatorname{St}(\mathbf{x} \mid \boldsymbol{\mu}, \boldsymbol{\Lambda}, \nu)=\int_{0}^{\infty} \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu},(\eta \boldsymbol{\Lambda})^{-1}\right) \operatorname{Gam}(\eta \mid \nu / 2, \nu / 2) \mathrm{d} \eta
St(x∣μ,Λ,ν)=∫0∞N(x∣μ,(ηΛ)−1)Gam(η∣ν/2,ν/2)dη
求出积分后得到:
St
(
x
∣
μ
,
Λ
,
ν
)
=
Γ
(
D
/
2
+
ν
/
2
)
Γ
(
ν
/
2
)
∣
Λ
∣
1
/
2
(
π
ν
)
D
/
2
[
1
+
Δ
2
ν
]
−
D
/
2
−
ν
/
2
Δ
2
=
(
x
−
μ
)
T
Λ
(
x
−
μ
)
\begin{aligned} &\operatorname{St}(\mathbf{x} \mid \boldsymbol{\mu}, \boldsymbol{\Lambda}, \nu)=\frac{\Gamma(D / 2+\nu / 2)}{\Gamma(\nu / 2)} \frac{|\boldsymbol{\Lambda}|^{1 / 2}}{(\pi \nu)^{D / 2}}\left[1+\frac{\Delta^{2}}{\nu}\right]^{-D / 2-\nu / 2}\\ &\Delta^{2}=(\mathbf{x}-\boldsymbol{\mu})^{\mathrm{T}} \boldsymbol{\Lambda}(\mathbf{x}-\boldsymbol{\mu}) \end{aligned}
St(x∣μ,Λ,ν)=Γ(ν/2)Γ(D/2+ν/2)(πν)D/2∣Λ∣1/2[1+νΔ2]−D/2−ν/2Δ2=(x−μ)TΛ(x−μ)
其中
D
D
D 为
x
\mathbf{x}
x 的维度,
Δ
2
\Delta^2
Δ2 为 squared Mahalanobis distance,具体推导过程如下:
最后,我们求解 t-分布的期望、方差以及最值点:
E
[
x
]
=
μ
,
if
ν
>
1
cov
[
x
]
=
ν
(
ν
−
2
)
Λ
−
1
,
if
ν
>
2
mode
[
x
]
=
μ
\begin{aligned} \mathbb{E}[\mathbf{x}] &=\boldsymbol{\mu}, & & \text { if } \quad \nu>1 \\ \operatorname{cov}[\mathbf{x}] &=\frac{\nu}{(\nu-2)} \mathbf{\Lambda}^{-1}, & & \text { if } \quad \nu>2 \\ \operatorname{mode}[\mathbf{x}] &=\boldsymbol{\mu} & & \end{aligned}
E[x]cov[x]mode[x]=μ,=(ν−2)νΛ−1,=μ if ν>1 if ν>2
2.3.8 周期变量 (Periodic variables)
本章以角度 ( θ \theta θ) 为例,讲解如何用高斯分布来度量周期变量。
我们用单位圆上的
x
=
(
x
1
,
x
2
)
\mathbf{x}=(x_1,x_2)
x=(x1,x2) 来代表
θ
\theta
θ,则假设
x
∼
N
(
μ
=
(
μ
1
,
μ
2
)
,
Σ
=
σ
2
I
)
\mathbf{x}\sim \mathcal{N}(\boldsymbol{\mu}=\left(\mu_{1}, \mu_{2}\right),\boldsymbol{\Sigma}=\sigma^{2} \mathbf{I})
x∼N(μ=(μ1,μ2),Σ=σ2I),即:
p
(
x
1
,
x
2
)
=
1
2
π
σ
2
exp
{
−
(
x
1
−
μ
1
)
2
+
(
x
2
−
μ
2
)
2
2
σ
2
}
p\left(x_{1}, x_{2}\right)=\frac{1}{2 \pi \sigma^{2}} \exp \left\{-\frac{\left(x_{1}-\mu_{1}\right)^{2}+\left(x_{2}-\mu_{2}\right)^{2}}{2 \sigma^{2}}\right\}
p(x1,x2)=2πσ21exp{−2σ2(x1−μ1)2+(x2−μ2)2}
接下来用极角坐标
(
r
,
θ
)
(r,\theta)
(r,θ) 来代替
x
\mathbf{x}
x,即令
x
1
=
r
cos
θ
,
x
2
=
r
sin
θ
,
μ
1
=
r
0
cos
θ
0
,
μ
2
=
r
0
sin
θ
0
x_1=r\cos\theta,x_2=r\sin\theta,\mu_1=r_0\cos\theta_0,\mu_2=r_0\sin\theta_0
x1=rcosθ,x2=rsinθ,μ1=r0cosθ0,μ2=r0sinθ0,则高斯分布 exp 中的内容变换为:
−
1
2
σ
2
{
(
r
cos
θ
−
r
0
cos
θ
0
)
2
+
(
r
sin
θ
−
r
0
sin
θ
0
)
2
}
=
−
1
2
σ
2
{
1
+
r
0
2
−
2
r
0
cos
θ
cos
θ
0
−
2
r
0
sin
θ
sin
θ
0
}
=
r
0
σ
2
cos
(
θ
−
θ
0
)
+
const
\begin{aligned} -& \frac{1}{2 \sigma^{2}}\left\{\left(r \cos \theta-r_{0} \cos \theta_{0}\right)^{2}+\left(r \sin \theta-r_{0} \sin \theta_{0}\right)^{2}\right\} \\ &=-\frac{1}{2 \sigma^{2}}\left\{1+r_{0}^{2}-2 r_{0} \cos \theta \cos \theta_{0}-2 r_{0} \sin \theta \sin \theta_{0}\right\} \\ &=\frac{r_{0}}{\sigma^{2}} \cos \left(\theta-\theta_{0}\right)+\text { const } \end{aligned}
−2σ21{(rcosθ−r0cosθ0)2+(rsinθ−r0sinθ0)2}=−2σ21{1+r02−2r0cosθcosθ0−2r0sinθsinθ0}=σ2r0cos(θ−θ0)+ const
令
m
=
r
0
/
σ
2
m=r_0/\sigma^2
m=r0/σ2,
p
(
θ
)
p(\theta)
p(θ) 转换为:
p
(
θ
∣
θ
0
,
m
)
=
1
2
π
I
0
(
m
)
exp
{
m
cos
(
θ
−
θ
0
)
}
I
0
(
m
)
=
1
2
π
∫
0
2
π
exp
{
m
cos
θ
}
d
θ
\begin{aligned} & p\left(\theta \mid \theta_{0}, m\right)=\frac{1}{2 \pi I_{0}(m)} \exp \left\{m \cos \left(\theta-\theta_{0}\right)\right\}\\ & I_{0}(m)=\frac{1}{2 \pi} \int_{0}^{2 \pi} \exp \{m \cos \theta\} \mathrm{d} \theta \end{aligned}
p(θ∣θ0,m)=2πI0(m)1exp{mcos(θ−θ0)}I0(m)=2π1∫02πexp{mcosθ}dθ
上述即为 von Mises distribution, or the circular normal,其中
θ
0
,
m
\theta_0,m
θ0,m 分别表示均值与浓度参数 (concentration parameter)。如下为分布图:
该分布有两条性质:
- 当 m 很大时,分布图集中在峰值 θ 0 \theta_0 θ0 处
- 当 m 很大时,将变为高斯分布
接下来对该分布进行极大似然估计:
ln
p
(
D
∣
θ
0
,
m
)
=
−
N
ln
(
2
π
)
−
N
ln
I
0
(
m
)
+
m
∑
n
=
1
N
cos
(
θ
n
−
θ
0
)
\ln p\left(\mathcal{D} \mid \theta_{0}, m\right)=-N \ln (2 \pi)-N \ln I_{0}(m)+m \sum_{n=1}^{N} \cos \left(\theta_{n}-\theta_{0}\right)
lnp(D∣θ0,m)=−Nln(2π)−NlnI0(m)+mn=1∑Ncos(θn−θ0)
分别对
θ
0
\theta_0
θ0 和
m
m
m 求导,可以得到如下估计值:
θ
0
M
L
=
tan
−
1
{
∑
n
sin
θ
n
∑
n
cos
θ
n
}
A
(
m
)
=
1
N
∑
n
=
1
N
cos
(
θ
n
−
θ
0
M
L
)
=
I
1
(
m
)
I
0
(
m
)
,
I
1
(
m
)
=
I
0
′
(
m
)
A
(
m
M
L
)
=
(
1
N
∑
n
=
1
N
cos
θ
n
)
cos
θ
0
M
L
−
(
1
N
∑
n
=
1
N
sin
θ
n
)
sin
θ
0
M
L
\begin{aligned} & \theta_{0}^{\mathrm{ML}}=\tan ^{-1}\left\{\frac{\sum_{n} \sin \theta_{n}}{\sum_{n} \cos \theta_{n}}\right\}\\ & A(m)=\frac{1}{N} \sum_{n=1}^{N} \cos \left(\theta_{n}-\theta_{0}^{\mathrm{ML}}\right)=\frac{I_{1}(m)}{I_{0}(m)},I_{1}(m)=I_{0}^{\prime}(m) \\ & A\left(m_{\mathrm{ML}}\right)=\left(\frac{1}{N} \sum_{n=1}^{N} \cos \theta_{n}\right) \cos \theta_{0}^{\mathrm{ML}}-\left(\frac{1}{N} \sum_{n=1}^{N} \sin \theta_{n}\right) \sin \theta_{0}^{\mathrm{ML}} \end{aligned}
θ0ML=tan−1{∑ncosθn∑nsinθn}A(m)=N1n=1∑Ncos(θn−θ0ML)=I0(m)I1(m),I1(m)=I0′(m)A(mML)=(N1n=1∑Ncosθn)cosθ0ML−(N1n=1∑Nsinθn)sinθ0ML
最后,von Mises distribution 有一个局限,即其是单峰分布;但我们可以通过将多个该分布组合 (mixtures),来得到能够处理多峰分布的周期变量。
2.3.9 高斯混合 (Mixtures of Gaussians)
高斯混合,即将多个高斯分布线性组合。理论上来说,如果高斯分布的数量足够多,且可以任意调节均值、方差以及线性组合系数,几乎任何连续分布图都可以在任意的精度上拟合出来。
高斯混合 (mixture of Gaussians):
- N ( x ∣ μ k , Σ k ) \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Sigma}_{k}\right) N(x∣μk,Σk) 为混合的一个部分 (component)
-
π
k
\pi_{k}
πk 为混合系数 (mixing coefficients),满足
∑
k
=
1
K
π
k
=
1
\sum_{k=1}^{K} \pi_{k}=1
∑k=1Kπk=1 与
0
⩽
π
k
⩽
1
0 \leqslant \pi_{k} \leqslant 1
0⩽πk⩽1
p ( x ) = ∑ k = 1 K π k N ( x ∣ μ k , Σ k ) p(\mathbf{x})=\sum_{k=1}^{K} \pi_{k} \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Sigma}_{k}\right) p(x)=k=1∑KπkN(x∣μk,Σk)
从贝叶斯推断视角出发:
- 其中
p
(
k
∣
x
)
p(k|\mathbf{x})
p(k∣x) 被称为 responsibilities
p ( x ) = ∑ k = 1 K p ( k ) p ( x ∣ k ) p ( k ) = π k p ( x ∣ k ) = N ( x ∣ μ k , Σ k ) γ k ( x ) ≡ p ( k ∣ x ) = p ( k ) p ( x ∣ k ) ∑ l p ( l ) p ( x ∣ l ) = π k N ( x ∣ μ k , Σ k ) ∑ l π l N ( x ∣ μ l , Σ l ) \begin{aligned} & p(\mathbf{x})=\sum_{k=1}^{K} p(k) p(\mathbf{x} \mid k) \\ & p(k)=\pi_k\\ & p(\mathbf{x} \mid k) = \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}_{k}, \mathbf{\Sigma}_{k}\right) \\ & \begin{aligned} \gamma_{k}(\mathbf{x}) & \equiv p(k \mid \mathbf{x}) \\ &=\frac{p(k) p(\mathbf{x} \mid k)}{\sum_{l} p(l) p(\mathbf{x} \mid l)} \\ &=\frac{\pi_{k} \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}_{k}, \mathbf{\Sigma}_{k}\right)}{\sum_{l} \pi_{l} \mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}_{l}, \mathbf{\Sigma}_{l}\right)} \end{aligned} \end{aligned} p(x)=k=1∑Kp(k)p(x∣k)p(k)=πkp(x∣k)=N(x∣μk,Σk)γk(x)≡p(k∣x)=∑lp(l)p(x∣l)p(k)p(x∣k)=∑lπlN(x∣μl,Σl)πkN(x∣μk,Σk)
从极大似然估计出发:
- 无法直接得到必式解,需要使用迭代求解的优化方法
ln p ( X ∣ π , μ , Σ ) = ∑ n = 1 N ln { ∑ k = 1 K π k N ( x n ∣ μ k , Σ k ) } \ln p(\mathbf{X} \mid \boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\Sigma})=\sum_{n=1}^{N} \ln \left\{\sum_{k=1}^{K} \pi_{k} \mathcal{N}\left(\mathbf{x}_{n} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Sigma}_{k}\right)\right\} lnp(X∣π,μ,Σ)=n=1∑Nln{k=1∑KπkN(xn∣μk,Σk)}
2.4 指数族 (The Exponential Family)
指数族形式:
- x \mathbf{x} x 可以是标量或向量
- η \eta η 是分布的自然参数 (natural parameters)
-
g
(
η
)
g(\eta)
g(η) 用来标准化分布
g
(
η
)
∫
h
(
x
)
exp
{
η
T
u
(
x
)
}
d
x
=
1
g(\boldsymbol{\eta}) \int h(\mathbf{x}) \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \mathbf{u}(\mathbf{x})\right\} \mathrm{d} \mathbf{x}=1
g(η)∫h(x)exp{ηTu(x)}dx=1
p ( x ∣ η ) = h ( x ) g ( η ) exp { η T u ( x ) } p(\mathbf{x} \mid \boldsymbol{\eta})=h(\mathbf{x}) g(\boldsymbol{\eta}) \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \mathbf{u}(\mathbf{x})\right\} p(x∣η)=h(x)g(η)exp{ηTu(x)}
接下来以三个指数族分布进行举例。
首先是伯努利分布 (Bernoulli distribution):
p
(
x
∣
μ
)
=
Bern
(
x
∣
μ
)
=
μ
x
(
1
−
μ
)
1
−
x
=
exp
{
x
ln
μ
+
(
1
−
x
)
ln
(
1
−
μ
)
}
=
(
1
−
μ
)
exp
{
ln
(
μ
1
−
μ
)
x
}
\begin{aligned} p(x \mid \mu) &=\operatorname{Bern}(x \mid \mu)\\ & =\mu^{x}(1-\mu)^{1-x}\\ &=\exp \{x \ln \mu+(1-x) \ln (1-\mu)\} \\ &=(1-\mu) \exp \left\{\ln \left(\frac{\mu}{1-\mu}\right) x\right\} \end{aligned}
p(x∣μ)=Bern(x∣μ)=μx(1−μ)1−x=exp{xlnμ+(1−x)ln(1−μ)}=(1−μ)exp{ln(1−μμ)x}
先关注 exp 中的形式,求出
η
\eta
η:
η
=
ln
(
μ
1
−
μ
)
μ
=
σ
(
η
)
=
1
1
+
exp
(
−
η
)
\begin{aligned} & \eta=\ln \left(\frac{\mu}{1-\mu}\right) \\ & \mu=\sigma(\eta)=\frac{1}{1+\exp (-\eta)}\\ \end{aligned}
η=ln(1−μμ)μ=σ(η)=1+exp(−η)1
最后得到指数族的标准形式:
p
(
x
∣
η
)
=
σ
(
−
η
)
exp
(
η
x
)
u
(
x
)
=
x
h
(
x
)
=
1
g
(
η
)
=
σ
(
−
η
)
\begin{aligned} & p(x \mid \eta)=\sigma(-\eta) \exp (\eta x) \\ & u(x)=x \\ & h(x)=1 \\ & g(\eta)=\sigma(-\eta) \end{aligned}
p(x∣η)=σ(−η)exp(ηx)u(x)=xh(x)=1g(η)=σ(−η)
接下来是多项式分布 (multinomial distribution):
p
(
x
∣
μ
)
=
∏
k
=
1
M
μ
k
x
k
=
exp
{
∑
k
=
1
M
x
k
ln
μ
k
}
p(\mathbf{x} \mid \boldsymbol{\mu})=\prod_{k=1}^{M} \mu_{k}^{x_{k}}=\exp \left\{\sum_{k=1}^{M} x_{k} \ln \mu_{k}\right\}
p(x∣μ)=k=1∏Mμkxk=exp{k=1∑Mxklnμk}
其指数族标准形式如下:
p
(
x
∣
η
)
=
exp
(
η
T
x
)
u
(
x
)
=
x
h
(
x
)
=
1
g
(
η
)
=
1
\begin{aligned} &p(\mathbf{x} \mid \boldsymbol{\eta})=\exp \left(\boldsymbol{\eta}^{\mathrm{T}} \mathbf{x}\right) \\ & \mathbf{u}(\mathbf{x})=\mathbf{x} \\ & h(\mathbf{x})=1 \\ & g(\boldsymbol{\eta})=1 \end{aligned}
p(x∣η)=exp(ηTx)u(x)=xh(x)=1g(η)=1
最后是高斯分布 (Gaussian distribution):
p
(
x
∣
μ
,
σ
2
)
=
1
(
2
π
σ
2
)
1
/
2
exp
{
−
1
2
σ
2
(
x
−
μ
)
2
}
=
1
(
2
π
σ
2
)
1
/
2
exp
{
−
1
2
σ
2
x
2
+
μ
σ
2
x
−
1
2
σ
2
μ
2
}
\begin{aligned} p\left(x \mid \mu, \sigma^{2}\right) &=\frac{1}{\left(2 \pi \sigma^{2}\right)^{1 / 2}} \exp \left\{-\frac{1}{2 \sigma^{2}}(x-\mu)^{2}\right\} \\ &=\frac{1}{\left(2 \pi \sigma^{2}\right)^{1 / 2}} \exp \left\{-\frac{1}{2 \sigma^{2}} x^{2}+\frac{\mu}{\sigma^{2}} x-\frac{1}{2 \sigma^{2}} \mu^{2}\right\} \end{aligned}
p(x∣μ,σ2)=(2πσ2)1/21exp{−2σ21(x−μ)2}=(2πσ2)1/21exp{−2σ21x2+σ2μx−2σ21μ2}
其对应参数为:
η
=
(
μ
/
σ
2
−
1
/
2
σ
2
)
u
(
x
)
=
(
x
x
2
)
h
(
x
)
=
(
2
π
)
−
1
/
2
g
(
η
)
=
(
−
2
η
2
)
1
/
2
exp
(
η
1
2
4
η
2
)
\begin{aligned} \boldsymbol{\eta} &=\left(\begin{array}{c} \mu / \sigma^{2} \\ -1 / 2 \sigma^{2} \end{array}\right) \\ \mathbf{u}(x) &=\left(\begin{array}{c} x \\ x^{2} \end{array}\right) \\ h(\mathbf{x}) &=(2 \pi)^{-1 / 2} \\ g(\boldsymbol{\eta}) &=\left(-2 \eta_{2}\right)^{1 / 2} \exp \left(\frac{\eta_{1}^{2}}{4 \eta_{2}}\right) \end{aligned}
ηu(x)h(x)g(η)=(μ/σ2−1/2σ2)=(xx2)=(2π)−1/2=(−2η2)1/2exp(4η2η12)
2.4.1 极大似然与充分统计量 (Maximum likelihood and sufficient statistics)
g ( η ) ∫ h ( x ) exp { η T u ( x ) } d x = 1 g(\boldsymbol{\eta}) \int h(\mathbf{x}) \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \mathbf{u}(\mathbf{x})\right\} \mathrm{d} \mathbf{x}=1 g(η)∫h(x)exp{ηTu(x)}dx=1
上式是指数族分布需要满足的基础条件,对上式等号两边同时对
η
\eta
η 求导,则可以得到下式:
∇
g
(
η
)
∫
h
(
x
)
exp
{
η
T
u
(
x
)
}
d
x
+
g
(
η
)
∫
h
(
x
)
exp
{
η
T
u
(
x
)
}
u
(
x
)
d
x
=
0
\begin{array}{l} \nabla g(\boldsymbol{\eta}) \int h(\mathbf{x}) \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \mathbf{u}(\mathbf{x})\right\} \mathrm{d} \mathbf{x} \\ +g(\boldsymbol{\eta}) \int h(\mathbf{x}) \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \mathbf{u}(\mathbf{x})\right\} \mathbf{u}(\mathbf{x}) \mathrm{d} \mathbf{x}=0 \end{array}
∇g(η)∫h(x)exp{ηTu(x)}dx+g(η)∫h(x)exp{ηTu(x)}u(x)dx=0
进一步化简可以得到
η
\eta
η 的表达式:
−
1
g
(
η
)
∇
g
(
η
)
=
g
(
η
)
∫
h
(
x
)
exp
{
η
T
u
(
x
)
}
u
(
x
)
d
x
=
E
[
u
(
x
)
]
−
∇
ln
g
(
η
)
=
E
[
u
(
x
)
]
\begin{aligned} & -\frac{1}{g(\boldsymbol{\eta})} \nabla g(\boldsymbol{\eta})=g(\boldsymbol{\eta}) \int h(\mathbf{x}) \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \mathbf{u}(\mathbf{x})\right\} \mathbf{u}(\mathbf{x}) \mathrm{d} \mathbf{x}=\mathbb{E}[\mathbf{u}(\mathbf{x})]\\ & -\nabla \ln g(\boldsymbol{\eta})=\mathbb{E}[\mathbf{u}(\mathbf{x})] \end{aligned}
−g(η)1∇g(η)=g(η)∫h(x)exp{ηTu(x)}u(x)dx=E[u(x)]−∇lng(η)=E[u(x)]
接下来对于数据
X
=
{
x
1
,
.
.
.
,
x
n
}
\mathbf{X}=\{\mathbf{x_1},...,\mathbf{x_n}\}
X={x1,...,xn} 进行极大似然估计:
p
(
X
∣
η
)
=
(
∏
n
=
1
N
h
(
x
n
)
)
g
(
η
)
N
exp
{
η
T
∑
n
=
1
N
u
(
x
n
)
}
p(\mathbf{X} \mid \boldsymbol{\eta})=\left(\prod_{n=1}^{N} h\left(\mathbf{x}_{n}\right)\right) g(\boldsymbol{\eta})^{N} \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \sum_{n=1}^{N} \mathbf{u}\left(\mathbf{x}_{n}\right)\right\}
p(X∣η)=(n=1∏Nh(xn))g(η)Nexp{ηTn=1∑Nu(xn)}
对
η
\eta
η 求导并令其等于 0 后,得到下式:
−
∇
ln
g
(
η
M
L
)
=
1
N
∑
n
=
1
N
u
(
x
n
)
-\nabla \ln g\left(\boldsymbol{\eta}_{\mathrm{ML}}\right)=\frac{1}{N} \sum_{n=1}^{N} \mathbf{u}\left(\mathbf{x}_{n}\right)
−∇lng(ηML)=N1n=1∑Nu(xn)
可以发现对于 η \eta η 的求解,我们只需要保存 ∑ n = 1 N u ( x n ) \sum_{n=1}^{N} \mathbf{u}\left(\mathbf{x}_{n}\right) ∑n=1Nu(xn) 的数值,而不用保存整个数据集,因此 ∑ n = 1 N u ( x n ) \sum_{n=1}^{N} \mathbf{u}\left(\mathbf{x}_{n}\right) ∑n=1Nu(xn) 是指数族分布的充分统计量 (sufficient statistic)。
2.4.2 共轭先验 (Conjugate priors)
In general, for a given probability distribution p ( x ∣ η ) p(x|\eta) p(x∣η), we can seek a prior p ( η ) p(\eta) p(η) that is conjugate to the likelihood function, so that the posterior distribution has the same functional form as the prior.
对于指数族分布,其共轭先验形式如下:
-
f
(
χ
,
ν
)
f(\boldsymbol{\chi}, \nu)
f(χ,ν) 为归一化参数
p ( η ∣ χ , ν ) = f ( χ , ν ) g ( η ) ν exp { ν η T χ } p(\boldsymbol{\eta} \mid \boldsymbol{\chi}, \nu)=f(\boldsymbol{\chi}, \nu) g(\boldsymbol{\eta})^{\nu} \exp \left\{\nu \boldsymbol{\eta}^{\mathrm{T}} \boldsymbol{\chi}\right\} p(η∣χ,ν)=f(χ,ν)g(η)νexp{νηTχ}
验证一下上述形式是否成立,似然函数形式:
p
(
X
∣
η
)
=
(
∏
n
=
1
N
h
(
x
n
)
)
g
(
η
)
N
exp
{
η
T
∑
n
=
1
N
u
(
x
n
)
}
p(\mathbf{X} \mid \boldsymbol{\eta})=\left(\prod_{n=1}^{N} h\left(\mathbf{x}_{n}\right)\right) g(\boldsymbol{\eta})^{N} \exp \left\{\boldsymbol{\eta}^{\mathrm{T}} \sum_{n=1}^{N} \mathbf{u}\left(\mathbf{x}_{n}\right)\right\}
p(X∣η)=(n=1∏Nh(xn))g(η)Nexp{ηTn=1∑Nu(xn)}
由此可以得到后验形式:
p
(
η
∣
X
,
χ
,
ν
)
∝
p
(
X
∣
η
)
p
(
η
∣
χ
,
ν
)
∝
g
(
η
)
ν
+
N
exp
{
η
T
(
∑
n
=
1
N
u
(
x
n
)
+
ν
χ
)
}
p(\boldsymbol{\eta} \mid \mathbf{X}, \boldsymbol{\chi}, \nu) \propto p(\mathbf{X} \mid \boldsymbol{\eta})p(\boldsymbol{\eta} \mid \boldsymbol{\chi}, \nu)\propto g(\boldsymbol{\eta})^{\nu+N} \exp \left\{\boldsymbol{\eta}^{\mathrm{T}}\left(\sum_{n=1}^{N} \mathbf{u}\left(\mathbf{x}_{n}\right)+\nu \boldsymbol{\chi}\right)\right\}
p(η∣X,χ,ν)∝p(X∣η)p(η∣χ,ν)∝g(η)ν+Nexp{ηT(n=1∑Nu(xn)+νχ)}
与先验形式一致,因此上述共轭形式成立。
2.4.3 无信息先验 (Noninformative priors)
在很多情况下,我们都对先验概率的形式一无所知,这时我们倾向于让先验形式对后验分布的影响越小越好,此时的先验就被称为无信息先验 (noninformative prior)。
当参数的取值有界时,可以根据均匀分布制定无信息先验;但当参数的取值*时,均匀分布则无法归一化,是 improper 的。
但在实际应用中,improper priors 也可以使用,只要对应的后验分布可以归一化。
书中举了两个例子,第一个是已知 σ 2 \sigma^2 σ2 求 μ \mu μ 分布的无信息先验 p ( μ ∣ μ 0 , σ 0 2 ) = N ( μ ∣ μ 0 , σ 0 2 ) p \left(\mu \mid \mu_{0}, \sigma_{0}^{2}\right)=\mathcal{N}\left(\mu \mid \mu_{0}, \sigma_{0}^{2}\right) p(μ∣μ0,σ02)=N(μ∣μ0,σ02),此时令 σ 0 2 → ∞ \sigma_0^2\rightarrow \infty σ02→∞ 则先验不对后验产生影响。
第二个是已知 μ \mu μ 求 σ 2 \sigma^2 σ2 的情况,此时先验符合 Gam ( λ ∣ a 0 , b 0 ) \text{Gam}(\lambda|a_0,b_0) Gam(λ∣a0,b0) 分布,令 a 0 = b 0 = 0 a_0=b_0=0 a0=b0=0,则先验不对后验产生影响,满足无信息先验的条件。
2.5 非参数化方法 (Nonparametric Methods)
参数化方法 (parametric methods):
- 预先确定好数据分布的函数形式,再根据数据来确定函数中的参数
- 缺点:预先选好的形式对数据的拟合能力可能很差
非参数化方法 (nonparametric methods):
- 对数据的分布形式不做或做很少的假设
本节从最直观的直方图统计方法引入,介绍非参数化方法:
- 将数据的定义域分为若干个等大的 bins
- 每个 bins 的概率用频率 p i = n i N Δ i p_{i}=\frac{n_{i}}{N \Delta_{i}} pi=NΔini 来定义
在上图中,绿线为原始数据分布, Δ \Delta Δ 为 bins 选取的大小。不难发现, Δ \Delta Δ 太大和太小都不能很好的拟合数据,大小适中的 Δ \Delta Δ 反而效果好一些。
总结一下,直方图统计法的优缺点:
- 【优点 1】直方图计算出来后,无需再保存原始数据
- 【优点 2】比较方便处理顺序到达的数据 (online learning)
- 【缺点 1】由于每个 bins 区间的左右开闭性,导致分布不连续
- 【缺点 2】如果数据是 D 维,每维开 M 个 bins,则 bins 的数量将达到 M D M^D MD,指数爆炸
在实际中,直方图统计法主要用于在某一两个维度上对数据分布有一个粗糙的直观理解,一般不用于概率分布估计。
但从该方法依然给了我们一些启发:
- 估计某点的概率密度,我们需要考虑该点附近的数据;直方图统计中,bins 的宽度就相当于一个自然平滑参数 (natural smoothing parameter),对局部区域的空间大小进行了描述
- 在选取自然平滑参数时,太大和太小都不合适
2.5.1 核密度估计 (Kernel density estimators)
根据之前得到的启发,假设一个小区域
R
\mathcal{R}
R 包含
x
\mathbf{x}
x,则
x
\mathbf{x}
x 落入
R
\mathcal{R}
R 的概率为:
P
=
∫
R
p
(
x
)
d
x
P=\int_{\mathcal{R}} p(\mathbf{x}) \mathrm{d} \mathbf{x}
P=∫Rp(x)dx
进一步地,我们可以求出
N
N
N 个数据,共
K
K
K 个数据落入
R
\mathcal{R}
R 中的概率:
Bin
(
K
∣
N
,
P
)
=
N
!
K
!
(
N
−
K
)
!
P
K
(
1
−
P
)
N
−
K
\operatorname{Bin}(K \mid N, P)=\frac{N !}{K !(N-K) !} P^{K}(1-P)^{N-K}
Bin(K∣N,P)=K!(N−K)!N!PK(1−P)N−K
不难发现这是一个二项分布,因此可以直接求出下述期望与方差:
E
[
K
/
N
]
=
E
[
K
]
N
=
P
var
[
K
/
N
]
=
var
[
K
]
N
2
=
P
(
1
−
P
)
/
N
\begin{aligned} & \mathbb{E}[K / N]=\displaystyle\frac{\mathbb{E}[K]}{N}=P\\ & \operatorname{var}[K / N]=\displaystyle\frac{\operatorname{var}[K]}{N^2}=P(1-P) / N \end{aligned}
E[K/N]=NE[K]=Pvar[K/N]=N2var[K]=P(1−P)/N
继续观察可以发现,当
N
→
∞
N\rightarrow\infty
N→∞ 时,
K
/
N
K/N
K/N 的分布将聚集在均值
P
P
P 的附近,即:
K
≃
N
P
K \simeq N P
K≃NP
另外,如果区域
R
\mathcal{R}
R 足够小,则假设该区域内概率处处相等,令其容量为
V
V
V,得到:
P
≃
p
(
x
)
V
P \simeq p(\mathbf{x}) V
P≃p(x)V
最终,依据两条看似有些矛盾的假设,得到
p
(
x
)
p(x)
p(x) 表达式:
- 【假设 1】 N N N 足够大
- 【假设 2】
R
\mathcal{R}
R 容量足够小,即
V
V
V 足够小
p ( x ) = K N V p(\mathbf{x})=\frac{K}{N V} p(x)=NVK
根据这两条假设,可以衍生出两种非参数估计方法:
- 【方法 1】固定 K K K,根据数据确定 V V V,即 K K K 近邻方法 (K-nearest-neighbour technique)
- 【方法 2】固定 V V V,根据数据确定 K K K,即核方法 (kernel approach)
It can be shown that both the K-nearest-neighbour density estimator and the kernel density estimator converge to the true probability density in the limit N → ∞ N\rightarrow \infty N→∞ provided V V V shrinks suitably with N N N, and K K K grows with N N N (Duda and Hart, 1973).
本小节主要介绍核方法,即固定
V
V
V,由数据确定
K
K
K,因此首先定义
k
(
u
)
k(\mathbf{u})
k(u) 函数 (kernel function / Parzen window),用于
K
K
K 的确定:
k
(
u
)
=
{
1
,
∣
u
i
∣
⩽
1
/
2
,
i
=
1
,
…
,
D
0
,
otherwise
k(\mathbf{u})=\left\{\begin{array}{ll} 1, & \left|u_{i}\right| \leqslant 1 / 2, \quad i=1, \ldots, D \\ 0, & \text { otherwise } \end{array}\right.
k(u)={1,0,∣ui∣⩽1/2,i=1,…,D otherwise
进一步地,我们可以根据直接确定
K
K
K:
K
=
∑
n
=
1
N
k
(
x
−
x
n
h
)
K=\sum_{n=1}^{N} k\left(\frac{\mathbf{x}-\mathbf{x}_{n}}{h}\right)
K=n=1∑Nk(hx−xn)
上述式子思想较简单,即根据有多少个数据落在了以
x
\mathbf{x}
x 为中心的各维度为
h
h
h 的超几何体内,来确定
K
K
K。其中
h
D
h^D
hD 可以当作
V
V
V 的度量,因此
p
(
x
)
p(\mathbf{x})
p(x) 如下:
p
(
x
)
=
1
N
∑
n
=
1
N
1
h
D
k
(
x
−
x
n
h
)
p(\mathbf{x})=\frac{1}{N} \sum_{n=1}^{N} \frac{1}{h^{D}} k\left(\frac{\mathbf{x}-\mathbf{x}_{n}}{h}\right)
p(x)=N1n=1∑NhD1k(hx−xn)
但这种直接的度量方法有一个缺点,即在每个超几何体的边界仍然是不连续的,因此可以将 kernel function 改为高斯分布,得到如下模型:
p
(
x
)
=
1
N
∑
n
=
1
N
1
(
2
π
h
2
)
1
/
2
exp
{
−
∥
x
−
x
n
∥
2
2
h
2
}
p(\mathbf{x})=\frac{1}{N} \sum_{n=1}^{N} \frac{1}{\left(2 \pi h^{2}\right)^{1 / 2}} \exp \left\{-\frac{\left\|\mathbf{x}-\mathbf{x}_{n}\right\|^{2}}{2 h^{2}}\right\}
p(x)=N1n=1∑N(2πh2)1/21exp{−2h2∥x−xn∥2}
上述模型可以理解为每一个数据点都代表了一个高斯分布,新数据点的概率等于,其它所有数据点对应的高斯分布对其概率的贡献均值。下图中蓝色为高斯 kernel function 的拟合结果,绿色为真实分布。
当然 kernel function 也可以进行修改,只需满足下述两个条件即可:
k
(
u
)
⩾
0
∫
k
(
u
)
d
u
=
1
\begin{array}{r} k(\mathbf{u}) \geqslant 0 \\ \int k(\mathbf{u}) \mathrm{d} \mathbf{u}=1 \end{array}
k(u)⩾0∫k(u)du=1
最后总结下该方法的优缺点:
- 【优点】无需训练集
- 【缺点】需要保存数据集的所有数据
2.5.2 最近邻法 (Nearest-neighbour methods)
固定 K,对于新的数据点不断扩大 V 使得包含 K 个数据,这种方法被称为 K 近邻 (K nearest neighbours),该方法拟合效果如下:
另外需要注意,K 近邻模型不是真正的概率密度模型,因为在所有空间上的积分是发散的,不满足归一化的性质。
接下来考虑用 KNN 进行分类,假设共
N
N
N 个数据,其中类别为
C
k
C_k
Ck 的数据有
N
k
N_k
Nk 个。先对于一个新数据,其包含
K
K
K 个数据的空间大小为
V
V
V,且其中
C
k
C_k
Ck 类别的数据共
K
k
K_k
Kk 个,则可根据贝叶斯进行类别估计:
p
(
x
∣
C
k
)
=
K
k
N
k
V
p
(
x
)
=
K
N
V
p
(
C
k
)
=
N
k
N
p
(
C
k
∣
x
)
=
p
(
x
∣
C
k
)
p
(
C
k
)
p
(
x
)
=
K
k
K
\begin{aligned} & p\left(\mathbf{x} \mid \mathcal{C}_{k}\right)=\frac{K_{k}}{N_{k} V}\\ & p(\mathbf{x})=\frac{K}{N V}\\ & p\left(\mathcal{C}_{k}\right)=\frac{N_{k}}{N}\\ & p\left(\mathcal{C}_{k} \mid \mathbf{x}\right)=\frac{p\left(\mathbf{x} \mid \mathcal{C}_{k}\right) p\left(\mathcal{C}_{k}\right)}{p(\mathbf{x})}=\frac{K_{k}}{K} \end{aligned}
p(x∣Ck)=NkVKkp(x)=NVKp(Ck)=NNkp(Ck∣x)=p(x)p(x∣Ck)p(Ck)=KKk
为了让分类错误率最低,该数据点的类别由最近的 K 个数据中,最多出现的类别决定。应用该模型,可以得到如下的示例结果图:
当 K = 1 时,且
N
→
∞
N\rightarrow\infty
N→∞,K 近邻方法的错误率不高于最优分类器错误率的两倍,证明如下:
K 近邻方法与 kernel 方法一样,当平滑参数较小时,会拟合噪音;太大时,则可能丢失结构。另外这两种方法都需要固定 K 或 V,但在实际数据分布中,数据空间中的不同区域,需要不同的 K、V 选择才能实现更准确的估计,而这也是非参数估计的一种局限。