简介这份资源面向空间计量经济学的研究者与高年级学生系统整理了截面数据下的主流空间回归估计方法帮助解决空间依赖建模中模型选择与参数估计的实操难题。包内覆盖空间滞后模型SLM、空间误差模型SEM、空间杜宾模型SDM及其误差形式并延伸至自变量空间滞后、Kelejian-Prucha模型、一般嵌套空间模型、空间扩展模型与地理加权回归等十余种估计程序同时提供拉格朗日乘子检验LM用于模型诊断。资源附带原始数据并集成jplv7与Elhorst_codes两套常用工具包便于直接复现与二次开发。压缩包约6.15MB以MATLAB脚本与数据文件为主脚本对应各模型的估计流程数据文件支撑实证演练。目前已有2073人学习下载适合希望从经典模型过渡到扩展设定、并借助现成代码快速完成截面空间计量实证的读者参考。1. 空间计量模型估计方法截面数据里那些“看着像回归其实差很远”的坑如果你手头有一份带地理坐标的截面数据比如各区县的人均产出、房价、污染浓度直接跑 OLS 大概率会翻车。原因不复杂相邻地区的变量往往互相影响误差项也可能在空间上抱团。空间计量模型估计方法就是来解决这个问题的核心思路是把“空间关系”显式写进模型——空间滞后模型SAR刻画因变量的溢出效应空间误差模型SEM处理误差项的空间自相关杜宾模型SDM则同时包含自变量和因变量的空间滞后项。截面数据是这套方法最经典的用武之地因为每个观测单元只有一个时间点空间权重矩阵直接决定了模型能不能用。这篇文章会从权重矩阵构造讲到三种模型的估计代码再到参数解释和踩坑记录适合手里有截面数据、想跑通空间回归的从业者。2. 空间权重矩阵与三种模型先搞清楚你的数据该用哪个2.1 空间权重矩阵不是随便填个邻接就行空间计量的第一步永远是构造空间权重矩阵 W。截面数据没有时间维度W 就是 n×n 的矩阵元素 w_ij 表示地区 i 和 j 的空间关系。最常见的三种构造方式邻接矩阵共享边界为 1否则为 0、距离阈值矩阵距离小于 d 为 1、反距离矩阵w_ij 1/d_ij。我一般会先做邻接矩阵因为它对边界敏感但解释直观如果研究区域边界碎片化严重再换反距离矩阵。构造完必须做行标准化让每行和为 1。这一步不做后面估计出来的空间滞后系数 ρ 会偏得离谱。行标准化的含义是“每个地区的邻居影响取平均”而不是简单加总。import numpy as np import pandas as pd from scipy.spatial.distance import cdist # 假设 df 里有 lon, lat 两列 coords df[[lon, lat]].values n len(coords) # 反距离矩阵对角线设为 0 D cdist(coords, coords, metriceuclidean) W 1.0 / (D 1e-12) np.fill_diagonal(W, 0) # 距离阈值超过 200 公里的邻居直接切断避免远距离伪相关 W[D 200] 0 # 行标准化 row_sum W.sum(axis1, keepdimsTrue) W_std W / row_sum这段代码的关键参数是距离阈值 200 公里你需要根据研究区域的实际尺度调整。阈值太小会导致孤立点某行全为 0太大则失去空间区分度。跑完检查np.any(row_sum 0)如果有孤立点要么放宽阈值要么手动指定最近邻。2.2 SAR、SEM、SDM 的公式差异决定了你的假设三个模型写出来差别不大但假设完全不同。SAR 是 y ρWy Xβ ε它假设因变量之间存在溢出效应比如相邻地区的房价会互相拉动。SEM 是 y Xβ uu λWu ε它认为溢出发生在误差项里可能是遗漏了某些空间相关的变量。SDM 是 y ρWy Xβ WXθ ε它把自变量也做了空间滞后适合你想同时检验“本地自变量”和“邻居自变量”的影响。选哪个模型不能靠拍脑袋。常见做法是先用 LM 检验Lagrange Multiplier判断残差是否存在空间自相关再用 LR 或 Wald 检验判断 SDM 能否退化成 SAR 或 SEM。如果 SDM 的 θ 和 ρ 都显著就别退化成 SAR否则会遗漏自变量的空间溢出。# 用 spreg 包做模型估计PySAL 生态 from spreg import ML_Lag, ML_Error, ML_Lag_Regimes # 假设 y 是因变量X 是自变量矩阵W_std 是行标准化权重 # SAR 模型 sar ML_Lag(y, X, wW_std, name_yy, name_x[x1,x2]) print(sar.summary) # SEM 模型 sem ML_Error(y, X, wW_std, name_yy, name_x[x1,x2]) print(sem.summary) # SDM 模型把 WX 也加入自变量 WX W_std X X_sdm np.hstack([X, WX]) sdm ML_Lag(y, X_sdm, wW_std, name_yy, name_x[x1,x2,Wx1,Wx2]) print(sdm.summary)ML_Lag和ML_Error分别对应 SAR 和 SEMML_Lag加上 WX 就是 SDM。注意name_x的顺序必须和矩阵列顺序一致否则输出表里的变量名会错位。估计方法默认是最大似然截面数据样本量小于 500 时建议用methodfull避免迭代不收敛。2.3 参数解释ρ、λ、θ 分别意味着什么ρ 是空间滞后系数衡量因变量的溢出强度。ρ 显著为正说明邻居的 y 越高本地的 y 也越高。λ 是误差项的空间自相关系数显著说明遗漏变量在空间上聚集。θ 是自变量的空间滞后系数显著意味着邻居的 X 会影响本地的 y。解释时要注意直接效应和间接效应的分解。SDM 里一个自变量对 y 的总效应 直接效应 间接效应直接效应不等于 β间接效应也不等于 θ。LeSage 和 Pace 提出的分解方法需要用 (I - ρW)^{-1} 计算spreg 包不直接输出需要自己写。# 计算 SDM 的直接效应和间接效应 I np.eye(n) rho sdm.betas[0] # 第一个系数是 rho beta sdm.betas[1:3] # x1, x2 的系数 theta sdm.betas[3:5] # Wx1, Wx2 的系数 # 总效应矩阵 S_total np.linalg.inv(I - rho * W_std) (beta.reshape(-1,1) theta.reshape(-1,1)) direct np.mean(np.diag(S_total)) total np.mean(S_total) indirect total - direct print(f直接效应: {direct:.4f}, 间接效应: {indirect:.4f})这段代码算的是平均效应更严谨的做法是逐观测计算再取平均。ρ 接近 1 时 (I - ρW) 接近奇异求逆会不稳定这时候要检查模型是否过度拟合。3. 截面数据实操从原始数据到模型输出的完整链路3.1 数据准备与缺失值处理截面数据的空间计量对缺失值极其敏感。如果某个地区缺了自变量整个 W 矩阵的那一行和那一列都得删掉否则估计结果会有偏。我一般会先做三件事检查 y 和 X 的缺失比例缺失超过 10% 的变量直接换掉检查空间权重矩阵的孤立点检查 y 的分布严重右偏就取对数。# 数据清洗 df df.dropna(subset[y, x1, x2, lon, lat]) # 检查 y 的偏度 from scipy.stats import skew if abs(skew(df[y])) 2: df[y] np.log1p(df[y]) # 重新构造 W确保和 df 的行顺序一致 coords df[[lon, lat]].values D cdist(coords, coords) W 1.0 / (D 1e-12) np.fill_diagonal(W, 0) W[D 200] 0 W_std W / W.sum(axis1, keepdimsTrue)注意df的行顺序必须和 W 的行顺序严格对应。如果中间做了排序或筛选W 必须重新构造。这个坑我踩过不止一次模型跑出来 ρ 是负的查了半天才发现是行错位。3.2 模型选择LM 检验与 LR 检验的代码实现LM 检验用来判断要不要做空间模型。原假设是残差不存在空间自相关。如果 LM-Lag 显著而 LM-Error 不显著选 SAR反之选 SEM两个都显著选 SDM。from spreg import OLS # 先跑 OLS 拿残差 ols OLS(y, X, name_yy, name_x[x1,x2]) print(ols.summary) # LM 检验需要手动计算或者用 GeoDa 等工具 # 这里给出简化版用残差和 W 计算 Morans I resid ols.u moran_I (n / W.sum()) * (resid.T W resid) / (resid.T resid) print(fMorans I: {moran_I:.4f})Morans I 显著为正说明残差有空间聚集需要上空间模型。更严格的 LM 检验需要计算统计量建议用 GeoDa 或 R 的 spdep 包做初步筛选再用 Python 估计。3.3 估计结果输出与稳健性检查跑完模型别急着写结论先做三件事检查 ρ 或 λ 的显著性检查对数似然值检查残差的空间自相关是否消除。# 检查 SAR 残差的空间自相关 sar_resid sar.u moran_sar (n / W.sum()) * (sar_resid.T W sar_resid) / (sar_resid.T sar_resid) print(fSAR 残差 Morans I: {moran_sar:.4f}) # 对比三个模型的 AIC print(fSAR AIC: {sar.aic}, SEM AIC: {sem.aic}, SDM AIC: {sdm.aic})如果 SAR 残差的 Morans I 仍然显著说明 SAR 没抓干净空间效应换 SEM 或 SDM。AIC 越小越好但别只看 AIC还要看系数是否符合经济含义。ρ 大于 1 或小于 -1 直接判定模型有问题。4. 避坑与排查空间计量截面数据最常见的 5 个翻车现场4.1 现象ρ 估计出来是负的且绝对值很大原因W 矩阵行标准化之前没有把对角线设为 0或者行顺序和 df 错位。对角线不为 0 会导致自己影响自己ρ 的符号完全乱掉。解决构造 W 之后立刻np.fill_diagonal(W, 0)然后检查df.index和 W 的行索引是否一一对应。如果 df 做过sort_valuesW 必须重新构造。4.2 现象模型不收敛报“Matrix is singular”原因W 矩阵有孤立点某一行全为 0导致 (I - ρW) 不可逆。或者自变量之间存在完全共线性。解决检查W.sum(axis1)是否有 0有就放宽距离阈值或手动指定最近邻。自变量共线性用 VIF 检查VIF 大于 10 的变量删掉。4.3 现象SDM 的 θ 和 ρ 都不显著但 SAR 的 ρ 显著原因SDM 引入了太多参数截面数据样本量不够自由度损失严重。解决如果 LR 检验接受 SDM 退化成 SAR就直接用 SAR。别硬上 SDM截面数据一般建议样本量至少 100 以上再考虑 SDM。4.4 现象直接效应和间接效应的分解结果和 β、θ 差很远原因ρ 较大时(I - ρW)^{-1} 会放大 β 和 θ 的效应直接效应不等于 β 是正常的。解决解释结果时以分解后的直接效应和间接效应为准不要直接拿 β 说事。如果 ρ 接近 1分解结果会很不稳定这时候要考虑模型是否过度拟合。4.5 现象换了距离阈值ρ 的显著性变了原因空间权重矩阵的构造对结果影响很大阈值太小导致邻居太少阈值太大导致空间效应被稀释。解决做敏感性分析试 3 到 5 个阈值看 ρ 的符号和显著性是否稳定。如果只在某个阈值下显著结论不可靠。5. 进阶技巧用空间计量做政策评估的边界与验证截面数据做空间计量最容易被人质疑的就是内生性。ρWy 里的 y 和误差项相关OLS 估计有偏所以必须用 ML 或 IV。ML 是默认选择但样本量小于 50 时 ML 的小样本性质不好这时候可以考虑 GMM。不过截面数据做 GMM 需要找工具变量一般用 WX 或 W 的高阶矩实操中不太好找。另一个进阶方向是空间断点回归。如果你有一个政策在某个边界上实施边界两侧的样本可以看作准自然实验这时候空间权重矩阵可以构造为“是否在边界同一侧”ρ 的解释就变成了政策溢出效应。这个做法在区域经济评估里越来越常见但对数据要求很高边界两侧的样本要足够多。验证模型是否靠谱我一般会做两个检查。第一把 W 矩阵随机置换 1000 次看 ρ 的分布是否集中在 0 附近。如果随机置换后的 ρ 仍然显著说明模型抓到的不是空间效应而是某种全局趋势。第二留出法交叉验证随机删掉 10% 的样本用剩下的估计模型预测删掉的样本看预测误差和 OLS 比是否更小。# 随机置换检验 rho_perm [] for _ in range(1000): perm np.random.permutation(n) W_perm W_std[perm][:, perm] sar_perm ML_Lag(y, X, wW_perm) rho_perm.append(sar_perm.betas[0]) # 看真实 rho 在置换分布中的位置 p_value np.mean(np.abs(rho_perm) np.abs(sar.betas[0])) print(f置换检验 p 值: {p_value:.4f})这个检验跑起来慢1000 次大概要十几分钟但能有效排除伪空间效应。我一般只在论文投稿前跑一次日常分析用 Morans I 就够了。最后说个血泪经验空间计量的结果对 W 矩阵极其敏感同一份数据换一种 W 构造方式ρ 可能从 0.3 跳到 0.7。所以别只报一个 W 的结果至少报邻接矩阵和反距离矩阵两套让读者自己判断稳健性。截面数据没有时间维度你没法用固定效应去吸收空间异质性W 的选择就是你的核心假设写论文时要把构造理由讲清楚。希望帮到你。本文还有配套的精品资源点击获取