Behaviour of gstat::variogram when using geographic coordinates: defining CRS and setting CRS to NA
Opened this issue · 0 comments
alexyshr commented
Greetings Dr. Edzer!
I have learned that I must use projected coordinates to do geostatistics (instead of using geographic coordinates). Nonetheless, I just want to be sure that the behavior of gstat::variogram
when using geographic coordinates is right and not a bug.
Why the distances of the semivariogram are between 0
and 150
considering that the maximum distance between points is 4.590449
?
If I use geographic coordinates without defining the coordinate system st_crs(st_grd_sample) <- NA
, then the semivariogram distances seem fine.
Thank you!
library(sf)
#> Linking to GEOS 3.9.3, GDAL 3.5.2, PROJ 8.2.1; sf_use_s2() is TRUE
library(ggplot2)
# Polygon of the state of Ohio, USA
ohio = sf::st_as_sf(maps::map("state", "ohio", plot = FALSE, fill = TRUE))
#Sample 200 points
st_grd_sample = sf::st_sf(sf::st_sample(ohio, 200))
#It has WGS84 geographic coordinates
sf::st_crs(st_grd_sample)$input
#> [1] "EPSG:4326"
#Review max distances between points
the_coords = sf::st_coordinates(st_grd_sample)
the_dist = dist(the_coords)
max(the_dist)
#> [1] 4.590449
#assign random values between 0 and 100
st_grd_sample$value = runif(dim(st_grd_sample)[1], 0, 100)
#
#Geographic Coordinates using CRS
#
y.v2 = gstat::variogram(value~1, loc=st_grd_sample)
#The distances for the semivariogram are between 0 and more than 150 (wrong?)
plot(y.v2$dist, y.v2$gamma,
xlab="Distance", ylab=expression(gamma(h)))
#
#Geographic Coordinates without defining CRS
#
st_crs(st_grd_sample) <- NA
y.v2 = gstat::variogram(value~1, loc=st_grd_sample)
#The distances for the semivariogram are ok
plot(y.v2$dist, y.v2$gamma,
xlab="Distance", ylab=expression(gamma(h)))
Created on 2023-04-11 with reprex v2.0.2