PLS系列003 单因变量线性PLS


单因变量线性偏最小二乘法

  • 1 单因变量线性PLS
    • 1.1 计算推导
    • 1.2 辅助分析技术
  • Reference

1 单因变量线性PLS 1.1 计算推导 由于在多因变量线性偏最小二乘法中,我们已经讨论了计算推导,在此,我们将但因变量进行简化计算推导过程:
①样本数据 X X X与 Y Y Y标准化预处理
②记 t 1 {{t}_{1}} t1?是 X X X的第1个成分有 t 1 = X w 1 {{t}_{1}}=X{{w}_{1}} t1?=Xw1?,其中 w 1 {{w}_{1}} w1?是 X X X的第1个轴(单位列向量即 ∥ w 1 ∥ = 1 \left\| {{w}_{1}} \right\|\text{=}1 ∥w1?∥=1)。
u 1 {{u}_{1}} u1?是 Y Y Y的第1个成分有 u 1 = Y v 1 {{u}_{1}}=Y{{v}_{1}} u1?=Yv1?,其中 v 1 {{v}_{1}} v1?是 X X X的第1个轴(单位列向量即 ∥ v 1 ∥ = 1 \left\| {{v}_{1}} \right\|\text{=}1 ∥v1?∥=1)。
t 1 {{t}_{1}} t1?、 u 1 {{u}_{1}} u1?为列向量,行数为 n n n,即正好是样本集合数。
w 1 {{w}_{1}} w1?为列向量,行数为 p p p,即正好是自变量个数
v 1 {{v}_{1}} v1?为列向量,行数为 q q q,即正好是因变量个数
由于 Y Y Y只是1个变量,故 v 1 {{v}_{1}} v1?是1个标量。 ∥ v 1 ∥ = 1 ? v 1 = 1 \left\| {{v}_{1}} \right\|\text{=}1\Rightarrow {{v}_{1}}=1 ∥v1?∥=1?v1?=1,即: u 1 = Y {{u}_{1}}=Y u1?=Y
t 1 {{t}_{1}} t1?和 u 1 {{u}_{1}} u1?满足(1)中两个条件则有:
变异信息最大: V a r ( t 1 ) → max ? , V a r ( u 1 ) → max ? Var({{t}_{1}})\to \max ,Var({{u}_{1}})\to \max Var(t1?)→max,Var(u1?)→max
相关程度最大: r ( t 1 , u 1 ) → max ? r({{t}_{1}},{{u}_{1}})\to \max r(t1?,u1?)→max 相关程度最大, r ( t 1 , u 1 ) r({{t}_{1}},{{u}_{1}}) r(t1?,u1?)指的就是线性相关了

综合可得协方差最大: C o v ( t 1 , u 1 ) = r ( t 1 , u 1 ) V a r ( t 1 ) V a r ( u 1 ) → max ? Cov({{t}_{1}},{{u}_{1}})=r({{t}_{1}},{{u}_{1}})\sqrt{Var({{t}_{1}})Var({{u}_{1}})}\to \max Cov(t1?,u1?)=r(t1?,u1?)Var(t1?)Var(u1?) ?→max
由于 1 n < X w 1 , Y v 1 > = C o v ( t 1 , u 1 ) \frac{1}{n}=Cov({{t}_{1}},{{u}_{1}}) n1?=Cov(t1?,u1?)且 n n n为常数,则:
max ? < X w 1 , Y v 1 > = ( X w 1 ) T Y v 1 = w 1 T X T Y v 1 s . t { w 1 T w 1 = ∥ w 1 ∥ 2 = 1 v 1 T v 1 = ∥ v 1 ∥ 2 = 1 \begin{aligned} & \max ={{(X{{w}_{1}})}^{T}}Y{{v}_{1}}=w_{_{1}}^{T}{{X}^{T}}Y{{v}_{1}} \\ & s.t\left\{ \begin{matrix} w_{_{1}}^{T}{{w}_{1}}={{\left\| {{w}_{1}} \right\|}^{2}}=1 \\ v_{_{1}}^{T}{{v}_{1}}={{\left\| {{v}_{1}} \right\|}^{2}}=1 \\ \end{matrix} \right. \\ \end{aligned} ?max=(Xw1?)TYv1?=w1?T?XTYv1?s.t{w1?T?w1?=∥w1?∥2=1v1?T?v1?=∥v1?∥2=1??
根据拉格朗日算法有:
f = w 1 T X T Y v 1 ? λ ( w 1 T w 1 ? 1 ) ? μ ( v 1 T v 1 ? 1 ) f=w_{_{1}}^{T}{{X}^{T}}Y{{v}_{1}}-\lambda (w_{_{1}}^{T}{{w}_{1}}-1)-\mu (v_{_{1}}^{T}{{v}_{1}}-1) f=w1?T?XTYv1??λ(w1?T?w1??1)?μ(v1?T?v1??1)
对 f f f分别求关于 w 1 , v 1 , λ , μ {{w}_{1}},{{v}_{1}},\lambda ,\mu w1?,v1?,λ,μ的偏导且置0(求),有:
{ ? f ? w 1 = X T Y v 1 ? 2 λ w 1 = 0 ? f ? v 1 = Y T X w 1 ? 2 μ v 1 = 0 ? f ? λ = ? ( w 1 T w 1 ? 1 ) = 0? f ? μ = ? ( v 1 T v 1 ? 1 ) = 0\left\{ \begin{matrix} \frac{\partial f}{\partial {{w}_{1}}}={{X}^{T}}Y{{v}_{1}}-2\lambda {{w}_{1}}=0 \\ \frac{\partial f}{\partial {{v}_{1}}}={{Y}^{T}}X{{w}_{1}}-2\mu {{v}_{1}}=0 \\ \frac{\partial f}{\partial \lambda }=-(w_{_{1}}^{T}{{w}_{1}}-1)=0\ \ \ \ \\ \frac{\partial f}{\partial \mu }=-(v_{_{1}}^{T}{{v}_{1}}-1)=0\ \ \ \ \ \\ \end{matrix} \right. ???????????w1??f?=XTYv1??2λw1?=0?v1??f?=YTXw1??2μv1?=0?λ?f?=?(w1?T?w1??1)=0?μ?f?=?(v1?T?v1??1)=0?
由上式可推出:
2 λ = 2 μ = w 1 T X T Y v 1 = ( X w 1 ) T Y v 1 = < X w 1 , Y v 1 > 2\lambda =2\mu =w_{_{1}}^{T}{{X}^{T}}Y{{v}_{1}}={{(X{{w}_{1}})}^{T}}Y{{v}_{1}}\text{=} 2λ=2μ=w1?T?XTYv1?=(Xw1?)TYv1?=

θ 1 = 2 λ = 2 μ = w 1 T X T Y v 1 {{\theta }_{1}}=2\lambda =2\mu =w_{_{1}}^{T}{{X}^{T}}Y{{v}_{1}} θ1?=2λ=2μ=w1?T?XTYv1?
则 θ 1 {{\theta }_{1}} θ1?是优化问题的目标函数且使是 θ 1 {{\theta }_{1}} θ1?达到最大必须有有:
{ X T Y v 1 = θ 1 w 1 Y T X w 1 = θ 1 v 1 \left\{ \begin{aligned} & {{X}^{T}}Y{{v}_{1}}={{\theta }_{1}}{{w}_{1}} \\ & {{Y}^{T}}X{{w}_{1}}\text{=}{{\theta }_{1}}{{v}_{1}} \\ \end{aligned} \right. {?XTYv1?=θ1?w1?YTXw1?=θ1?v1??
将上面组合式结合得:
X T Y ( 1 θ 1 Y T X w 1 ) = θ 1 w 1 ? X T Y Y T X w 1 = θ 1 2 w 1 {{X}^{T}}Y(\frac{1}{{{\theta }_{1}}}{{Y}^{T}}X{{w}_{1}})={{\theta }_{1}}{{w}_{1}}\Rightarrow {{X}^{T}}Y{{Y}^{T}}X{{w}_{1}}=\theta _{_{1}}^{2}{{w}_{1}} XTY(θ1?1?YTXw1?)=θ1?w1??XTYYTXw1?=θ1?2?w1?
同理可得:
Y T X X T Y v 1 = θ 1 2 v 1 {{Y}^{T}}X{{X}^{T}}Y{{v}_{1}}=\theta _{_{1}}^{2}{{v}_{1}} YTXXTYv1?=θ1?2?v1?

v 1 = 1 ? Y T X X T Y = θ 1 2 ? θ 1 2 = ∥ X T Y ∥ 2 ? θ 1 = ∥ X T Y ∥ {{v}_{1}}=1\Rightarrow {{Y}^{T}}X{{X}^{T}}Y=\theta _{_{1}}^{2}\Rightarrow \theta _{_{1}}^{2}\text{=}{{\left\| {{X}^{T}}Y \right\|}^{2}}\Rightarrow {{\theta }_{1}}\text{=}\left\| {{X}^{T}}Y \right\| v1?=1?YTXXTY=θ1?2??θ1?2?=∥∥?XTY∥∥?2?θ1?=∥∥?XTY∥∥?
.根据联合式(1)我们可以得到:
w 1 = 1 θ 1 X T Y v 1 = 1 θ 1 X T Y {{w}_{1}}=\frac{1}{{{\theta }_{1}}}{{X}^{T}}Y{{v}_{1}}=\frac{1}{{{\theta }_{1}}}{{X}^{T}}Y w1?=θ1?1?XTYv1?=θ1?1?XTY
结合上面两个表达式,我们可以推出:
w 1 = 1 θ 1 X T Y = X T Y ∥ X T Y ∥ {{w}_{1}}=\frac{1}{{{\theta }_{1}}}{{X}^{T}}Y=\frac{{{X}^{T}}Y}{\left\| {{X}^{T}}Y \right\|} w1?=θ1?1?XTY=∥XTY∥XTY?
由于 X X X、 Y Y Y是标准化矩阵,有:
X T Y = ( x 1 , x 2 , ? ? , x p ) T Y = ( x 1 , x 2 , ? ? , x p ) T y= ( x 1 T x 2 T ? x p T ) y = ( x 1 T y x 2 T y ? x p T y ) = ( r ( x 1 , y ) r ( x 2 , y ) ? r ( x p , y ) ) \begin{aligned} & {{X}^{T}}Y={{({{x}_{1}},{{x}_{2}},\cdots ,{{x}_{p}})}^{T}}Y={{({{x}_{1}},{{x}_{2}},\cdots ,{{x}_{p}})}^{T}}y \\ & \ \ \ \ \ \ \ =\left( \begin{matrix} x_{1}^{T} \\ x_{2}^{T} \\ \vdots \\ x_{p}^{T} \\ \end{matrix} \right)y=\left( \begin{matrix} x_{1}^{T}y \\ x_{2}^{T}y \\ \vdots \\ x_{p}^{T}y \\ \end{matrix} \right)=\left( \begin{matrix} r({{x}_{1}},y) \\ r({{x}_{2}},y) \\ \vdots \\ r({{x}_{p}},y) \\ \end{matrix} \right) \\ \end{aligned} ?XTY=(x1?,x2?,?,xp?)TY=(x1?,x2?,?,xp?)Ty=??????x1T?x2T??xpT????????y=??????x1T?yx2T?y?xpT?y???????=??????r(x1?,y)r(x2?,y)?r(xp?,y)????????
进一步,有:
w 1 = X T Y ∥ X T Y ∥ = 1 ∑ i = 1 p r 2 ( x i , y ) ( r ( x 1 , y ) r ( x 2 , y ) ? r ( x p , y ) ) {{w}_{1}}=\frac{{{X}^{T}}Y}{\left\| {{X}^{T}}Y \right\|}=\frac{1}{\sqrt{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},y)}}}\left( \begin{matrix} r({{x}_{1}},y) \\ r({{x}_{2}},y) \\ \vdots \\ r({{x}_{p}},y) \\ \end{matrix} \right) w1?=∥XTY∥XTY?=i=1∑p?r2(xi?,y) ?1???????r(x1?,y)r(x2?,y)?r(xp?,y)???????
t 1 = X w 1 = 1 ∑ i = 1 p r 2 ( x i , y ) X ( r ( x 1 , y ) r ( x 2 , y ) ? r ( x p , y ) ) = 1 ∑ i = 1 p r 2 ( x i , y ) ( x 1 r ( x 1 , y ) + x 2 r ( x 2 , y ) + ? + x p r ( x p , y ) ) {{t}_{1}}=X{{w}_{1}}=\frac{1}{\sqrt{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},y)}}}X\left( \begin{matrix} r({{x}_{1}},y) \\ r({{x}_{2}},y) \\ \vdots \\ r({{x}_{p}},y) \\ \end{matrix} \right)=\frac{1}{\sqrt{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},y)}}}({{x}_{1}}r({{x}_{1}},y)+{{x}_{2}}r({{x}_{2}},y)+\cdots +{{x}_{p}}r({{x}_{p}},y)) t1?=Xw1?=i=1∑p?r2(xi?,y) ?1?X??????r(x1?,y)r(x2?,y)?r(xp?,y)???????=i=1∑p?r2(xi?,y) ?1?(x1?r(x1?,y)+x2?r(x2?,y)+?+xp?r(xp?,y))
在 t 1 {{t}_{1}} t1?中,关于 x i {{x}_{i}} xi?的线性组合系数是:
r ( x i , y ) ∑ i = 1 p r 2 ( x i , y ) \frac{r({{x}_{i}},y)}{\sqrt{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},y)}}} i=1∑p?r2(xi?,y) ?r(xi?,y)?
在 t 1 {{t}_{1}} t1?中,关于 x i {{x}_{i}} xi?的线性组合系数中,若 x i {{x}_{i}} xi?与 y y y的相关程度越强,则在 t 1 {{t}_{1}} t1?成分中的组合系数越大。
而此时目标函数的优化值: θ 1 = ∥ X T Y ∥ = ∑ i = 1 p r 2 ( x i , y ) {{\theta }_{1}}\text{=}\left\| {{X}^{T}}Y \right\|\text{=}\sqrt{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},y)}} θ1?=∥∥?XTY∥∥?=i=1∑p?r2(xi?,y) ?
这里和多因变量不同在于不用求特征向量了,因为特征向量就是1,而特征值可以直接求出来
我们通过求得 w 1 {{w}_{1}} w1?和 v 1 = 1 {{v}_{1}}\text{=}1 v1?=1之后即可得到第1成分:
{ t 1 = X w 1 u 1 = Y v 1 = Y \left\{ \begin{aligned} & {{t}_{1}}=X{{w}_{1}} \\ & {{u}_{1}}=Y{{v}_{1}}\text{=}Y \\ \end{aligned} \right. {?t1?=Xw1?u1?=Yv1?=Y?
由(1)式我们可以进一步推导出:
θ 1 = < t 1 , u 1 > = w 1 T X T Y v 1 = w 1 T X T Y {{\theta }_{1}}\text{=}<{{t}_{1}},{{u}_{1}}>=w_{1}^{T}{{X}^{T}}Y{{v}_{1}}=w_{1}^{T}{{X}^{T}}Y θ1?==w1T?XTYv1?=w1T?XTY
然后分别进行 X X X、 Y Y Y对 t 1 {{t}_{1}} t1?的回归(这里 Y Y Y对 t 1 {{t}_{1}} t1?的回归):
{ X = t 1 p 1 T + X 1 Y = t 1 r 1 T + Y 1 \left\{ \begin{aligned} & X={{t}_{1}}p_{1}^{T}+{{X}_{1}} \\ & Y={{t}_{1}}r_{1}^{T}+{{Y}_{1}} \\ \end{aligned} \right. {?X=t1?p1T?+X1?Y=t1?r1T?+Y1??
其中,回归系数向量:
{ p 1 = X T t 1 ∥ t 1 ∥ 2 r 1 = Y T t 1 ∥ t 1 ∥ 2 \left\{ \begin{aligned} & {{p}_{1}}=\frac{{{X}^{T}}{{t}_{1}}}{{{\left\| {{t}_{1}} \right\|}^{2}}} \\ & {{r}_{1}}=\frac{{{Y}^{T}}{{t}_{1}}}{{{\left\| {{t}_{1}} \right\|}^{2}}} \\ \end{aligned} \right. ???????????p1?=∥t1?∥2XTt1??r1?=∥t1?∥2YTt1???
另外, X 1 {{X}_{1}} X1?、 Y 1 {{Y}_{1}} Y1?则为 X X X、 Y Y Y的残差信息矩阵。(回归系数向量可利用PLS回归性质推导?)
③用残差信息矩阵 X 1 {{X}_{1}} X1?、 Y 1 {{Y}_{1}} Y1?取代 X X X、 Y Y Y,用同样的方法重复第②步,得到:
w 2 = X 1 T Y 1 ∥ X 1 T Y 1 ∥ = 1 ∑ i = 1 p r 2 ( x i , y ) ( r ( x 1 , y ) r ( x 2 , y ) ? r ( x p , y ) ) {{w}_{2}}=\frac{{{X}_{1}}^{T}{{Y}_{1}}}{\left\| {{X}_{1}}^{T}{{Y}_{1}} \right\|}=\frac{1}{\sqrt{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},y)}}}\left( \begin{matrix} r({{x}_{1}},y) \\ r({{x}_{2}},y) \\ \vdots \\ r({{x}_{p}},y) \\ \end{matrix} \right) w2?=∥∥?X1?TY1?∥∥?X1?TY1??=i=1∑p?r2(xi?,y) ?1???????r(x1?,y)r(x2?,y)?r(xp?,y)???????
这里注意 x i ( i = 1 , 2 , ? ? , p ) {{x}_{i}}(i=1,2,\cdots ,p) xi?(i=1,2,?,p)是指 X 1 {{X}_{1}} X1?中的列向量,而 y y y是指 Y 1 {{Y}_{1}} Y1?中的列向量。
求第2个成分 t 2 {{t}_{2}} t2?、 u 2 {{u}_{2}} u2?和第2个轴 w 2 {{w}_{2}} w2?、 v 2 =1 {{v}_{2}}\text{=1} v2?=1,即:
{ t 2 = X 1 w 2 u 2 = Y 1 v 2 = Y 1 \left\{ \begin{aligned} & {{t}_{2}}={{X}_{1}}{{w}_{2}} \\ & {{u}_{2}}={{Y}_{1}}{{v}_{2}}\text{=}{{Y}_{1}} \\ \end{aligned} \right. {?t2?=X1?w2?u2?=Y1?v2?=Y1??
θ 2 = < t 2 , u 2 > = w 2 T X 1 T Y 1 v 2 = w 2 T X 1 T Y 1 {{\theta }_{2}}=<{{t}_{2}},{{u}_{2}}>=w_{2}^{T}X_{1}^{T}{{Y}_{1}}{{v}_{2}}=w_{2}^{T}X_{1}^{T}{{Y}_{1}} θ2?==w2T?X1T?Y1?v2?=w2T?X1T?Y1?
得到回归方程:
{ X 1 = t 2 p 2 T + X 2 Y 1 = t 2 r 2 T + Y 2 \left\{ \begin{aligned} & {{X}_{1}}={{t}_{2}}p_{2}^{T}+{{X}_{2}} \\ & {{Y}_{1}}={{t}_{2}}r_{2}^{T}+{{Y}_{2}} \\ \end{aligned} \right. {?X1?=t2?p2T?+X2?Y1?=t2?r2T?+Y2??
其中,回归系数向量:
{ p 2 = X 1 T t 2 ∥ t 2 ∥ 2 r 2 = Y 1 T t 2 ∥ t 2 ∥ 2 \left\{ \begin{aligned} & {{p}_{2}}=\frac{X_{1}^{T}{{t}_{2}}}{{{\left\| {{t}_{2}} \right\|}^{2}}} \\ & {{r}_{2}}=\frac{{{Y}_{1}}^{T}{{t}_{2}}}{{{\left\| {{t}_{2}} \right\|}^{2}}} \\ \end{aligned} \right. ???????????p2?=∥t2?∥2X1T?t2??r2?=∥t2?∥2Y1?Tt2???
④如此利用剩下的残差信息矩阵不断迭代计算,我们假设 X X X的秩为 m m m(即可以有A个成分):
{ X = t 1 p 1 T + t 2 p 2 T + ? + t m p m T Y = t 1 r 1 T + t 2 r 2 T + ? + t m r m T + Y m \left\{ \begin{aligned} & X={{t}_{1}}p_{1}^{T}+{{t}_{2}}p_{2}^{T}+\cdots +{{t}_{m}}p_{m}^{T} \\ & Y={{t}_{1}}r_{1}^{T}+{{t}_{2}}r_{2}^{T}+\cdots +{{t}_{m}}r_{m}^{T}\text{+}{{Y}_{m}} \\ \end{aligned} \right. {?X=t1?p1T?+t2?p2T?+?+tm?pmT?Y=t1?r1T?+t2?r2T?+?+tm?rmT?+Ym??
而 t 1 , t 2 , ? ? , t m {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{m}} t1?,t2?,?,tm?可表示成 X =???? { ????x 1 , x 2 , ? ? , x p???? } ????X\text{= }\!\!\{\!\!\text{ }{{x}_{1}},{{x}_{2}},\cdots ,{{x}_{p}}\text{ }\!\!\}\!\!\text{ } X= { x1?,x2?,?,xp? } 的线性组合
其中 Y m {{Y}_{m}} Ym?为第 m m m个残差矩阵
由于 w h ? = ∏ k = 1 h ? 1 ( E ? w k p k T ) w h&t h = X w h ? w_{h}^{*}=\prod\limits_{k=1}^{h-1}{(E-{{w}_{k}}p_{k}^{T})}{{w}_{h}}\ \And \ \ {{t}_{h}}=Xw_{h}^{*} wh??=k=1∏h?1?(E?wk?pkT?)wh? &th?=Xwh?? (在多因变量线性偏最小二乘法性质中)则有:
Y = t 1 r 1 T + t 2 r 2 T + ? + t m r m T + Y m= ( X w 1 ? ) r 1 T + ( X w 2 ? ) r 2 T + ? + ( X w m ? ) r m T + Y m= X ( ∑ i = 1 m w i ? r i T ) + Y m \begin{aligned} & Y={{t}_{1}}r_{1}^{T}+{{t}_{2}}r_{2}^{T}+\cdots +{{t}_{m}}r_{m}^{T}+{{Y}_{m}} \\ & \ \ \ =(Xw_{1}^{*})r_{1}^{T}+(Xw_{2}^{*})r_{2}^{T}+\cdots +(Xw_{m}^{*})r_{m}^{T}+{{Y}_{m}} \\ & \ \ =X\left( \sum\limits_{i=1}^{m}{w_{i}^{*}r_{i}^{T}} \right)+{{Y}_{m}} \\ \end{aligned} ?Y=t1?r1T?+t2?r2T?+?+tm?rmT?+Ym?=(Xw1??)r1T?+(Xw2??)r2T?+?+(Xwm??)rmT?+Ym?=X(i=1∑m?wi??riT?)+Ym??
令 B = ∑ i = 1 m w i r i T B=\sum\limits_{i=1}^{m}{{{w}_{i}}r_{i}^{T}} B=i=1∑m?wi?riT?即为PLS回归方程的回归系数向量,有:
Y = X B + F m Y=XB\text{+}{{F}_{m}} Y=XB+Fm?
1.2 辅助分析技术 ①与典型相关分析对应的分析技术
ⅰ.精度分析
在PLS计算推导中,在 X X X提取的自变量成分 t h {{t}_{h}} th?不仅要尽可能多的携带 X X X中的变异信息,而且要尽可能与 Y Y Y相关程度达到最大来解释 Y Y Y的信息。我们模仿典型相关分析中的精度分析,为了测量 t h {{t}_{h}} th?对 X X X和 Y Y Y的解释能力,定义 t h {{t}_{h}} th?的各种解释能力,有:
t h {{t}_{h}} th?对某自变量 x i {{x}_{i}} xi?的解释能力:
R d ( x i ; t h ) = r 2 ( x i ; t h ) Rd({{x}_{i}}; {{t}_{h}})={{r}^{2}}({{x}_{i}}; {{t}_{h}}) Rd(xi?; th?)=r2(xi?; th?)
t h {{t}_{h}} th?对 X X X的解释能力:
R d ( X ; t h ) = 1 p ∑ i = 1 p R d ( x i ; t h ) = 1 p ∑ i = 1 p r 2 ( x i , t h ) Rd(X; {{t}_{h}})=\frac{1}{p}\sum\limits_{i=1}^{p}{Rd({{x}_{i}}; {{t}_{h}})}=\frac{1}{p}\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},{{t}_{h}})} Rd(X; th?)=p1?i=1∑p?Rd(xi?; th?)=p1?i=1∑p?r2(xi?,th?)
t 1 , t 2 , ? ? , t h {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}} t1?,t2?,?,th?对 X X X的累计解释能力:
R d ( X ; t 1 , t 2 , ? ? , t h ) = ∑ h = 1 m R d ( X ; t h ) = 1 p ∑ h = 1 m ∑ i = 1 p r 2 ( x i , t h ) Rd(X; {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}})=\sum\limits_{h=1}^{m}{Rd(X; {{t}_{h}})}=\frac{1}{p}\sum\limits_{h=1}^{m}{\sum\limits_{i=1}^{p}{{{r}^{2}}({{x}_{i}},{{t}_{h}})}} Rd(X; t1?,t2?,?,th?)=h=1∑m?Rd(X; th?)=p1?h=1∑m?i=1∑p?r2(xi?,th?)
t h {{t}_{h}} th?对某因变量 y j {{y}_{j}} yj?的解释能力:
R d ( y j ; t h ) = r 2 ( y j ; t h ) Rd({{y}_{j}}; {{t}_{h}})={{r}^{2}}({{y}_{j}}; {{t}_{h}}) Rd(yj?; th?)=r2(yj?; th?)
t h {{t}_{h}} th?对 Y Y Y的解释能力:
R d ( Y ; t h ) = 1 q ∑ j = 1 q R d ( y j ; t h ) = 1 q ∑ j = 1 q r 2 ( y j , t h ) Rd(Y; {{t}_{h}})=\frac{1}{q}\sum\limits_{j=1}^{q}{Rd({{y}_{j}}; {{t}_{h}})}=\frac{1}{q}\sum\limits_{j=1}^{q}{{{r}^{2}}({{y}_{j}},{{t}_{h}})} Rd(Y; th?)=q1?j=1∑q?Rd(yj?; th?)=q1?j=1∑q?r2(yj?,th?)
t 1 , t 2 , ? ? , t h {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}} t1?,t2?,?,th?对 Y Y Y的累计解释能力:
R d ( Y ; t 1 , t 2 , ? ? , t h ) = ∑ h = 1 m R d ( Y ; t h ) = 1 q ∑ h = 1 m ∑ j = 1 q r 2 ( y j , t h ) Rd(Y; {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}})=\sum\limits_{h=1}^{m}{Rd(Y; {{t}_{h}})}=\frac{1}{q}\sum\limits_{h=1}^{m}{\sum\limits_{j=1}^{q}{{{r}^{2}}({{y}_{j}},{{t}_{h}})}} Rd(Y; t1?,t2?,?,th?)=h=1∑m?Rd(Y; th?)=q1?h=1∑m?j=1∑q?r2(yj?,th?)
ⅱ.测量自变量 x i {{x}_{i}} xi?对因变量集合 Y Y Y的解释能力
x i {{x}_{i}} xi?在解释 Y Y Y时作用的重要性,我们可以通过变量投影重要性指标( V I P i VI{{P}_{i}} VIPi?)来测量(Variable Importance in Projection),有:
V I P i = p R d ( Y ; t 1 , t 2 , ? ? , t h ) ∑ h = 1 m R d ( Y ; t h ) w h i 2 VI{{P}_{i}}=\sqrt{\frac{p}{Rd(Y; {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}})}\sum\limits_{h=1}^{m}{Rd(Y; {{t}_{h}})w_{hi}^{2}}} VIPi?=Rd(Y; t1?,t2?,?,th?)p?h=1∑m?Rd(Y; th?)whi2? ?
这里Y可看成单个因变量,也可看成因变量集合。
其中 w h i {{w}_{hi}} whi?是轴 w h {{w}_{h}} wh?第 i i i个分量(就是一个标量,其有 p p p个分量, w h {{w}_{h}} wh?是一个列向量,行数 p p p),由于针对 x i {{x}_{i}} xi?,在 t h = X h ? 1 w h {{t}_{h}}={{X}_{h-1}}{{w}_{h}} th?=Xh?1?wh?中, w h {{w}_{h}} wh?的第 i i i个分量(标量)对应解释 X h ? 1 {{X}_{h-1}} Xh?1?中的 x i {{x}_{i}} xi?,则 V I P i VI{{P}_{i}} VIPi?对应于 x i {{x}_{i}} xi?对 Y Y Y的解释时起到的作用程度,有:
∑ i = 1 p w h i 2 = w h T w h = 1 \sum\limits_{i=1}^{p}{w_{hi}^{2}}=w_{h}^{T}{{w}_{h}}=1 i=1∑p?whi2?=whT?wh?=1
上面可以如此解释: x i {{x}_{i}} xi?对 Y Y Y的解释是通过 t h {{t}_{h}} th?来实现的,则若 R d ( Y ; t h ) Rd(Y; {{t}_{h}}) Rd(Y; th?)值很大即 t h {{t}_{h}} th?对 Y Y Y的解释能力很强,由于 x i {{x}_{i}} xi?在构造 t h {{t}_{h}} th?起到非常重要作用,则 x i {{x}_{i}} xi?对 Y Y Y的解释能力就被视为很大。另外, x i {{x}_{i}} xi?是通过 w h {{w}_{h}} wh?来构造 t h {{t}_{h}} th?的,当 w h i {{w}_{hi}} whi?取很大值时,则 x i {{x}_{i}} xi?对 Y Y Y的解释能力就被视为很大,有:
V I P i 2 = p ∑ h = 1 m R d ( Y ; t h ) w h i 2 R d ( Y ; t 1 , t 2 , ? ? , t h ) VIP_{i}^{2}=\frac{p\sum\limits_{h=1}^{m}{Rd(Y; {{t}_{h}})w_{hi}^{2}}}{Rd(Y; {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}})} VIPi2?=Rd(Y; t1?,t2?,?,th?)ph=1∑m?Rd(Y; th?)whi2??
通过上面分析,当 R d ( Y ; t h ) Rd(Y; {{t}_{h}}) Rd(Y; th?)很大时,则有 w h i 2 w_{hi}^{2} whi2?很大,进一步有 V I P i 2 VIP_{i}^{2} VIPi2?很大。
∑ i p V I P i 2 = ∑ i p p ∑ h = 1 m R d ( Y ; t h ) w h i 2 R d ( Y ; t 1 , t 2 , ? ? , t h ) = p ∑ h = 1 m R d ( Y ; t h ) ∑ i p w h i 2 R d ( Y ; t 1 , t 2 , ? ? , t h ) = p \sum\limits_{i}^{p}{VIP_{i}^{2}}=\sum\limits_{i}^{p}{\frac{p\sum\limits_{h=1}^{m}{Rd(Y; {{t}_{h}})w_{hi}^{2}}}{Rd(Y; {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}})}}=\frac{p\sum\limits_{h=1}^{m}{Rd(Y; {{t}_{h}})\sum\limits_{i}^{p}{w_{hi}^{2}}}}{Rd(Y; {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{h}})}=p i∑p?VIPi2?=i∑p?Rd(Y; t1?,t2?,?,th?)ph=1∑m?Rd(Y; th?)whi2??=Rd(Y; t1?,t2?,?,th?)ph=1∑m?Rd(Y; th?)i∑p?whi2??=p
从上面分析我们可以知道,若针对所有的 x i {{x}_{i}} xi?与之对应的 V I P i ( i = 1 , 2 , ? ? , p ) VI{{P}_{i}}(i=1,2,\cdots ,p) VIPi?(i=1,2,?,p)均相等即在解释 Y Y Y时的作用相同,则所有的 V I P i VI{{P}_{i}} VIPi?均为1,否则对于 V I P i > 1 VI{{P}_{i}}>1 VIPi?>1的 x i {{x}_{i}} xi?在解释 Y Y Y时起到更加重要的作用。上面我们定义了 V I P i VI{{P}_{i}} VIPi?指标,均定性的能够分析出哪些自变量的起到的作用更大。
②与主成分分析对应的分析技术
ⅰ.特异点分析
我们可以模仿主成分分析定义第 i i i个样本点对地 h h h成分 t h {{t}_{h}} th?的贡献率 T h i 2 T_{hi}^{2} Thi2?以此来发现样本点集合中的特异点,有:
T h i 2 = t h i 2 ( n ? 1 ) s h 2 T_{hi}^{2}=\frac{t_{hi}^{2}}{(n-1)s_{h}^{2}} Thi2?=(n?1)sh2?thi2??
其中: t h i {{t}_{hi}} thi?是列向量 t h {{t}_{h}} th?(行数 n n n)Xscores的第 i i i个样本点对应的值, s h 2 s_{h}^{2} sh2?是成分 T H {{T}_{H}} TH?的方差
则样本点 I I I对成分 T 1 , T 2 … … T M {{T}_{1}},{{T}_{2}}……{{T}_{M}} T1?,T2?……TM?的累计贡献率:
T i 2 = 1 n ? 1 ∑ h = 1 m t h i 2 s h 2 T_{i}^{2}=\frac{1}{n-1}\sum\limits_{h=1}^{m}{\frac{t_{hi}^{2}}{s_{h}^{2}}} Ti2?=n?11?h=1∑m?sh2?thi2??
我们模仿主成分分析,由于一个样本点如果对成分构成贡献很大,则其存在会使分析造成比较大的误差,所以一个样本点对成分构成的贡献不可以很大,在SIMCA-P软件中利用特雷西等人证明的统计量:
n 2 ( n ? m ) m ( n 2 ? 1 ) T i 2 ~ F ( m , n ? m ) \frac{{{n}^{2}}(n-m)}{m({{n}^{2}}-1)}T_{i}^{2}\sim F(m,n-m) m(n2?1)n2(n?m)?Ti2?~F(m,n?m)
根据 F F F统计量检验,当 T i 2 ≥ m ( n 2 ? 1 ) n 2 ( n ? m ) F 0.05 ( m , n ? m ) T_{i}^{2}\ge \frac{m({{n}^{2}}-1)}{{{n}^{2}}(n-m)}{{F}_{0.05}}(m,n-m) Ti2?≥n2(n?m)m(n2?1)?F0.05?(m,n?m)我们认为在 95 95% 95的检验水平上,样本点 i i i对成分 t 1 , t 2 , ? ? , t m {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{m}} t1?,t2?,?,tm?的贡献过大,我们称之为样本点 I I I为一个特异点。
我们一般如果选择 M = 2 M=2 M=2即PLS回归中只采用了2个主成分或者 ( X ) = 2 (X)=2 (X)=2,此时有:
T i 2 = 1 n ? 1 ( t 1 i 2 s 1 2 + t 2 i 2 s 2 2 ) ≥ 2 ( n 2 ? 1 ) n 2 ( n ? 2 ) F 0.05 ( 2 , n ? 2 ) T_{i}^{2}\text{=}\frac{1}{n-1}\left( \frac{t_{1i}^{2}}{s_{1}^{2}}+\frac{t_{2i}^{2}}{s_{2}^{2}} \right)\ge \frac{2({{n}^{2}}-1)}{{{n}^{2}}(n-2)}{{F}_{0.05}}(2,n-2) Ti2?=n?11?(s12?t1i2??+s22?t2i2??)≥n2(n?2)2(n2?1)?F0.05?(2,n?2)
最后我们得到:
t 1 i 2 s 1 2 + t 2 i 2 s 2 2 ≥ 2 ( n 2 ? 1 ) ( n ? 1 ) n 2 ( n ? 2 ) F 0.05 ( 2 , n ? 2 ) \frac{t_{1i}^{2}}{s_{1}^{2}}+\frac{t_{2i}^{2}}{s_{2}^{2}}\ge \frac{2({{n}^{2}}-1)(n-1)}{{{n}^{2}}(n-2)}{{F}_{0.05}}(2,n-2) s12?t1i2??+s22?t2i2??≥n2(n?2)2(n2?1)(n?1)?F0.05?(2,n?2)
令 c = 2 ( n 2 ? 1 ) ( n ? 1 ) n 2 ( n ? 2 ) F 0.05 ( 2 , n ? 2 ) c=\frac{2({{n}^{2}}-1)(n-1)}{{{n}^{2}}(n-2)}{{F}_{0.05}}(2,n-2) c=n2(n?2)2(n2?1)(n?1)?F0.05?(2,n?2),有:
t 1 i 2 s 1 2 + t 2 i 2 s 2 2 = c \frac{t_{1i}^{2}}{s_{1}^{2}}+\frac{t_{2i}^{2}}{s_{2}^{2}}\text{=}c s12?t1i2??+s22?t2i2??=c
判断提取多个主成分是否在椭圆内外关系可通过:
t 1 I 2 s 1 2 + t 2 I 2 s 2 2 + ? t m i 2 s m 2 \frac{t_{1I}^{2}}{s_{1}^{2}}+\frac{t_{2I}^{2}}{s_{2}^{2}}+\cdots \frac{t_{mi}^{2}}{s_{m}^{2}} s12?t1I2??+s22?t2I2??+?sm2?tmi2??

m ( n 2 ? 1 ) ( n ? 1 ) n 2 ( n ? m ) F 0.05 ( m , n ? m ) \frac{m({{n}^{2}}-1)(n-1)}{{{n}^{2}}(n-m)}{{F}_{0.05}}(m,n-m) n2(n?m)m(n2?1)(n?1)?F0.05?(m,n?m)
计算方法:
m ( n 2 ? 1 ) ( n ? 1 ) n 2 ( n ? m ) f 0.05 ( m , n ? m ) = ( n 2 ? 1 ) ( n ? 1 ) n 2 ? m n ? m ? f 0.05 ( m , n ? m ) \frac{m({{n}^{2}}-1)(n-1)}{{{n}^{2}}(n-m)}{{f}_{0.05}}(m,n-m)\text{=}\frac{({{n}^{2}}-1)(n-1)}{{{n}^{2}}}\centerdot \frac{m}{n-m}\centerdot {{f}_{0.05}}(m,n-m) n2(n?m)m(n2?1)(n?1)?f0.05?(m,n?m)=n2(n2?1)(n?1)??n?mm??f0.05?(m,n?m)
MATLAB计算式:(n-1)* (n^2-1)/( n^2) * j*finv(0.95,j , n-j)/(n-j) j从1开始
三维:
t 1 i 2 s 1 2 + t 2 i 2 s 2 2 + t 3 i 2 s 3 2 = c ? t 1 i 2 ( s 1 c ) 2 + t 2 i 2 ( s 2 c ) 2 + t 3 i 2 ( s 3 c ) 2 =1 \frac{t_{1i}^{2}}{s_{1}^{2}}+\frac{t_{2i}^{2}}{s_{2}^{2}}+\frac{t_{3i}^{2}}{s_{3}^{2}}\text{=}c\Leftrightarrow \frac{t_{1i}^{2}}{{{\left( {{s}_{1}}\sqrt{c} \right)}^{2}}}+\frac{t_{2i}^{2}}{{{\left( {{s}_{2}}\sqrt{c} \right)}^{2}}}+\frac{t_{3i}^{2}}{{{\left( {{s}_{3}}\sqrt{c} \right)}^{2}}}\text{=1} s12?t1i2??+s22?t2i2??+s32?t3i2??=c?(s1?c ?)2t1i2??+(s2?c ?)2t2i2??+(s3?c ?)2t3i2??=1
上式是一个椭圆,所以,我们以 t 1 i {{t}_{1i}} t1i?和 t 2 i {{t}_{2i}} t2i?作为坐标轴,在 t 1 / t 2 {{t}_{1}}/{{t}_{2}} t1?/t2?平面图上,可以得到这个 t 2 {{t}^{2}} t2椭圆图,若所有样本点都落在这个椭圆内部,则认为所有样本点分布是均匀的,否则落在外部,则称这些点为特异点,即这个样本点远离所有样本集合的平均水平。
ⅱ.PLS后的数据质量分析
我们通过主成分分析可以知道,在PLS回归中有以下同样情况产生:由于特异点的存在或者仍然有一些样本点在PLS模型分析中得不到很好地表示,对于此类样本点,就无法根据PLS回归的表现来判断其特征,对于这类样本点分析必须十分小心。
由于在PLS模型分析中去除了一部分原始信息( m < ( A ) m<(A) m<(A))而使得一些样本点在 y j {{y}_{j}} yj?上的拟合值与原始值差异比较大。
由PLS模型计算推导我们可以知道,当提取了 m m m个成分 t 1 , t 2 , ? ? , t m {{t}_{1}},{{t}_{2}},\cdots ,{{t}_{m}} t1?,t2?,?,tm?后,有:
{ X ^ = t 1 p 1 T + t 2 p 2 T + ? + t m p m T Y ^ = t 1 r 1 T + t 2 r 2 T + ? + t m r m T + Y m \left\{ \begin{aligned} & \hat{X}={{t}_{1}}p_{1}^{T}+{{t}_{2}}p_{2}^{T}+\cdots +{{t}_{m}}p_{m}^{T} \\ & \hat{Y}={{t}_{1}}r_{1}^{T}+{{t}_{2}}r_{2}^{T}+\cdots +{{t}_{m}}r_{m}^{T}\text{+}{{Y}_{m}} \\ \end{aligned} \right. {?X^=t1?p1T?+t2?p2T?+?+tm?pmT?Y^=t1?r1T?+t2?r2T?+?+tm?rmT?+Ym??
我们定义样本点 i ( i = 1 , 2 , ? ? , n ) i(i=1,2,\cdots ,n) i(i=1,2,?,n)在 X X X空间与PLS模型的距离 D M o d X i ( s i ) DMod{{X}_{i}}({{s}_{i}}) DModXi?(si?):
s i = D M o d X i = ∑ j = 1 p e i j 2 p ? m ? n n ? m ? 1 {{s}_{i}}=DMod{{X}_{i}}=\sqrt{\frac{\sum\limits_{j=1}^{p}{e_{ij}^{2}}}{p-m}}\cdot \sqrt{\frac{n}{n-m-1}} si?=DModXi?=p?mj=1∑p?eij2?? ??n?m?1n? ?
其中 e i j 2 = ( x i j ? x ^ i j ) 2 e_{ij}^{2}={{({{x}_{ij}}-{{\hat{x}}_{ij}})}^{2}} eij2?=(xij??x^ij?)2, x ^ i j {{\hat{x}}_{ij}} x^ij?是重构矩阵 X ^ \hat{X} X^中样本点 i i i在变量 x j {{x}_{j}} xj?上的取值。
从上式我们可以知道,参入PLS模型的成分个数越多( m m m越大), s i {{s}_{i}} si?就越小即数据重构的误差就越小。可是,有时候 m m m过大,PLS模型的预测能力反而会降低,这和多元回归分析中一样,使用成分个数过多即使用变量个数过多,模型拟合效果看起来非常完美,但是模型却不能够识别系统信息与噪声,有时候如果我们把噪声加在了模型中,那这样的拟合效果反而更差
为此我们模型多元回归分析,定义一个调整复测定系数 R ˉ 2 {{\bar{R}}^{2}} Rˉ2,则由此我们我们这里定义模型距离的概念。
所有样本点重构的平均质量: S X 2 = 1 n ∑ i = 1 n s i 2 S_{X}^{2}\text{=}\frac{1}{n}\sum\limits_{i=1}^{n}{s_{i}^{2}} SX2?=n1?i=1∑n?si2?,进一步所有样本点的重构平均距离 S X {{S}_{X}} SX?
S X = 1 n ∑ i = 1 n ∑ j = 1 p e i j 2 p ? m ? n n ? m ? 1 = ∑ i = 1 n ∑ j = 1 p e i j 2 ( p ? m ) ( n ? m ? 1 ) {{S}_{X}}=\sqrt{\frac{1}{n}\sum\limits_{i=1}^{n}{\frac{\sum\limits_{j=1}^{p}{e_{ij}^{2}}}{p-m}\cdot \frac{n}{n-m-1}}}\text{=}\sqrt{\frac{\sum\limits_{i=1}^{n}{\sum\limits_{j=1}^{p}{e_{ij}^{2}}}}{(p-m)(n-m-1)}} SX?=n1?i=1∑n?p?mj=1∑p?eij2???n?m?1n? ?=(p?m)(n?m?1)i=1∑n?j=1∑p?eij2?? ?
则根据上面的 s i {{s}_{i}} si?和 S X {{S}_{X}} SX?我们可以推出样本点 i i i的标准化模型距离:
( D M o d X , N ) i = s i S X = D M o d X i S X {{(DModX,N)}_{i}}=\frac{{{s}_{i}}}{{{S}_{X}}}=\frac{DMod{{X}_{i}}}{{{S}_{X}}} (DModX,N)i?=SX?si??=SX?DModXi??
上式表明同所有样本点的重构平均质量相比,样本点 i i i是否偏大。
同理我们可以得到:
样本点 i ( i = 1 , 2 , ? ? , n ) i(i=1,2,\cdots ,n) i(i=1,2,?,n)在 Y Y Y空间与PLS模型的距离 D M o d Y i DMod{{Y}_{i}} DModYi?:
D M o d Y i = ∑ k = 1 q f i k 2 q ? m ? n n ? m ? 1 DMod{{Y}_{i}}=\sqrt{\frac{\sum\limits_{k=1}^{q}{f_{ik}^{2}}}{q-m}}\cdot \sqrt{\frac{n}{n-m-1}} DModYi?=q?mk=1∑q?fik2?? ??n?m?1n? ?
其中, f i j 2 = ( y i k ? y ^ i k ) 2 f_{ij}^{2}={{({{y}_{ik}}-{{\hat{y}}_{ik}})}^{2}} fij2?=(yik??y^?ik?)2, y ^ i k {{\hat{y}}_{ik}} y^?ik?是重构矩阵 Y ^ \hat{Y} Y^中样本点 i i i在变量 y j {{y}_{j}} yj?上的取值。
所有样本点重构的平均质量:
S X 2 = 1 n ∑ i = 1 n s i 2 S_{X}^{2}\text{=}\frac{1}{n}\sum\limits_{i=1}^{n}{s_{i}^{2}} SX2?=n1?i=1∑n?si2?
进一步所有样本点的重构平均距离 S X {{S}_{X}} SX?:
S Y = 1 n ∑ i = 1 n ∑ k = 1 q f i k 2 q ? m ? n n ? m ? 1 = ∑ i = 1 n ∑ k = 1 q e i k 2 ( q ? m ) ( n ? m ? 1 ) {{S}_{Y}}=\sqrt{\frac{1}{n}\sum\limits_{i=1}^{n}{\frac{\sum\limits_{k=1}^{q}{f_{ik}^{2}}}{q-m}\cdot \frac{n}{n-m-1}}}\text{=}\sqrt{\frac{\sum\limits_{i=1}^{n}{\sum\limits_{k=1}^{q}{e_{ik}^{2}}}}{(q-m)(n-m-1)}} SY?=n1?i=1∑n?q?mk=1∑q?fik2???n?m?1n? ?=(q?m)(n?m?1)i=1∑n?k=1∑q?eik2?? ?
则根据上面的 s i {{s}_{i}} si?和 S X {{S}_{X}} SX?我们可以推出样本点 i i i的标准化模型距离:
( D M o d Y , N ) i = D M o d Y i S Y {{(DModY,N)}_{i}}=\frac{DMod{{Y}_{i}}}{{{S}_{Y}}} (DModY,N)i?=SY?DModYi??
上式表明同所有样本点的重构平均质量相比,样本点 i i i是否偏大,若偏大,则说明数据重构质量不理想即PLS模型不好或者说 m m m的取值不理想即成分个数选取不适当。
Reference
【PLS系列003 单因变量线性PLS】王惠文.偏最小二乘方法原理及其应用
郭建校. 改进的高维非线性PLS回归方法及应用研究[D]. 天津大学, 2010.

    推荐阅读