中内存解决我的大标准方程使用`outer`时最小二乘估计

Question

考虑R中下面的例子：

x1 <- rnorm(100000) 
x2 <- rnorm(100000) 
g <- cbind(x1, x2, x1^2, x2^2) 
gg <- t(g) %*% g 
gginv <- solve(gg) 
bigmatrix <- outer(x1, x2, "<=") 
Gw <- t(g) %*% bigmatrix 
beta <- gginv %*% Gw 
w1 <- bigmatrix - g %*% beta 

如果我试图在我的电脑上运行这样的事情，它会抛出内存错误（因为bigmatrix太大）。

你知道我怎么能达到相同的目的，而不会遇到这个问题？

Answer 1

这是一个具有100,000个响应的最小二乘问题。您的bigmatrix是响应（矩阵），beta是系数（矩阵），而w1是残差（矩阵）。

bigmatrix，以及w1，如果形成明确，将每个成本

(100,000 * 100,000 * 8)/(1024^3) = 74.5 GB 

这是过于庞大。

由于对每个响应的估计是独立的，所以实际上不需要一次性形成bigmatrix并尝试将其存储在RAM中。我们可以形成它瓷砖通过瓷砖，并使用迭代过程：形成一个瓷砖，使用瓷砖，然后丢弃它。例如，下面考虑尺寸100,000 * 2,000的瓦，与存储器大小：

(100,000 * 2,000 * 8)/(1024^3) = 1.5 GB 

通过这样的迭代过程，存储器使用实际上是在控制之下。

x1 <- rnorm(100000) 
x2 <- rnorm(100000) 
g <- cbind(x1, x2, x1^2, x2^2) 
gg <- crossprod(g) ## don't use `t(g) %*% g` 
## we also don't explicitly form `gg` inverse 

## initialize `beta` matrix (4 coefficients for each of 100,000 responses) 
beta <- matrix(0, 4, 100000) 

## we split 100,000 columns into 50 tiles, each with 2000 columns 
for (i in 1:50) { 
    start <- 2000 * (i-1) + 1 ## chunk start 
    end <- 2000 * i ## chunk end 
    bigmatrix <- outer(x1, x2[start:end], "<=") 
    Gw <- crossprod(g, bigmatrix) ## don't use `t(g) %*% bigmatrix` 
    beta[, start:end] <- solve(gg, Gw) 
    } 

注意，不要试图计算剩余矩阵w1，因为这将花费74.5 GB。如果您在以后的工作中需要残余矩阵，您仍然应该尝试将其分解为多个瓷砖并逐一工作。

你不需要担心这里的循环。每次迭代中的计算足够昂贵，以分摊循环开销。