以R爲基礎減少線性迴歸

我正在R for Hadoop中進行分佈式線性迴歸計算，但在實現它之前，我想驗證一下我的計算是否與lm函數的結果一致。以R爲基礎減少線性迴歸

我有以下功能，試圖實施Andrew Ng等人討論的通用「求和」框架。在論文Map-Reduce for Machine Learning on Multicore中。

對於線性迴歸，這涉及到映射每一行Y_I和X_I到P_I和Q_I使得：

P_i = x_i * transpose(x_i) 
Q_i = x_i * y_i

然後減少來求解係數，θ-： theta = (sum(P_i))^-1 * sum(Q_i)

R的功能，以做到這一點是：

calculate_p <- function(dat_row) { 
    dat_row %*% t(dat_row) 
} 

calculate_q <- function(dat_row) { 
    dat_row[1,1] * dat_row[, -1] 
} 

calculate_pq <- function(dat_row) { 
    c(calculate_p(matrix(dat_row[-1], nrow=1)), calculate_q(matrix(dat_row, nrow=1))) 
} 

map_pq <- function(dat) { 
    t(apply(dat, 1, calculate_pq)) 
} 

reduce_pq <- function(pq) { 
    (1/sum(pq[, 1])) * apply(pq[, -1], 2, sum) 
}

您可以通過運行實現它的一些合成數據：

X <- matrix(rnorm(20*5), ncol = 5) 
y <- as.matrix(rnorm(20)) 
reduce_pq(map_pq(cbind(y, X))) 
[1] 0.010755882 -0.006339951 -0.034797768 0.067438662 -0.033557351 
coef(lm.fit(X, y)) 
      x1   x2   x3   x4   x5 
-0.038556283 -0.002963991 -0.195897701 0.422552974 -0.029823962

不幸的是，輸出不匹配，所以顯然我做錯了什麼。任何想法如何解決它？

來源

2012-10-10 Erik Shilts

你看過'biglm'包是否能爲你提供一些方法嗎？ – ako

我知道我可以將塊發送到'biglm'並使用'update'函數來更新結果，但它要求每個塊中的每個因子級別都可用。我將在一大組數據上運行這個工具，因此爲了確保足夠的覆蓋率而對塊進行隨機化的開銷將是太大的實施難題。 –

你在reduce_pq中的逆需要是一個矩陣逆。我也改變了一些功能。

calculate_p <- function(dat_row) { 
    dat_row %*% t(dat_row) 
} 

calculate_q <- function(dat_row) { 
    dat_row[1] * dat_row[-1] 
} 

calculate_pq <- function(dat_row) { 
    c(calculate_p(dat_row[-1]), calculate_q(dat_row)) 
} 

map_pq <- function(dat) { 
    t(apply(dat, 1, calculate_pq)) 
} 

reduce_pq <- function(pq) { 
    solve(matrix(apply(pq[, 1:(ncol(X) * ncol(X))], 2, sum), nrow=ncol(X))) %*% apply(pq[, 1:ncol(X) + ncol(X)*ncol(X)], 2, sum) 
} 


set.seed(1) 
X <- matrix(rnorm(20*5), ncol = 5) 
y <- as.matrix(rnorm(20)) 

t(reduce_pq(map_pq(cbind(y, X)))) 
      [,1]  [,2]  [,3]  [,4]  [,5] 
[1,] 0.1236914 0.2482445 0.5120975 -0.1104451 -0.04080922 

coef(lm.fit(X,y)) 
     x1   x2   x3   x4   x5 
0.12369137 0.24824449 0.51209753 -0.11044507 -0.04080922 

> all.equal(as.numeric(t(reduce_pq(map_pq(cbind(y, X))))), as.numeric(coef(lm.fit(X,y)))) 
[1] TRUE

來源

2012-10-11 04:20:53 Sameer

你們都可能想看'crossprod'和'tcrossprod'，因爲這些功能都集中精力優化。 –

你能解釋一下reduce_pq步驟嗎？ P的總和應該是一個標量，所以我不明白爲什麼減少步驟的前半部分有'1（ncol（X）* ncol（X））'。 –

由於你有5個變量，P是一個5X5矩陣。 P的總和也是一個5X5的矩陣。它們存儲在前25列的pq中。這就是爲什麼你需要從pq中選擇1：（ncol（X）* ncol（x）），以便獲得可以塑造成5x5矩陣的東西。 – Sameer

以R爲基礎減少線性迴歸

回答

相關問題