studio - superponer graficas en r

Cálculo de la matriz de distancia por pares dispersos en R

Question

Tengo una matriz NxM y quiero calcular la matriz NxN de distancias euclidianas entre los puntos M En mi problema, N es alrededor de 100.000. Como planeo usar esta matriz para un algoritmo vecino k más cercano, solo necesito mantener las k distancias más pequeñas, por lo que la matriz NxN resultante es muy escasa. Esto contrasta con lo que sale de dist() , por ejemplo, lo que resultaría en una matriz densa (y probablemente en problemas de almacenamiento para mi tamaño N ).

Los paquetes para kNN que he encontrado hasta ahora ( knnflex , kknn , etc.) parecen usar matrices densas. Además, el paquete Matrix no ofrece una función de distancia por pares.

Más cerca de mi objetivo, veo que el paquete de spam tiene una función más nearest.dist() que le permite a uno solo considerar distancias menores que algún umbral, delta . En mi caso, sin embargo, un valor particular de delta puede producir demasiadas distancias (por lo que tengo que almacenar la matriz NxN densamente) o muy pocas distancias (por lo que no puedo usar kNN).

He visto una discusión previa sobre bigmemory/biganalytics intentar realizar clústeres de k-means utilizando los paquetes bigmemory/biganalytics , pero no parece que pueda aprovechar estos métodos en este caso.

¿Alguien sabe una función / implementación que calculará una matriz de distancia de forma dispersa en R? Mi (temido) plan de copia de seguridad es tener dos for bucles y guardar resultados en un objeto Matrix .

Answer 1

Bueno, no podemos hacer que recurras a los bucles for, ahora podemos :)

Por supuesto, está la cuestión de cómo representar la matriz dispersa. Una forma simple es que solo contenga los índices de los puntos más cercanos (y vuelva a calcularlos según sea necesario). Pero en la solución a continuación, coloco la distancia (''d1'', etc.) y el índice (''i1'', etc.) en una sola matriz:

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
}

Probando en 9 2d-puntos:

> 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

Y probándolo en un problema mayor (10k puntos). Aún así, en 100k puntos y más dimensiones tomará mucho tiempo (como 15-30 minutos).

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

PS acaba de señalar que publicaste una respuesta mientras escribía esto: la solución aquí es aproximadamente el doble de rápida porque no calcula la misma distancia dos veces (la distancia entre los puntos 1 y 13 es la misma que entre los puntos 13 y 1) .

Answer 2

Por ahora estoy usando lo siguiente, inspirado en esta respuesta . La salida es una matriz nxk donde el elemento (i,k) es el índice del punto de datos que es el k th más cercano a 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

Si uno está preocupado por la forma en que se manejan los vínculos y otras cosas, tal vez debería incorporarse el rank() .

El código anterior parece algo rápido, pero estoy seguro de que podría mejorarse (aunque no tengo tiempo para ir por la ruta C o fortran ). Así que todavía estoy abierto a implementaciones rápidas y dispersas de lo anterior.

A continuación incluyo una versión paralelizada que terminé usando:

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))
}

Answer 3

Si desea mantener la lógica de su función min.k.dist y devolver distancias duplicadas, puede considerar modificarla un poco. Parece inútil devolver la primera línea con 0 distancia, ¿verdad? ... y al incorporar algunos de los trucos en mi otra respuesta, puedes acelerar tu versión en un 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