Calculate Moran's I for a large df (12k+ columns, 3.5k rows) and store results in a separate df

Viewed 89

I'd like to run Moran's I test on a data frame with 12044 columns (and 3400 rows) and store the results in a df or a list of sorts. The first three columns are ID, Lat, and Long respectively. The rest are variables that I am interested in.

I understand that lapply is designed to do what I want, but I am not very good at storing the results in a separate df. The result of the test has four variable: observed, expected, sd, and p.value.

Here is a sample of the df and the function itself.

set.seed(1)
df <- data.frame(
     ID = 1:15,
     LATITUDE = c(42.6, 42.5, 42.3, 42.8, 42.4, 42.4, 42.4, 42.3, 42.4, 42.4, 41.4, 41.6, 41.8, 43.7, 47.3),
     LONGITUDE = c(-71.5, -71.6, -71.9, -71.0, -71.1, -71.1, -71.1, -71.1, -71.2, -71.2, -70.5, -70.3, -71.2, -70.3, -68.3),
     x1 = runif(15, min=0, max=1000),
     x2 = runif(15, min=0, max=1000),
     x3 = runif(15, min=0, max=1000),
     x4 = runif(15, min=0, max=1000),
     x5 = runif(15, min=0, max=1000),   
     x6 = runif(15, min=0, max=1000),
     x7 = runif(15, min=0, max=1000),
     x8 = runif(15, min=0, max=1000)  )
require(ape)

dists <- as.matrix(dist(cbind(df$LONGITUDE, df$LATITUDE)))
dists.inv <- 1/dists
diag(dists.inv) <- 0
#check
dists.inv[1:5, 1:5]
#deal with the infinite values in the matrix
dists.inv[is.infinite(dists.inv)] <- 0
#calculate Moran's I
Moran.I(df$x1, dists.inv)

Thank you all

2 Answers

Consider using tidyverse. Select only the columns that starts_with 'x' or matches("^x\\d+$"), loop over those columns with map, apply the Moran.I with the already created 'dists.inv' and the looped column, and return a tibble by row binding the list elements (_dfr)

library(purrr)
library(dplyr)
df %>% 
    select(starts_with('x')) %>%
     map_dfr(~ ape::Moran.I(.x, dists.inv))

-output

# A tibble: 8 x 4
  observed expected     sd p.value
     <dbl>    <dbl>  <dbl>   <dbl>
1 -0.0305   -0.0714 0.0745  0.583 
2 -0.0854   -0.0714 0.0739  0.850 
3 -0.185    -0.0714 0.0712  0.111 
4 -0.237    -0.0714 0.0737  0.0250
5 -0.109    -0.0714 0.0736  0.612 
6 -0.0280   -0.0714 0.0749  0.562 
7  0.00361  -0.0714 0.0731  0.305 
8 -0.177    -0.0714 0.0737  0.152 

Another way using a base R could be

output <- apply(df[-(1:3)], 2, function(x) Moran.I(x, dists.inv))

and then bind this list into a data.frame:

do.call("rbind.data.frame", output)

This returns

       observed    expected         sd    p.value
x1 -0.030529141 -0.07142857 0.07452502 0.58314178
x2 -0.085369231 -0.07142857 0.07390247 0.85037818
x3 -0.184828111 -0.07142857 0.07123184 0.11138959
x4 -0.236554103 -0.07142857 0.07367464 0.02500791
x5 -0.108772142 -0.07142857 0.07359794 0.61187441
x6 -0.028012329 -0.07142857 0.07485506 0.56191185
x7  0.003612685 -0.07142857 0.07309663 0.30460722
x8 -0.177143267 -0.07142857 0.07372009 0.15157193
Related