0 引言

主成分分析^[1]是考察多个变量间相关性的一种多元统计方法, 通过几个主成分来揭示多个变量间的内部结构, 即从原始变量中导出几个主成分, 使其尽可能多地保留原始变量的信息, 且彼此间互不相关.该方法被广泛应用于金融、经济、管理等领域^[2-5].但通过主成分分析所得的每个主成分为所有原始变量的线性组合, 使所得主成分难于解释, 而解决实际问题时, 有时只需考虑与主成分关系比较密切的一些原始变量.为了凸显主成分和原始变量的关系, 一些学者将稀疏性引入主成分分析. 2003年, Jolliffe受Lasso^[6]的启发, 将L₁罚引入主成分^[7], 提出了模型

$ \mathit{\boldsymbol{\hat \beta }}{\rm{ = }}\mathop {\arg \min }\limits_\mathit{\boldsymbol{\beta }} {\left( {{\mathit{\boldsymbol{Y}}_i} - \mathit{\boldsymbol{X}}{\mathit{\boldsymbol{\beta }}_j}} \right)^2},{\rm{s}}{\rm{.t}}\;\;\sum\limits_{j = 1}^p {\left| {{\mathit{\boldsymbol{\beta }}_j}} \right|} < t. $

(1)

该模型实现主成分对原始变量的自动选择, 保留与主成分关系密切的原始变量, 剔除与主成分关系不密切的原始变量. 2005年, Zou^[8]等将模型(1) 中主成分的求解问题直接转化为Lasso回归问题, 有效地把主成分的求解转化为线形模型的变量选择问题.但当观测变量的个数远远大于样本的个数时, 通过L₁罚所得的解过于稀疏, 导致大量信息被损失.为了克服该缺点, Zou将“elastic net”惩罚结构引入主成分, 提出了稀疏主成分分析^[9].其模型为

$ \left( {\mathit{\boldsymbol{A}},\mathit{\boldsymbol{B}}} \right) = \mathop {\arg \min }\limits_{\mathit{\boldsymbol{A}},\mathit{\boldsymbol{B}}} \sum\limits_{i = 1}^n {{{\left\| {{\mathit{\boldsymbol{X}}_i} - \mathit{\boldsymbol{A}}{\mathit{\boldsymbol{B}}^{\rm{T}}}{\mathit{\boldsymbol{X}}_i}} \right\|}^2}} + \lambda \sum\limits_{j = 1}^k {{{\left\| {{\mathit{\boldsymbol{\beta }}_j}} \right\|}^2}} + \sum\limits_{j = 1}^k {{\lambda _{1,j}}{{\left\| {{\mathit{\boldsymbol{\beta }}_j}} \right\|}_1}} ,{\rm{s}}{\rm{.t}}\;{\mathit{\boldsymbol{A}}^{\rm{T}}}\mathit{\boldsymbol{A = }}{\mathit{\boldsymbol{I}}_{k \times k}}. $

(2)

其中, A是通过主成分分析得到的前k个主成分的系数组成的矩阵, β_j为矩阵B中第j列.

稀疏主成分分析具有L₁罚和L₂罚的优点, 能凸显所得主成分与原始变量的关系, 并有效解决变量个数大于样本个数的优化问题.在稀疏主成分分析中, L₁罚是保证所得主成分系数稀疏, L₂罚是为了克服当样本量小于变量个数时, 所得的主成分系数过稀疏的缺点.一些学者将稀疏主成分分析应用到综合评价、股票研究及其他方面^[10-12], 取得较好的效果.然而, 基于L₁罚的优化问题仍存在至少两方面的不足:第一, 数据之间可能存在很大的冗余难以去除; 第二, 无法区分稀疏尺度的位置, 即会出现低尺度的能量转移到高尺度的现象, 因而易出现高频震荡现象.

2000年, Nikolova在研究变量选择时提出了变换L₁罚^[13](transformed L₁ penalty, 简称TL₁罚), 并指出TL₁罚函数$\frac{{\left( {a + 1} \right)\left| \mathit{\boldsymbol{\beta }} \right|}}{{a + \left| \mathit{\boldsymbol{\beta }} \right|}}$是介于L₀罚和L₁罚的一种罚函数.事实上, 当a→0时, TL₁罚逼近L₀罚; 当a→1时, TL₁罚逼近L_1/2罚; 当a→∞时, TL₁罚逼近L₁罚. 2009年, Lv和Fan^[14]研究了TL₁罚的无偏性、稀疏性和连续性.Xin等在2014年把TL₁罚应用到压缩感知^[15-16], 并提出了TL₁罚函数的阈值迭代函数等.从上述的研究中可以发现, 基于TL₁罚优化问题的解要优于基于L₁罚优化问题.

本文将TL₁罚应用于稀疏主成分分析, 即用TL₁罚替代稀疏主成分分析中的L₁罚, 对稀疏主成分分析进行改进, 以克服基于L₁罚优化问题的不足, 并给出优化模型.然后通过2-范数的性质把TL₁罚稀疏主成分分析的求解进行转化.最后利用阈值迭代算法对优化问题(4) 进行求解, 并将该方法应用到蔬菜选择实例中, 结果表明TL₁罚稀疏主成分分析具有效果更优.

1 TL₁罚稀疏主成分分析

用TL₁罚$\frac{{\left( {a + 1} \right)\left| \mathit{\boldsymbol{\beta }} \right|}}{{a + \left| \mathit{\boldsymbol{\beta }} \right|}}$替代稀疏主成分分析中的L₁罚函数, 将用TL₁罚$\frac{{\left( {a + 1} \right)\left| \mathit{\boldsymbol{\beta }} \right|}}{{a + \left| \mathit{\boldsymbol{\beta }} \right|}}$替换后所得的稀疏主成分分析称为TL₁罚稀疏主成分分析, 则TL₁罚稀疏主成分分析模型为

$ \mathit{\boldsymbol{\hat \beta }}{\rm{ = }}\mathop {\arg \min }\limits_\mathit{\boldsymbol{\beta }} \left\| {{\mathit{\boldsymbol{Y}}_i} - \mathit{\boldsymbol{X}}{\mathit{\boldsymbol{\beta }}_j}} \right\|_2^2 + {\lambda _1}\left\| {{\mathit{\boldsymbol{\beta }}_j}} \right\|_2^2 + {\lambda _2}\frac{{\left( {a + 1} \right)\left| {{\mathit{\boldsymbol{\beta }}_j}} \right|}}{{a + \left| {{\mathit{\boldsymbol{\beta }}_j}} \right|}}. $

(3)

其中, Y_i是第i个主成分, λ₁, λ₂, a是参数.

利用2-范数的性质, TL₁罚稀疏主成分分析可转换为

$ \mathit{\boldsymbol{\hat \beta }}_i^* = \mathop {\arg \min }\limits_\mathit{\boldsymbol{\beta }} \left\| {\mathit{\boldsymbol{Y}}_i^* - {\mathit{\boldsymbol{X}}^*}\mathit{\boldsymbol{\beta }}_i^*} \right\|_2^2 + {\lambda _2}\frac{{\left( {a + 1} \right)\left| {\mathit{\boldsymbol{\beta }}_i^*} \right|}}{{a + \left| {\mathit{\boldsymbol{\beta }}_i^*} \right|}}. $

(4)

其中, $\mathit{\boldsymbol{Y}}_i^* = \left( \begin{gathered} {\mathit{\boldsymbol{Y}}_i} \hfill \\ 0 \hfill \\ \end{gathered} \right),{\mathit{\boldsymbol{X}}^*} = \left( \begin{gathered} \mathit{\boldsymbol{X}} \hfill \\ \sqrt {{\lambda _1}\mathit{\boldsymbol{I}}} \hfill \\ \end{gathered} \right),\mathit{\boldsymbol{\beta }}_i^* = {\mathit{\boldsymbol{\beta }}_i}$.

为了求解优化问题(4), 利用文献[15-16]中提出的阈值迭代算法, 给出该问题的阈值迭代函数, 即

$ {\mathit{\boldsymbol{\beta }}^ * } = \left\{ {\begin{array}{*{20}{c}} {0,}&{\left| \mathit{\boldsymbol{\beta }} \right| < t,}\\ {{g_{{\lambda _2}}}\left( \mathit{\boldsymbol{\beta }} \right),}&{\left| \mathit{\boldsymbol{\beta }} \right| > t.} \end{array}} \right. $

其中,

$ \begin{array}{*{20}{c}} {t = \left\{ {\begin{array}{*{20}{c}} {{\lambda _2}\mu \frac{{a + 1}}{a},}&{{\lambda _2} \le \frac{{{a^2}}}{{2\left( {a + 1} \right)\mu }},}\\ {\sqrt {2{\lambda _2}\mu \left( {a + 1} \right)} - \frac{a}{2},}&{{\lambda _2} > \frac{{{a^2}}}{{2\left( {a + 1} \right)\mu }},} \end{array}} \right.}\\ {{g_{{\lambda _2}}}\left( \mathit{\boldsymbol{\beta }} \right) = {\rm{sign}}\left( \mathit{\boldsymbol{\beta }} \right)\left\{ {\frac{2}{3}\left( {a + \left| \mathit{\boldsymbol{\beta }} \right|} \right)\cos \left( {\frac{{\varphi \left( \mathit{\boldsymbol{\beta }} \right)}}{3}} \right) - \frac{{2a}}{3} + \frac{{\left| \mathit{\boldsymbol{\beta }} \right|}}{3}} \right\},}\\ {\varphi \left( \mathit{\boldsymbol{\beta }} \right) = \arccos \left( {1 - \frac{{27{\lambda _2}a\left( {a + 1} \right)}}{{2{{\left( {a + \left| \mathit{\boldsymbol{\beta }} \right|} \right)}^3}}}} \right).} \end{array} $

下面, 给出TL₁罚稀疏主成分分析的阈值迭代算法, 具体步骤如下:

(1) 对原始数据进行主成分分析, 并按方差累计贡献率提取k个主成分;

(2) 初始化:将第i主成分的系数初始为x⁰, 给一个合适的a, ε, μ₀;

(3) 计算zⁿ=B_μ(xⁿ)=xⁿ+μA^T(y-Axⁿ), 令λ₂=λ₀, μ=μ₀;

(4) 判断λ₂与$\frac{{{a}^{2}}}{2\left( a+1 \right)\mu }$的大小:如果${{\lambda }_{2}}\le \frac{{{a}^{2}}}{2\left( a+1 \right)\mu }$, 令$t={{\lambda }_{2}}\mu \frac{a+1}{a}$, 判断|z_n(i)|与t的大小, 若|z_n(i)|>t, 则xⁿ⁺¹(i)=g_λμ(zⁿ(i)), 若|z_n(i)|≤t, 则xⁿ⁺¹(i)=0;如果${{\lambda }_{2}}>\frac{{{a}^{2}}}{2\left( a+1 \right)\mu }$, 令$t=\sqrt{2{{\lambda }_{2}}\mu \left( a+1 \right)}-\frac{a}{2}$, 判断|z_n(i)|与t的大小, 若|z_n(i)|>t, 则xⁿ⁺¹(i)=g_λμ(zⁿ(i)), 若|z_n(i)|≤t, 则xⁿ⁺¹(i)=0;

(5) 重复步骤(3)~(4), 当|xⁿ⁺¹(i)-xⁿ(i)| < ε, 或n>3 000时, 停止迭代, 输出xⁿ⁺¹.

2 数值模拟

对研究数据(2014年数学建模D题中常见蔬菜营养成分表中数据) 进行处理, 将蔬菜的种类作为样本, 蔬菜中所含的各种膳食纤维的含量作为变量, 分别进行主成分分析、稀疏主成分分析和TL₁罚稀疏主成分分析.分析时，按方差累计贡献率的80%来提取主成分, 提取主成分的个数为4.在此基础上, 利用阈值迭代算法对所提取的主成分进行稀疏主成分分析, 得到稀疏主成分对应的系数矩阵及相应的方差贡献率, 如表 1所示.

稀疏主成分形式为

$ {F_i} = \mathit{\boldsymbol{\alpha }}_i^{\rm{T}}\mathit{\boldsymbol{x}}, $

(5)

其中, F_i为第i个主成分, α_i为第i个主成分的系数, i=1, 2, 3, 4, x=(x₁, x₂, …, x₁₆)′.

由表 1和式(5) 可知, 稀疏主成分分析的稀疏性表现在稀疏主成分的系数上, 系数中零的个数越多, 稀疏主成分越稀疏.

由得到的稀疏主成分和方差贡献率, 可以得到稀疏主成分分析的模型为

$ F = 0.531{F_1} + 0.209{F_2} + 0.135{F_3} + 0.125{F_4}. $

利用TL₁罚稀疏主成分分析的阈值迭代函数对主成分分析的主成分进行处理, 得到TL₁罚稀疏主成分对应的系数及相应的方差贡献率, 如表 2所示.

TL₁罚稀疏主成分的形式为

$ {{\hat F}_i} = \mathit{\boldsymbol{\beta }}_i^{\rm{T}}\mathit{\boldsymbol{x}}. $

其中, ${{\hat{F}}_{i}}$为第i个稀疏主成分, β_i为第i个稀疏主成分的系数, i=1, 2, 3, 4, x=(x₁, x₂, …，x₁₆)′.

与表 1相同, 表 2中主成分具有稀疏性, 且零的个数越多, 稀疏主成分越稀疏.由表 1和表 2可知, 稀疏主成分分析和TL₁罚稀疏主成分分析都具有稀疏性, 且两者稀疏性基本相同, 但在方差累计贡献率方面, TL₁罚稀疏主成分分析略高于稀疏主成分分析.

根据所得的稀疏主成分和方差贡献率, 可以得到TL₁罚稀疏主成分分析的模型为

$ \hat F = 0.651{{\hat F}_1} + 0.164{{\hat F}_2} + 0.086{{\hat F}_3} + 0.099{{\hat F}_4}. $

由上述所得的两种模型, 根据常见蔬菜各种膳食纤维营养的含量, 计算每种蔬菜的各主成分得分, 再利用主成分分析(PCA)、稀疏主成分分析(SPCA) 和TL₁罚稀疏主成分分析(TLPCA) 的主成分模型对已知的蔬菜主成分得分进行排序.结果如表 3所示.

由表 3可知, 主成分分析、稀疏主成分分析和TL₁罚稀疏主成分分析对蔬菜的排序结果相差不大.其中, 排名前五的蔬菜均为蘑菇、榨菜、木耳、香菇和茄子.

3 结束语

用TL₁罚函数$\frac{{\left( {a + 1} \right)\left| \mathit{\boldsymbol{\beta }} \right|}}{{a + \left| \mathit{\boldsymbol{\beta }} \right|}}$来替代稀疏主成分分析中的L₁罚函数, 利用阈值迭代算法对稀疏主成分分析和TL₁罚稀疏主成分分析进行数值实验, 根据所得稀疏主成分和TL₁罚稀疏主成分的系数可知, 所得主成分都具有稀疏性, 且两者的稀疏度基本相同, 但TL₁罚稀疏主成分分析的方差累计贡献率略高于稀疏主成分分析, 即TL₁罚稀疏主成分分析提取原始数据的信息量略多于稀疏主成分分析.利用稀疏主成分分析和TL₁罚稀疏主成分分析对常见蔬菜进行营养含量排序, 排序结果基本相同, 即在主要蔬菜选择方面, TL₁罚稀疏主成分分析是有效的.综上所述, TL₁罚稀疏主成分分析比稀疏主成分分析更优.