基于Haversine函数统计指定半径内POI的最优方案咨询

Hey Jan, great question! When dealing with 5k stores and 200k POIs, brute-force nested loops (that's 1 billion distance calculations!) are going to be painfully slow. Below are two optimized approaches, ordered from fast to blazingly fast, to solve this problem efficiently:

Approach 1: Numpy Vectorization (100x Faster Than Pure Python Loops) #

This leverages numpy's broadcasting to avoid slow Python-level loops. We'll compute distances between a single store and all POIs in one go, then count how many fall within your target radius. No extra libraries needed—just build on the code you already have:

import pandas as pd
import numpy as np

def haversine_np(lon1, lat1, lon2, lat2):
    lon1, lat1, lon2, lat2 = map(np.radians, [lon1, lat1, lon2, lat2])
    dlon = lon2 - lon1
    dlat = lat2 - lat1
    a = np.sin(dlat/2.0)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2.0)**2
    c = 2 * np.arcsin(np.sqrt(a))
    km = 6367 * c
    return km

# Load your data
df1 = pd.read_excel(r'Stores.xlsx')
df2 = pd.read_excel(r'POI.xlsx')

# Extract POI coordinates as numpy arrays (avoids repeated dataframe lookups)
poi_lons = df2['longitude'].values  # Replace with your actual column name
poi_lats = df2['latitude'].values   # Replace with your actual column name

# Set your target radius (e.g., 5 kilometers)
TARGET_RADIUS_KM = 5

# Vectorized function to count nearby POIs for a single store
def count_nearby_pois(store_lon, store_lat):
    # Compute distances from this store to ALL POIs in one vectorized operation
    distances = haversine_np(store_lon, store_lat, poi_lons, poi_lats)
    # Count how many distances are within the target radius
    return np.sum(distances <= TARGET_RADIUS_KM)

# Apply to all stores and add the result as a new column
df1['nearby_poi_count'] = df1.apply(
    lambda row: count_nearby_pois(row['longitude'], row['latitude']),
    axis=1
)

Why this works so well: #

• Numpy offloads calculations to optimized C code, skipping the overhead of Python loops
• Broadcasting lets us compute distances for all POIs at once without writing nested loops

Approach 2: Spatial Index (R-tree) + Local Queries (10-100x Faster Than Vectorization) #

For 200k POIs, even vectorized calculations mean 200k distance checks per store. With an R-tree spatial index, we first filter down to only POIs that are potentially within the radius (usually just a few hundred), then compute precise distances for that small subset. This cuts computation time drastically. We'll use geopandas (which wraps R-tree under the hood) for this:

Step 1: Install dependencies #

pip install geopandas rtree

Step 2: Implementation Code #

import pandas as pd
import geopandas as gpd
from shapely.geometry import Point

# Load your data
df1 = pd.read_excel(r'Stores.xlsx')
df2 = pd.read_excel(r'POI.xlsx')

# Set target radius (convert to meters, since we'll use a projected CRS)
TARGET_RADIUS_M = 5 * 1000  # 5 km = 5000 m

# Convert data to GeoDataFrames with WGS84 coordinate system (EPSG:4326)
gdf_stores = gpd.GeoDataFrame(
    df1,
    geometry=gpd.points_from_xy(df1['longitude'], df1['latitude']),
    crs="EPSG:4326"
)
gdf_pois = gpd.GeoDataFrame(
    df2,
    geometry=gpd.points_from_xy(df2['longitude'], df2['latitude']),
    crs="EPSG:4326"
)

# Convert to a local projected CRS (UTM) for accurate distance calculations
# This automatically picks the right UTM zone for your data
utm_crs = gdf_stores.estimate_utm_crs()
gdf_stores = gdf_stores.to_crs(utm_crs)
gdf_pois = gdf_pois.to_crs(utm_crs)

# Build R-tree spatial index for POIs (this is a one-time cost)
poi_sindex = gdf_pois.sindex

# Function to count nearby POIs using the spatial index
def count_nearby_pois_with_index(row):
    store_point = row.geometry
    # Create a buffer around the store and get the bounds for the index query
    search_buffer = store_point.buffer(TARGET_RADIUS_M)
    # Get indices of POIs that fall within the buffer's bounds (candidates)
    candidate_indices = list(poi_sindex.intersection(search_buffer.bounds))
    candidate_pois = gdf_pois.iloc[candidate_indices]
    # Filter candidates to those actually within the radius (precise check)
    valid_pois = candidate_pois[candidate_pois.distance(store_point) <= TARGET_RADIUS_M]
    return len(valid_pois)

# Apply to all stores
df1['nearby_poi_count'] = gdf_stores.apply(count_nearby_pois_with_index, axis=1)

Why this is the best option for your data: #

• The R-tree index eliminates 99% of irrelevant POIs before any distance calculations
• Using a projected CRS gives more accurate distance measurements than the Haversine spherical approximation
• Scales incredibly well—this will stay fast even if your POI dataset grows to millions of records

Quick Comparison of Approaches #

内容的提问来源于stack exchange，提问作者Jan