I'm trying to do some landscape analyses using R, and am having issues where using the pointDistance function from the raster package and the costDistance function from the gDistance package can give significantly different results depending on the raster properties. Here is a trivial example:
library(raster) library(gdistance) x <- c(5, 370 , 900) y <- c(5, 370 , 900) coords <- cbind(x,y) pd <- pointDistance(coords, lonlat = FALSE) #works; dist_matrix matches pd r <- raster(ncol=1000,nrow=1000, crs=NA, xmn=0, xmx=1000, ymn=0, ymx=1000) r <- setValues(r,1) tr <- transition(r, mean, directions=8) tr1 <- geoCorrection(tr, type="c", multpl=FALSE) dist = costDistance(tr1,coords) dist_matrix = as.matrix(dist) #doesn't work r <- raster(ncol=100,nrow=100, crs=NA, xmn=0, xmx=1000, ymn=0, ymx=1000) r <- setValues(r,1) tr <- transition(r, mean, directions=8) tr1 <- geoCorrection(tr, type="c", multpl=FALSE) dist = costDistance(tr1,coords) dist_matrix = as.matrix(dist)
When the raster has spatial resolution of 1 for each pixel, it works, but anything else, it starts to fail. In my case, it's resulting in a 1-10% difference between each site, which is a lot. I'm adding barriers (lakes) to my rasters which have less of an impact...
So how do I deal with this? Depending on the extent I'm working with, my rasters range from tens of megabytes to hundreds of megabytes in size, with each pixel representing a 50x50 meter area. Eventually, when I have the entire process ironed out, it's going to drop to a 10x10 meter area for each pixel, which is going to be the limit of what I can work with computationally. I won't be able to get it down to a 1x1 resolution like the trivial example that works.
Edit: Turns out that the trivial example only works when the x component and y component of a single point are the same value. Also, changing spatial resolutions of the raster doesn't seem to affect the costDistance results, so at least that function seems consistent regardless of the size of the pixels in the raster.