大家好,又见面了,我是你们的朋友全栈君。
因子分析是一种数据简化技术,是一种数据的降维方法。
因子分子可以从原始高维数据中,挖掘出仍然能表现众多原始变量主要信息的低维数据。此低维数据可以通过高斯分布、线性变换、误差扰动生成原始数据。
因子分析基于一种概率模型,使用EM算法来估计参数。
主成分分析(PCA)也是一种特征降维的方法。
学习理论中,特征选择是要剔除与标签无关的特征,比如“汽车的颜色”与“汽车的速度”无关;
PCA中要处理与标签有关、但是存在噪声或者冗余的特征,比如在一个汽车样本中,“千米/小时”与“英里/小时”中有一个冗余了。
PCA的方法比较直接,只要计算特征向量就可以降维了。
独立成分分析(ICA)是一种主元分解的方法。
其基本思想是从一组混合的观测信号中分离出独立信号。比如在一个大房间里,很多人同时在说话,样本是这个房间里各个位置的一段录音,ICA可以从这些混合的录音中分离出每个人独立的说话的声音。
ICA认为观测信号是若干个统计独立的分量的线性组合,ICA要做的是一个解混过程。
因为因子分析、PCA、ICA都是对数据的处理方法,就放在这同一份总结里了。
1、因子分析(Factor analysis)
1.1、因子分析的直观理解
因子分析认为高维样本点实际上是由低维样本点经过高斯分布、线性变换、误差扰动生成的。让我们来看一个简单例子,对低维数据如何生成高维数据有一个直观理解。
假设我们有m=5个2维原始样本点如下:
图一
那么按照因子分析的做法,原始数据可以由以下过程生成:
①在一个低维空间(此处是1维)中,存在着由高斯分布生成的 m 个点
z(i)
图二
②使用某个
Λ=(a,bT)
将1维的
z
映射到2维的空间中:
图三
③加上
μ(μ1,μ2)T
,让直线过
μ
——实际上是将样本点横坐标加
μ1
,纵坐标加
μ2
:
图四
④对直线上的点做一定的扰动,其扰动为
ε
~
N(0,ψ)
:
图五
黑点就是图一中的原始数据。
1.2、因子分析的一般过程
因子分析认为m个n维特征的训练样例 (x(1),x(2),⋯,x(m)) 的产生过程如下:
①在一个 k 维空间中,按照多元高斯分布生成m个
z(i)
k<n
z(i) ~ N(0,I)
②存在一个变换矩阵 Λ∈Rn∗k ,将 z(i) 映射到 n 维空间中,即
Λz(i)
③将 Λz(i) ( n 维)加上一个均值
μ
μ+Λz(i)
④对每个点加上符合多元高斯分布的扰动 ε ~ N(0,ψ) ( n 维向量),即
x(i)=μ+Λz(i)+ε
1.3、因子分析模型
模型与参数概述
由上面的分析,我们定义因子分析的模型为:
z
~
N(0,I)
ε
~
N(0,ψ)
x=μ+Λz+ε(1)
其中
z
和
ε
是相互独立的。并且由上面的分析过程,我们可以直观地感受到我们的
参数是 μ∈Rn 、 Λ∈Rn∗k 、 ψ∈Rn∗n 。
另一个等价的假设是, (x,z) 联合分布如下,其中 z∈Rk 是一个隐藏随机变量:
x∣z
~
N(μ+Λz,ψ)
(2)
这个假设会在使用EM算法求解因子分析参数,E步中迭代
Q
分布的时候用到。
接下来的课程,是使用高斯模型的矩阵表示法来对模型进行分析。矩阵表示法认为
z
[zx]
~
N(μzx,Σ)
多元高斯分布的原始模型是:
其中
x
是
k
维向量,
μ
是
k
维向量,
Σ
是
k∗k
协方差矩阵。
很明显在多元高斯分布模型下,参数是
μzx,Σ
——它们是由
x,z
的联合分布生成的,所以我们可以用我们的原始参数
μ,Λ,ψ
来表示
μzx,Σ
,求得
x
的边缘分布,再把相关参数带入式(3),这就得到了关于我们参数的概率分布,然后就可以通过最大似然估计来求取我们的参数。
求取
μzx
μzx 是 x,z 联合分布的期望值(期望的定义:所有结果*相应概率的总和):
μzx=E[zx]=[E(z)E(x)](4)
由 z ~
N(0,I)
我们可以简单获得 E(z)=0 。
类似地由 ε ~ N(0,ψ) , x=μ+Λz+ε , μ 是一个常数,我们有:
E[x]=E[μ+Λz+ε]=E[μ]+ΛE[z]+E[ε]=μ+0+0=μ(5)
所以:
μzx=[0⃗ μ](6)
求取 Σ
Σ 是 x,z 联合分布的协方差矩阵。
方差,度量随机变量与期望之间的偏离程度,定义如下:
Var(X)=E((X−E(X))2)=E(X2)−(E(X)2)(7)
协方差,两个变量总体误差的期望,定义如下:
Cov(X,Y)=E((X−E(X))(Y−E(Y)))(8)
协方差、方差、期望之间的一些相互关系如下:
Cov(X,X)=Cov(X)=Var(X)=E(XXT)=σ2(9)
下面开始求取 Σ 。
Σ=Cov[zx]=[ΣzzΣxzΣzxΣxx]=E[(z−E(z))(z−E(z))T(x−E(x))(z−E(z))T(z−E(z))(x−E(x))T(x−E(x))(x−E(x))T](10)
由 z ~
N(0,I)
,可以简单得到:
Σzz=Cov(z)=σ2=I(11)
由 ε ~ N(0,ψ) , x=μ+Λz+ε , E(x)=μ ,并且 z 和
ε
是相互独立,有:
Σzx=E[(z−E(z))(x−E(x))T]=E[(z−0)(μ+Λz+ε−μ)T]=E[zzT]ΛT+E[zεT]=IΛT+0=ΛT(12)
类似地,我们可以得到:
Σxx=E[(x−E(x))(x−E(x))T]=E[(μ+Λz+ε−μ)(μ+Λz+ε−μ)T]=ΛE[zzT]ΛT+E[εεT]=ΛIΛT+ψ=ΛΛT+ψ(13)
用最大似然估计法求解参数
经过上面的步骤,我们就把 μzx,Σ 用我们的参数 μ,Λ,ψ 表示出来了:
[zx]
~
N(μzx,Σ)
~
N([0⃗ μ],[IΛΛTΛΛT+ψ])
然后我们可以求得 x 的边缘分布:
x
~ N(μ,ΛΛT+ψ)
因此,给定一个训练集 {
x(i);i=1,2,⋯,m} ,把参数带入式(3),我们可以写出下面的似然函数:
l(μ,Λ,ψ)=log∏i=1m12πn∣∣ΛΛT+ψ∣∣−−−−−−−−−−−√exp(−12(x(i)−μ)T(ΛΛT+ψ)−1(x(i)−μ))(14)
对此似然函数做最大似然估计,就能求得我们的参数。
1.4、因子分析的EM求解
可以感受到,直接对这个似然函数求解是很困难的,在这个情况下,用EM算法就登场了——当一个似然函数难以直接求解其最大值的时候,可以通过EM算法不断建立下界、最大化下界的方式不断逼近该似然函数真实的最大值,当EM算法收敛,我们就认为已经求得了此最大值。
E-step
对于EM算法的E-step,我们有:
Qi(z(i)):=p(z(i)∣x(i);μ,Λ,ψ)(15)
进一步地:
Qi(z(i))=12πk∣∣Σz(i)∣x(i)∣∣−−−−−−−−−−√exp(−12(z(i)−μz(i)∣x(i))TΣ−1z(i)∣x(i)(z(i)−μz(i)∣x(i)))(16)
其中:
μz(i)∣x(i)Σz(i)∣x(i)=ΛT(ΛΛT+ψ)−1(x(i)−μ)=I−ΛT(ΛΛT+ψ)−1Λ(17)
μz(i)∣x(i),Σz(i)∣x(i) 是讲义与课上直接给出的,这里也不进行推导。
M-step
在M-step中,我们需要最大化如下公式来求取参数 μ,Λ,ψ :
∑i=1m∫z(i)Qi(z(i))logp(x(i),z(i);μ,Λ,ψ)Qi(z(i))dz(i)(18)
视为期望,打开log
在这里,因为 z 是连续的,所以使用积分;如果是离散的,则使用累加。
并且,积分部分可以当成
z
服从 Q 分布时,函数
logp(x(i),z(i);μ,Λ,ψ)Qi(z(i))
的期望,这里将会用E表示,省略 z(i) ~ Qi 的下标;对于函数中 x,z 的联合分布,我们可以用贝叶斯公式把它打开 p(x,z)=p(x∣x)p(z) ;为了方便计算我们还要把log函数打开——经过这些分析,我们有如下推导:
∑i=1m∫z(i)Qi(z(i))logp(x(i),z(i);μ,Λ,ψ)Qi(z(i))dz(i)=∑i=1mE[logp(x(i),z(i);μ,Λ,ψ)Qi(z(i))]=∑i=1mE[logp(x(i)∣z(i);μ,Λ,ψ)p(z(i))Qi(z(i))]=∑i=1mE[logp(x(i)∣z(i);μ,Λ,ψ)+logp(z(i))−logQi(z(i))](19)
去掉无关项后带入具体分布
这就比较清爽了,然后,记住我们的目标是求得参数 μ,Λ,ψ ,但是它们不能一起求解,所以下面以参数 Λ 为例,对公式进行求解——在式(19)中,对参数 Λ 求偏导。另外式(19)中的 p(z(i) 与 Qi(z(i)) 与 Λ 无关,可以忽略掉,所以实际上就是对下式求偏导:
∑i=1mE[logp(x(i)∣z(i);μ,Λ,ψ)](20)
在对式(20)求偏导之前,还可以对其进行一些处理——由式(2),并且 x∣z 服从多元高斯分布,所以有:
∑i=1mE[logp(x(i)∣z(i);μ,Λ,ψ)]=∑i=1mE[log12πn|ψ|−−−−−−√exp(−12(x(i)−(μ+Λz(i)))Tψ−1(x(i)−(μ+Λz(i))))]=∑i=1mE[−12log|ψ|−n2log(2π)−12(x(i)−μ−Λz(i))Tψ−1(x(i)−μ−Λz(i)))](21)
去掉无关项后求偏导
同样地,我们的目标是与 Λ 有关的项,所以忽略掉前面的无关项之后,我们实际上是对下式求偏导并求解:
∇Λ∑i=1mE[−12(x(i)−μ−Λz(i))Tψ−1(x(i)−μ−Λz(i)))]=∑i=1m∇Λ−E[12((x(i)T−μT)ψ−1A−z(i)TΛTψ−1B)((x(i)−μ)C−Λz(i)D)]=∑i=1m∇Λ−E[12(AC−AD−BC+BD)](22)
打开后我们可以发现, AC 这一项是与 Λ 无关的,把这一项忽略掉,所以由式(22)继续推导有:
∑i=1m∇Λ−E[12(−AD−BC+BD)]=∑i=1m∇Λ−E[12(−(x(i)T−μT)ψ−1Λz(i)E−z(i)TΛTψ−1(x(i)−μ)F+z(i)TΛTψ−1Λz(i))](23)
因为期望是一个常数,又因为 a=tr(a) ,所以可以直接对上式求迹;
因为 trA=trAT ,可以对E求转置,又因为对角矩阵的转置是它本身—— (ψ−1)T=ψ−1 ,所以有 trE=trET=trF ,对式(23)继续推导有:
∑i=1m∇Λ−E[12(−(x(i)T−μT)ψ−1Λz(i)E−z(i)TΛTψ−1(x(i)−μ)F+z(i)TΛTψ−1Λz(i))]=∑i=1m∇Λ−E[tr12(−z(i)TΛTψ−1(x(i)−μ)ET−z(i)TΛTψ−1(x(i)−μ)F+z(i)TΛTψ−1Λz(i))]=∑i=1m∇ΛE[−tr12z(i)TΛTψ−1Λz(i)+trz(i)TΛTψ−1(x(i)−μ)](24)
然后利用 trAB=trBA ,把式(25)中的 z(i)T 放到它们自己的后面,再把求导切换到期望里面——求导是针对 Λ ,期望是针对 z(i) ,所以是可以切换的:
∑i=1m∇ΛE[−tr12z(i)TΛTψ−1Λz(i)+trz(i)TΛTψ−1(x(i)−μ)]=∑i=1m∇ΛE[−tr12ΛTψ−1Λz(i)z(i)T+trΛTψ−1(x(i)−μ)z(i)T]=∑i=1m(E[−∇Λtr12ΛTψ−1Λz(i)z(i)T]+E[∇ΛtrΛTψ−1(x(i)−μ)z(i)T])(25)
对于式(25),先用矩阵的迹的性质 ∇ATf(A)=(∇Af(A))T 处理一下:
∑i=1m(E[−(∇ΛTtr12ΛTψ−1Λz(i)z(i)T)T]+E[(∇ΛTtrΛTψ−1(x(i)−μ)z(i)T)T])(26)
对式(26)的第一项使用 ∇AtrABATC=CAB+CTABT 的性质,对第二项使用 ∇AtrAB=BT 的性质:
∑i=1m⎛⎝⎜⎜E⎡⎣⎢⎢−⎛⎝⎜∇ΛTAtr12ΛTAψ−1BΛATz(i)z(i)TC⎞⎠⎟T⎤⎦⎥⎥+E⎡⎣⎢⎢⎛⎝⎜∇ΛTAtrΛTAψ−1(x(i)−μ)z(i)TB⎞⎠⎟T⎤⎦⎥⎥⎞⎠⎟⎟=∑i=1mE((−12z(i)z(i)TΛTψ−1−12(z(i)z(i)T)TΛT(ψ−1)T)T+((ψ−1(x(i)−μ)z(i)T)T)T)=∑i=1mE((−z(i)z(i)TΛTψ−1)T+ψ−1(x(i)−μ)z(i)T)=∑i=1mE(−ψ−1Λz(i)z(i)T+ψ−1(x(i)−μ)z(i)T)(27)
令式(27)=0,并化简,就可以求得参数 Λ :
⟹⟹⟹∑i=1mE(−ψ−1Λz(i)z(i)T+ψ−1(x(i)−μ)z(i)T)=0∑i=1m−ψ−1ΛE[z(i)z(i)T]+∑i=1m−ψ−1(x(i)−μ)E[z(i)T]=0∑i=1mΛE[z(i)z(i)T]=∑i=1m(x(i)−μ)E[z(i)T]Λ=(∑i=1m(x(i)−μ)E[z(i)T])(∑i=1mE[z(i)z(i)T])−1(28)
我们发现,这里的公式与线性回归中最小二乘法的矩阵形式相似。
相似原因:在因子分析中, x 是
z
的线性函数,在E-step中给出 z 的
Q
分布之后,在M-tep中寻找 x 与
z
的映射关系 Λ ;在线性回归的最小二乘中,也是寻找 x 与
y
的线性关系。
不同之处:最小二乘只用到了 z 的最优估计,因子分析还用到了
z(i)z(i)T
的估计。
对于参数 Λ ,这里还有未知数 E[z(i)T] 与 E[z(i)z(i)T] ,并且此处的期望是在 z(i) 服从 Qi 前提下计算的,所以对于前者,通过式(17)我们有:
E[z(i)T]=μTz(i)∣x(i)(29)
对于后者,由式(7)~式(9)方差与协方差的性质,我们有:
Cov(z)⟹E[z(i)z(i)T]=E(zzT)−E(z)E(zT)=E(z(i))E(z(i)T)+Cov(z)=μz(i)∣x(i)μTz(i)∣x(i)+Σz(i)∣x(i)(30)
注意这里的 E[z(i)z(i)T] 不仅仅等于 E(z(i))E(z(i)T) ,后面还有加上一个后验概率 p(z∣x) 协方差,要特别注意。
到这里,我们就可以把参数 Λ 的最终形式给出来了:
Λ=(∑i=1m(x(i)−μ)μTz(i)∣x(i))(∑i=1mμz(i)∣x(i)μTz(i)∣x(i)+Σz(i)∣x(i))−1(31)
另外对于其他两个参数 Λ,ψ ,使用相同的方法可以求得,这里直接给出结果
:
μ=1m∑i=1mx(i)(32)
ϕ=1m∑i=1mx(i)x(i)T−x(i)μTz(i)∣x(i)ΛT−Λμz(i)∣x(i)x(i)T+Λμz(i)∣x(i)μTz(i)∣x(i)+Σz(i)∣x(i)ΛT(33)
注意这里的参数是 ϕ 不是 ψ ,得到 ϕ 之后还需要将 ψii=ϕii ,因为 ϕ 不是对角矩阵,所以只需要取 ϕ 对角线上的值即可。
2、主成分分析(Principal component analysis,简称PCA
PCA的意义
PCA技术的一大好处是对数据进行降维的处理。我们可以对新求出的“主元”向量的重要性进行排序,根据需要取前面最重要的部分,将后面的维数省去,可以达到降维从而简化模型或是对数据进行压缩的效果。同时最大程度的保持了原有数据的信息。
PCA将n个特征降维到k个,可以用来进行数据压缩,如果100维的向量最后可以用10维来表示,那么压缩率为90%。同样图像处理领域的KL变换使用PCA做图像压缩。但PCA要保证降维后,还要保证数据的特性损失最小。
预处理
运行PCA算法之前,数据一般要进行预处理。预处理步骤如下:
①令 μ=1m∑mi=1x(i)
②用 (x(i)−μ) 代替 x(i)
③令 σ2=1m∑mi=1(x(i)j)2
④用 x(i)j/σ2 代替 x(i)j
步骤①与②将数据的均值变为0;
步骤③与④将数据每个维度的方差变为1,使每个维度都在同一个维度下被度量。
如果已知数据均值为0,或者数据已在同样一个维度下,就无需进行上面的步骤。
最大方差理论解释PCA
Ng在课上说,PCA有9到10种解释方法,这里仅用其中的最大方差理论来解释PCA
在信号处理中认为信号具有较大的方差,噪声有较小的方差,信噪比就是信号与噪声的方差比,越大越好。如前面的图,样本在横轴上的投影方差较大,在纵轴上的投影方差较小,那么认为纵轴上的投影是由噪声引起的。
因此我们认为,最好的k维特征是将n维样本点转换为k维后,每一维上的样本方差都很大。
下面举例来说明如何寻找主方向。
图六
图六左图是经过预处理的样本点,均值为0,特征方差归一;
图六的中图和右图是将样本投影到某个维度上的表示, 该维度用一条过原点的直线表示,前处理的过程实质是将原点移到样本点的中心点。
图七
图七中,原点到样本在直线上的投影的距离 x(i)Tu 既是方差;样本点到直线上的距离是平方误差。
我们可以直观感受到,图六中间的图方差和是最大的,平方误差是最小的。下面给出用最大方差理论寻找该最佳方向的公式定义。
形式化
设 x(i) 是数据集中的点, u 是要求解的单位向量,那么方差最大化可以形式化为最大化下式:
1m∑i=1m(x(i)Tu)2=1m∑i=1muTx(i)x(i)Tu=uT(1m∑i=1mx(i)x(i)T)u(34)
因为归一化后的数据,投影值的均值也为0,所以在方差计算中直接平方。
同时这个式子还有一个约束条件,即 ∥u∥2=1 。熟悉的最大化某个带约束的函数,引入拉格朗日乘子来求解:
l=uT(1m∑i=1mx(i)x(i)T)u−λ(∥u∥2−1)=uTΣu−λ(uTu−1)(35)
对 u 求导:
∇ul=∇u(uTΣu−λ(uTu−1))=∇utr(uTΣu)−λ∇utr(uTu)=(∇uTtr(uTΣu))T−(λ∇uTtr(uTu))T=((Σu)T)T−(λuT)T=Σu−λu(36)
这里主要运用了 ∇ATf(A)=(∇Af(A))T 与 ∇AtrAB=BT 的性质,处理方式与上面因子分析中式(26)与(27)的方式类似,不再赘述。
令导数为0,可知 u 是
Σ
的特征向量—— (A−λ)x=0 看起来会眼熟一些。
因为 Σ=1m∑mi=1x(i)x(i)T 是对称矩阵,因而可得到 n 个相正交的特征向量:
⎡⎣⎢⎢⎢⎢u1u2⋮un⎤⎦⎥⎥⎥⎥
此时如何达到降维的效果?选取 最大的 k 个(
k<n
)特征值对应的特征向量即可。 u∈Rn , x(i)∈Rn ,降维后的数据为:
y(i)=⎡⎣⎢⎢⎢⎢⎢⎢uT1x(i)uT2x(i)⋮uTkx(i)⎤⎦⎥⎥⎥⎥⎥⎥∈Rk(37)
这样,对于 n 维的原始数据
x(i)
,我们有了更低维的 k 维
y(i)
去表示它。
以上就是PCA的数学表达式。
3、独立成分分析(Independent component analysis,简称ICA)
解决的问题
ICA解决的是原始数据分解的问题,更具体一些的例子就是将声音信号分离的问题。
经典的鸡尾酒宴会问题(cocktail party problem)。假设在party中有n个人,他们可以同时说话,我们也在房间中一些角落里共放置了个声音接收器用来记录声音。宴会过后,我们从 n 个麦克风中得到了一组数据
{x(i);i=1,2,⋯,m}
,i表示采样的时间顺序,也就是说共得到了m组采样,每一组采样都是n维的。我们的目标是单单从这m组采样数据中分辨出每个人说话的信号。
同时我们假设 n 个信号源为
s=(s1,s2,⋯,sn)T
, s∈Rn ,每一维都是一个人的声音信号,每个人发出的声音信号独立。A是一个未知的混合矩阵(mixing matrix),用来组合叠加信号s,那么:
x=⎡⎣⎢|x(1)||x(2)|⋯|x(m)|⎤⎦⎥=⎡⎣⎢|As(1)||As(2)|⋯|As(m)|⎤⎦⎥=A⎡⎣⎢⎢⎢⎢⎢⎢s(1)1s(1)2⋮s(1)ns(2)1s(2)2⋮s(2)n⋯⋯⋱⋯s(m)1s(m)2⋮s(m)n⎤⎦⎥⎥⎥⎥⎥⎥=As(38)
x 不是一个向量,是一个矩阵。其中每个列向量是
x(i)
, x(i)=As(i) 。 A 和
s
都是未知的, x 是观测收集到的声音信号,是已知的,我们要想办法从
x
推出 s 。这个过程也称为盲信号分离。
令
W=A−1
,有 s(i)=A−1x(i)=Wx(i) ,所以 W 可以表示为:
W=⎡⎣⎢⎢⎢⎢⎢———WT1WT2⋮WTn———⎤⎦⎥⎥⎥⎥⎥(39)
其中 Wi∈Rn ,其实就是将 Wi 写成行向量形式。那么得到: s(i)j=Wx(i)j 。
ICA不确定性(Ambiguities)
下面是两个原信号不确定性:
第一个,由于 w 和
s
都不确定,那么在没有先验知识的情况下,无法同时确定这两个相关参数。比如上面的公式 s=wx 。当 w 扩大两倍时,
s
只需要同时扩大两倍即可,等式仍然满足,因此无法得到唯一的 s 。
第二个,如果将
s
的顺序打乱,变成另外一个顺序,那么只需要调换 A 的列向量顺序即可,因此也无法单独确定
s
。
另外,还有一种ICA不适用的情况,那就是信号不能是高斯分布的。当源信号是高斯分布的时候,可以由不同的混合矩阵 A 乘上
s
,得到相同分布的 x ,一样无法确定源信号。
密度函数和线性变换
假设我们的随机变量
s
有概率密度函数 ps(s) (连续值是概率密度函数,离散值是概率)。为了简单,我们再假设 s 是实数,还有一个随机变量
x=As
, A 和
x
都是实数。令 px(x) 是x的概率密度,那么怎么求 px(x) ?
求解之前先看一个密度函数与分布函数之间的关系。
设一个概率密度函数(Probability density function,简称pdf)为 f(x) ,一个累积分布函数(Cumulative distribution function,简称cdf)为 F(x) ,它们之间的关系为:
F(x)f(x)=∫xf(x)dx=F′(x)(40)
以正态分布为例,它的概率密度函数与累积分布函数的图像如下图:
图八
一个很直观的感受就是,后者是由前者积分得到的,这也是符合式(40)的公式表示的。
关于概率密度的直接定义是:
Fx(a)=P(X≤a)=∫a−∞f(x)dx(41)
关于 px(x) 的推导如下:
Fx(x)px(x)=P(X≤x)=P(As≤x)=P(s≤Wx)=Fs(Wx)=F′x(x)=F′s(Wx)=ps(Wx)|W|(42)
ICA算法
这里使用最大似然估计来解释算法, 我们假定每个 si 有概率密度 ps ,那么给定时刻原信号的联合分布就是:
p(s)=∏i=1nps(si)(43)
这个式子有一个假设前提:每个人发出的声音信号各自独立。由式(42),我们有:
p(x)=ps(Wx)|W|=|W|∏i=1nps(WTix)(44)
前面提到过,如果没有先验知识,我们无法求得W和s。因此我们需要知道 ps(sj) ,我们打算选取一个概率密度函数赋给s,但是我们不能选取高斯分布的密度函数。在概率论里我们知道密度函数 p(x) 由累计分布函数(cdf) F(x) 求导得到。 F(x) 要满足两个性质是:单调递增和在[0,1]。我们发现sigmoid函数很适合,定义域负无穷到正无穷,值域0到1,缓慢递增。我们假定s的累积分布函数符合sigmoid函数:
g(s)=11+e−s(45)
求导后有:
ps(s)=g′(s)=es(1+es)2(46)
知道 ps(s) 之后,就剩下参数 W 了,给定训练样本
{x(i);i=1,2,⋯,m}
后,样本对数似然估计如下:
l(W)=log∏i=1mp(x(i))=∑i=1m⎛⎝∑j=1nlogg′(WTix(i))+log|W|⎞⎠(47)
接下来就是对W求导了,这里牵涉一个问题是对行列式|W|进行求导的方法,属于矩阵微积分,这里直接给出结果:
∇W|W|=|W|(W−1)T(48)
W 最后的迭代公式如下,
logg′(s)
的导数为 1−2g(s) , α 是梯度上升速率,人为指定。:
W:=W+α⎛⎝⎜⎜⎜⎜⎜⎡⎣⎢⎢⎢⎢⎢1−2g(WT1x(i))1−2g(WT2x(i))⋮1−2g(WTnx(i))⎤⎦⎥⎥⎥⎥⎥x(i)T+(WT)−1⎞⎠⎟⎟⎟⎟⎟(49)
4、后记
从11月初到3月底学习这份课程,再从4月到现在5月初把这6份总结写完,前后也是用了半年时间——这都是在业余时间完成的,鬼知道我经历了什么。最后写到这里,竟然是不知道说什么了。清明三天早10晚11把第一份总结写完,在假期最后一天9点多躺在床上的时候,没有什么快感,只有一种身体和脑子被掏空的感觉,现在也是。而且一直压制着想要直接做项目的冲动,我知道一旦直接去操作了,我应该就没有办法写下这样的东西了。
想要学习决策树,Boosting,神经网络,深度学习,还有这个课程中缺失的最后一部分关于强化学习的完整内容,还有其他的很多很多东西,到时候会边学边写,一个个来,像这次积累了这么多东西然后一次写完真是太难过了。但是回头一看,这些东西也还是太少了。接下来主要是做点项目,没学到又要用到的算法和知识到时候先求会用。
时间是最大的敌人,世界在变好,弱是原罪。
告一段落。
μzx
N(0,I)
类似地由 ε ~ N(0,ψ) , x=μ+Λz+ε , μ 是一个常数,我们有:
所以:
方差,度量随机变量与期望之间的偏离程度,定义如下:
协方差,两个变量总体误差的期望,定义如下:
协方差、方差、期望之间的一些相互关系如下:
N(0,I)
由 ε ~ N(0,ψ) , x=μ+Λz+ε , E(x)=μ ,并且 z 和
ε
类似地,我们可以得到:
[zx]
~
N(μzx,Σ)
~
N([0⃗ μ],[IΛΛTΛΛT+ψ])
x
因此,给定一个训练集 {
x(i);i=1,2,⋯,m} ,把参数带入式(3),我们可以写出下面的似然函数:
对此似然函数做最大似然估计,就能求得我们的参数。
1.4、因子分析的EM求解
可以感受到,直接对这个似然函数求解是很困难的,在这个情况下,用EM算法就登场了——当一个似然函数难以直接求解其最大值的时候,可以通过EM算法不断建立下界、最大化下界的方式不断逼近该似然函数真实的最大值,当EM算法收敛,我们就认为已经求得了此最大值。
E-step
对于EM算法的E-step,我们有:
进一步地:
其中:
μz(i)∣x(i),Σz(i)∣x(i) 是讲义与课上直接给出的,这里也不进行推导。
M-step
在M-step中,我们需要最大化如下公式来求取参数 μ,Λ,ψ :
视为期望,打开log
在这里,因为 z 是连续的,所以使用积分;如果是离散的,则使用累加。
并且,积分部分可以当成
z
logp(x(i),z(i);μ,Λ,ψ)Qi(z(i))
去掉无关项后带入具体分布
这就比较清爽了,然后,记住我们的目标是求得参数 μ,Λ,ψ ,但是它们不能一起求解,所以下面以参数 Λ 为例,对公式进行求解——在式(19)中,对参数 Λ 求偏导。另外式(19)中的 p(z(i) 与 Qi(z(i)) 与 Λ 无关,可以忽略掉,所以实际上就是对下式求偏导:
在对式(20)求偏导之前,还可以对其进行一些处理——由式(2),并且 x∣z 服从多元高斯分布,所以有:
去掉无关项后求偏导
同样地,我们的目标是与 Λ 有关的项,所以忽略掉前面的无关项之后,我们实际上是对下式求偏导并求解:
打开后我们可以发现, AC 这一项是与 Λ 无关的,把这一项忽略掉,所以由式(22)继续推导有:
因为期望是一个常数,又因为 a=tr(a) ,所以可以直接对上式求迹;
因为 trA=trAT ,可以对E求转置,又因为对角矩阵的转置是它本身—— (ψ−1)T=ψ−1 ,所以有 trE=trET=trF ,对式(23)继续推导有:
然后利用 trAB=trBA ,把式(25)中的 z(i)T 放到它们自己的后面,再把求导切换到期望里面——求导是针对 Λ ,期望是针对 z(i) ,所以是可以切换的:
对于式(25),先用矩阵的迹的性质 ∇ATf(A)=(∇Af(A))T 处理一下:
对式(26)的第一项使用 ∇AtrABATC=CAB+CTABT 的性质,对第二项使用 ∇AtrAB=BT 的性质:
令式(27)=0,并化简,就可以求得参数 Λ :
我们发现,这里的公式与线性回归中最小二乘法的矩阵形式相似。
相似原因:在因子分析中, x 是
z
Q
z
y
不同之处:最小二乘只用到了 z 的最优估计,因子分析还用到了
z(i)z(i)T
对于参数 Λ ,这里还有未知数 E[z(i)T] 与 E[z(i)z(i)T] ,并且此处的期望是在 z(i) 服从 Qi 前提下计算的,所以对于前者,通过式(17)我们有:
对于后者,由式(7)~式(9)方差与协方差的性质,我们有:
注意这里的 E[z(i)z(i)T] 不仅仅等于 E(z(i))E(z(i)T) ,后面还有加上一个后验概率 p(z∣x) 协方差,要特别注意。
到这里,我们就可以把参数 Λ 的最终形式给出来了:
另外对于其他两个参数 Λ,ψ ,使用相同的方法可以求得,这里直接给出结果
:
注意这里的参数是 ϕ 不是 ψ ,得到 ϕ 之后还需要将 ψii=ϕii ,因为 ϕ 不是对角矩阵,所以只需要取 ϕ 对角线上的值即可。
2、主成分分析(Principal component analysis,简称PCA
PCA的意义
PCA技术的一大好处是对数据进行降维的处理。我们可以对新求出的“主元”向量的重要性进行排序,根据需要取前面最重要的部分,将后面的维数省去,可以达到降维从而简化模型或是对数据进行压缩的效果。同时最大程度的保持了原有数据的信息。
PCA将n个特征降维到k个,可以用来进行数据压缩,如果100维的向量最后可以用10维来表示,那么压缩率为90%。同样图像处理领域的KL变换使用PCA做图像压缩。但PCA要保证降维后,还要保证数据的特性损失最小。
预处理
运行PCA算法之前,数据一般要进行预处理。预处理步骤如下:
①令 μ=1m∑mi=1x(i)
②用 (x(i)−μ) 代替 x(i)
③令 σ2=1m∑mi=1(x(i)j)2
④用 x(i)j/σ2 代替 x(i)j
步骤①与②将数据的均值变为0;
步骤③与④将数据每个维度的方差变为1,使每个维度都在同一个维度下被度量。
如果已知数据均值为0,或者数据已在同样一个维度下,就无需进行上面的步骤。
最大方差理论解释PCA
Ng在课上说,PCA有9到10种解释方法,这里仅用其中的最大方差理论来解释PCA
在信号处理中认为信号具有较大的方差,噪声有较小的方差,信噪比就是信号与噪声的方差比,越大越好。如前面的图,样本在横轴上的投影方差较大,在纵轴上的投影方差较小,那么认为纵轴上的投影是由噪声引起的。
因此我们认为,最好的k维特征是将n维样本点转换为k维后,每一维上的样本方差都很大。
下面举例来说明如何寻找主方向。
图六
图六左图是经过预处理的样本点,均值为0,特征方差归一;
图六的中图和右图是将样本投影到某个维度上的表示, 该维度用一条过原点的直线表示,前处理的过程实质是将原点移到样本点的中心点。
图七
图七中,原点到样本在直线上的投影的距离 x(i)Tu 既是方差;样本点到直线上的距离是平方误差。
我们可以直观感受到,图六中间的图方差和是最大的,平方误差是最小的。下面给出用最大方差理论寻找该最佳方向的公式定义。
形式化
设 x(i) 是数据集中的点, u 是要求解的单位向量,那么方差最大化可以形式化为最大化下式:
1m∑i=1m(x(i)Tu)2=1m∑i=1muTx(i)x(i)Tu=uT(1m∑i=1mx(i)x(i)T)u(34)
因为归一化后的数据,投影值的均值也为0,所以在方差计算中直接平方。
同时这个式子还有一个约束条件,即 ∥u∥2=1 。熟悉的最大化某个带约束的函数,引入拉格朗日乘子来求解:
对 u 求导:
∇ul=∇u(uTΣu−λ(uTu−1))=∇utr(uTΣu)−λ∇utr(uTu)=(∇uTtr(uTΣu))T−(λ∇uTtr(uTu))T=((Σu)T)T−(λuT)T=Σu−λu(36)
这里主要运用了 ∇ATf(A)=(∇Af(A))T 与 ∇AtrAB=BT 的性质,处理方式与上面因子分析中式(26)与(27)的方式类似,不再赘述。
令导数为0,可知 u 是
Σ
因为 Σ=1m∑mi=1x(i)x(i)T 是对称矩阵,因而可得到 n 个相正交的特征向量:
⎡⎣⎢⎢⎢⎢u1u2⋮un⎤⎦⎥⎥⎥⎥
此时如何达到降维的效果?选取 最大的 k 个(
k<n
这样,对于 n 维的原始数据
x(i)
y(i)
以上就是PCA的数学表达式。
3、独立成分分析(Independent component analysis,简称ICA)
解决的问题
ICA解决的是原始数据分解的问题,更具体一些的例子就是将声音信号分离的问题。
经典的鸡尾酒宴会问题(cocktail party problem)。假设在party中有n个人,他们可以同时说话,我们也在房间中一些角落里共放置了个声音接收器用来记录声音。宴会过后,我们从 n 个麦克风中得到了一组数据
{x(i);i=1,2,⋯,m}
同时我们假设 n 个信号源为
s=(s1,s2,⋯,sn)T
x 不是一个向量,是一个矩阵。其中每个列向量是
x(i)
s
x
令
W=A−1
W=⎡⎣⎢⎢⎢⎢⎢———WT1WT2⋮WTn———⎤⎦⎥⎥⎥⎥⎥(39)
其中 Wi∈Rn ,其实就是将 Wi 写成行向量形式。那么得到: s(i)j=Wx(i)j 。
ICA不确定性(Ambiguities)
下面是两个原信号不确定性:
第一个,由于 w 和
s
s
第二个,如果将
s
s
另外,还有一种ICA不适用的情况,那就是信号不能是高斯分布的。当源信号是高斯分布的时候,可以由不同的混合矩阵 A 乘上
s
密度函数和线性变换
假设我们的随机变量
s
x=As
x
求解之前先看一个密度函数与分布函数之间的关系。
设一个概率密度函数(Probability density function,简称pdf)为 f(x) ,一个累积分布函数(Cumulative distribution function,简称cdf)为 F(x) ,它们之间的关系为:
以正态分布为例,它的概率密度函数与累积分布函数的图像如下图:
图八
一个很直观的感受就是,后者是由前者积分得到的,这也是符合式(40)的公式表示的。
关于概率密度的直接定义是:
关于 px(x) 的推导如下:
ICA算法
这里使用最大似然估计来解释算法, 我们假定每个 si 有概率密度 ps ,那么给定时刻原信号的联合分布就是:
这个式子有一个假设前提:每个人发出的声音信号各自独立。由式(42),我们有:
前面提到过,如果没有先验知识,我们无法求得W和s。因此我们需要知道 ps(sj) ,我们打算选取一个概率密度函数赋给s,但是我们不能选取高斯分布的密度函数。在概率论里我们知道密度函数 p(x) 由累计分布函数(cdf) F(x) 求导得到。 F(x) 要满足两个性质是:单调递增和在[0,1]。我们发现sigmoid函数很适合,定义域负无穷到正无穷,值域0到1,缓慢递增。我们假定s的累积分布函数符合sigmoid函数:
求导后有:
知道 ps(s) 之后,就剩下参数 W 了,给定训练样本
{x(i);i=1,2,⋯,m}
接下来就是对W求导了,这里牵涉一个问题是对行列式|W|进行求导的方法,属于矩阵微积分,这里直接给出结果:
W 最后的迭代公式如下,
logg′(s)
4、后记
从11月初到3月底学习这份课程,再从4月到现在5月初把这6份总结写完,前后也是用了半年时间——这都是在业余时间完成的,鬼知道我经历了什么。最后写到这里,竟然是不知道说什么了。清明三天早10晚11把第一份总结写完,在假期最后一天9点多躺在床上的时候,没有什么快感,只有一种身体和脑子被掏空的感觉,现在也是。而且一直压制着想要直接做项目的冲动,我知道一旦直接去操作了,我应该就没有办法写下这样的东西了。
想要学习决策树,Boosting,神经网络,深度学习,还有这个课程中缺失的最后一部分关于强化学习的完整内容,还有其他的很多很多东西,到时候会边学边写,一个个来,像这次积累了这么多东西然后一次写完真是太难过了。但是回头一看,这些东西也还是太少了。接下来主要是做点项目,没学到又要用到的算法和知识到时候先求会用。
时间是最大的敌人,世界在变好,弱是原罪。
告一段落。
发布者:全栈程序员-用户IM,转载请注明出处:https://javaforall.cn/145363.html原文链接:https://javaforall.cn
【正版授权,激活自己账号】: Jetbrains全家桶Ide使用,1年售后保障,每天仅需1毛
【官方授权 正版激活】: 官方授权 正版激活 支持Jetbrains家族下所有IDE 使用个人JB账号...