Warning message when doing spatial join nearest neighbor on GeoPandas

Question

I am trying a nearest neighbor spatial join with two point Geodataframes but a warning message is appearing.

My code:

import geopandas as gpd

yield_2013=gpd.read_file('https://github.com/kevinkuranyi/archive/raw/main/Brazil_Maize_yield_2013.shp')
temperature=gpd.read_file('https://github.com/kevinkuranyi/archive/raw/main/Temperature_2013.shp')
print(yield_2013.crs)
print(temperature.crs)
df=yield_2013.sjoin_nearest(temperature)

The following warning message is shown:

UserWarning: Geometry is in a geographic CRS. Results from 'sjoin_nearest' are likely incorrect. Use 'GeoSeries.to_crs()' to re-project geometries to a projected CRS before this operation.

When I print both CRS the result is the same: EPSG:4674. If I follow the warning message instructions and do:

df=yield_2013.to_crs(epsg=4674).sjoin_nearest(temperature.to_crs(epsg=4674))

The same warning message appears. I manually checked a dozen of points from the resulting file and it seems to be ok, but I have thousands of points and this message (likely incorrect) is making me afraid.

Am I doing the spatial join in the wright way?

Answer 1

A geographic CRS like EPSG:4674 wraps your coordinates around the surface of a sphere or ellipsoid, but that really complicates things when quantifying distance. Now that you're on a sphere, your distances are going to be represented by an angle (degrees) which are of little practical use. And then you have to use a set of trigonometric functions (geodesy) just to convert them to a human-interpretable unit of measurement. Spatial projections (projected CRS) have been developed to solve this problem. They simplify the measurement of distance, shape and/or area and save you from having to perform mathematical transformations.

As for geopandas, if we refer to the docs about sjoin_nearest, we understand that

Every operation in GeoPandas is planar, i.e. the potential third dimension is not taken into account.

In other words, we should project our coordinates onto a 2D plane if we want the spatial join to be accurate. A good choice for Brazil would be EPSG:32723 , because UTM projections are conformal and the units of measurement are meters.

yield_2013_projected = yield_2013.to_crs(epsg=32723)
temperature_projected = temperature.to_crs(epsg=32723)
df=yield_2013_projected.sjoin_nearest(temperature_projected)