4

所以我有这些巨大的矩阵 X 和 Y。X 和 Y 都有 1 亿行,X 有 10 列。我正在尝试用这些矩阵实现线性回归,我需要数量(X^T*X)^-1 * X^T * Y。我怎样才能尽可能节省空间地计算它?

现在我有

X = readMatrix("fileX.txt")
Y = readMatrix("fileY.txt")
return (X.getT() * X).getI() * X.getT() * Y

这里有多少矩阵存储在内存中?是否同时存储两个以上的矩阵?有更好的方法吗?

我有大约 1.5 GB 的内存用于这个项目。如果我关闭所有其他程序,我可能可以将其拉伸到 2 或 2.5。理想情况下,该过程也会在很短的时间内运行,但内存限制更加严格。

我尝试过的另一种方法是将计算的中间步骤保存为文本文件,并在每一步后重新加载它们。但这很慢。

4

3 回答 3

3

X 的大小是 100e6 x 10 Y 的大小是 100e6 x 1

所以最终大小(X^T*X)^-1 * X^T * Y是 10 x 1

您可以通过以下步骤计算它:

  1. 计算a = X^T*X-> 10 x 10
  2. 计算b = X^T*Y-> 10 x 1
  3. 计算a^-1 * b

第 3 步中的矩阵非常小,因此您只需要执行一些中间步骤即可计算 1 和 2。

例如,您可以读取 X 和 Y 的第 0 列,并通过numpy.dot(X0, Y).

对于 float64 dtype,X0 和 Y 的大小约为 1600M,如果无法容纳内存,您可以分别调用 numpy.dot 两次 X0 & Y 的前半部分和后半部分。

所以计算X^T*Y你需要调用 numpy.dot 20 次,计算X^T*X你需要调用 numpy.dot 200 次。

于 2012-10-31T00:42:50.963 回答
2

普通最小二乘回归的一个简洁属性是,如果您有两个数据集 X1、Y1 和 X2、Y2,并且您已经计算了所有

  • X1' * X1
  • X1' * Y1
  • X2' * X2
  • X2' * Y2

你现在想要对组合数据集 X = [X1; 进行回归。X2] 和 Y = [Y1; Y2],您实际上不必重新计算太多。关系

  • X' * X = X1' * X1 + X2' * X2
  • X' * Y = X1' * Y1 + X2' * Y2

持有,所以有了这些计算,你只需计算

  • β = inv(X' * X) * (X' * Y)

你就完成了。这导致了一个在非常大的数据集上用于 OLS 的简单算法:

  • 加载部分数据集(例如,前一百万行)并计算 X' * X 和 X' * Y(它们是非常小的矩阵)并存储它们。
  • 继续对接下来的一百万行执行此操作,直到处理完整个数据集。
  • 将您存储的所有 X' * Xs 和 X' * Ys 加在一起
  • 计算 beta = inv(X' * X) \ (X' * Y)

这并不比一次加载整个数据集慢得多,而且它使用的内存要少得多。

最后一点:你不应该通过首先计算 (X' * X) 并找到它的逆来计算 beta(有两个原因 - 1. 它很慢,2. 它容易出现数值错误)。

相反,您应该解决线性系统 -

  • (X' * X) * 贝塔 = X' * Y

在 MATLAB 中,这是一个简单的单行

beta = (X' * X) \ (X' * Y);

我希望 numpy 具有类似的解决线性系统的方法,而无需反转矩阵。

于 2015-01-12T21:07:53.127 回答
1

RAM 相当便宜——你应该考虑投资。具有 24 Gig RAM 的系统不再需要花费一臂之力 - 戴尔的低端服务器之一可以容纳这么多。

如果矩阵是稀疏的(很多零),请使用稀疏矩阵类来节省大量 RAM。

如果矩阵不是稀疏的,您要么需要更多 RAM(或至少更多虚拟内存),要么使用磁盘文件进行矩阵运算。

磁盘文件当然比 RAM 慢一个数量级,并且根据您的访问模式,颠簸您的虚拟内存系统实际上可能比这更糟。

于 2012-10-30T22:11:17.520 回答