Just Rookie
Just Rookie

Reputation: 209

Efficient way to compute standard deviation of nearest neighbours of each element in matrix

I would like to compute the standard deviation of the nearest neighbors (3*3 moving window) of each element in a matrix. I wrote some code in R to implement it:

library(FNN)    
df <- matrix(1:10000, nrow = 100, ncol = 100, byrow = TRUE)

df_ <- reshape2::melt(df)
df_index <- df_[, c(1,2)]

df_query <- df_index
neighbor_index <- knnx.index(df_index, df_query, k = 9, algorithm = 'kd_tree')

neighbor_coor<- apply(neighbor_index, 1, function(x) df_query[x, ])

neighbor_sd <- lapply(neighbor_coor, function(x) sd(df[x[, 1], x[, 2]]))

sd <- do.call(rbind, neighbor_sd)

But the speed is too slow. Would you give me some advice to speed up? Are there other ways to implement it?

Upvotes: 1

Views: 425

Answers (1)

loki
loki

Reputation: 10350

As @romanlustrik proposed in his comment, we can use a raster::focal() for this problem.

library(raster)

df <- matrix(1:10000, nrow = 100, ncol = 100, byrow = TRUE)
dfR <- raster(df)

dfSD <- as.matrix(focal(dfR, w = matrix(1,3,3), fun = sd))

where, w is the a matrix representing the nearest neighbors and their weighting within fun (in this case 3x3 which is the cell itself and it 8 neighbors). Thus, any neighborhood pattern is imaginable as long as it it can be represented by a matrix.

matrix(1,3,3)
#      [,1] [,2] [,3]
# [1,]    1    1    1
# [2,]    1    1    1
# [3,]    1    1    1

An example with only the 4 neighbors (excluding diagonals and the cell itself):

matrix(c(0,1,0,1,0,1,0,1,0), 3, 3)
#      [,1] [,2] [,3]
# [1,]    0    1    0
# [2,]    1    0    1
# [3,]    0    1    0

Upvotes: 1

Related Questions