在地图上看见数据，空间中的故事
地理空间分析（Geospatial Analysis）是处理带有地理位置信息的数据的科学，从城市规划到流行病学，从物流优化到气候变化，地理信息无处不在。
| 数据类型 | 格式 | 示例 | 工具 |
|---|---|---|---|
| 矢量-点 | CSV/GeoJSON | 城市位置 | geopandas |
| 矢量-线 | Shapefile | 道路网络 | osmnx |
| 矢量-面 | GeoJSON | 行政区划 | geopandas |
| 栅格 | GeoTIFF | 卫星影像 | rasterio |
| 网格 | NetCDF | 气象数据 | xarray |
地理空间分析的第一步是理解坐标系统和距离计算方法。
import numpy as np
def haversine(lat1, lon1, lat2, lon2):
"""计算两点间球面距离(km)"""
R = 6371 # 地球半径(km)
dlat = np.radians(lat2 - lat1)
dlon = np.radians(lon2 - lon1)
a = np.sin(dlat/2)**2 + np.cos(np.radians(lat1)) * np.cos(np.radians(lat2)) * np.sin(dlon/2)**2
c = 2 * np.arcsin(np.sqrt(a))
return R * c
# 北京→上海
dist = haversine(39.9042, 116.4074, 31.2304, 121.4737)
print(f"北京→上海: {dist:.1f} km")
# 中国城市间距离矩阵
cities = [
("北京", 39.90, 116.41), ("上海", 31.23, 121.47),
("广州", 23.13, 113.26), ("深圳", 22.54, 114.06),
("成都", 30.57, 104.07),
]
for name, lat, lon in cities:
for name2, lat2, lon2 in cities:
if name < name2:
d = haversine(lat, lon, lat2, lon2)
print(f" {name}→{name2}: {d:.0f} km")folium是Python最流行的交互式地图库，基于Leaflet.js，可在浏览器中查看。
import folium
import numpy as np
# 中国主要城市数据
cities = [
("北京", 39.90, 116.41, 2171, 4027),
("上海", 31.23, 121.47, 2489, 4321),
("广州", 23.13, 113.26, 1868, 2882),
("深圳", 22.54, 114.06, 1756, 3246),
("成都", 30.57, 104.07, 2094, 1992),
("杭州", 30.27, 120.16, 1194, 1810),
("武汉", 30.59, 114.31, 1233, 1890),
("重庆", 29.43, 106.91, 3205, 2789),
("西安", 34.34, 108.94, 1295, 1141),
("南京", 32.06, 118.80, 931, 1681),
]
# 创建暗色主题地图
m = folium.Map(
location=[35.0, 110.0], # 中国中心
zoom_start=4,
tiles='CartoDB dark_matter' # 暗色主题
)
# 添加城市气泡标记
for name, lat, lon, pop, gdp in cities:
folium.CircleMarker(
location=[lat, lon],
radius=pop / 400, # 人口越大圆越大
popup=f"{name}: 人口{pop}万, GDP{gdp}亿",
color='#3b82f6',
fill=True,
fill_color='#3b82f6',
fill_opacity=0.6
).add_to(m)
# 添加城市连线
for i, (n1, lat1, lon1, _, _) in enumerate(cities[:5]):
for j, (n2, lat2, lon2, _, _) in enumerate(cities[5:], 5):
dist = haversine(lat1, lon1, lat2, lon2)
if dist < 1500: # 只画1500km以内的连线
folium.PolyLine(
locations=[[lat1, lon1], [lat2, lon2]],
color='#3b82f6',
weight=1,
opacity=0.3
).add_to(m)
m.save('china_cities.html')
print("地图已保存! 用浏览器打开china_cities.html查看")from scipy.stats import pearsonr
from sklearn.cluster import KMeans
# 空间相关性分析
pops = np.array([c[3] for c in cities])
gdps = np.array([c[4] for c in cities])
lats = np.array([c[1] for c in cities])
r, p = pearsonr(pops, gdps)
print(f"人口-GDP相关: r={r:.3f}, p={p:.4f}")
r2, p2 = pearsonr(lats, gdps)
print(f"纬度-GDP相关: r={r2:.3f}, p={p2:.4f}")
# 人均GDP分析
gdp_per_capita = gdps / pops * 10000
sorted_idx = np.argsort(gdp_per_capita)[::-1]
print(f"\n人均GDP排名:")
for i in sorted_idx:
print(f" {cities[i][0]:>4}: {gdp_per_capita[i]:.2f}万元/人")
# 空间聚类
coords = np.array([[c[1], c[2]] for c in cities])
kmeans = KMeans(n_clusters=3, random_state=42, n_init=10)
clusters = kmeans.fit_predict(coords)
for cid in range(3):
members = [cities[i][0] for i in range(len(cities)) if clusters[i] == cid]
center = kmeans.cluster_centers_[cid]
print(f"\n聚类{cid}: {', '.join(members)}")
print(f" 中心: {center[0]:.1f}°N, {center[1]:.1f}°E")# GeoPandas核心概念 (伪代码，需安装geopandas)
# import geopandas as gpd
# from shapely.geometry import Point, Polygon
# 1. 创建GeoDataFrame
# geometry = [Point(lon, lat) for lat, lon in coords]
# gdf = gpd.GeoDataFrame(cities_df, geometry=geometry, crs='EPSG:4326')
# 2. 空间连接
# result = gpd.sjoin(cities_gdf, provinces_gdf, how='left', op='within')
# 3. 缓冲区分析
# buffer = gdf.buffer(distance=10000) # 10km缓冲区
# 4. 面积计算（需转投影坐标系）
# gdf_proj = gdf.to_crs('EPSG:3857') # 转Web Mercator
# areas = gdf_proj.area / 1e6 # km²
# 5. 可视化
# gdf.plot(column='population', cmap='Blues', legend=True)
print("GeoPandas核心操作: GeoDataFrame创建、空间连接、缓冲区、面积计算、choropleth可视化")