I'm using the correlog
function from the ncf
package to figure out the distance between points at which spatial autocorrelation isn't detected. I was wondering what unit was being returned in the correlogram or even in the object's mean.of.class
element.
RNGkind('Mersenne-Twister')
set.seed(42)
library(sp)
library(raster)
x <- seq_len(1000) #Zonal coordinates
y <- seq_len(1000) #Meridionial coordinates
z <- runif(1000, -5, 5) #Simulated residuals
my.coords <- data.frame(x, y)
coordinates(my.coords) <- c('x', 'y')
crs(my.coords) <- CRS('+proj=merc +datum=WGS84 +units=m +ellps=WGS84')
library(ncf)
cor.obj <- correlog(x = x, y = y, z = z, latlon = FALSE, increment = 20, resamp = 10)
Where x
are the zonal coordinates, y
are the meridionial coordinates, z
are the model residuals I'm running the Moran's I test on, increment
is the "increment for the uniformly distributed distance classes" (from ?ncf::correlog
) and latlon
being a logical value (TRUE
or FALSE
) depending on what coordinate systems x
and y
are.
The output object produces a few things, but I'm interested in correlation
and x.intercept
outputs:
Value
An object of class "correlog" is returned, consisting of the following components:
correlation
the value for the moran (or Mantel) similarity....
x.intercept
the interpolate x.intercept of Epperson (1993).
Calling ?ncf::correlog
doesn't return much in regards to what units are being returned, but I'm convinced that it returns meters if I'm using World Mercator and set the latlon
argument to FALSE
. But when I plot the object and call the for the mean.of.class
, it doesn't make sense:
plot(cor.obj)
car.obj$mean.of.class
1 2 3 4 5 6 7 8 9 10 11
10.58345 30.38211 50.18075 69.97939 89.77802 109.57664 129.37525 149.87711 170.38660 190.18517 209.98373
12 13 14 15 16 17 18 19 20 21 22
229.78226 249.58079 269.37929 289.88057 310.39043 330.18889 349.98732 369.78573 389.58412 409.38248 429.88299
23 24 25 26 27 28 29 30 31 32 33
450.39330 470.19158 489.98983 509.78805 529.58624 549.38439 569.88385 590.39471 610.19273 629.99070 649.78861
34 35 36 37 38 39 40 41 42 43 44
669.58647 689.38427 709.88220 730.39377 750.19135 769.98883 789.78622 809.58351 829.38068 849.87619 870.38868
45 46 47 48 49 50 51 52 53 54 55
890.18541 909.98196 929.77832 949.57445 969.37033 989.86150 1010.37510 1030.16994 1049.96432 1069.75818 1089.55141
56 57 58 59 60 61 62 63 64 65 66
1109.34390 1129.82530 1150.33951 1170.12850 1189.91572 1209.70066 1229.48259 1249.26042 1269.70497 1290.20724 1309.95620
67 68 69 70 71
1329.68162 1349.36177 1368.93628 1388.15163 1404.31407
I've only had 1000 meters in the original x
and y
vectors, so why is it plotting it all the way to 1400?
Changing latlon
to TRUE
produces this:
cor.obj <- correlog(x = my.coords$x, y = my.coords$y, z = z, latlon = TRUE, increment = 20, resamp = 10)
plot(cor.obj)
Converting my.coords
to latitude/longitude produces this monstrosity:
my.coords <- spTransform(my.coords, CRS('+proj=longlat +ellps=WGS84 +datum=WGS84'))
cor.obj <- correlog(x = x, y = y, z = z, latlon = TRUE, increment = 20, resamp = 10)
plot(cor.obj)
Changing the increment
to 0.02
produces something a bit more recognizable but the distance (mean-of-class) has changed to what I think are kilometers:
cor.obj <- correlog(x = my.coords$x, y = my.coords$y, z = z, latlon = T, increment = 0.02, resamp = 10)
plot(cor.obj)
But again, I'm a bit thrown by the x-axis (e.g. call cor.obj$mean.of.class
), because I only started off with 1000 values in my.coords
, so why does it go up to 1400? Am I missing something? What unit is being returned?