計算稀疏成對距離矩陣R

我有一個NxM矩陣，我想計算M點之間的歐幾里得距離矩陣。在我的問題中，N約爲100,000。由於我打算將這個矩陣用於k-最近鄰算法，我只需要保持最小距離，所以得到的矩陣非常稀疏。這與例如dist()的結果相反，這將導致密集的矩陣（並且對於我的尺寸N可能存在存儲問題）。計算稀疏成對距離矩陣R

迄今爲止我發現的kNN包（knnflex，kknn等）都顯示爲使用密集矩陣。另外，Matrix包不提供成對距離功能。

接近我的目標，我看到spam包有一個nearest.dist()函數，允許人們只考慮小於某個閾值的距離，delta。然而，在我的情況下，特定值delta可能會產生太多距離（因此我必須密集存儲密度矩陣）或距離太短（以至於我不能使用kNN）。

我已經見過以前的討論，試圖使用bigmemory/biganalytics軟件包來執行k-means clustering，但在這種情況下似乎沒有可以利用這些方法。

有沒有人知道一個函數/實現將在R中以稀疏方式計算距離矩陣？我的（可怕的）備份計劃是有兩個for循環並將結果保存在Matrix對象中。

來源

2011-04-06 Christopher DuBois

只要確保...你知道'dist' http：// stat。 ethz.ch/R-manual/R-patched/library/stats/html/dist.html，對嗎？ – Benjamin 2011-04-06 17:08:21

對不起，我不清楚爲什麼dist（）不適合我的情況。它導致了一個稠密的矩陣，並且存儲NxN矩陣有點煩人。 – 2011-04-07 01:10:41

您應該或者接受其中一個答案，您認爲它實際上回答了問題（如果您認爲它最合適，那麼是您自己的問題），或者編輯您的問題以澄清問題的原因。 – Tommy 2011-07-25 16:45:58

好了，我們不能讓你訴諸for循環，我們現在可以:)

當然有如何表示稀疏矩陣的問題。一個簡單的方法是讓它只包含最近點的索引（並根據需要重新計算）。但在下面的溶液中，我把二者的距離（「D1」等）和索引（「I1」等）在一個單一的矩陣：

sparseDist <- function(m, k) { 
    m <- t(m) 
    n <- ncol(m) 
    d <- vapply(seq_len(n-1L), function(i) { 
     d<-colSums((m[, seq(i+1L, n), drop=FALSE]-m[,i])^2) 
     o<-sort.list(d, na.last=NA, method='quick')[seq_len(k)] 
     c(sqrt(d[o]), o+i) 
     }, numeric(2*k) 
    ) 
    dimnames(d) <- list(c(paste('d', seq_len(k), sep=''), 
     paste('i', seq_len(k), sep='')), colnames(m)[-n]) 
    d 
}

嘗試出來9 2D點：

> m <- matrix(c(0,0, 1.1,0, 2,0, 0,1.2, 1.1,1.2, 2,1.2, 0,2, 1.1,2, 2,2), 
       9, byrow=TRUE, dimnames=list(letters[1:9], letters[24:25])) 
> print(dist(m), digits=2) 
    a b c d e f g h 
b 1.1        
c 2.0 0.9       
d 1.2 1.6 2.3      
e 1.6 1.2 1.5 1.1     
f 2.3 1.5 1.2 2.0 0.9    
g 2.0 2.3 2.8 0.8 1.4 2.2   
h 2.3 2.0 2.2 1.4 0.8 1.2 1.1  
i 2.8 2.2 2.0 2.2 1.2 0.8 2.0 0.9 
> print(sparseDist(m, 3), digits=2) 
    a b c d e f g h 
d1 1.1 0.9 1.2 0.8 0.8 0.8 1.1 0.9 
d2 1.2 1.2 1.5 1.1 0.9 1.2 2.0 NA 
d3 1.6 1.5 2.0 1.4 1.2 2.2 NA NA 
i1 2.0 3.0 6.0 7.0 8.0 9.0 8.0 9.0 
i2 4.0 5.0 5.0 5.0 6.0 8.0 9.0 NA 
i3 5.0 6.0 9.0 8.0 9.0 7.0 NA NA

並試圖解決更大的問題（10k點）。儘管如此，在100k點和更多維度上需要很長時間（比如15-30分鐘）。

n<-1e4; m<-3; m=matrix(runif(n*m), n) 
system.time(d <- sparseDist(m, 3)) # 9 seconds on my machine...

P.S.剛纔注意到，當我寫這篇文章時，你發佈了一個答案：這裏的解決方案大概是兩倍的速度，因爲它不會計算兩次相同的距離（點1和點13之間的距離與點13和點1之間的距離相同）。

來源

2011-04-06 16:19:07 Tommy

感謝您的回答。我同意它快兩倍。然而，對於我的應用（kNN），我認爲只有距離矩陣的上三角實際上稍微不方便。我想我可以堅持我提交的代碼的並行版本。不過謝謝你！ – 2011-04-07 01:06:09

現在我使用以下內容，靈感來自this answer。輸出是一個n x k矩陣，其中元素(i,k)是數據點的索引，即k最接近i。

n <- 10 
d <- 3 
x <- matrix(rnorm(n * d), ncol = n) 

min.k.dists <- function(x,k=5) { 
    apply(x,2,function(r) { 
    b <- colSums((x - r)^2) 
    o <- order(b) 
    o[1:k] 
    }) 
} 

min.k.dists(x) # first row should be 1:ncol(x); these points have distance 0 
dist(t(x))  # can check answer against this

如果是因爲擔心關係是如何被處理和諸如此類的東西，也許rank()應納入。

上面的代碼似乎有點快，但我相信它可以改進（雖然我沒有時間去C或fortran路線）。所以我仍然願意快速和稀疏地實施上述。

下面我包括我最終使用並行版本：

min.k.dists <- function(x,k=5,cores=1) { 
    require(multicore) 
    xx <- as.list(as.data.frame(x)) 
    names(xx) <- c() 
    m <- mclapply(xx,function(r) { 
    b <- colSums((x - r)^2) 
    o <- order(b) 
    o[1:k] 
    },mc.cores=cores) 
    t(do.call(rbind,m)) 
}

來源

2011-04-06 16:03:00

你需要做dist（t（x））來獲得可比的答案。 – Tommy 2011-04-06 16:44:57

如果您想保留min.k.dist函數的邏輯並返回重複的距離，您可能需要考慮修改它。用0距離返回第一條線似乎毫無意義，對吧？ ...並通過在我的其他答案中加入一些技巧，你可以加快你的版本30％：

min.k.dists2 <- function(x, k=4L) { 
    k <- max(2L, k + 1L) 
    apply(x, 2, function(r) { 
    sort.list(colSums((x - r)^2), na.last=NA, method='quick')[2:k] 
    }) 
} 

> n<-1e4; m<-3; m=matrix(runif(n*m), n) 
> system.time(d <- min.k.dists(t(m), 4)) #To get 3 nearest neighbours and itself 
    user system elapsed 
    17.26 0.00 17.30 
> system.time(d <- min.k.dists2(t(m), 3)) #To get 3 nearest neighbours 
    user system elapsed 
    12.7  0.0 12.7

來源

2011-04-07 15:54:59 Tommy

計算稀疏成對距離矩陣R

回答

相關問題