Calculate distance between a coordinate and a county in GeoPandas

One of the packages that Geopandas depends on is Shapely. Shapely uses 2D Cartesian coordinates (x, y) and calculations are performed in 2D plane.

If you use "projected coordinate system", no problem. The distance that you get is the distance on the map (not on the spherical earth). When using "geographic coordinate system - GCS", the distance that you get will be the shortest distance in 3D space.

So you should use a formula to calculate distance on the sphere, and that is Haversine formula. But it calculates great-circle distance between two points on a sphere given their longitudes and latitudes. However, you have a point and a polygon. So you need a point on the polygon which has nearest distance to the point.

I assume that the coordinates of county polygon are on GCS, not projected. If projected, you need to transform to GCS using .to_crs() method. Then use that code to find nearest point coordinates on the polygon. (Code Reference: Find Coordinate of Closest Point on Polygon Shapely)

from shapely.geometry import LinearRing

pol_ext = LinearRing(county.geometry.exterior.coords)
d = pol_ext.project(point.geometry)
p = pol_ext.interpolate(d)
closest_point_coords = list(p.coords)[0]

And then, use the function below to calculate spherical distance that you need. (Code Reference: Haversine Formula in Python)

from math import radians, cos, sin, asin, sqrt

def haversine(lon1, lat1, lon2, lat2):
    """
    Calculate the great circle distance between two points 
    on the earth (specified in decimal degrees)
    """
    # convert decimal degrees to radians 
    lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])

    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
    c = 2 * asin(sqrt(a)) 
    r = 3956 # Radius of earth in miles. Use 6371 for kilometers
    return c * r

You need to reproject your geodataframes into a CRS that adequately preserves distances; UTM is a common choice, and the results would be in meters. You must determine which UTM zone to use, find the appropriate EPSG, then use the .to_crs(epsg=XXXX) method on your geodataframes to convert them.

Choosing the correct CRS in many cases is nontrivial; I have just provided a starting point that may work for you. UTM projections do not preserve distances, but over small areas, the distortion is often acceptable. The Plate Caree apparently preserves distances (epsg 32662).