CN104699660A - 概述非线性回归问题中的结构化矩阵 - Google Patents

Abstract

公开了一种用于快速并近似解决结构化回归问题的系统、方法和计算机程序产品。在一方面，该系统、方法和计算机程序产品适用于在各种统计建模设置中自然出现的问题(当设计矩阵是范德蒙德矩阵或是一系列这样的矩阵时)。利用范德蒙德矩阵结构还加速解决了回归问题，实现比“输入稀疏”更快的运行时间。建模框架加速有益于用于解决结构化回归问题的随机回归。

Description

概述非线性回归问题中的结构化矩阵

政府合同

本发明是在政府支持下完成的，合同号为FA8750-12-C-0323，由国防部高级研究计划局授予。政府对本发明享有某些权利。

技术领域

本发明一般涉及预测分析领域，特别地，涉及基于概述(sketch)非线性回归问题中的结构化矩阵的高效和可伸缩的预测分析系统。

背景技术

数据分析应用中普遍出现的一个基本统计问题是非线性回归。对于大的数据集，解决这样的问题从计算上来说很有挑战性。

尽管存在解决线性最小二乘回归和最小绝对偏差问题的随机技术的发展，它们并不会利用在回归问题中经常发现的结构。例如，在多项式拟合问题和某些核回归问题中，设计矩阵的结构化很强，可想象其能导致解决回归问题的更快方法。

需要一种方法来利用这种结构。

发明内容

本发明提供了一种运行对于数字线性代数问题的基于概述的算法的系统和方法，包括多项式拟合和结构化回归。

在一个实施例中，当设计矩阵是范德蒙德矩阵或这种矩阵的序列时，该方法快速并近似地解决结构化回归问题，其中问题以各种统计建模设置自然出现，包括经典多项式拟合问题、加法模型(additive model)和近似，以及用于可伸缩核方法(kernel method)的最近开发的随机技术。

该范德蒙德矩阵结构可被开发，以进一步加速解决回归问题，实现比“输入稀疏”更快的运行时间(即，其写下描述回归问题的矩阵所需的时间)。

在一方面，提供了一种计算机实施的用于解决快速非线性回归和分类问题的方法。该方法包括：利用稀疏嵌入矩阵和结构化随机矩阵来概述输入和输出数据；利用输入和输出数据中的结构来加速概述，所述结构是范德蒙德矩阵形式；并解决概述的数据的回归问题。

在又一个方面，提供了一种用于解决快速非线性回归和分类问题的系统，包括：存储器；以及耦合到存储器的硬件处理器设备，被配置为执行包括以下步骤的方法：利用稀疏嵌入矩阵和结构化随机矩阵来概述输入和输出数据；利用输入和输出数据中的结构来加速概述，所述结果是范德蒙德矩阵形式；并解决概述的数据的回归问题。

在又一个方面，提供了一种用于执行操作的计算机程序产品。计算机程序产品包括存储介质，其可被处理电路读取，并存储由处理电路运行的指令以运行方法。该方法与如上描述的一样。

附图说明

现在将参考以下附图作为非限制性的例子描述实施例。

图1示出了根据在此描述的方法解决的回归问题的概述；

图2A到2C示出了本发明的预测分析方法的方面：图2A示出了应用概述矩阵S的一般方法的第一步骤20；图2B示出了在应用概述矩阵S后改良的回归问题30；且图2C示出了概述矩阵S的应用的描述以及S和A都被结构化以允许增加的处理速度的指示。

图3示出了在一个实施例中，实施处理器“p”55a，...55n的描述50，处理器被配置为将相同的概述矩阵S局部应用到其数据，并将相应的结果传播到中央处理器。

图4示出了在范德蒙德矩阵中使用的数据结构70以有助于实现更快的处理方法，并指示了几个应用，诸如多项式拟合和多变量多项式拟合；

图5示出了在一个实施例中被命名为算法1(StructRegression-2)100的第一概述矩阵乘法算法；

图6示出了在一个实施例中被命名为算法2(StructRegression-1)150的第二概述矩阵乘法算法；

图7示出了概括(generalization)性能表175形式的不同方法的测试误差和训练时间的比较；

图8A-8C示出了概述的示例性性能曲线，其报告了概述在运行时间方面的好处，以及在准确率方面的权衡；且

图9示出了其中可运行本方法的计算系统架构200的示例性硬件配置。

具体实施方式

系统和方法使用某种随机概述和采样转换来大大压缩大数据集，同时保持其关键特性。这允许对在较小的数据概述上很快地运行分析，但达到相同或类似的输出质量，就好象在整个数据集上运行那样。

系统和方法通过应用用于非线性回归的概述/采样概念，并针对特定结构化矩阵的问题，执行预测分析。

图1示出了以代数形式写出的回归问题的概念描述。在图1中，回归问题10是发现矩阵“x”，使得给定回归问题类，矩阵“Ax”与矢量“b”(值列表)尽可能接近：min_x|Ax-b|_p。

在此描述的实施例中，矩阵A是块范德蒙德(block-Vandermonde)结构化矩阵。接近性是根据Lp度量而度量的，它是曲线有多好的一种度量的回归。“lp”回归中的p的不同值对应于不同的误差，例如p＝1是每个数据点对曲线的距离和，p＝2是每个数据点到曲线的距离的平方和。

图2A-2C示出了在此讨论的实施例的方法。

图2A示出了将概述矩阵S应用到矩阵A和矢量b的一般方法的第一步骤20。即，由编程计算机系统实施的方法首先计算以下乘积：S*A和S*b，其分别是矩阵-矩阵乘积和矩阵-矢量乘积。概述矩阵S在每列具有单个非零项。矩阵A具有与其相关的快速矢量-矩阵方法，因为它是范德蒙德结构化矩阵。矩阵SA被缩小尺寸为比“A”还小的矩阵，如同矢量“sb”被缩小尺寸为比“b”还小的矢量。

图2B示出了应用概述矩阵S后的改良的即新的回归问题30，即，解决min_x|SAx-Sb|_p。为了解决该缩减的问题，用于“lp”回归的任何算法可被用作“黑盒”处理。如果“S”是一种概述矩阵，且如果多项式矩阵A具有范德蒙德形式，“SA”矩阵可被更快得多地计算。“Sb”矢量也可被更快得多地计算，这样回归问题整体都被处理地更快得多。

图2C示出了概述矩阵S的应用的描述40，以及S和A都被结构化以允许增加的处理速度的指示45。即，使用下列事实，S*A矩阵-矩阵计算被很快速地执行：(1)S在每列具有单个非零项；以及(2)A具有与其相关的快速矢量-矩阵方法，因为它是范德蒙德矩阵。

更具体地，如以下将详细描述的，考虑以下等式1)中的类型的回归问题类：

$\arg \underset{x\∈C}{\min }{\left|\right|\mathrm{Zx}-b\left|\right|}_p,$ 其中p＝1,2(1)

其中C是凸约束集，Z∈R^n×k是特征样本(sample-by-feature)设计矩阵，而b∈Rⁿ是目标矢量。此后假设相对于数据维度，样本数量很大(n＞＞k)。设置p＝2对应于经典最小二乘回归，而p＝1导致了最小绝对偏差拟合，由于其鲁棒性，其非常重要。约束集合C可纳入正则化。当C＝R^k且p＝2时，可以以时间O(nk log k)+(kε^-1)使用随机化来得到ε-最佳解，当ε不是很小时，这比O(nk²)确定性求解器快很多(对ε的依赖可被改善为O(1/ε)，如果需要更高精确度)。类似地，以时间O(nk log n)+(kε^-1)来运行用于l₁回归的随机求解器。

在许多设置中，使得这样的加速变成可能的是因为存在合适的不经意子空间嵌入(Oblivious Subspace Embedding，OSE)。OSE可被认为是独立于数据的随机“概述”矩阵S∈R^t×n，其在子空间上(例如，在Z,b的列空间上)的近似等距特性意味着：

||S(Zx-b)||_p≈||Zx-b||_p对于所有x∈C,

这转而允许x相对于更小得多尺寸的“概述的”数据集被优化，而不失去解的质量。概述矩阵包括高斯随机矩阵、结构化随机矩阵(其允许经由类似于FFT操作的快速矩阵乘法)等。

在一方面，本方法还不平凡地(non-trivially)利用了Z中的额外结构，以进一步加速运行时间。例如，已知当Z是稀疏的且具有nnz(Z)＝nk个非零项时，使用基于散列的概述矩阵，有可能实现更快得多的“输入稀疏性”运行时间。而且，该方法还提供了更准确的概述技术，其被设计用于非线性和非参数回归。例如，考虑将多项式模型拟合到一组样本(z_i,b_i)∈R×R,i＝1,...,n的基本问题。于是，设计矩阵Z具有可在回归求解器中利用的范德蒙德结构。在大型数据集上估计非参数模型特别有吸引力。

为了能精确描述Z上的结构，这里的描述参考以下定义:

(范德蒙德矩阵)，假设x₀,x₁,...,x_n-1是实数。由V_q,n(x₀,x₁,...,x_n-1)表示的范德蒙德矩阵具有如下形式：

$V_{q,n}(x_{1,}{x}_1,...,x_{n-1})=\left(\begin{array}{cccc}1& 1& ...& 1\\ x_0& x_1& ...& x_{n-1}\\ ...& ...& ...& ...\\ x_0^{q-1}& x_1^{q-1}& ...& x_{n-1}^{q-1}\end{array}\right)$

维数q×n的范德蒙德矩阵仅需要O(n)隐式存储，并允许O(n log²q)矩阵-矢量乘积时间。现在定义以下矩阵运算符T_q，其将矩阵A映射到块范德蒙德结构化矩阵：

(矩阵运算符)给定矩阵A∈R^n×d，定义以下矩阵，T_q(A)＝[V_q,n(A_1,1,...,A_n,1)^T V_q,n(A_1,2,...,A_n,2)^T … V_q,n(A_1,d,...,A_n,d)^T]。在此考虑等式(1)的类型的回归问题，其中Z可被写为：

Z＝T_q(A)   (2)

对于n×d矩阵A，从而k＝dq。运算符T_q通过应用直到q-1阶的单项式转换，将原始数据集A的每个特征(列)扩展为Z的q列。这将块范德蒙德结构给于Z。这样的结构自然出现在多项式插入问题中，但也更广泛地适用于如下讨论的非参数加法模型和核心方法。考虑到这些，一个实施例是解决以下“结构化回归”问题：

结构化回归：给定A和b，以恒定概率输出矢量x′∈C，对于该矢量，针对准确度参数>0，

||T_q(A)x′-b||_p≤(1+ε)||T_q(A)x^*-b||_p

其中x^*＝argmin_x∈C||T_q(A)x-b||_p。

在一个实施例中，系统和方法提供了以下：对于p＝2，以时间O(nnz(A)log²q)+(dqε^-1)或O(nnz(A)log²q)+poly(dq/ε)来解决以上结构化回归问题的算法。通过将概述方法与预调整(pre-conditioned)迭代求解器进行组合，可获取关于ε的对数依赖性。对于p＝1，提供了具有运行时间O(nnz(A)log n log²q)+(dqε^-1log n)的算法。这意味着从线性(即Z＝A)移动到非线性回归(Z＝T_q(A)))仅引起轻度附加的log²q运行时成本，同时不要求额外存储。由于nnz(T_q(A))＝qnnz(A)，这提供了概述方法，其比“输入稀疏”时间运行得更快。对于最小二乘回归及其鲁棒的l₁回归对等物两者，这些算法适用于广泛的线性模型类。尽管多项式插入和具有单项式基函数的加法模型被这里的方法立即覆盖，也可示出在约束集合C的合适选择下，具有用于高斯随机矩阵G的Z＝T_q(AG)结构化回归问题近似使用高斯核的非参数回归。当与用于核心方法的随机傅立叶映射相比时，本系统和方法提供了更灵活的建模框架。实验结果确认了建模框架的实际值，以及概述的加速好处。

在一个实施例中，现在将考虑凸约束集C描述用于块范德蒙德结构l_p回归问题的随机求解器的开发。

给定矩阵M∈R^n×d，令M₁,...,M_d是M的列，且M¹,...,M_n是M的行。将||M||₁定义为M的逐元素l₁范数。即，||M||₁＝∑_i∈[d]||M_i||₁。令是M的弗罗宾尼斯范数。令[n]＝{1,...,n}。

A矩阵的良好调整(well-conditioning)和采样

给定矩阵M∈R^n×d，称M是(α,β,1)良好调整的，如果(1)对于任何x∈R^d，||x||_∞≤β||Mx||₁，以及(2)||M||₁≤α。

假设S是r×n矩阵，从而对于所有的x∈R^d，||Mx||₁≤||SMx||₁≤κ||Mx||₁。令Q·R是SM的QR分解，从而QR＝SM，且Q具有标准正交列。于是MR^-1是()良好调整的。

假设U是n×d矩阵A的(α,β,1)良好调整的基。对于每个i∈[n]，令 $p_i\≥\min (1,\frac{{\left|\right|U_i\left|\right|}_1}{t{\left|\right|U\left|\right|}_1}),$ 其中 $t\≥32\mathrm{\α\β}(d\ln \left(\frac{12}{\mathrm{\ϵ}}\right)+\ln \left(\frac{2}{\mathrm{\δ}}\right))/\left({\mathrm{\ϵ}}^2\right).$ 假设每行以概率p_i来独立采样，并创建对角矩阵S，其中如果i没有被采样，则S_i,i＝0，且如果i被采样，则S_i,i＝1/p_i。于是，以至少1-δ的概率，同时对于所有的x∈R^d，这导致：

|||SAx||₁-||Ax||₁|≤ε||Ax||₁.

可使用如下方法来快速获取对p_i的近似。使得U∈R^n×d是n×d矩阵A的(α,β,1)良好调整的基。假设G是i.i.d.高斯的d×O(log n)矩阵。令对于所有的i，其中t如上定义。然后以概率1-1/n，选择G，出现如下情况：如果以p_i概率来采样每行，且S被创建，则在选择采样的行时以至少1-δ的概率，同时用于所有的x∈R^d，存在：

|||SAx||₁-||Ax||₁|≤||Ax||₁.

这在此被称为“定理6”，如在图6中示出的算法所使用的。

不经意子空间嵌入

此外，令A∈R^n×d并假设n>d，nnz(A)表示A的非零项的数目。还假设nnz(A)≥n，且在A中没有全零行或列，则针对l₂范数：对于参数t，定义了根据如下矩阵族的随机线性映射ΦD:Rⁿ→R^t：

是随机映射，从而对于每个i∈[n]，对于t′∈[t]以概率1/t，h(i)＝t′。

Φ∈{0,1}^t×n是t×n二元矩阵，其中Φ_h(i),i＝1，且所有剩余项为0。

D是n×n随机对角矩阵，其中每个对角项以相同概率被独立地选择为+1或-1。

量Π＝ΦD被称为稀疏嵌入矩阵。

对于特定t，且以对于选择Φ和D至少.99的概率，对于任何固定A∈R^n×d，对于所有的x∈R^d都同时存在：

(1-ε)·||Ax||₂≤||ΠAx||₂≤(1+ε)·||Ax||₂,

即，A的整个列空间被保留。t的最佳已知值是t＝O(d²/²)。

还使用了不经意子空间嵌入，被称为子采样随机阿达玛变换(subsampled randomized Hadamard transform)或SRHT。

对于线性映Π′射存在分布，使得对于Π′的选择以.99的概率，针对任何固定A∈R^n×d，对于所有x∈R^d都同时存在：

(1-ε)·||Ax||₂≤||Π′Ax||₂≤(1+ε)·||Ax||₂,

其中Π′的行数是且计算Π′A的时间是O(ndlogt′)。这在此被称为“定理7”，如在图5中示出的算法所使用的。

对于l₁范数：结果可被概括为针对l₁范数的子空间嵌入。在一个实施例中，目标是设计矩阵Ψ上的分布，使得以至少.99的概率，针对任何固定A∈R^n×d，对于所有x∈R^d都同时存在：

||Ax||₁≤||ΨAx||₁≤κ||Ax||₁,

其中κ>1是失真参数。跟d无关的κ的最佳已知值是κ＝O(d²log²d)。其矩阵Ψ族被选为是Π·E形式，其中Π如上具有参数t＝d^1+γ，对于任意小常数γ>0，且E是具有E_i,i＝1/u_i的对角矩阵，其中u₁,...,u_n是独立标准指数分布的随机变量。

指数分布具有支持x∈[0,∞)、概率密度函数(PDF)f(x)＝e^-x以及累积分布函数(CDE)F(x)＝1-e^-x。如果X是从指数分布中选择的，随机变量X是指数的。

快速范德蒙德乘法

令x₀,...,x_n-1∈R且V＝V_q,n(x₀,...,x_n-1)。对于任何y∈Rⁿ和z∈R^q，可以以O((n+q)log²q)时间来计算矩阵-向量乘积Vy和V^Tz。

主要引理：范数l₂和l₁被分开处理。图5中示出了被命名为算法1(StructRegression-2)100的第一概述矩阵相乘算法。图6示出了被命名为算法2(StructRegression-1)150的第二概述矩阵相乘算法。本方法和系统根据以下讨论的引理使用这些子例程：

对于图5中涉及的引理(引理9)：稀疏概述和T_q(A)的高效相乘：令A∈R^n×d且令Π＝ΦD为用于l₂范数(其具有对于任意值t的相关散列函数h:[n]→[t])的稀疏嵌入矩阵，并且令E是任意对角矩阵。存在图5的确定性算法100以O((nnz(A)+dtq)log²q)时间来计算乘积Φ·D·E·T_q(A)。

图5示出了用于lp-2回归的被命名为算法1(StructRegression-2)100的第一算法(例如，对应于每个数据到曲线的距离的平方和)。

对于图6中涉及的引理(引理10)：右侧的T_q(A)的高效相乘：令A∈R^n×d。对于任意矢量z，存在确定性算法以O((nnz(A)+dq)log²q)时间来计算矩阵矢量乘积T_q(A)·z。

对于引理(引理11)：左侧的T_q(A)的高效乘积：假设A∈R^n×d。对于任意矢量z，存在确定性算法以O((nnz(A)+dq)log²q)时间来计算矩阵矢量乘积z·T_q(A)。

如上所述，对于快速l₂回归，考虑可以由图5的算法100解决的p＝2情形下的结构化回归问题。即，算法StuctRegression-2以下列时间解决了p＝2的结构化回归：

O(nnz(A)log²q)+(dq/ε)或O(nnz(A)log²q)+poly(dq/ε).

证明：根据稀疏嵌入矩阵的性质(参见不经意子空间嵌入)，以至少0.99的概率，对于t＝O((dq)²/ε²)，对于与b结合的T_q(A)的列范围内所有y，同时存在：

(1-ε)||y||₂≤||Πy||₂≤(1+ε)||y||₂,

因为该空间的范围具有至少dq+1的维度。通过先前引用定理7，以0.99的概率，对于的列范围内的所有矢量z，同时存在：

(1-ε)||z||₂≤||Π'z||₂≤(1+ε)||z||₂.

于是对于所有的矢量x∈R^d,

(1-O(ε))||T_q(A)x-b||₂≤||Pi'Π(T_q(A)x-B)||₂≤(1+O(ε))||T_q(A)x-b||₂.

于是，通过至少.98概率的联合界(union bound)，StructRegression-2算法100的输出是(1+ε)近似。

对于时间复杂度，可以通过引理9以O((nnz(A)+dtq)log²q)来计算ΠT_q(A)，而可以以O(n)时间来计算Πb。可以以(dq/ε)时间来执行剩余步骤，于是总时间为O(nnz(A)log²q)+(dq/ε)。或者，可以以poly(dq/ε)时间来执行这些步骤，因此总时间为O(nnz(A)log²q)+poly(dq/ε)。

对1/ε的对数依赖

StructRegression-2算法100可被修改，以通过将基于概述的方法与迭代方法相组合来得到对ε有对数依赖的运行时间。

算法I以时间O((nnz(A)+dq)log(1/ε))+(dq)解决了p＝2的结构化回归问题。或者，可以以O((nnz(A)+dq)log(1/ε))+poly(dq)的总时间来解决该问题。

快速l₁回归

现在考虑p＝1情形下的结构化回归。该情形下的算法150(图6)比p＝2的情形更复杂。

算法StructRegression-1以时间O(nnz(A)log n log²q)+(dqε^-1log n).解决了p＝1的问题中的结构化回归。

需要注意，当存在凸约束集合C时，上述算法中的唯一改变是在x′∈C上进行优化。

多项式拟合、加法模型和随机傅立叶映射

在一个实施例中，用于各类块范德蒙德结构化回归问题的各种概述方法通过表明这些问题自然出现在各种统计建模设置中而得到促进。

最基本的应用是一维(d＝1)多项式回归。

在多变量加法回归模型中，连续的目标变量y∈R和输入变量z∈R^d通过模型来关联，其中μ是截距项，ε_i为均值为零的高斯误差项，且f_i是平滑单变量函数。使用基函数h_i,t(·)，每个函数可被扩展为并且典型地通过受限或受罚最小二乘模型来估算未知的参数矢量其中，在训练样本(z_i,y_i),i＝1...n上，b＝(y₁...y_n)^T，并且对于(H_i)_j,t＝h_i,t(z_j)，Z＝[H₁...H_q]∈R^n×dq。约束集合C典型地带来平滑、稀疏和组稀疏约束。容易看到，选择单项式基h_i,s(u)＝u^s立即将设计矩阵Z映射到等式2)的结构化回归形式。对于p＝1，这里展示的算法还提供了用于鲁棒多项式加法模型的快速求解器。

加法模型带来了单变量非线性的受限形式，其忽略了协变量之间的交互。将交互项表示为其中∑_iα_i＝q，α_i∈{0...q}，通过{z^α,α∈{0,...q}^d,∑_iα_i≤q}来扩展q阶多变量多项式函数空间P_q。P_q允许所有可能的q阶交互，但具有维度d^q，这在计算上是不可能显式地操作的，除非对于低阶和低维度或稀疏数据集。具有多项式核心k(z,z′)＝(z^Tz′)^q＝∑_αz^αz′^α的核心方法提供了隐含的机制来计算与P_q关联的特征空间中的内积。但是，它们需要O(n³)的计算来解决相关的核化(脊状)回归问题以及稠密n×n格拉姆矩阵K(给定K_ij＝k(z_i,z_j))的O(n²)存储，且因此不能很好地伸缩。

对于d×D矩阵G，令S_G是由下式来跨越的子空间

$\{{\left(\underset{i=1}{\overset{d}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{z}_i\right)}^t,t=1...q,j=1...D\}.$

假设D＝d^q并且G是i.i.d高斯变量的随机矩阵，则几乎确切地存在S_G＝P_q。于是直观地吸引人的显式可伸缩方法是使用D＝d^q。在该情形下，S_G基本上跨越P_q的随机空间。用于解决受限于S_G的多变量多项式回归的设计矩阵具有Z＝T_q(AG)的形式，其中

该方案可涉及近似平移不变核心函数的上下文中的随机傅立叶特征的概念，以高斯核作为主要例子。根据博赫纳定理(Bochner’s Theorem)，表明高斯核是具有协方差矩阵σ^-1I_d的均值为零的多变量高斯分布的傅立叶变换，其中，I_d表示d-维单位矩阵，

$k(z,z^{\′})=\exp (-{\left|\right|z-z^{\′}\left|\right|}_2^2/2{\mathrm{\σ}}^2)=E_{\mathrm{\ω}:N(0_d,{\mathrm{\σ}}^{-1}{I}_d)}\left[{\mathrm{\φ}}_{\mathrm{\ω}}\left(z\right){\mathrm{\φ}}_{\mathrm{\ω}}{\left(z^{\′}\right)}^{*}\right]$

其中，φ_ω(z)＝e^iω′z。可以通过采样D频率ω:N(0_d,σ^-1I_d)并设置来获取该期望的经验近似。这意味着可以用K≈RR^T的高强度来近似高斯核的格拉姆矩阵其中，R＝[cos(AG)sin(AG)]∈R^n×2D(作为标量函数来逐元素地应用正弦和余弦)。该高斯核的随机显式特征映射意味着以R作为设计矩阵的标准线性回归于是可被用于获取时间O(nD²)的解决方案。通过对正弦和余弦进行直到q阶的迈克劳林级数(Maclaurin series)展开，可以看到argmin_x∈range(Q)||T_q(AG)x-b||_p(其中，矩阵Q∈R^2Dq×2D包含迈克劳林级数的适当系数)形式的受限结构化的回归问题将非常近似于随机傅立叶特征构造。通过放弃或修改约束集合x∈range(Q)，上述设置原则上可以定义更丰富的模型类。

现在描述示例性分类和回归数据集的两组试验。第一组试验将结构化非线性最小二乘回归模型的概括性能与标准线性回归进行比较，并将非线性回归与随机傅立叶特征进行比较。第二组试验集中在概述的可伸缩性收益。

图7示出了形式为概括性能表175的不同方法的测试误差和训练时间的比较。

在表中，n是训练实例的数量，d是每个实例的特征数量，且k是测试集中的实例的数量。“普通回归”表示普通的`2回归。“加法多项式回归”表示加法多项式`2回归。针对分类任务，报告了不正确预测的测试点的百分比。对于回归任务，报告了||y_p-y||₂/||y||，其中，y_p是预测值，且y是地面真值(ground truth)。

在列176所示的普通l₂线性回归中非常快，特别是如果矩阵是稀疏的。但是，它仅提供了平庸的结果。通过列177的加法多项式回归，该结果有所改进。列177的加法多项式回归保持稀疏结构，从而它可以利用快速的稀疏求解器。一旦在列178引入随机特征(D)，由此引入交互项，结果明显改进。在与179的如图7所示的随机傅立叶特征进行比较时，针对相同数量的随机特征D，具有随机特征的加法多项式回归比具有随机傅立叶特征的回归得到更好的结果。如果随机特征的数量不同，则如果D_Fourier＝D_Poly·q(其中D_Fourier是傅立叶特征的数量，且D_Poly是加法多项式回归中的随机特征的数量)，则与具有随机特征的加法多项式回归相比，具有随机傅立叶特征的回归表现出改进。但是，计算随机特征是最昂贵的步骤之一，因此需要以更少的随机特征来计算更好的近似。

图8A-8C示出了概述的示例性性能图，并报告了概述在运行时间方面的好处，以及在准确率方面的权衡。在该试验中，使用具有300,000个样本的mnist数据集的更大采样。计算了1,500随机特征的样本数量，然后通过精确地且概述到不同数量的行、以q＝4来解决相应的加法多项式回归问题。另外测试了基于采样的方法，其简单随机地采样了行的子集(没有概述)。图8A是相对于精确解的概述方法的加速的图180。在这些试验中，使用非最优直接实现，其不会利用快速范德蒙德乘法或并行处理。因此，使用顺序执行来测量运行时间。在一个实施例中，测量解决回归问题所需的时间。为了此试验，使用了具有两个四核Intel E54102.33GHz和32GBDDR2800MHz RAM的机器。图8B以图185来探索解决回归问题的次优性。更具体而言，绘制了其中，Y是标签矩阵，是最佳近似(精确解)，且Y_p是概述的解。可以看到误差随着采样的矩阵大小增加而减少，并且使用不太大的概述大小，存在大约大10％的目标。在图8C中，看到在图190中这变成了误差率的增加。令人鼓舞的是，小到15％的例子数量的概述足以具有非常小的误差率增长，而仍然能更快超过5倍地解决回归问题(对于更大的数据集，期望加速会增加)。

于是，在一个实施例中，提供了减少浮点数运算并减少并行和分布式系统中的通信成本的方法。

此外，所使用的方法通过以这样的方式降低大小或“概述”数据从而近似地保持下游算法的统计性能，来做到这一点。大小的降低带来FLOPS的减少。概述的数据集可以适合于单个节点或集群中的多个节点的存储器。

图3示出了包含多个处理器55a,…55n的处理系统50的示例性实现，其被配置将相同的概述矩阵S局部应用到其数据，并将相应的结果发送到中央处理器，该中央处理器被配置为计算并获取局部计算机器上的整个矩阵60的概述。这涉及到每个处理器针对矩阵A的其本地份额Aⁱ来计算S*(Aⁱ)，其中，

该方法同时不再需要进行内核外部的计算以及数据的I/O扫描成本。例如，可以通过概述数据以使它适合于单个节点或聚合存储器，来避免通信或内核外部解决方案。

这里描述的方法可以直接加速很多种实际预测分析系统，例如，语音识别、图像分类、文本分类等。

即，在数据矩阵具有特定结构的涉及快速非线性回归和分类的预测分析系统例如语音识别、图像分类、文本分类等系统中，处理方法包括：使用稀疏嵌入矩阵(散列)和结构化随机矩阵来概述数据(输入、输出)；利用数据中的结构来加速概述。

图4示出了在范德蒙德矩阵中使用的示例性数据结构70，其有助于实现更快的处理方法，并指示几种应用，例如多项式拟合和多变量多项式拟合。

这里描述的方法减少了FLOPS、通信成本和/或内核外部通信的数据扫描成本，并在被概述的数据上解决了简化的回归问题。在一个实施例中，从原始数据特征(加法模型)来隐式地生成结构化矩阵。在另一实施例中，从原始数据特征(交互项)的随机组合来隐式地生成结构化矩阵，其中，建模问题是最小二乘回归(l2)或鲁棒回归(l1)。对于鲁棒回归l1，适用修改步骤。例如，这种处理解决下列类型的问题：(1)例如形式的多项式拟合；(2)加法模型(多变量多项式拟合)，例如Y＝以及(3)交互项，例如X_i＝(ZG)_i。这些技术可以扩展到其他结构化矩阵以及其他类的问题(例如低秩近似)。

图9示出了计算或移动设备计算系统基础设施200的示例性硬件配置，在其中运行本方法。一方面，计算系统200从映射应用中的用户输入搜索查询接收搜索结果，并被编程为执行图2A-2C、3等的方法处理步骤。硬件配置优选地具有至少一个处理器或中央处理单元(CPU)211。CPU 211通过系统总线212互联到随机存取存储器(RAM)214、只读存储器(ROM)216、输入/输出(I/O)适配器218(用于将外围设备例如磁盘单元221和磁带驱动器240连接到总线212)、用户接口适配器222(用于将键盘224、鼠标226、扬声器228、磁盘驱动器设备232和/或其他用户接口设备连接到总线212)、通信适配器234(用于将系统200连接到数据处理网络、互联网、内联网、局域网(LAN)等、以及显示适配器236(用于将总线212连接到显示设备238和/或打印机239(例如数字打印机等))。

本发明可以是一种系统、方法和/或计算机程序产品。计算机程序产品可包括计算机可读存储介质(或媒介)，其上有计算机可读程序指令以使得处理器能执行本发明的方面。

计算机可读存储介质可以是有形设备，其可保持并存储供指令执行设备使用的指令。计算机可读存储介质可以是包括但不限于电存储设备、磁存储设备、光存储设备、电磁存储设备、半导体存储设备或前述任意设备的合适组合。计算机可读存储介质的更具体的例子的非穷尽列表包括以下：便携式计算机软盘、硬盘、随机存取存储器(RAM)、只读存储器(ROM)、可擦除可编程只读存储器(EPROM或闪存存储器)、静态随机存取存储器(SRAM)、便携式紧致盘只读存储器(CD-ROM)、数字多功能盘(DVD)、存储棒、软盘、诸如冲卡的机械编程设备或其上记录有指令的凹槽中的凸起结构，以及前述设备的任意合适的组合。在此使用的计算机可读存储介质不应被认为是暂时信号本身，诸如无线波或其他自由传播的电磁波，通过波导或其他传输介质(例如穿过光纤电缆的光脉冲)传播的电磁波，或通过电线被传输的电信号。

在此描述的计算机可读程序指令可经由网络，例如因特网、局域网、广域网和/或无线网，从计算机可读存储介质被下载到各个计算/处理设备，或到外部计算机或外部存储设备。网路可包括铜传输线缆、光传输光纤、无线传输、路由器、防火墙、交换机、网关计算机和/或边缘服务器。每个计算/处理设备中的网络适配卡或网络接口从网络接收计算机可读程序指令，并转发计算机可读程序指令以存储在各个计算/处理设备中的计算机可读存储介质中。

用于执行本发明的操作的计算机可读程序指令可以是汇编指令、指令系统结构(ISA)指令、机器指令、机器相关指令、伪代码、固件指令、状态设定数据，或以一种或多种编程语言的任意组合编程的其他源代码或目标代码，所述编程语言包括面向对象的编程语言—诸如Java、Smalltalk、C++，还包括常规的过程式编程语言—诸如”C”语言或类似的编程语言。程序代码可以完全地在用户计算机上执行、部分地在用户计算机上执行、作为一个独立的软件包执行、部分在用户计算机上部分在远程计算机上执行、或者完全在远程计算机或服务器上执行。在后一种情形中，远程计算机可通过任意类型的网络(包括局域网(LAN)或广域网(WAN))连接到用户计算机，或可连接到外部计算机(例如，使用因特网服务提供商通过因特网)。在一些实施例中，电路包括例如可编程逻辑电路、场可编程门阵列(FPGA)或可编程逻辑阵列(PLA)，可通过使用计算机可读程序指令的状态信息以个性化电路来执行计算机可读程序指令，以便执行本发明的方面。

在此可参考根据本发明实施例的方法、设备(系统)和计算机程序产品的流程图描述和/或框图来描述本发明的方面。将理解，流程图描述和/或框图的每个框，以及流程图描述和/或框图的结合，可被计算机可读程序指令实施。

这些计算机可读程序指令可被提供给通用计算机、专用计算机或其他可编程数据处理设备的处理器，以生成机器，这样经由计算机的处理器或其他可编程数据处理设备的指令可创建执行在流程图和/或框图的框中执行的功能/动作的装置。这些计算机可读程序指令也可被存储在计算机可读存储介质中，其可指示计算机、可编程数据处理设备和/或其他设备以特定方式起作用，这样其中存储有指令的计算机可读存储介质包括制造品，其包括实现在流程图和/或框图的框中指定的功能/行为的方面的指令。

也可以把计算机程序指令加载到计算机、其它可编程数据处理装置、或其它设备上，使得在计算机、其它可编程数据处理装置或其它设备上执行一系列操作步骤，以产生计算机实现的过程，从而使得在计算机或其它可编程装置上执行的指令能够实现流程图和/或框图中的方框中规定的功能/操作。

附图中的流程图和框图显示了根据本发明的多个实施例的系统、方法和计算机程序产品的可能实现的体系架构、功能和操作。在这点上，流程图或框图中的每个方框可以代表一个模块、程序段或代码的一部分，所述模块、程序段或代码的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令。也应当注意，在有些作为替换的实现中，方框中所标注的功能也可以以不同于附图中所标注的顺序发生。例如，两个连续的方框实际上可以基本并行地执行，它们有时也可以按相反的顺序执行，这依所涉及的功能而定。也要注意的是，框图和/或流程图中的每个方框、以及框图和/或流程图中的方框的组合，可以用执行规定的功能或操作的专用的基于硬件的系统来实现，或者可以用专用硬件与计算机指令的组合来实现。

以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。因此本技术领域的普通技术人员可进行许多修改和变更，而不偏离所附权利要求所限定的本发明的范围和精神。

Claims (16)

1.一种计算机实施的用于解决快速非线性回归和分类问题的方法，包括：

使用稀疏嵌入矩阵和结构化随机矩阵来概述输入数据和输出数据；

利用数据中的结构来加速概述，所述结构是范德蒙德矩阵形式；以及

在概述的数据上解决回归问题，硬件处理器执行一个或多个所述概述、结构利用和问题解决。

2.如权利要求1所述的计算机实施的方法，其中回归问题是min_x|Ax-b|_p的形式，其中x是矩阵，A是块范德蒙德结构化矩阵，b是矢量，而p是1范数或欧几里得范数。

3.如权利要求1所述的计算机实施的方法，还包括：

生成对于所述回归问题的输出x’，其中输出x’满足|Ax'-b|_p≤(1+eps)min_x|Ax-b|_p，其中eps>0是用户指定的精确度参数。

4.如权利要求1所述的计算机实施的方法，其中回归问题是min_x|T_q(A)x-b|_p的形式，其中x是矩阵，A是任意(n x d)矩阵，且T_q(A)通过用q-元组(1,A_i,j,A_i,j ²,…,A_i,j ^q-1)替换每个项A_i,j将A扩展为(n x(dq))矩阵。

5.如权利要求4所述的计算机实施的方法，还包括：

生成对于回归问题的输出x’，其中输出x’满足|Ax'-b|_p≤(1+eps)min_x|Ax-b|_p，其中eps(ε)>0是用户指定的精确度参数。

6.如权利要求5所述的计算机实施的方法，还包括：

在根据O(nnz(A)log²q)+poly(dq/ε)的时间中以p＝2来解决所述回归问题，其中nnz(A)表示矩阵A的非零项的数量。

7.如权利要求5所述的计算机实施的方法，还包括：

在根据O((nnz(A)+dq)log(1/ε))+poly(dq)的时间中以p＝2来解决所述回归问题，其中nnz(A)表示矩阵A的非零项的数量。

8.如权利要求5所述的计算机实施的方法，还包括：

在根据O(nnz(A)lognlog²q)+(dqε^-1logn)的时间中以p＝1来解决所述回归问题，其中nnz(A)表示矩阵A的非零项的数量。

9.一种用于解决快速非线性回归和分类问题的系统，包括：

存储器；

硬件处理器，其耦合到存储器并被配置为执行包括如下步骤的方法：

使用稀疏嵌入矩阵和结构随机矩阵来概述输入数据和输出数据；

利用输入数据和输出数据中的结构来加速概述，所述结构是范德蒙德矩阵形式；以及

在概述的数据上解决回归问题。

10.如权利要求9所述的系统，其中回归问题是minx|Ax-b|_p的形式，其中x是矩阵，A是块范德蒙德结构化矩阵，b是矢量，而p是1范数或欧几里得范数。

11.如权利要求9所述的系统，其中，所述硬件处理器还被配置为：

生成对于所述回归问题的输出x’，其中输出x’满足|Ax'-b|_p≤(1+eps)min_x|Ax-b|_p，其中eps>0是用户指定的精确度参数。

12.如权利要求9所述的系统，其中回归问题是min_x|T_q(A)x-b|_p的形式，其中x是矩阵，A是任意(n x d)矩阵，且T_q(A)通过用q-元组(1,A_i,j,A_i,j ²,…,A_i,j ^q-1)替换每个项A_i,j将A扩展为(n x(dq))矩阵。

13.如权利要求9所述的系统，其中，所述硬件处理器还被配置为：

生成对于回归问题的输出x’，其中输出x’满足|Ax'-b|_p≤(1+eps)min_x|Ax-b|_p，其中eps(ε)>0是用户指定的精确度参数。

14.如权利要求9所述的系统，其中，所述硬件处理器还被配置为：

在根据O(nnz(A)log²q)+poly(dq/ε)的时间中以p＝2来解决所述回归问题，其中nnz(A)表示矩阵A的非零项的数量。

15.如权利要求9所述的系统，其中，所述硬件处理器还被配置为：

在根据O((nnz(A)+dq)log(1/ε))+poly(dq)的时间中以p＝2来解决所述回归问题，其中nnz(A)表示矩阵A的非零项的数量。

16.如权利要求9所述的系统，其中，所述硬件处理器还被配置为：

在根据O(nnz(A)lognlog²q)+(dqε^-1logn)的时间中以p＝1来解决所述回归问题，其中nnz(A)表示矩阵A的非零项的数量。