做空间计量这几年我踩过最多的坑几乎都集中在“不知道用哪个命令、权重矩阵到底怎么构建”这些最基础的地方。很多人打开Stata第一反应是去搜“stata下载”“stata安装包”装完就急着跑回归结果连莫兰指数都没算就开始做空间杜宾模型最后审稿人一问“你的空间依赖性是真实存在的吗”就直接卡住。今天这篇就把莫兰指数和Stata里16个sp系列空间计量命令一次性讲透从数据准备到结果解读给你一条能直接照做的路径。先说结论在Stata里做空间计量你不需要学一堆杂七杂八的外部命令核心就围绕官方sp系列命令展开再配合几个经典的第三方spat命令做诊断。官方命令负责建权重矩阵、估计模型、分解效应spat系列负责算莫兰指数和做空间依赖诊断。两者搭配起来横截面数据也好、面板数据也好基本都能覆盖日常研究的全部需求。1. 莫兰指数动手建模前必须先回答的问题1.1 莫兰指数到底在检验什么空间计量和普通计量最大的区别就是承认“邻居会影响你”。一个地区的房价高往往不是因为那个地区本身条件好而是因为隔壁地区房价高、经济活跃这种空间溢出效应在普通OLS里是被忽略的。莫兰指数Morans I就是用来回答“这种空间依赖性到底存不存在”的统计量。它的计算逻辑很直白。你先要有一个空间权重矩阵WW的第i行第j列表示地区j是不是地区i的邻居或者地区i和地区j之间的距离有多近。莫兰指数的公式是[ I \frac{n}{\sum_{i}\sum_{j}w_{ij}} \cdot \frac{\sum_{i}\sum_{j}w_{ij}(x_i-\bar{x})(x_j-\bar{x})}{\sum_{i}(x_i-\bar{x})^2} ]看着复杂其实核心思想就是把某个变量x的取值和“邻居的x取值”做相关分析。如果I为正且显著说明高值地区旁边也是高值低值旁边也是低值存在正向空间聚集如果I为负说明高值旁边是低值空间分布存在排斥效应如果I接近0说明变量在空间上是随机分布的做不做空间计量就值得商榷。用生活场景打个比方你在一个菜市场里观察摊位的价格如果贵的摊位总是扎堆在一起便宜的也总是聚在另一个角落那这个菜市场的价格就有正向空间自相关。如果贵的和便宜的交替分布那就是负相关。如果价格高低完全随机那就没有空间自相关。莫兰指数做的正是这件事只不过把“摊位”换成了“地区”。1.2 全局莫兰指数的Stata实现步骤先提醒一句计算莫兰指数之前你必须有空间权重矩阵。很多人忽略了这个前置步骤直接拿变量跑spatgsa结果报错“variable not found”或者干脆算出一个无法解释的数值。在Stata里计算全局莫兰指数我推荐直接用第三方命令spatgsa。使用方法如下ssc install spatwmat, replace ssc install spatgsa, replace先安装命令然后准备一个包含地区ID、变量值和权重矩阵的数据。权重矩阵可以用spatwmat生成它支持从shapefile文件读取邻接关系也支持根据经纬度计算距离权重。use yourdata.dta, clear spatwmat using yourdata_shape.dta, name(W) xcoord(x) ycoord(y) spatgsa var1 var2, weights(W) moran这里的var1、var2是你想检验的变量weights(W)指定权重矩阵moran选项输出莫兰指数。运行结果会给出Morans I的值、期望值、标准差和z统计量。z值大于1.65单侧5%显著性水平基本可以认为存在显著的空间自相关。如果你用的是官方sp矩阵也可以从spreg回归后的estat moran得到残差的莫兰指数这个后面第4节会细讲。1.3 局部莫兰指数与LISA找出“热点”在哪里全局莫兰指数只能告诉你“整体上有没有空间聚集”但它不会告诉你具体是哪个地区在带动这个聚集。这时候就要用局部莫兰指数Local Morans I也就是LISA。局部莫兰的计算逻辑类似于把全局莫兰拆解到每个地区对每个地区i计算它和邻居之间的局部相关性。Stata里实现这个的思路比较简单用spatlsa命令ssc install spatlsa, replace spatlsa var1, weights(W) moran这个命令不仅会输出每个地区的局部莫兰指数值还能生成莫兰散点图。散点图分成四个象限第一象限是高-高聚集区也就是所谓的热点第三象限是低-低聚集区是冷点第二象限是低-高区域第四象限是高-低区域这两种都是空间异质性比较明显的地方。搞研究时如果审稿人问你“空间聚集主要发生在哪些区域”你就可以把LISA图贴出来配合显著性地图说明。实操时要注意spatlsa默认的显著性检验是基于随机置换的所以跑的Bootstrap次数越多越稳定。如果数据量比较大建议把permutation次数设成999以上得到的结果才经得住推敲。2. 为什么建议直接使用Stata官方sp系列命令2.1 官方sp命令 vs 第三方spat命令很多人一搜“Stata空间计量”搜出来的全是spatwmat、spatreg这类第三方命令。这些命令确实老牌在Stata 15之前几乎是唯一选择。但Stata 15之后出了一套官方空间计量命令也就是sp开头的那一批包括spset、spmatrix、spgenerate、spreg、spxtreg等。这套命令和第三方命令最大的区别是它把空间权重矩阵、空间滞后变量、回归估计整合到了一个统一的框架里数据管理更规范输出也更友好。我个人的建议是能用官方命令做估计的尽量用官方命令第三方spat系列命令主要用于做莫兰指数和空间依赖诊断这些官方命令覆盖得不够细的部分。两者不是对立关系而是互补关系。2.2 理解官方sp系列的操作主线官方sp命令的设计逻辑其实就三步定义空间数据、构建空间权重矩阵、估计空间模型。这三个步骤对应三组命令。第一步是用spset命令告诉Stata“这份数据是空间数据”需要指定唯一ID变量、坐标变量经度纬度或者链接外部shapefile。如果不先spset后面的spmatrix和spreg根本没法用。第二步是用spmatrix命令创建和管理空间权重矩阵。spmatrix create contiguity会基于邻接关系生成0-1权重矩阵spmatrix create idistance会基于逆距离生成权重矩阵。创建后可以用spmatrix save保存spmatrix dir查看有哪些矩阵spmatrix export导出。第三步是估计模型。横截面数据用spreg面板数据用spxtreg工具变量场景用spivreg。回归完之后还能用estat moran检验残差是否还有空间自相关用estat impact分解直接效应和间接效应。这套流程的好处是权重矩阵和估计命令高度集成不容易出错。很多第三方命令之间依赖关系混乱一个矩阵格式不匹配就会报各种莫名其妙的错误。2.3 牢记这3个基础命令其他都顺理成章开始学习sp系列时不要想着一次把16个命令全记住先掌握三个核心命令即可。第一个是spset这是空间数据声明命令。第二个是spmatrix create这是权重矩阵构建命令。第三个是spreg这是模型估计命令。把这三个命令串起来一个最简空间自回归模型就出来了use data.dta, clear spset id, modify spmatrix create contiguity W spreg x1 x2 x3, dv(y) ml weights(W)运行结果会报告空间自回归系数rho如果rho显著为正说明存在正向空间溢出效应。学到这里空间计量的核心操作就已经完成了。剩下的命令更多是在这个基础上做扩展比如处理面板数据、加入工具变量、做更复杂的权重矩阵等。3. 16个sp系列命令逐一拆解3.1 数据准备类spset、spbalance、spgeneratespset是进入空间计量世界的第一道门。这个命令的作用是声明数据的空间属性也就是告诉Stata“哪个变量是唯一ID哪个变量是经度哪个变量是纬度”。spset id, coords(xcoord ycoord) modify如果你有shapefile可以写成spset id, modify replaceStata会自动匹配坐标信息。spset执行后可以用spset summarize查看空间数据概况比如有多少个唯一多边形、坐标变量的范围。spbalance是用来处理空间面板数据平衡性的。面板数据经常存在某些年份缺值的情况缺值会导致权重矩阵和模型估计出现问题。spbalance的作用就是检查并转换数据让每个时空单元在时间维度上保持平衡。它有几个选项fill表示填充缺失的观测base()表示指定基准年份generate()可以把是否平衡的信息存成新变量。spgenerate是一个被低估的神器。它的功能是生成空间滞后变量。所谓空间滞后就是邻居变量的加权平均。比如你想创建一个“邻居人均GDP”的变量可以先创建权重矩阵W然后spgenerate lag_gdp W * gdp这里的spgenerate会把W矩阵和gdp变量做矩阵乘法生成一个lag_gdp变量。这个变量在做空间滞后模型的稳健性检验、画莫兰散点图、甚至做空间异质性分析时都很有用。3.2 权重矩阵类spmatrix、spdistancespmatrix是整个官方sp系列的地基。它负责创建、导入、导出、管理空间权重矩阵。常用的创建方式有* 邻接矩阵根据共同边界或顶点相连 spmatrix create contiguity W1 * 逆距离矩阵距离越近权重越大可设定阈值 spmatrix create idistance W2, cutoff(100) power(1)邻接矩阵适合处理行政区域数据比如省、市、县这类有明显边界的数据逆距离矩阵适合处理像企业、银行网点这类没有行政边界但有经纬度坐标的数据。spmatrix还支持spmatrix import从外部文件导入自定义矩阵比如从GeoDa导出的gal文件或者从ArcGIS导出的权重矩阵。这个功能非常实用因为很多审稿人要求你用不同权重矩阵做稳健性检验你就需要反复创建和切换矩阵。spdistance是计算两两地区之间距离的命令。它会生成一个距离矩阵你可以把它当作构建自定义权重矩阵的基础。通常情况下如果你用spmatrix create idistance就不需要单独用spdistance但当你要做距离衰减阈值分析、或者要检验距离的某种非线性影响时spdistance能给你更多灵活性。它生成的距离矩阵可以用spmatrix save保存也可以导出成数据文件供其他软件使用。3.3 横截面估计类spreg、spregcs、spivregspreg是官方横截面空间回归命令支持三种主要模型空间自回归模型SAR、空间误差模型SEM、空间杜宾模型SDM。它的语法很直观* 空间自回归模型 spreg y x1 x2, dv(y) ml weights(W) * 空间误差模型 spreg y x1 x2, error(1) dv(y) ml weights(W) * 空间杜宾模型 spreg y x1 x2, dv(y) ml weights(W) durbin(x1 x2)选择哪个模型主要看研究假设。如果你认为被解释变量存在空间溢出效应比如一个地区的房价会影响邻近地区房价用SAR如果你认为误差项存在空间相关可能是遗漏了某个空间相关的变量用SEM如果你认为解释变量的空间滞后也影响被解释变量用SDM。常用的模型选择策略是先用spatdiag做LM检验再结合经济学理论决定。spregcs是spreg的扩展它把几种常见空间模型纳入一个更灵活的框架可以通过选项组合空间滞后和空间误差项还能输出直接效应、间接效应和总效应的分解。这个命令特别适合做政策评估类研究因为你不仅要关心某个变量的回归系数还要关心它对本地区的影响直接效应和对邻近地区的影响间接效应。spivreg是空间工具变量回归命令。当你的解释变量存在内生性时普通spreg会得到有偏估计。spivreg允许你指定工具变量同时处理空间滞后项和内生解释变量。语法类似ivregress但多了一个weights(W)选项来指定空间权重矩阵。这个命令用得相对少但一旦遇到内生性问题它就是救场的存在。3.4 面板估计类spxtreg、spxtregsspxtreg是面板数据的空间回归命令。语法和xtreg非常相似区别在于多了权重矩阵和空间效应选项。它支持固定效应和随机效应模型* 固定效应空间自回归模型 spxtreg y x1 x2, fe dv(y) weights(W) * 随机效应空间自回归模型 spxtreg y x1 x2, re dv(y) weights(W)面板数据的核心优势是能控制个体异质性。固定效应模型主要消除不随时间变化的地区特征影响随机效应模型则假设个体效应与解释变量不相关。选择fe还是re可以用传统的Hausman检验但要注意空间滞后项的存在可能影响检验统计量的渐进性质建议配合理论判断。spxtregs是spxtreg的扩展版本它支持更复杂的空间效应设定包括空间固定效应、时间固定效应以及空间误差项。如果你的数据时间跨度比较长面板结构复杂spxtregs会是更好的选择。比如同时控制地区和年份固定效应的空间杜宾模型spxtregs就可以通过选项组合实现。3.5 经典诊断类spatwmat、spatgsa、spatlsa、spatdiag、spatcorr、spatreg这六个命令来自第三方但它们在空间计量诊断环节价值极高。spatwmat用于生成空间权重矩阵支持从shapefile读取邻接关系也支持基于坐标计算距离矩阵。虽然官方spmatrix功能更全但spatwmat生成的矩阵格式是spat系列命令通用的所以在做莫兰检验时我常常用spatwmat生成权重矩阵再喂给spatgsa。spatgsa是全局莫兰和Geary指数的计算命令输出结果简洁清晰是做空间自相关初步检验的首选。spatlsa用于局部莫兰指数计算和莫兰散点图绘制能识别出高-高聚集、低-低聚集、高-低离群等不同类型。spatdiag是OLS回归后的空间依赖诊断命令对OLS残差做LM检验和稳健LM检验帮你判断该用SAR还是SEM。spatcorr可以计算空间相关图在不同距离范围内考察空间自相关随距离变化的趋势。spatreg是第三方空间回归估计命令虽然估计方法相对传统但对于某些特殊模型设定仍有自己的优势适合作为稳健性检验的补充。这个夜间诊断组合拳的思路是先用spatgsa确认存在空间依赖再用spatlsa观察聚集位置然后跑OLS用spatdiag判断模型形式最后再用spreg或spatreg正式估计。4. 完整实操案例从莫兰指数到空间回归模型4.1 数据准备和权重矩阵构建为了讲清楚操作路径我用一个模拟的横截面案例来演示。假设有100个地区每个地区有GDP、教育支出等变量以及经纬度坐标。研究目的是检验GDP是否存在空间聚集以及某解释变量对GDP的影响是否具有空间外溢。第一步先把数据读入Stata并用spset声明空间属性。由于没有shapefile我使用经纬度坐标来声明use spatial_data.dta, clear spset id, coords(xcoord ycoord) modify这一步至关重要。如果不做spset后续一切sp开头命令都会报错。spset执行后用spdescribe查看数据状态确认ID无重复、坐标变量有效。第二步创建空间权重矩阵。我选择基于距离的逆距离矩阵因为100个地区分布比较分散单纯用邻接矩阵会导致很多地区没有邻居权重矩阵过于稀疏spmatrix create idistance W, cutoff(50) power(1)cutoff(50)表示只有距离在50单位以内的地区才算邻居power(1)表示权重取距离的倒数即距离越近权重越大。4.2 全局莫兰指数检验与结果解读直接用官方spmatrix创建的矩阵计算莫兰指数一个快捷方法是先用spgenerate生成空间滞后变量再用corr命令计算Morans I。更标准的做法是使用spatgsa。两者各有利弊spatgsa更正式能给出期望值和标准化统计量。* 方式一用spatgsa spatwmat using spatial_data.dta, name(W) xcoord(xcoord) ycoord(ycoord) spatgsa gdp, weights(W) moran结果可以看到Morans I约为0.35z值约为6.8p值小于0.001说明GDP在空间上存在显著的正向聚集效应。这意味着GDP高的地区倾向于和GDP高的地区相邻空间计量建模的合理性得到了初步验证。如果你想使用官方spmatrix创建的权重矩阵来完成检验可以通过spgenerate Wgdp W * gdp corr gdp Wgdp这个相关系数虽然不等于莫兰指数但也能反映空间滞后的相关程度。要获得比较严格的检验结果还是建议使用spatgsa。4.3 从OLS残差诊断到模型选择在莫兰指数显著之后下一步不是直接套模型而是先做一个普通OLS回归然后对残差做空间依赖诊断。这一步的作用是判断应该采用SAR、SEM还是SDM。reg gdp edu invest predict e, resid spatdiag, weights(W)spatdiag会输出多种检验统计量包括Morans I、LM-error、LM-lag以及对应的稳健版本。如果LM-lag显著而LM-error不显著优先考虑空间滞后模型SAR如果LM-error显著而LM-lag不显著优先考虑空间误差模型SEM如果两者都显著再看稳健版本稳健lm-lag显著则选SAR稳健lm-error显著则选SEM。如果两个稳健统计量都显著就要考虑SDM或者更复杂的SAC模型。本例中LM-lag和稳健LM-lag均显著LM-error不显著因此选择空间自回归模型SAR是合理的。4.4 估计空间回归模型并解读直接与间接效应选定SAR模型后使用官方spreg命令估计spreg gdp edu invest, dv(gdp) ml weights(W) estat impact, nose结果中的rho为空间自回归系数显著为正说明GDP存在正向的空间溢出效应。edu和invest的回归系数是直接的边际效应但空间模型里解释变量对被解释变量的总影响不仅限于本地区还有对邻居地区的间接影响。estat impact会输出直接效应、间接效应和总效应。假设edu的直接效应是0.42间接效应是0.18总效应是0.60。你可以在论文里这样描述教育支出每提高1%不仅会显著促进本地区GDP增长0.42%还会通过空间溢出效应促进邻近地区GDP增长0.18%。这就是空间计量区别于普通回归的核心价值。如果你的数据能支撑面板模型可以把流程换成spset配合spxtreg其他诊断思路完全一致。面板数据还能控制地区固定效应缓解遗漏变量偏误。5. 常见问题与排查技巧5.1 shapefile导入报错spset无法识别坐标很多人拿到一份shapefile直接在Stata里use然后spset结果报错“variable not found”。原因很可能是Stata没有直接读取shapefile你需要先把shapefile转换为Stata格式。办法是用shp2dta命令shp2dta using province.shp, data(province_data.dta) coords(province_coords.dta) gencentroids(centroids)这个命令会把shapefile的属性数据和坐标数据分开导出。然后再用spset去连接。转换时要注意shapefile的坐标系如果是经纬度spset里不用额外设置如果投影坐标建议先处理坐标单位保证距离计算符合实际。5.2 权重矩阵行标准化警告在做空间计量时Stata经常会提示“row-standardization recommended”提醒你把权重矩阵标准化。这是一个关键问题。行标准化的目的是让每行权重之和为1这样空间滞后变量的含义就是“邻居变量的加权平均”解释起来更直观。如果你用的是spmatrix创建的矩阵可以在创建后使用spmatrix normalize W, row进行行标准化。如果你用的是spatwmat可以给spatgsa或spatdiag指定标准化选项比如standardize。忽视这个步骤可能导致莫兰指数计算出来的值偏大或偏小显著性检验也会失真。5.3 面板数据spxtreg报错面板ID变量不匹配spxtreg对数据结构很挑剔。它的要求是经spset声明的空间ID必须和面板ID完全对应且不能在时间维度上有缺失。常见的报错是“panel data must be strongly balanced”或者“spatial panel data are not balanced”。解决办法是先用xtset设置面板结构然后用xtdescribe查看是否有缺失年份再用spbalance调整。还有一些情况你的地区ID和形状数据里的ID变量名字不同导致spset时无法自动匹配。解决办法是在spset前先看一眼shapefile数据的ID变量比如用describe查看再用rename改成一致的名字或者用spset的id选项显式指定。5.4 局部莫兰指数的结果无法保存spatlsa输出的是每个地区的局部莫兰值、p值和散点图但它不像普通回归那样自动把所有结果存入e()或r()。你如果要导出LISA表格需要手动把结果保存下来。一个技巧是使用spatlsa的graph选项生成莫兰散点图然后用Stata的graph save保存图片。如果要导出数值可以在spatlsa运行的临时数据区里手动复制或者用log记录来截取结果。另一个更高效的方式是直接用spgenerate计算空间滞后变量然后用常见的统计命令画散点图spgenerate Wgdp W * gdp scatter Wgdp gdp, mlabel(id)这个散点图本质上就是莫兰散点图横轴是本地区值纵轴是邻居平均值。斜率为正且越陡空间正向相关越强。这个方法的好处是可以用到官方命令生成的权重矩阵不需要来回切换工具。5.5 一个容易忽视的细节莫兰指数对权重矩阵高度敏感最后分享一个容易被忽视的细节莫兰指数的结果对权重矩阵的选择非常敏感。同一份数据用邻接矩阵算出来可能显著用逆距离矩阵算出来可能不显著。这不是软件问题而是空间计量的固有特性。所以在论文里一定要报告你所用的权重矩阵类型和构建方式并且建议做多种权重矩阵的稳健性检验让结论不那么依赖单一阵设定。我在实际研究中通常会同时使用邻接矩阵和不同的距离阈值逆距离矩阵如果莫兰指数和空间回归系数在多种设定下都保持一致的方向和显著性结果才算是真正站得住脚。如果结果在不同矩阵下差异很大就需要仔细分析是什么原因造成的是矩阵太稀疏还是距离阈值选得不合理而不是直接挑一个“好看”的结果发出去。这个习惯帮我避开了不少审稿质疑。