2013-03-18 91 views
2

我大約六個月前開始使用R,並且我在R中獲得了一些經驗。最近,我遇到了有關矩陣內子集的問題,希望能夠幫助您制定解決方案我有更高的效率。'R'沒有循環的矩陣子集

我想要做的是以下幾點。假設我有一個矩陣和兩個向量如下:

# matrix 
a <- matrix(seq(1,100,by=1),10,10) 
# vector (first column of matrix a) 
b <- c(2,4,5,6,7,8) 
# vector (column numbers of matrix a) 
c <- c(5,3,1,4,6,2) 

只是重申,

  • 矢量b指矩陣a的第一列。
  • 向量c是指矩陣的列號a

我想獲得tmp99 <- a[b,c:8]。但是,當我這樣做時,我收到以下警告消息。

Warning message: 
In c:8 : numerical expression has 6 elements: only the 
     first used (index has to be scalar and not vector) 

所以,我試着解決問題,使用循環和列表,我得到我想要的解決方案。我假設有一個比這更有效的解決方案。該解決方案是我到目前爲止是這樣的:

a <- matrix(seq(1,100,by=1),10,10) 
b <- c(2,4,5,6,7,8) 
c <- c(5,3,1,4,6,2) 
tmp <- list() 
for (i in 1:length(b)) tmp[[i]] <- c(a[b[i],(c[i]:8)]) 
tmp99 <- t(sapply(tmp, '[', 1:max(sapply(tmp, length)))) 
tmp99[is.na(tmp99)] <- 0 

我想知道什麼是如果有辦法避免使用循環實現上述,因爲我的矩陣尺寸爲200000 x 200,因爲我有做這個很多(在我的問題中,bc被確定爲代碼的另一部分的一部分,所以我不能使用絕對索引號),我想減少相同的時間。任何幫助將不勝感激。謝謝。

+0

這是爲什麼標有'html',只有是什麼? – CBroe 2013-03-18 11:04:04

+0

作爲一般的良好實踐,您可能希望避免通過函數名稱調用變量(如'c') – ds440 2013-03-18 15:09:40

回答

1

以下是使用base程序包執行此操作的一種方法。有可能是更好的解決方案使用data.table但以下工作:)

a <- matrix(seq(1, 100, by = 1), 10, 10) 
b <- c(2, 4, 5, 6, 7, 8) 
c <- c(5, 3, 1, 4, 6, 2) 

res <- t(sapply(X = mapply(FUN = function(b, c) expand.grid(b, seq(from = c, to = 8)), b, c, SIMPLIFY = FALSE), FUN = function(x) { 
    c(a[as.matrix(x)], rep(0, 8 - nrow(x))) 
})) 

res 
##  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] 
## [1,] 42 52 62 72 0 0 0 0 
## [2,] 24 34 44 54 64 74 0 0 
## [3,] 5 15 25 35 45 55 65 75 
## [4,] 36 46 56 66 76 0 0 0 
## [5,] 57 67 77 0 0 0 0 0 
## [6,] 18 28 38 48 58 68 78 0 



# Let's break it down in multiple steps. 

coordinates <- mapply(FUN = function(b, c) expand.grid(b, seq(from = c, to = 8)), b, c, SIMPLIFY = FALSE) 

# below sapply subsets c using each element in coordinates and pads result with additional 0s such that total 8 elements are returned. 

res <- sapply(X = coordinates, FUN = function(x) { 
    c(a[as.matrix(x)], rep(0, 8 - nrow(x))) 
}) 
res 
##  [,1] [,2] [,3] [,4] [,5] [,6] 
## [1,] 42 24 5 36 57 18 
## [2,] 52 34 15 46 67 28 
## [3,] 62 44 25 56 77 38 
## [4,] 72 54 35 66 0 48 
## [5,] 0 64 45 76 0 58 
## [6,] 0 74 55 0 0 68 
## [7,] 0 0 65 0 0 78 
## [8,] 0 0 75 0 0 0 


# you probably need result as traspose 
res <- t(res) 

res 
##  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] 
## [1,] 42 52 62 72 0 0 0 0 
## [2,] 24 34 44 54 64 74 0 0 
## [3,] 5 15 25 35 45 55 65 75 
## [4,] 36 46 56 66 76 0 0 0 
## [5,] 57 67 77 0 0 0 0 0 
## [6,] 18 28 38 48 58 68 78 0 
2

你可以嘗試某種矩陣索引的解決方案,是這樣的。目前尚不清楚實際上是否會更快;在小的情況下,我認爲它肯定會是,但是在大的情況下,創建矩陣到索引的開銷可能需要比遍歷for循環更長的時間。爲了得到更好的答案,編制一個類似於我們可以測試的數據集。

idx.in <- cbind(rep(b, 8-c+1), unlist(lapply(c, function(x) x:8))) 
idx.out <- cbind(rep(seq_along(b), 8-c+1), unlist(lapply(c, function(x) 1:(8-x+1)))) 
tmp99 <- array(0, dim=apply(idx.out, 2, max)) 
tmp99[idx.out] <- a[idx.in] 

這是一個帶有矩陣索引的版本,但是它爲每一行分別進行。這可能會更快,具體取決於要替換的行數和列數。你想避免的是內存不足,for循環可以提供幫助,因爲它不會同時在內存中保存每一步的所有細節。

out <- array(0, dim=c(length(b), 8-min(c)+1)) 
for(idx in seq_along(b)) { 
    out[cbind(idx, 1:(8-c[idx]+1))] <- a[cbind(b[idx], c[idx]:8)] 
} 
out 
+0

非常感謝Aaron,@geektrader,Roland和Arun向我展示瞭如何加速解決方案。我嘗試了4種方法中的3種(還沒有嘗試過Arun的方法),並且它們比當前的'for循環'解決方案更慢和/或需要更多內存。爲了完整性,我有16GB RAM i7系統。與此同時,我將嘗試構建一個數據集,像Aaron建議的那樣,並且發佈相同的內容,看看是否有幫助。謝謝大家花時間爲我提供幫助。我已經明確瞭解解決這個問題的不同方法。 – Ram 2013-03-19 08:10:41

0
tmp <- lapply(seq_len(length(b)),function(i) { 
    res <- a[b[i],c[i]:8] 
    res <- c(res,rep(0,c[i]-1)) 
    res 
               }) 
tmp99 <- do.call("rbind",tmp) 
#  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] 
# [1,] 42 52 62 72 0 0 0 0 
# [2,] 24 34 44 54 64 74 0 0 
# [3,] 5 15 25 35 45 55 65 75 
# [4,] 36 46 56 66 76 0 0 0 
# [5,] 57 67 77 0 0 0 0 0 
# [6,] 18 28 38 48 58 68 78 0