当前位置：首页 > news >正文

别再只用GWR了！用Python的mgtwr包搞定时空地理加权回归(GTWR)实战

news 2026/6/7 5:42:28

时空地理加权回归实战：用Python的mgtwr包突破传统GWR局限

当我们需要分析房价波动、污染物扩散或流行病传播这类同时具有空间和时间演变特征的数据时，传统的地理加权回归(GWR)往往力不从心。我曾在一个城市空气质量分析项目中深有体会——当试图用GWR模型解释PM2.5浓度的空间分布时，发现模型在跨季节数据上的预测效果波动很大，这才意识到忽略了时间维度的动态影响。这正是时空地理加权回归(GTWR)大显身手的场景。

1. GWR与GTWR的核心差异解析

GWR模型通过地理加权的方式，让回归系数随着空间位置变化而调整，这在处理纯空间异质性数据时表现优异。但在现实世界中，很少有现象是纯粹空间性的。以房价为例，同一地段在不同经济周期下的增值规律可能完全不同，这就是典型的时空双重异质性。

GTWR的核心创新在于其权重矩阵的构造。与传统GWR仅考虑空间距离不同，GTWR的权重计算同时融合了空间距离dᵢⱼ和时间距离dₜⱼ：

Wᵢ = exp[-(dᵢⱼ²/bw² + dₜⱼ²/τ²)]

其中bw是空间带宽，τ是时间带宽参数。这种设计带来了几个关键优势：

时空耦合效应捕捉：能够识别出某些影响因素在特定时期对特定区域的作用强度
动态过程建模：比如交通建设对房价的影响可能随时间衰减，这种动态可以被量化
跨期数据整合：避免将不同时期数据简单堆砌导致的空间自相关误判

注意：实际计算时，时空距离需要标准化处理，因为空间单位(如公里)和时间单位(如月)的量纲不同。

2. mgtwr包实战全流程

Python的mgtwr包(1.0.4版本)是目前实现GTWR较为成熟的工具。下面通过一个模拟数据集演示完整工作流。

2.1 数据准备的关键细节

首先导入必要库并生成模拟数据：

import numpy as np import pandas as pd from mgtwr.gtwr import GTWR from mgtwr.sel_bws import Sel_bws # 生成时空数据 np.random.seed(42) n_samples = 1000 u = np.random.uniform(0, 10, n_samples) # 空间x坐标 v = np.random.uniform(0, 10, n_samples) # 空间y坐标 t = np.random.randint(0, 24, n_samples) # 时间点(0-23小时) coords = np.column_stack([u, v]) # 生成解释变量和响应变量 X1 = np.random.normal(5, 2, n_samples) X2 = np.random.poisson(3, n_samples) X = np.column_stack([X1, X2]) # 构造时空变化的回归系数 beta1 = 2 + 0.5*u + 0.1*v + 0.05*t beta2 = 3 - 0.3*u + 0.2*v - 0.02*t y = 5 + beta1*X1 + beta2*X2 + np.random.normal(0, 1, n_samples)

关键点：坐标数据应当使用投影坐标系(如UTM)的直角坐标，而非经纬度。这是因为包内计算时空距离时直接使用欧式距离公式。若原始数据是经纬度，需要先进行投影转换。

2.2 带宽参数优化实战

GTWR有两个关键超参数需要确定：

空间带宽(bw)：控制空间影响的衰减速度
时间带宽(τ)：控制时间影响的衰减速度

使用黄金分割搜索法自动优化：

selector = Sel_bws(coords, t, y, X, kernel='gaussian', fixed=True) bw, tau = selector.search(bw_max=20, tau_max=10, verbose=True) print(f"最优参数 - bw: {bw:.2f}, tau: {tau:.2f}")

输出示例：

bw, tau, score: 12.4 6.2 4821.65 bw, tau, score: 7.7 6.2 4789.32 ... 最优参数 - bw: 5.8, tau: 3.1

参数选择对结果影响显著。过大的带宽会导致模型退化为全局回归，过小则可能过拟合。建议：

先使用自动搜索获得初始值
在最优值附近手动微调
通过交叉验证评估稳定性

2.3 模型拟合与结果解读

使用优化后的参数拟合模型：

model = GTWR(coords, t, y, X, bw=5.8, tau=3.1, kernel='gaussian', fixed=True).fit()

关键输出解析：

R²：模型整体解释力，GTWR通常比GWR提高0.1-0.3
betas：各点的回归系数矩阵，形状为(n_samples, n_features+1)
tvalues：系数的t统计量，评估显著性

可视化系数空间分布：

import matplotlib.pyplot as plt plt.scatter(coords[:,0], coords[:,1], c=model.betas[:,1], cmap='coolwarm', s=10) plt.colorbar(label='X1系数') plt.title('X1系数的空间分布') plt.show()

3. 进阶技巧与避坑指南

3.1 时空权重函数的选型

mgtwr支持多种核函数，通过kernel参数指定：

核函数类型	公式	适用场景
Gaussian	exp(-d²/h²)	平滑连续变化
Bisquare	(1-d²/h²)² if d<h else 0	局部突变明显
Exponential	exp(-d/h)	长尾衰减

# 比较不同核函数 kernels = ['gaussian', 'bisquare', 'exponential'] for kern in kernels: model = GTWR(..., kernel=kern).fit() print(f"{kern}: R²={model.R2:.3f}")

3.2 大规模数据优化策略

当数据量>10,000时，可采取以下优化：

空间分块：将研究区域划分为子区域分别建模
时间分段：按季节/时段建立分模型
采样策略：使用空间随机森林等方法筛选重要点位

# 空间分块示例 from sklearn.cluster import KMeans # 将空间划分为5个子区域 kmeans = KMeans(n_clusters=5).fit(coords) for i in range(5): mask = kmeans.labels_ == i sub_model = GTWR(coords[mask], t[mask], y[mask], X[mask]).fit()

3.3 模型诊断与验证

推荐进行以下验证：

残差空间自相关检验：使用Moran's I指数
时间自相关检验：Ljung-Box检验
多重共线性诊断：VIF分析

from mgtwr.diagnostics import morans_i residuals = model.resid moran = morans_i(coords, residuals) print(f"Moran's I: {moran.I:.3f}, p-value: {moran.p:.3f}")