EN中文

地球信息科学学报 >

2024 , Vol. 26 >Issue 3: 567 - 590

DOI: https://doi.org/10.12082/dqxxkx.2024.230413

地球信息科学理论与方法

基于多带宽局部多项式的时空地理加权分位数回归分析

  • 王守芬 , 1, 3, * ,
  • 王守霞 2 ,
  • 顾建祥 1, 3
展开
  • 1.上海市测绘院，上海 200063
  • 2.北京大学数学科学学院，北京 100871
  • 3.自然资源部超大城市自然资源时空大数据分析应用重点实验室，上海 200063

王守芬（1990— ），女，山东临沂人，硕士，工程师，注册测绘师，主要研究方向为测绘地理信息技术应用。E-mail:

Copy editor: 蒋树芳

收稿日期: 2023-07-19

  修回日期: 2023-09-20

  网络出版日期: 2024-03-31

基金资助

上海市2021年度“科技创新行动计划”社会发展科技攻关项目(21DZ1204100)

收起

Geographically and Temporally Weighted Quantile Regression Analysis Based on Multi-bandwidth Local Polynomial

  • WANG Shoufen , 1, 3, * ,
  • WANG Shouxia 2 ,
  • GU Jianxiang 1, 3
Expand
  • 1. Shanghai Surveying and Mapping Institute, Shanghai 200063, China
  • 2. School of Mathematical Sciences, Peking University, Beijing 100871, China
  • 3. Key Laboratory of Spatial-temporal Big Data Analysis and Application of Natural Resources in Megacities, Ministry of Natural Resources, Shanghai 200063, China
*WANG Shoufen, E-mail:

Received date: 2023-07-19

  Revised date: 2023-09-20

  Online published: 2024-03-31

Supported by

2021 Shanghai "Science and Technology Innovation Action Plan" Social Development Science and Technology Research Project(21DZ1204100)

Fold

摘要

在高斯-马尔可夫假设条件下，最小二乘估计是具有最小方差的最优线性无偏估计量，因此基于最小二乘估计的时空地理加权回归方法在满足此假设时可以获得最优估计，但现实中这些条件有时得不到满足。如果样本数据中存在异常值或者呈厚尾分布，最小二乘回归模型估计值可能会存在较大偏误。而分位数回归受异常值影响较小，相比最小二乘回归更为稳健且应用条件相对更为宽松。更为重要的是最小二乘回归模型只能探索解释变量对响应变量条件均值的影响，而分位数回归可以探索解释变量对响应变量分布的影响（如响应变量的多个分位数），可以挖掘到更为丰富的信息。本文在局部多项式估计原理基础上，提出了基于多带宽局部多项式的时空地理加权分位数回归模型，利用两步迭代估计方法得到系数估计，并且允许不同自变量（影响因素)的最优带宽可以不同。本文通过数值模拟，将该模型与时空地理加权最小二乘回归进行对比，基于分位数回归的系数估计的均方误差和平均绝对误差均比最小二乘估计量小（例如，在0.75分位数，基于最小二乘回归得到的系数估计的均方误差和平均绝对误差分别是基于分位数回归的10倍和4倍），说明本文的分位数回归具有稳健性且可以探索影响响应变量分布的因素。最后以上海市2017—2021年商品房住宅小区为案例对象，应用该方法，探究不同影响因素对不同分位数的住宅价格（如高位房价、中等房价、低位房价）的影响，说明了本文方法的实用性。实际数据研究表明同一个影响因素对不同水平房价的影响效果不同，即同一影响因素系数的时间分布和空间分布在高位房价、中等房价和低位房价存在明显差异，并且不同影响因素的最优带宽也存在差异；与基于最小二乘回归的MGTWR相比，本文的分位数回归模型对于异常值的存在更为稳健（删除1%极端值后基于分位数回归模型拟合的平均绝对误差的变化比基于最小二乘回归模型小1%），并且分位数回归模型可以探究多个水平房价如高位房价、中等房价和低位房价的影响因素。

关键词： 时空地理加权; 多带宽; 局部多项式; 分位数回归; 稳健; 全局分布; 异常值; 异方差

本文引用格式

王守芬 , 王守霞 , 顾建祥 . 基于多带宽局部多项式的时空地理加权分位数回归分析[J]. 地球信息科学学报, 2024 , 26(3) : 567 -590 . DOI: 10.12082/dqxxkx.2024.230413

Abstract

The geographically and temporally weighted regression method based on weighted least squares estimation achieves optimal estimates under the assumption of Gauss-Markov independent identical distributions. However, these conditions cannot be always satisfied. If there are outliers or heavy-tailed distributions in the data, the least squares estimates may be significantly biased. On the other hand, quantile regression is less affected by outliers and is more robust than least squares regression, which can be applied in a broader range of applications under more relaxed conditions. More importantly, the least squares regression model only focuses on the mean of the response, while quantile regression explores the global distribution of the response variable (e.g., quantiles of the response variable) and can obtain richer information. In this paper, we propose the geographically and temporally weighted quantile regression model based on the local polynomial estimation. This model allows for different optimal bandwidths for different explanatory variables and use a two-step estimation method to obtain the estimates of the coefficients. To illustrate the superiority of the proposed method, we compare the proposed method with the geographically and temporally weighted least squares regression through numerical simulations. The simulation results show that the mean square error and the mean absolute error of the coefficient estimates for the proposed quantile regression model are both smaller than those of the least squares regression model. For example, at the 0.75 quantile, the mean square error and mean absolute error of the coefficient estimates based on the least squares regression are 10 times and 4 times those based on the quantile regression, respectively. This indicates that our proposed method is robust and can explore the global distribution of the response variable compared to the least squares regression model. Finally, to illustrate the practical ability of the method, we apply it to the data of Shanghai's commercial residential neighborhoods from 2017 to 2021 to investigate the effects of different factors on residential prices at different quantiles (e.g., high house prices, medium house prices, and low house prices). The results show that the explanatory variables have different effects on house prices at different quantiles. The spatial and temporal distributions of the coefficients of the variables differ significantly among high house prices, medium house prices, and low house prices, and the optimal bandwidths for different explanatory variables also differ. Compared to the MGTWR based on least squares regression, the quantile regression model proposed in this paper is more robust with the presence of outliers. After removing 1% of extreme values, the change in the mean absolute error of the fitting based on the quantile regression model is 1% smaller than that based on the least squares regression model. Additionally, the quantile regression model can explore the factors affecting the different price levels of the housing such as the high house prices, medium house prices, and low house prices.

Key words： geographically and temporally weighted regression; multi-bandwidths; local polynomials estimate; quantile regression; robustness; global distribution; outlier; heteroscedasticity

1 引言

回归分析作为常用的空间统计方法，可以对空间关系进行建模、检查和探究。早期空间统计分析方法多从全局的角度假设空间变量关系是固定的，忽略了现实地理世界空间关系的异质性，因而并不能反映回归参数的真实空间特征。1996年，Brunsdon等[1]在总结前人关于局部回归和变参数研究的基础上，提出了基于空间变化关系建模的局部线性回归方法——地理加权回归模型(Geographically Weighted Regression，GWR)。地理加权回归通过建立空间范围内每个点处的局部回归方程，来探索研究对象在某一尺度下的空间变化及相关驱动因素，很好地解释了变量的局部空间关系与空间异质性，因而在相关领域得到了广泛应用。空间位置是影响城市住宅价格的关键因素，近年来地理加权回归分析在住宅价格影响因子建模分析方面受到了广泛关注。地理加权回归模型虽然考虑了事物的空间变化，却未考虑时间变化。而现实世界中变量间的关系或者结构往往还会随着时间而改变，例如市场趋势、通货膨胀等驱动因子的时间效应对于住宅价格的影响也不可忽略[2]。为了考虑时间动态性， Huang等[3]在地理加权回归模型中加入了时间因素，把时间和空间信息整合到权矩阵中，并提出了在建立的椭球坐标系中以观测位置的时空距离函数为权值的估计方法，称之为时空地理加权回归方法(Geographically and Temporally Weighted Regression)简称GTWR。近年来，GTWR一直被不断地创新和演化[4]，时空地理加权自回归模型[5]和基于局部多项式求解的GTWR[6]等分析方法相继被提出。
GTWR的加权最小二乘估计方法在高斯-马尔可夫假设条件下可以取得最优估计，现实中这些条件有时得不到满足。如现实中因变量的方差或者说模型的随机误差项的方差是不相同的，即异方差。异方差在时空分析中普遍存在。例如一线城市城区和郊区的房价变化存在明显差异。另一方面GTWR估计方法本质上属于NW(Nadaraya-Watson)核估计方法，即局部常数估计。我们知道NW方法存在边界效应，即系数函数在边界区域的估计偏误较大.由于GTWR模型的系数函数是时空坐标的三元函数，GTWR方法的边界效应会更加严重。局部多项式估计的偏误比NW估计小，且边界效应更小[7,8]，另一方面局部多项式估计可以消除异方差的影响，减小回归系数估计值偏误，提升拟合精度，因此，Wang等[9]将局部多项式回归与GWR方法相结合，提出了局部线性地理加权回归方法。对于时间异质性和空间异质性，最初考虑利用混合GWR(mixed GWR)模型，其参数分为全局和局部两种尺度特征[10]。尽管混合GWR能够对GWR模型参数估计尺度进行差异化呈现，但尺度特征区分仍然相对单一，缺乏对细节差异的精准呈现[11]。Fotheringham等[12]提出了MGWR(Multiscale Geographically Weighted Regression)，为模型中的每个变量产生单独的优化带宽，从而得到更精确的局部参数估计。Wu等[13]在GTWR的基础上提出了多尺度MGTWR(Multiscale Geographically and Temporally Weighted Regression)，该方法不仅考虑时间和空间因素，而且充分考虑了不同变量的不同尺度。
研究数据在收集过程中可能存在失误从而导致异常值的出现。另一方面，所研究数据也不一定为高斯分布，可能呈现厚尾分布。最小二乘估计对此类数据稳健性较差，因此Zhang等[14]在局部线性地理加权回归方法的基础上研究了回归方法的稳健性，并证明在采用局部多项式改进后，能有效消除异方差，提升估计精度。普通最小二乘 (Ordinary Least Squares，OLS) 回归模型关注的是均值，研究在给定某些解释变量的条件下响应变量的期望均值，模型估计方法为最小二乘法。有时候我们不仅需要研究解释变量对响应变量条件均值的影响，而且希望能够探索解释变量对响应变量分布的影响(如响应变量的多个分位数)，这种情况下则需要分位数回归。分位数回归在1978 年由 Koenker和 Bassett[15]提出，利用响应变量的条件分位数来建模，采用加权的最小绝对离差和(Weighted Least Absolute Deviation，WLAD)法进行回归系数估计，应用条件相对更为宽松，可以描述响应变量的全局特征，可以挖掘到更为丰富的信息。另一方面，最小二乘回归的高斯-马尔可夫条件在现实中有时不能满足。如果样本数据中存在异常值、异方差、或厚尾情况，最小二乘回归模型估计值可能会存在较大偏误。而分位数回归估计采用加权的最小绝对离差和法进行估计，通常不受离群点的影响，结果更为稳健[16]。
考虑到分位数回归、局部多项式估计和多带宽的优点，本文结合局部多项式估计和分位数回归，提出了基于多带宽局部多项式的时空地理加权分位数回归模型，并通过数值模拟和上海市住宅价格实际数据分析，说明所提出方法的优越性和实用性。

2 研究方法

本研究的基本框架如图1所示，包含建立模型及模型估计、蒙特卡洛仿真模拟及实际数据分析等部分。
图1 研究基本框架

Fig. 1 Framework of this study

2.1 时空地理加权回归模型

时空地理加权回归模型具有如下形式。
      Y i = β 0 u i , v i , t i + ∑ j = 1 p β j u i , v i , t i X i j + ε i                   （ i = 1,2 , ⋯ , n ）
式中： Y i和 X i 1 , X i 2 , ⋯ , X i p是因变量 Y和自变量 X 1 , X 2 , ⋯ , X p在观测位置 u i , v i , t i处的观测值； p表示自变量的个数, 在本文实证分析中因变量Y代表上海市住宅价格, 自变量 X 1 , X 2 , ⋯ , X p是影响上海市住宅价格 Y的 p个影响因素, 如楼龄、容积率、绿化率、地铁站距离和数量、学校距离、商场距离和三甲医院距离； β 0 , β 1 , ⋯ , β p为自变量 X 1 , X 2 , ⋯ , X p在时空观测位置 u i , v i , t i的未知效应函数（系数函数）, 即每个影响因素X如何影响上海市住宅价格Y； ε i , i = 1,2 , ⋯ , n为独立同分布的随机误差项, 其均值为零, 方差为 σ 2, n为观测样本量。上述模型本质上是一个半参数变系数模型[8,17]。下面具体介绍时空地理加权回归方法。对于研究区域内的任意一点 u , v , t, 首先需构造 u , v , t处的一组权重 w i u , v , t , i = 1,2 , ⋯ , n, 其中第 i个权重值表示第 i组观测值对拟合 Y值的贡献。令 X = 1 , X 1 , ⋯ , X p为 n × p + 1维自变量观测矩阵, W u , v , t = d i a g w 1 u , v , t , w 2 u , v , t , ⋯ , w n u , v , t为权重矩阵, Y = Y 1 , ⋯ , Y n T为观测因变量向量。由最小二乘法可得 u , v , t处的系数函数的估计为
β ^ u , v , t = β ^ 0 u , v , t , β ^ 1 u , v , t , ⋯ , β ^ p u , v , t T X T W u , v , t X - 1 X T W u , v , t Y
在估计系数函数的过程中，距离点 u i , v i , t i近的观测数据对估计 β k u i , v i , t i的影响比距离点 u i , v i , t i远的数据大，因此常用的权重函数 w u i , v i , t i的计算方法有距离阈值法、距离反比法、高斯(Gauss)核函数和截尾型(Bi-square)核函数等方法[18]。本文考虑高斯核函数法，即( u i, v i , t i)处的权重取为高斯核函数，如式（3）所示。
w i u i , v i , t i = e x p - d i j 2 h 2 = e x p - u i - u j 2 + v i - v j 2 h s 2 + t i - t j 2 h T 2 = w i u i , v i w i t i
式中： d i j为 u j , v j , t j到第 i组观测所在时空位置 u i , v i , t i之间的时空距离； h > 0为时空带宽参数， h S = h 2 / λ , h T = h 2 / μ分别为待定的空间和时间带宽参数， λ和 μ为比例因子,代表空间和时间距离对时空距离不同权重的影响，且 λ + μ = 1。 h S和 h T可以用交叉验证(Cross Validation, CV)、广义交叉验证(Generalized Cross Validation, GCV)、或AIC(Akaike Information Criterion)、BIC (Bayesian Information Criterion)准则来确定。

2.2 基于局部多项式估计的GTWR模型

时空地理加权回归模型(1)的上述估计方法本质上属于NW核估计方法，该方法存在边界效应[7-8]，因此在这里考虑局部多项式估计[7,8]。假设式(1)中的每个系数函数 β j u , v , t j = 1,2 , ⋯ , p关于空间位置坐标 u , v和时间坐标 t均有连续的 q + 1阶偏导数。对于研究区域内给定的任一点 u 0 , v 0 , t 0，由泰勒公式知，在 u 0 , v 0 , t 0的邻域内 β j u , v , t可用 u , v和 t的 q次多项式函数逼近，式(4)给出 q = 2的泰勒展开形式，即对 j = 0,1 , 2 , ⋯ , p，有：
β j u , v , t ≈ β j u 0 , v 0 , t 0 + β j u u 0 , v 0 , t 0 u - u 0 + β j v u 0 , v 0 , t 0 v - v 0 + β j t u 0 , v 0 , t 0 t - t 0 + 1 2 β j u u u 0 , v 0 , t 0 ( u - u 0 ) 2 + 1 2 β j v v u 0 , v 0 , t 0 ( v - v 0 ) 2 + 1 2 β j t t u 0 , v 0 , t 0 ( t - t 0 ) 2 + β j u v u 0 , v 0 , t 0 u - u 0 v - v 0 + β j v t u 0 , v 0 , t 0 v - v 0 t - t 0 + β j t u u 0 , v 0 , t 0 t - t 0 u - u 0
式中: β j u u 0 , v 0 , t 0 ,   β j v u 0 , v 0 , t 0和 β j t u 0 , v 0 , t 0分别表示 β j u , v , t关于 u , v和 t的一阶偏导数在 u 0 , v 0 , t 0处的值，为方便展示，后文取 q = 1，也就是常用的局部线性估计（也可利用高阶多项式函数逼近）。前面提到时空地理加权模型（式(1)）是一个半参数变系数模型，因此根据变系数模型的局部线性拟合[8,17]及GTWR方法，需最小化以下目标函数来得到系数函数 β j的估计：
∑ i = 1 n { Y i - [ β j u 0 , v 0 , t 0 + β j u u 0 , v 0 , t 0 u i - u 0 + β j v u 0 , v 0 , t 0 v i - v 0 + β j t u 0 , v 0 , t 0 t i - t 0 ] X i j } 2 K d 0 i h s t
式中: X i 0 ≡ 1 ， d 0 i为 u 0 , v 0 , t 0到 u i , v i , t i的时空距离， i = 1,2 , ⋯ , n ;   K ∙为高斯权重函数。为了方便展示，令 X T u 0 , v 0 , t 0 = X 1 u 0 , v 0 , t 0 , ⋯ , X n u 0 , v 0 , t 0， 其中 X i T u 0 , v 0 , t 0 = X 0 i , ⋯ , X p i ⊗ 1 , u i - u 0 , v i - v 0 , t i - t 0 , i = 1,2 , ⋯ , n，则上述最小二乘估计问题的解可用矩阵表示为：
P ^ u 0 , v 0 , t 0 = X T u 0 , v 0 , t 0 W u 0 , v 0 , t 0 X u 0 , v 0 , t 0 - 1 X T u 0 , v 0 , t 0 W u 0 , v 0 , t 0 Y
式中: P u 0 , v 0 , t 0 = β T u 0 , v 0 , t 0 , β u u 0 , v 0 , t 0 T , β v u 0 , v 0 , t 0 T , β t u 0 , v 0 , t 0 T T ， β u 0 , v 0 , t 0为各回归系数函数 β ^ j u , v , t j = 0,1 , 2 , ⋯ , p在 u 0 , v 0 , t 0处的值所构成的列向量。由式(6)可得系数估计：
β ^ u 0 , v 0 , t 0 = I p + 1 , 0 p + 1 , 0 p + 1 , 0 p + 1 P ^ u 0 , v 0 , t 0

2.3 多带宽局部多项式GTWR分位数回归模型

假定误差项 ϵ i { i = 1 , … , n }为独立同分布的 序列，且分布情况未知，则 Y的 τ条件分位数 为 m τ x = a r g m i n θ ∈ R E   ρ τ Y - θ ∣   X = x  。其中， ρ τ x = x τ I x ≥ 0 - 1 - τ I x < 0为分位数回归领域的分位数损失函数，有些地方也称之为对号函数(check function)。 I ∙为示性函数。分位数损失函数也可以等价地表示为：
ρ τ x = τ x                                             x ≥ 0 τ - 1 x                           x < 0
或 ρ τ x = x τ - I x < 0 ，希望 E ρ τ x - X ^尽可能小，最好等于0。多带宽局部多项式GTWR分位数回归模型由式(9)给出。
Y i = β τ 0 , h s 0 t 0 u i , v i , t i + ∑ j = 1 p β τ j , h s j t j u i , v i , t i X i j + ε τ i               （ i = 1,2 , ⋯ , n ）
其中 h s j t j > 0 ,   j = 0,1 , … , p为可变带宽，在此允许不同自变量 X 0 , ⋯ , X p的时空带宽可以不相同，这样更符合常理。根据变系数模型的局部多项式拟合及GTWR方法，分位数回归需最小化以下目标函数 Q n X , Y来得到系数函数 β τ j的估计：
Q n X , Y = ∑ i = 1 n ρ τ { Y i - ∑ j = 0 p [ β τ j u 0 , v 0 , t 0 + β τ j u u 0 , v 0 , t 0 u i - u 0 + β τ j v u 0 , v 0 , t 0 v i - v 0 + β τ j t u 0 , v 0 , t 0 t - t 0 ] X i j } K d 0 i h s j t j                
式中： K ∙ , d 0 i与最小二乘回归目标函数(5)中相同，最小二乘回归目标函数（式(5)）中的损失函数为平方损失，而分位数回归目标函数 Q n X , Y中的损失函数为分位数损失函数 ρ τ x，该损失函数在0处不可导，因此分位数回归不像最小二乘回归一样有显示解（式(6)）。解决以上最小化（式(10)）的问题一般通过迭代算法，如:单纯形算法（Simplex Method），内点算法（Interior Point Method），平滑算法（Smoothing Method）等。
首先假设 h s j t j = h，则可较容易通过最小化目标函数(10)得到分位数回归的系数估计和偏导数估计，记为 β ˜ τ u 0 , v 0 , t 0， β ˜ τ u u 0 , v 0 , t 0， β ˜ τ v u 0 , v 0 , t 0 β ˜ τ t u 0 , v 0 , t 0，此时带宽 h可以根据交叉验证或广义交叉验证法选取最优带宽。上文提到，允许不同自变量 X 0 , ⋯ , X p的时空带宽可不相同 h s j t j , j = 0,1 , … , p，对于多带宽的估计过程将在下一节给出。

2.4 多带宽局部多项式GTWR分位数回归模型的估计

对于多带宽分位数回归模型的估计，本文通过两步迭代估计法来实现。估计过程如下所示：
（1）首先假设 h s j t j = h，通过最小化目标函数(10)得到系数的初始估计：
β ˜ τ u 0 , v 0 , t 0 = β ˜ τ , 0 u 0 , v 0 , t 0 ,   β ˜ τ , 1 u 0 , v 0 , t 0 , ⋯ ,     β ˜ τ , p u 0 , v 0 , t 0 T
（2）对于自变量 X 0 ≡ 1 ，令 Y ˜ i 0 = Y i - ∑ j = 1 p β ˜ τ , j u i , v i , t i X i j。然后利用 { Y ˜ i 0 } i = 1 n对 { X i 0 } i = 1 n进行分位数回归，最小化(12)式目标函数 Q n 0 X 0 , Y ˜ 0得到系数 β τ , 0的估计 β ^ τ , 0。
Q n 0 X 0 , Y ˜ 0 = ∑ i = 1 n ρ τ { Y ˜ i 0 - [ β τ 0 u 0 , v 0 , t 0 + β τ 0 u u 0 , v 0 , t 0 u i - u 0 + β τ 0 v u 0 , v 0 , t 0 v i - v 0 + β τ 0 t u 0 , v 0 , t 0 t i - t 0 ] X i 0 } K d 0 i h s 0 t 0
通过交叉验证或广义交叉验证法选取自变量 X 0对应的最优带宽 h s 0 t 0(时间带宽 h t 0和空间带宽 h s 0)。
（3）对于自变量 X k ,   k = 1,2 , 3 , ⋯ ,   p，令 Y ˜ i k = Y i - ∑ j = 0 k - 1 β ^ τ , j u i , v i , t i X i j - ∑ j = k + 1 p β ˜ τ , j u i , v i , t i X i j。然后利用 { Y ˜ i k } i = 1 n对 { X i k } i = 1 n进行分位数回归，最小化式(13)目标函数 Q n k X k , Y ˜ k得到系数 β τ , k的估计 β ^ τ , k。
Q n k X k , Y ˜ k = ∑ i = 1 n ρ τ { Y ˜ i k - [ β τ k u 0 , v 0 , t 0 + β τ k u u 0 , v 0 , t 0 u i - u 0 + β τ k v u 0 , v 0 , t 0 v i - v 0 + β τ k t u 0 , v 0 , t 0 t i - t 0 ] X i k } K d 0 i h s k t k
通过交叉验证或广义交叉验证法选取自变量 X k对应的最优带宽 h s k t k（时间带宽 h t k和空间带宽 h s k）。
（4）计算 S O C f，若 S O C f大于给定的精度 δ （ δ在本文中取为 10 - 5），则重复步骤（2）和步 骤（3）；否则停止迭代。其中 S O C f表示各个可加项的相对变化，与Fotheringham等[12]的定义相同，由以下公式给出：
S O C f = ∑ j = 1 p ∑ i = 1 n ( β ^ τ , j n e w u i , v i , t i X i j - β ^ τ , j o l d u i , v i , t i X i j ) 2 / n ∑ i = 1 n ( ∑ j = 1 p β ^ τ , j n e w u i , v i , t i X i j ) 2
原始MGTWR的算法需要迭代计算得到最终不同变量的系数估计结果，由于本文的多带宽GTWR分位数回归模型不像最小二乘回归有显式解，因此在每次求解系数估计时需要迭代求解，这在一定程度上会增加计算负担。但分位数求解算法已经较为成熟，此外，对于分位数回归损失函数不可导造成的迭代负担也可通过Horowitz光滑方法[19]解决。通过一些方法加快上述多带宽迭代算法的收敛。影响算法收敛快慢的因素主要有：① 初值 β ˜ τ的选择； ② 最优带宽的选取方法； ③ 精度 δ的选取； ④ 指标 S O C f的选择。给定的评价指标 S O C f，若对精度需求较低，则可选择一个较大的 δ（如 10 - 3）；对于初值的选择，本文通过拟合GTWR给出，此外还可通过直接利用Y对每个变量做GTWR分位数回归并选择最优带宽给出；对于最优带宽的选择，在每次迭代中是通过对每个变量回归，利用广义交叉验证法选取，广义交叉验证相比交叉验证计算更快，另外，某些变量的最优带宽在几步迭代之后会趋于平稳，那么在之后的迭代过程中可不再为该变量选取最优带宽。
表1总结了现有的GWR及其变种以及本文方法的原理、公式、适用条件及参考文献。
表1 GWR及其变种以及本文方法的原理、公式、适用条件及参考文献.

Tab. 1 Principles, formulas, applicable conditions and references of GWR and its variants as well as the method in this paper

分类 原理 公式 编号 变量说明 适用条件 参考文献
地理加权回归（GWR） 基于局部最优思想
(局部常数估计）的
空间变系数回归模型		 $\begin{aligned}Y_{i}= & \beta_{0}\left(u_{i}, v_{i}\right)+\\ \sum_{j=1}^{p} \beta_{j}\left(u_{i}, v_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (15) X i j和 β j分别是位置( u i , v i)的第j个自变量及其系数； ε i是误差项； Y i是因变量 针对存在空间
非平稳性的空
间数据		 文献[1]
混合地理加权回归(Mixed-GWR) 在GWR模型基础上，结合OLS回归模型，
同时考虑数据的局部和全局变化情况，将
模型参数分为全局和局部两种尺度特征，
利用后向迭代算法
对模型进行估计		 $\begin{aligned}Y_{i}= & \sum_{k=1}^{P_{a}} \beta_{k}^{(a)} X_{i k}^{(a)}+\\ \sum_{l=1}^{P_{b}} \beta_{l}^{(b)}\left(u_{i}, v_{i}\right) X_{i l}^{(b)}+\varepsilon_{i}\end{aligned} $
$P_{a}+P_{b}=P$		 (16) Y i表示因变量， β ( a )表示常系数， ( X i 1 a , X i 2 a , ⋯ , X i P a a )表示全局自变量， ( X i 1 b , X i 2 b , ⋯ , X P b b )表示位置 u i , v i的局部自变量， β ( b )表示变系数 针对同时存在全
局平稳性和局部
非平稳性的空间
数据		 文献[20]
时空地理加权回归（GTWR） GWR模型基础上增
加时间维度		 $\begin{aligned}Y_{i}= & \beta_{0}\left(u_{i}, v_{i}, t_{i}\right)+ \\ \sum_{j=1}^{p} \beta_{j}\left(u_{i}, v_{i}, t_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (17) X i j 和 β j分别是位置( u i , v i , t i)的第j个自变量及其系数； ε i是误差项； Y i是因变量 针对存在时空
非平稳性的时
空数据		 文献[5]
多尺度地理加权回归（MGWR） 为模型中的每个关系产生单独的优化带宽，从而指示不同关系如何在不同空间尺度上操作，并产生更精确的局部参数估计 $\begin{aligned}Y_{i}= & \beta_{b w s_{0}}\left(u_{i}, v_{i}\right)+ \\ \sum_{j=1}^{p} \beta_{b w s_{j}}\left(u_{i}, v_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (18) X i j是位置( u i , v i)的第j个自变量； β b w s j u i , v i为位置 u i , v i处基于空间带宽 b w s j的第j个自变量的系数； ε i是误差项； Y i是因变量 针对多元空间数
据关系对应不同
变化尺度的情况		 文献[12]
距离度量-参数对应地理加权回归
(PSDM GWR)		 为模型中的每个空间关系产生不同的最优距离度量和空间带宽，从而更好地处理数据关系中的空间异质性，利用后向迭代算法对模型进行估计 $\begin{array}{l}Y_{i}=\beta_{0 i}\left(D M_{0}, b w_{0}\right)_{\left(u_{i}, v_{i}\right)+}\\ \sum_{j=1}^{p} \beta_{j i}\left(D M_{j}, b w_{j}\right)_{\times}\quad\left(u_{i}, v_{i}\right) X_{i j i j}+\varepsilon_{i}\end{array}$ (19) X i j为位置( u i , v i)的第j个自变量； β j i ( D M j , b w j )   u i , v i为第j个自变量在位置 u i , v i处基于距离度量 D M j和空间带宽 b w j的系数； ε i是误差项； Y i是因变量 针对多元空间数
据关系对应复杂
不同变化尺度的
情况		 文献[21]
文献[22]
局部线性地理加权回归（LGWR） 局部线性回归与地理加权回归方法相结合 $\begin{aligned}Y_{i}= & \beta_{0}\left(u_{i}, v_{i}\right)+ \\ \sum_{j=1}^{p} \beta_{j}\left(u_{i}, v_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (20) X i j 和 β j分别是位置( u i , v i)的第j个自变量及其系数； ε i是误差项； Y i是因变量 针对存在异方差
的空间数据		 文献[9]
混合时空地理加权回归(Mixed-GTWR) 结合OLS回归和GTWR，将特征变量
分成全局特征变量
和时空非平稳特征
变量，前者用于反映
全局平稳特性，可用OLS模型表示；后者
反映时空非平稳特性，可用GTWR模型表示		 $\begin{array}{c}Y_{i}=\sum_{k=1}^{P_{a}} \beta_{k}^{(a)} X_{i k}^{(a)}+ \\\sum_{i=1}^{P_{b}} \beta_{l}^{(b)}\left(u_{i}, v_{i}, t_{i}\right) X_{i l}^{(b)}+\varepsilon_{i} \\P_{a}+P_{b}=P\end{array}$ (21) Y i表示因变量， β ( a )表示常系数， β ( b )表示变系数， ε i是误差项； ( X i 1 a , X i 2 a , … , X i P a a ) , Y i是因变量，表示全局自变量 ( X i 1 b , X i 2 b … , X P b b )表示时空位置 u i , v i , t i的局部自变量， 针对全局时空平
稳特征和局部时
空非平稳特征同
时存在的现象		 文献[23]
多尺度时空地理加权回归(MGTWR) 通过后向拟合算法，
为模型中的每个自
变量与因变量的映
射关系产生独有的
优化时空带宽		 $\begin{aligned}Y_{i}= & \beta_{b w t_{0} s_{0}}\left(u_{i}, v_{i}, t_{i}\right)+ \\& \sum_{j=1}^{p} \beta_{b w t_{j} s_{j}}\left(u_{i}, v_{i}, t_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (22) X i j 为 第 j 个 自 变 量 β b w t j s j u i , v i , t i是位置 u i , v i , t i处基于特定空间带宽 b w s j和时间带宽 b w t j的第j个自变量 X i j的系数； ε i是误差项； Y i是因变量 GTWR限制了模型
内的每个局部关系
在单一尺度下进行
变化，MGTWR针对多元时空数据关系
对应不同变化尺度
的情况，解决空间
尺度效应问题		 文献[13]
局部多项式时空地理加权回归(LPGTWR) 局部多项式回归与
时空地理加权回归
方法相结合		 $\begin{aligned}Y_{i}= & \beta_{0}\left(u_{i}, v_{i}, t_{i}\right)+ \\& \sum_{j=1}^{p} \beta_{j}\left(u_{i}, v_{i}, t_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (23) X i j 和 β j分别是位置( u i , v i , t i)的第j个自变量及其系数； ε i是误差项； Y i是因变量 针对存在异方差的
时空数据		 文献[6]
基于多带宽局部多项式的时空地理加权分位数回归模型(本文方法) 分位数回归、局部多
项式估计、多尺度多
带宽与时空地理加
权回归方法相结合		 $\begin{aligned}Y_{i}= & \beta_{\tau 0, h_{s_{0} t_{0}}}\left(u_{i}, v_{i}, t_{i}\right)+ \\& \sum_{j=1}^{p} \beta_{\tau j, h_{s j_{j}}}\left(u_{i}, v_{i}, t_{i}\right) X_{i j}+\varepsilon_{i}\end{aligned}$ (24) X i j 为 第 j 个 自 变 量 β τ j , h s j t j表
示位置 u i , v i , t i处基于空
间带宽 h s j和时间带宽 h t j
的第j个自变量 X i j的 τ
分位数系数； ε i是误差
项； Y i是因变量		 针对不满足高斯-马尔可夫假设的时空数据，如存在异常值、异方差、或厚尾分布。针对多元时空数据关系对应不同解释变量的变化尺度的情况，解决时间、空间尺度效应问题。研究解释变量对响应变量分布的影响，研究解释变量对响应变量多个分位数的影响 本文

3 模拟研究

本节采用蒙特卡洛仿真实验来探索本文方法的有限样本表现，通过比较本文的基于多带宽局部多项式的时空地理加权分位数回归模型(Q)和多带宽时空地理加权最小二乘回归模型(LS)在不同仿真实验设定下的表现来展现本文方法的优越性。设置自变量 X的个数为 p = 3，通过以下式(25)产生数据Y。
Y i = β 0 u i , v i , t i + X i 1 β 1 u i , v i , t i + X i 2 β 2 u i , v i , t i + X i 3 β 3 u i , v i , t i + ε i τ
式中： u i , v i , t i从 0,1上的均匀分布中生成，自变量 X i 1, X i 2和 X i 3分别从标准正态分布，均匀分布和正态分布中产生，即 X i 1 ∽ N 0,1, X i 2 ∽ U 0,1, X i 3 ∽ N 1,4。误差项 ε i τ = ε i - F - 1 τ,其中 F为 ε i的条件分布函数， ε i从自由度为3的 t分布 t 3，也就是厚尾分布中生成。假定数据中存在4个异常值，即有4个 Y值与真实值相差较多。系数函数 β p ∙ , ∙ , ∙的值分别由式(26)函数产生。
β 0 u , v , t = 4 u v β 1 u , v , t = 3 β 2 u , v , t = 2 v t + 3 v 2 t 2 β 3 u , v , t = 4 u 2 v 2 t 2
在本小节中，设置样本量n=1 000，考虑分位数τ=0.5，τ=0.75和τ=0.25，重复模拟m=200次来估计系数函数。为评估本章估计方法的有限样本表现，使用均方误差(MSE)和平均绝对误差(MAE)评估系数函数的估计效果。
M S E β ^ p = 1 N t ∑ k = 1 N t 1 N v ∑ j = 1 N v 1 N u                                             ∑ i = 1 N u β ^ p u i , v j , t k - β p u i , v j , t k 2                     （ p = 0,1 , 2,3 ）
M A E β ^ p = 1 N t ∑ k = 1 N t 1 N v ∑ j = 1 N v 1 N u                                               ∑ i = 1 N u β ^ p u i , v j , t k - β p u i , v j , t k                   （ p = 0,1 , 2,3 ）
式中： u i ,   i = 1 , … , N u、 v j ,   j = 1 , … , N v和 t k , k = 1 , … , N t为 0.1 ,   0.9之间均匀的格子点，格子点个数选取为 N u = N v = N t =20。
表2给出了200次重复模拟下，系数函数分位数回归估计和最小二乘估计的平均MSE(AMSE)、平均MAE(AMAE)以及对应的标准差，括号里为200次模拟得到的系数函数估计的MSE(MAE)的标准差，其中MSEQ和MAEQ分别表示基于本文的时空地理加权分位数回归模型(Q)得到的系数函数估计的MSE和MAE，MSELS和MAELS表示基于时空地理加权最小二乘回归模型(LS)得到的系数函数估计的MSE和MAE。表2表明在0.5分位数( τ = 0.5)，0.75分位数( τ = 0.75)和0.25分位数( τ = 0.25)下，所有系数的分位数回归估计量的AMSE和AMAE均比最小二乘估计量小，且MSE和MAE的标准差也较小。当分位数( τ = 0.75)时，分位数回归与最小二乘回归结果差距较大，尤其是系数β0和β2的表现，基于最小二乘回归系数估计的AMSE是基于分位数回归系数估计的10倍，而基于最小二乘回归系数估计的AMAE是基于分位数回归系数估计的4倍，说明数据存在异常值，数据呈厚尾分布时，时空地理加权分位数回归模型比时空地理加权最小二乘回归模型效果表现更好，估计结果更为准确也更为稳健。
表2 通过2种方法得到的系数函数估计的 MSE 和 MAE 的均值和标准差(n=1000, m=200)

Tab. 2 The mean and standard deviation of MSE and MAE of the estimated coefficient function by two methods, Q and LS(n=1000, m=200)

τ β p MSEQ MSELS MAEQ MAELS
τ = 0.5 β 0 0.046 8 (0.022 5) 0.053 0 (0.026 0) 0.176 1 (0.044 1) 0.185 1 (0.046 8)
β 1 0.006 2 (0.003 6) 0.009 4 (0.006 5) 0.063 6 (0.019 5) 0.076 1 (0.026 6)
β 2 0.111 2 (0.049 9) 0.121 6 (0.061 3) 0.256 9 (0.062 3) 0.271 0 (0.064 5)
β 3 0.005 9 (0.002 0) 0.007 8 (0.002 8) 0.060 5 (0.012 0) 0.069 7 (0.012 9)
τ = 0.75 β 0 0.056 8 (0.022 2) 0.785 3 (0.170 1) 0.191 9 (0.044 6) 0.858 1 (0.099 1)
β 1 0.012 1 (0.009 0) 0.014 6 (0.009 5) 0.086 0 (0.032 9) 0.095 9 (0.031 3)
β 2 0.163 7 (0.088 5) 1.189 9 (0.390 5) 0.316 7 (0.086 4) 1.011 1 (0.178 5)
β 3 0.009 3 (0.003 2) 0.017 6 (0.007 3) 0.077 5 (0.014 8) 0.111 3 (0.027 4)
τ = 0.25 β 0 0.077 3 (0.048 0) 0.529 6 (0.158 3) 0.225 5 (0.073 5) 0.690 4 (0.112 8)
β 1 0.009 1 (0.005 8) 0.009 3 (0.007 2) 0.074 0 (0.024 8) 0.076 3 (0.027 9)
β 2 0.175 3 (0.101 4) 2.126 5 (0.598 8) 0.329 4 (0.105 8) 1.395 6 (0.210 1)
β 3 0.007 2 (0.002 7) 0.055 2 (0.013 0) 0.063 3 (0.014 8) 0.220 1 (0.028 5)
表3展示了当样本量增大到n=1 500时，系数函数0.5分位数( τ = 0.5)回归估计和最小二乘估计的平均MSE(AMSE)、平均MAE(AMAE)以及对应的标准差，可以看到，与表2中 τ = 0.5的结果相比，当样本量增大时，各个系数估计的平均MSE和平均MAE减小，说明样本量越大，系数估计结果越准确。表4展示了模拟次数增大到300次时，系数函数0.5分位数( τ = 0.5)回归估计和最小二乘估计的平均MSE(AMSE)、平均MAE(AMAE)以及对应的标准差，与表2中 τ = 0.5的结果相比，增大模拟次数，各个系数估计的平均MSE和平均MAE略有增大，但增加幅度十分小，这是随机重复次数增多造成的。
表3 通过2种方法得到的系数函数估计的MSE和MAE的均值和标准差 (n=1500, m=200)

Tab. 3 The mean and standard deviation of MSE and MAE of the estimated coefficient function by two methods, Q and LS (n=1500, m=200)

τ β p MSEQ MSELS MAEQ MAELS
τ = 0.5 β 0 0.040 4 (0.020 2) 0.041 0 (0.020 9) 0.165 1 (0.047 6) 0.167 2 (0.047 9)
β 1 0.004 5 (0.003 2) 0.009 3 (0.006 7) 0.052 5 (0.020 7) 0.076 1 (0.026 6)
β 2 0.094 1 (0.040 2) 0.115 5 (0.054 6) 0.245 4 (0.060 2) 0.274 2 (0.071 8)
β 3 0.006 3 (0.001 6) 0.008 9 (0.002 6) 0.058 3 (0.010 3) 0.074 0 (0.013 8)
表4 通过2种方法得到的系数函数估计的MSE和MAE的均值和标准差(n=1000, m=300)

Tab. 4 The mean and standard deviation of MSE and MAE of the estimated coefficient function by two methods, Q and LS(n=1000, m=300)

τ β p MSEQ MSELS MAEQ MAELS
τ = 0.5 β 0 0.048 7 (0.029 1) 0.054 9 (0.031 1) 0.177 0 (0.054 3) 0.199 0 (0.057 2)
β 1 0.006 6 (0.004 2) 0.010 9 (0.008 6) 0.064 5 (0.022 0) 0.081 2 (0.031 5)
β 2 0.128 4 (0.071 5) 0.194 3 (0.112 2) 0.280 7 (0.087 0) 0.353 4 (0.110 5)
β 3 0.006 3 (0.002 3) 0.009 8 (0.004 1) 0.060 8 (0.014 0) 0.080 5 (0.018 0)

4 实证分析

本文以上海市住宅价格时空分布格局与影响因素关系作为研究对象，通过时空地理加权回归的分位数回归模型研究分析上海市住宅价格的影响因素、每个因素如何影响住宅价格、不同分位数的住宅价格（高房价、中等房价、低房价）的影响因素及影响效果是否相同。

4.1 研究区域与数据源

作为我国的典型超大城市，上海的住宅市场一直备受关注，本文以上海市中心城区及商品房住宅交易活跃的郊区为主要研究区域，包括黄浦区、静安区、徐汇区、杨浦区、虹口区、长宁区、普陀区、宝山区、嘉定区、青浦区、松江区、闵行区及浦东新区13个区，如图2所示。以住宅小区为基本研究单元，选取时间区域为2017年1月1日—2021年12月31日，对研究区域内各住宅小区进行住宅价格信息采集。研究区域内相关住宅价格的数据均来自安居客（https://shanghai.anjuke.com/），以上海市出售的新房或二手房成交价作为住宅价格数据源，共计采集9 021个住宅小区57 189条住宅价格数据，即住宅价格为Y，样本量n=57 189，利用百度地图对各住宅小区进行空间定位，与住宅价格关联。
图2 研究区域

Fig. 2 Study area

区域住宅价格受到众多因素的影响，许多学者对住宅价格影响因素做了大量研究，这些研究文献主要考虑了建筑所在小区的基本属性、交通设施、医疗设施、生活设施等方面的影响指标[24⇓⇓⇓-28]。综合以往的研究文献，本文从住宅基本属性、小区基本属性、交通设施、医疗设施、教育设施、生活设施6个方面细化选取以下9个特征变量 X:住宅面积、住宅楼龄、小区容积率、小区绿化率、距离最近地铁站的距离、2 km以内的地铁站数量、距离最近综合商场的距离、距离最近三甲医院的距离以及距离最近小学的距离作为影响房价时空变动的主要因素，见表5。
表5 上海市住宅价格及影响因素的描述性统计分析

Tab.5 Variables influencing residential prices in Shanghai and their descriptive statistics

变量
类型		 特征变量 标签 最小值
(Min)		 第一四分位数(Q1) 中位数
(Median)		 均值(Mean) 第三四分位数(Q3) 最大值(Max)
因变量 住宅价格(price) Y 1.0 285 396 503 597 6 783
基本属性 面积(area) X 1 12 52 71 79.37 94 1 000
楼龄(old) X 2 -4 15 21 20.77 26 97
交通设施 距离最近地铁站的距离(subway) X 5 35 410 620 918.50 892 9 999
2 km以内的地铁站的数量(subway num) X 9 0 1 2 2.38 3 11
教育设施 距离最近小学的距离(school) X 8 0 468 741 868.9 1134 6 579
小区属性 容积率(volume rate) X 3 0.27 1.50 1.80 1.98 2.25 7.20
绿化率(greening) X 4 0.10 0.30 0.35 0.35 0.40 0.80
生活设施 距离最近综合商场的距离(mall) X 6 0 830 1 311 1 711 2 041 15 309
医疗设施 距离最近三甲医院的距离(hospital) X 7 226 3 971 6 163 8 515 11 193 38 684
表5给出了因变量Y和各个自变量的描述性统计量，即最小值，第一四分位数，中位数，第三四分位数，最大值。从表5中因变量Y的描述性统计量以及图3给出的2017—2021年上海市住宅价格Y的箱线图可以看出，上海市住宅价格数据存在明显的异常值。表5中住宅和小区基本属性因素数据（面积、楼龄、容积率、绿化率）来源于安居客（https://shanghai.anjuke.com/），交通、医疗、教育、生活设施特征变量（距离最近地铁站的距离、2 km以内的地铁站的数量、距离最近三甲医院的距离、距离最近小学的距离、距离最近综合商场的距离）通过各类POI点数据与住宅价格点数据计算所得，相关POI点数据由百度地图获取。地理加权回归分析中常采用欧氏距离作为度量空间距离的标准，欧氏距离未考虑空间可达性，本文采用路网距离，基于交通网络（来源于百度地图路网）通过路径 规划算法计算最短路径距离，更能科学地表达实际城市中的空间距离。本文所用数据汇总如表6所示。
图3 2017—2021年上海市住宅价格Y箱线图

Fig. 3 Boxplot of the residential price in Shanghai, 2017-2021

表6 研究数据汇总

Tab. 6 Summary table of the research data

数据名称 来源 详情
住宅价格、住宅楼龄
住宅面积、小区容积率
小区绿化率		 安居客
（https://shanghai.anjuke.com）		 通过2017—2021年上海各住宅小区成交记录采集成交价格、
面积、楼龄、小区容积率、小区绿化率、小区名称、所在区县7个属性（共9 021个小区，57 189条成交记录）
小区空间位置 百度地图 基于小区名称及所在区县获取小区空间位置坐标
综合商场、三甲医院
小学、地铁站、路网		 百度地图

4.2 结果及分析

2017—2021年上海市住宅价格的空间分布如图4所示。从房价的空间分布图可以看出，高位房价主要集中在市中心，郊区房价相对较低，且房价的空间分布也随着时间的改变有所改变。为进一步研究影响房价的因素、每个因素如何影响房价、不同分位数房价的影响因素及影响效果是否相同，以及多带宽影响等，利用第2节介绍的多带宽时空地理加权分位数回归模型，第2.4节介绍的两步估计方法得到不同分位数的系数函数估计。在此我们考虑3个分位数 τ = 0.25 、 0.50 、 0.75，分别对应低位房价，中等房价，高位房价。其中初始估计的最优时间、空间带宽和每个自变量 X j , j = 0,1 , … , p对应的最优时间、空间带宽均通过广义交叉验证选取。
图4 2017—2021年上海市住宅价格分布

Fig. 4 Residential price distribution in Shanghai, 2017-2021

（1）不同影响因素在低位房价、中等房价、高位房价的时间及空间带宽
表7展示了不同分位数房价回归中每个自变量 X j , j = 0,1 , … , p对应的最优时间带宽 h t j、空间带宽 h s j、时空带宽 h s j t j及比例因子 λ j。从表7可以看出常数项( X 0)、面积( X 1)、容积率( X 3)、绿化率( X 4)、距离最近三甲医院的距离( X 7)、距离小学距离( X 8)以及2 km以内的地铁站数( X 9)的时间带宽与空间带宽都相对较小，而楼龄( X 2)、距离最近地铁的距离( X 5)、商场距离( X 6)的时间带宽与空间带宽都相对较大。时间带宽与空间带宽较大表明变量没有明显的时间与空间异质性，说明这些变量对住宅价格的影响在时间与空间上是全局的。不同变量的时间带宽与空间带宽不相同则说明不同的变量对房价影响的时间异质性与空间异质性水平也不相同。
表7 不同自变量$X_{j}$在不同分位数的时间带宽$h_{t_{j}}$、空间带宽$h_{s_{j}}$、时空带宽$h_{s_{j}t_{j}}$及比例因子$\lambda_{j}$

Tab. 7 Temporal bandwidth$h_{t_{j}}$, spatial bandwidth $h_{s_{j}}$, spatio-temporal bandwidth $h_{s_{j}t_{j}}$ and scaling factor$\lambda_{j}$. for different explanatory variables $X_{j}$ at different quantiles

分位数 τ 自变量 X j 时间带宽 h t j 空间带宽 h s j 时空带宽 h s j t j 比例因子 λ j
常数 ( X 0) 0.174 0 0.039 9 0.038 9 0.95
面积 ( X 1) 0.174 0 0.039 9 0.038 9 0.95
楼龄 ( X 2) 4.347 0 0.997 3 0.972 0 0.95
容积率 ( X 3) 0.696 4 0.159 5 0.155 5 0.95
τ = 0.25 绿化率 ( X 4) 0.869 4 0.199 5 0.194 4 0.95
地铁距离 ( X 5) 1.205 6 1.643 0 0.972 0 0.35
商场距离 ( X 6) 4.347 0 0.997 3 0.972 0 0.95
医院距离 ( X 7) 0.289 8 0.262 1 0.194 4 0.55
小学距离 ( X 8) 2.260 2 0.518 5 0.505 4 0.95
地铁站数 ( X 9) 1.564 8 0.359 0 0.349 9 0.95
常数 ( X 0) 0.174 0 0.039 9 0.038 9 0.95
面积 ( X 1) 0.077 8 0.044 9 0.038 9 0.75
楼龄 ( X 2) 1.122 0 1.944 0 0.972 0 0.25
容积率 ( X 3) 0.869 4 0.199 5 0.194 4 0.95
τ = 0.5 绿化率 ( X 4) 0.869 4 0.199 5 0.194 4 0.95
地铁距离 ( X 5) 4.347 0 0.997 3 0.972 0 0.95
商场距离 ( X 6) 4.347 0 0.997 3 0.972 0 0.95
医院距离 ( X 7) 0.311 0 0.179 6 0.155 5 0.75
小学距离 ( X 8) 2.434 2 0.558 4 0.544 3 0.95
地铁站数 ( X 9) 1.390 8 0.319 0 0.311 0 0.95
常数 ( X 0) 0.174 0 0.039 9 0.038 9 0.95
面积 ( X 1) 0.100 4 0.042 2 0.038 9 0.85
楼龄 ( X 2) 1.122 0 1.944 0 0.972 0 0.25
容积率 ( X 3) 0.231 8 0.209 7 0.155 5 0.55
τ = 0.75 绿化率 ( X 4) 0.232 8 0.134 4 0.116 4 0.75
地铁距离 ( X 5) 1.205 6 1.643 0 0.972 0 0.35
商场距离 ( X 6) 4.347 0 0.997 3 0.972 0 0.95
医院距离 ( X 7) 0.520 6 0.119 4 0.116 4 0.95
小学距离 ( X 8) 2.260 2 0.518 5 0.505 4 0.95
地铁站数 ( X 9) 0.869 4 0.199 5 0.194 4 0.95
另一方面，不同分位数 τ下的时间带宽与空间带宽也有所差别，这说明同一个变量对低位房价、中等房价、高位房价的影响的时间异质性与空间异质性水平也不相同。如楼龄 X 2在低位房价 τ = 0.25下时间带宽远大于空间带宽，而在中等房价和高位房价 τ = 0.50,0.75正好相反；距离最近地铁的距离( X 5)在中等房价 τ = 0.50下时间带宽远大于空间带宽，而在低位房价和高位房价 τ = 0.50 ,   0.75正好相反；距离最近的三甲医院距离( X 7)的时间带宽在高位房价 τ = 0.75最大，在低位房价 τ = 0.25最小，而空间带宽正好相反，说明三甲医院距离( X 7)在高位房价的空间异质性更大、时间异质性更小，在低位房价的时间异质性更大、空间异质性更小。
表8展示了不同分位数下的同带宽和多带宽拟合的平均绝对误差(MAE)，其中同带宽对应第2.4节介绍的两步估计法中的初始估计，多带宽对应最终估计。很明显在3个分位数水平下，多带宽拟合的平均绝对误差均小于同带宽相应值，说明多带宽拟合效果更好，不同变量对房价的影响存在不同的时间异质性和空间异质性。表8最后一行为基于最小二乘回归的MGTWR的结果，结果显示中位数回归（即 τ = 0.50）的结果优于最小二乘回归的结果。为了进一步说明分位数回归的稳健性，删除原始数据1%的极端值（即删除原始数据0.5%的最小值和0.5%的最大值）后分别利用本文的分位数回归模型（ τ = 0.50）和基于最小二乘回归的MGTWR（在表9中用LSE表示）拟合，结果展示在表9。在同带宽和多带宽情况下，分位数回归模型的平均绝对误差的变化较小（相比原始数据的MAE的下降比例小），说明分位数回归对于极端值或者异常值比最小二乘回归更具稳健性。
表8 不同分位数下同带宽和多带宽拟合的平均绝对误差

Tab. 8 MAE under equal bandwidth and muti-bandwidth at different quantiles

分位数 τ 带宽 平均绝对误差(MAE)
τ = 0.25 多带宽 100
同带宽 104
τ = 0.5 多带宽 81
同带宽 85
τ = 0.75 多带宽 97
同带宽 99
LSE 多带宽 87
同带宽 90
表9 删除1%极端值后中位数回归和最小二乘回归拟合的平均绝对误差及下降比例

Tab. 9 The MAE and drop ratio under median regression and least squares regression after removing 1% of extreme values

分位数 τ 窗宽 平均绝对误差(MAE) 下降比例/%
τ = 0.5 多带宽 76.4 5.7
同带宽 80.8 5.0
LSE 多带宽 81.1 6.8
同带宽 84.5 6.1
（2）不同影响因素的系数在低位房价、中等房价、高位房价的时间分布
再看不同分位数房价对不同自变量系数函数随时间的变化。图5—图6分别描述了在一个固定的空间位置点 u 0 , v 0，不同分位数下的面积、楼龄、容积率、绿化率、三甲医院距离和学校距离等各系数估计值随时间的变化情况。从图5—图6看出，所有考虑的自变量的系数估计 β ^ j u 0 , v 0 , t均随时间有变化，即所考虑的9个自变量对不同分位数的房价影响均随时间变化。例如，从图5(a)—图5(c)可以看出从2017—2021年，3个分位数的面积系数估计值均先下降后上升，说明面积对房价的正向影响效果先减少后增加。图5(d)—图5(f)可看出从2017—2021年楼龄对低位、中等、高位房价影响随时间变化趋势存在明显不同。对于低位房价( τ = 0.25)，楼龄系数估计的绝对值随时间而减小，说明从2017—2021年楼龄对低位房价的影响效果减弱，而对于中等和高位房价( τ = 0.50 ,   0.75)，楼龄的影响效果越来越强，商场距离对房价影响的时间变化趋势与楼龄相同（图片省略）。图5(g)—图5(i)可看出容积率对不同水平房价的负向影响效果先增大后减弱。
图5 不同分位数下的系数估计的随时间的变化

Fig. 5 The changes of coefficient estimates along time at different quantiles

图6 不同分位数下的系数估计的随时间的变化

Fig. 6 The changes of coefficient estimates along time at different quantiles

从图6(a)—图6(c)可看出，从2017—2021年绿化率对不同分位数房价的正向影响效果均随着时间不断增大，地铁站数量的系数图与绿化率相同（图片省略），也就意味着人们在买房时越来越重视小区的绿化率和附近的交通便利程度。这也与2017—2021年几条新的地铁线的开通有关：2017年12月30日，17号线一期东方绿舟—虹桥火车站开通；2018年3月31日浦江线开通；2020年12月26日18号线长江南路—御桥开通；2021年1月23日15号线开通；2021年12月30日14号线开通。图6(d)—图6(i)可看出，从2017—2021年，三甲医院距离和学校距离对不同价格水平的房价的影响变化均存在差异，对于低位房价和中等房价( τ = 0.25,0.50)，三甲医院距离对低位房价和中等房价的影响效果先增加后减弱，学校距离的影响效果随时间减弱，而对于高位房价( τ = 0.75)，三甲医院距离和学校距离的影响效果越来越强。需要特别说明，本文研究的时间跨度为2017—2021年，在此期间会有通货膨胀，但目前在研究中并未考虑通货膨胀的效应，在未来的研究中可以将通货膨胀因素纳入研究。
（3）不同影响因素的系数在低位房价、中等房价、高位房价的空间分布
最后来看不同自变量系数函数在不同分位数房价的空间分布。图7,图8,图9,图10,图11,图12分别描述了每一年，不同分位数下的面积、楼龄、容积率、绿化率、商场距离、三甲医院距离、学校距离和地铁站数量的各系数估计 值的空间分布，很明显各个系数在不同水平房价的空间分布均存在差异，这也说明不同影响因素对不同水平房价的影响是不同的，体现了分位数回归可以探索y的全局分布的优点。
图7 2017年和2021年不同分位数下的面积( X 1)系数估计的空间分布

Fig. 7 Spatial distribution of estimated area ( X 1) coefficients at different quantiles for 2017and 2021

图8 2017年和2021年不同分位数下的楼龄( X 2)系数估计的空间分布

Fig. 8 Spatial distribution of estimated coefficients for X 2at different quantiles for 2017 and 2021

图9 2017年和2021年不同分位数下的容积率( X 3)系数估计的空间分布

Fig. 9 Spatial distribution of estimated coefficients for X 3at different quantiles for 2017 and 2021

图10 2017年和2021年不同分位数下的绿化率( X 4)系数估计的空间分布

Fig. 10 Spatial distribution of estimated coefficients for X 4at different quantiles for 2017 and 2021

图11 2021年不同分位数下的商场距离( X 6)、三甲医院距离( X 7)和学校距离( X 8 )系数估计的空间分布

Fig. 11 Spatial distribution of estimated coefficients for   X 6   , X 7   , X 8at different quantiles for 2021

图12 2017、2019和2021年不同分位数下的地铁站数量( X 8)系数估计的空间分布

Fig. 12 Spatial distribution of estimated coefficients for   X 8  at different quantiles for 2017, 2019 and 2021

图7描述了2017年和2021年不同分位数下面积系数估计 β ^ 1的空间分布。很明显，面积对不同水平房价均具有正向影响，且系数估计随着分位数增大而增大，意味着面积对高位房价的影响大于低位房价。从面积系数估计的空间分布看，不同分位数的系数估计的绝对值在郊区如浦东新区、青浦区、嘉定区、松江区、宝山区、闵行区均比较小，在市中心如黄浦、普陀等区较大，说明面积对市中心区域房价的影响大于郊区。
从图8描述的楼龄系数估计 β ^ 2的空间分布可以看出，不同分位数的楼龄系数估计均为负数，即楼龄越高房价越低。注意到，楼龄系数估计 β ^ 2的空间分布在不同分位数水平有明显不同，楼龄对低位房价( τ = 0.25)的影响从西北(青浦、嘉定)到东南增大(浦东新区)；对于中等房价( τ = 0.5)的影响从北(嘉定、宝山)向南增大(浦东新区、闵行区)；对于高位房价( τ = 0.75)的影响从西南(松江、青浦)到东北(宝山、杨浦)增大。
从图9的容积率系数的空间分布看出，整体上，容积率系数在市中心和松江区为正数，在其他区基本为负数。此外，容积率系数的空间分布在不同分位数水平存在明显的不同，这种差异在浦东新区南部尤为明显。在低位房价( τ = 0.25)，容积率系数在浦东新区南部为正数，而在中等房价和高位房价为负数。
图10展示了绿化率系数的空间分布在不同分位数水平也存在差异。在2021年，绿化率系数为负值的区域明显减小，说明人们越来越重视绿化率。整体上来说，绿化率对于市中心地区和浦东新区房价具有较强的正向影响。另一方面，绿化率对高位房价的影响大于低位房价。
图11(a)—图11(i)分别描述了2021年不同分位数下商场距离、三甲医院距离和学校距离系数估计的空间分布。很明显，相比自变量面积、楼龄、容积率、绿化率，这3个变量的系数估计的绝对值较小且基本为负值，其空间分布在不同水平房价差异较小，说明商场距离、三甲医院距离和学校距离对于房价有较小的负向影响。整体上，商场距离对于西北郊区如嘉定、宝山房价的影响较大，三甲医院距离对东北区域如杨浦、虹口的影响较大，而学校距离对于东南区如浦东新区及市中心区域的房价的影响较大。
图12描述了2017、2019和2021年每一年不同分位数下地铁站数量系数估计 β ^ 9的空间分布。相比地铁距离，地铁站数量对房价的影响较大,故省略地铁距离系数估计 β ^ 5的空间分布图。整体来说，地铁站数量的系数估计在除浦东新区之外的地区均为正数，说明地铁站数量对于市中心地区及西部郊区的房价有较大的正向影响，且对西部郊区（如青浦、松江、嘉定）房价的影响比市中心更大。
本节对上海市住宅价格进行了建模分析，分别对高位房价、中等房价、低位房价建立了基于多带宽局部多项式的时空地理加权分位数回归模型。从结果来看，基于多带宽的模型拟合结果更好，且不同影响因素的最优时间带宽和空间带宽存在差异，说明不同影响因素的时间异质性与空间异质性水平不同，这表明了考虑多带宽的必要性。为了说明分位数回归模型的稳健性，删除原始数据1%的极端值后分别利用本文的分位数回归模型和基于最小二乘回归的MGTWR拟合，结果显示分位数回归模型的平均绝对误差变化较小（相比原始数据的MAE下降比例小），说明分位数回归对于极端值或者异常值比最小二乘回归更具稳健性。此外，高位房价、中等房价、低位房价的影响因素存在差异，同一个影响因素对不同水平房价如高位房价、中等房价、低位房价的影响效应也明显不同，具体体现在其时间分布和空间分布存在差异，这说明了分位数回归模型的必要性，体现了分位数回归可以用来研究影响房价分布的因素，而最小二乘回归模型只能探索均值房价的影响因素，无法研究高位房价、低位房价等不同水平价格的影响因素。

5 结论与讨论

5.1 结论

首先，本文基于局部多项式估计和分位数回归，提出了基于多带宽局部多项式的时空地理加权分位数回归模型，且允许不同自变量拥有不同的最优带宽。局部多项式估计相比核估计边界估计更准确，且多带宽可以处理不同自变量的时间异质性和空间异质性。基于加权最小二乘估计的时空地理加权回归方法对异常值不稳健，而分位数回归受异常值影响较小，相比最小二乘回归更为稳健且应用条件相对更为宽松。此外，分位数回归最重要的优点在于可以探索解释变量与响应变量分布的影响关系(如响应变量的多个分位数)，可以挖掘到更为丰富的信息，而最小二乘回归模型只能研究解释变量与响应变量条件均值的关系。本文基于局部多项式估计，利用两步迭代估计法给出模型的系数估计。
其次，本文通过数值模拟，将本文的分位数回归模型与时空地理加权最小二乘回归进行对比，结果显示，基于分位数回归的系数估计的均方误差和平均绝对误差均比最小二乘估计量小，比如，在0.75分位数，基于最小二乘回归得到的系数估计的均方误差和平均绝对误差分别是基于分位数回归的10倍和4倍，说明了分位数回归具有稳健性且可以研究影响响应变量分布的因素。
最后，本文以上海市2017—2021年商品房住宅小区为案例对象，应用该模型，探究不同影响因素对不同分位数的住宅价格(如高位房价、中等房价、低位房价)的影响。研究表明影响因素对不同水平房价的影响效果不同，具体体现在同一个影响因素的时间分布和空间分布在高位房价、中等房价、低位房价存在明显差异，这说明了分位数回归的优点和必要性，因为基于最小二乘回归的时空地理加权回归只能研究房价的平均水平，无法对低位房价或高位房价的影响因素进行探究。另一方面，上海市住宅价格数据建模分析也说明不同影响因素的最优带宽存在差异，体现了多带宽的必要性。此外，与基于最小二乘回归的MGTWR相比，本文的分位数回归模型对于异常值的存在更为稳健（删除1%极端值后基于分位数回归模型拟合的平均绝对误差的变化比基于最小二乘回归模型小1%）。对不同价格水平的房价进行建模分析可以帮助我们进一步理解影响房价（尤其是高位房价和低位房价）分布的因素及其影响机制。

5.2 讨论

尽管分位数回归相比最小二乘回归更为稳健且应用条件相对更为宽松，但在数据满足高斯-马尔可夫条件时，最小二乘估计是最优估计，若此时仅仅想研究响应变量y的平均水平，不关注其全局分布（不同分位数），那么应该使用基于最小二乘估计的MGTWR。如果想探索响应变量不同分位数水平的影响因素，则需要利用本文的多带宽局部多项式GTWR分位数回归模型。此外，Lu等[21,22]考虑了距离度量-参数对应的GWR模型，为不同变量选择不同距离度量和空间带宽，在之后的研究中也可以考虑在本文多带宽GTWR分位数回归模型中加入不同距离度量，以便更好地处理不同分位数下不同变量的时空相关性。

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
Brunsdon C, Fotheringham A S, Charlton M E. Geographically weighted regression: A method for exploring spatial nonstationarity[J]. Geographical Analysis, 1996, 28(4):281-298. DOI:10.1111/j.1538-4632.1996.tb00936.x

[2]
叶健, 胡鑫, 徐鸿蒙, 等. 多尺度GTWR城市住宅价格建模与分析[J]. 测绘学报, 2021, 50(9):1266-1274.

DOI

[ Ye J, Hu X, Xu H M, et al. Modeling and analysis of urban housing price models based on multiscale geographically and temporally weighted regression[J]. Acta Geodaetica et Cartographica Sinica, 2021, 50(9):1266-1274. ]. DOI:10.13203/j.whugis20190346

[3]
Huang B, Wu B, Barry M. Geographically and temporally weighted regression for modeling spatio-temporal variation in house prices[J]. International Journal of Geographical Information Science, 2010, 24(3):383-401. DOI:10.1080/13658810802672469

[4]
Fotheringham A S, Crespo R, Yao J. Geographical and temporal weighted regression (GTWR)[J]. Geographical Analysis, 2015, 47(4):431-452. DOI:10.1111/gean.12071

[5]
Wu B, Li R R, Huang B. A geographically and temporally weighted autoregressive model with application to housing prices[J]. International Journal of Geographical Information Science, 2014, 28(5):1186-1204. DOI:10.1080/13658816.2013.878463

[6]
赵阳阳, 张小璐, 张福浩, 等. 一种局部多项式时空地理加权回归方法[J]. 测绘学报, 2018, 47(5):663-671.

DOI

[ Zhao Y Y, Zhang X L, Zhang F H, et al. A local polynomial geographically and temporally weight regression[J]. Acta Geodaetica et Cartographica Sinica, 2018, 47(5): 663-671. ]. DOI:10.11947/j.AGCS.2015.20150005

[7]
Fan J Q, Gijbels I. Local polynomial modelling and its applications[M]. London: Chapman & Hall, 1996

[8]
Li Q, Racine J S. Nonparametric econometrics: theory and practice[M]. Princeton University Press, 2010.

[9]
Wang N, Mei C L, Yan X D. Local linear estimation of spatially varying coefficient models: An improvement on the geographically weighted regression technique[J]. Environment and Planning A: Economy and Space, 2008, 40(4):986-1005. DOI:10.1068/a3941

[10]
Mei C L, Wang N, Zhang W X. Testing the importance of the explanatory variables in a mixed geographically weighted regression model[J]. Environment and Planning A: Economy and Space, 2006, 38(3):587-598. DOI: 10.1068/a3768

[11]
卢宾宾, 葛咏, 秦昆, 等. 地理加权回归分析技术综述[J]. 武汉大学学报·信息科学版, 2020, 45(9):1356-1366.

[ Lu B B, Ge Y, Qin K, et al. A review on geographically weighted regression[J]. Geomatics and Information Science of Wuhan University, 2020, 45(9):1356-1366. ] DOI:10.13203/j.whugis20190346

[12]
Fotheringham A S, Yang W B, Kang W. Multiscale geographically weighted regression (MGWR)[J]. Annals of the American Association of Geographers, 2017, 107(6):1247-1265. DOI:10.1080/24694452.2017.1352480

[13]
Wu C, Ren F, Hu W, et al. Multiscale geographically and temporally weighted regression: Exploring the spatiotemporal determinants of housing prices[J]. International Journal of Geographical Information Science, 2019, 33(3):489-511. DOI:10.1080/13658816.2018.1545158

[14]
Zhang H G, Mei C L. Local least absolute deviation estimation of spatially varying coefficient models: Robust geographically weighted regression approaches[J]. International Journal of Geographical Information Science, 2011, 25(9):1467-1489. DOI:10.1080/13658816.2010.528420

[15]
Koenker R, Bassett G. Regression quantiles[J]. Econometrica, 1978, 46(1):33-50. DOI: 10.2307/1913643.

[16]
田茂再. 分位回归与复杂分层结构数据分析[M]. 2015.

[ Tian M Z. Quantile regression and complex hierarchical data analysis[M]. 2015.]

[17]
Hastie T, Tibshirani R. Varying-coefficient models[J]. Journal of the Royal Statistical Society: Series B (Methodological), 1993, 55(4):757-779. DOI:10.1111/j.2517-6161.1993.tb01939.x.

[18]
杨毅. 顾及时空非平稳性的地理加权回归方法研究[D]. 武汉: 武汉大学, 2016. [ Yang Y. Research on geographically weighted regression method considering spatio-temporal nonstationarity[D]. Wuhan: Wuhan University, 2016. ]

[19]
Horowitz J L. Bootstrap methods for Median regression models[J]. Econometrica, 1998, 66(6):1327-1351. DOI: 10.2307/2999619

[20]
Brunsdon C, Fotheringham A S, Charlton M. Some notes on parametric significance tests for geographically weighted regression[J]. Journal of Regional Science, 1999, 39(3):497-524. DOI:10.1111/0022-4146.00146

[21]
Lu B B, Brunsdon C, Charlton M, et al. Geographically weighted regression with parameter-specific distance metrics[J]. International Journal of Geographical Information Science, 2017, 31(5):982-998. DOI:10.1080/13658816.2016.1263731

[22]
Lu B B, Yang W B, Ge Y, et al. Improvements to the calibration of a geographically weighted regression with parameter-specific distance metrics and bandwidths[J]. Computers, Environment and Urban Systems, 2018, 71:41-57. DOI:10.1016/j.compenvurbsys.2018.03.012

[23]
赵阳阳, 刘纪平, 杨毅, 等. 混合时空地理加权回归及参数的两步估计[J]. 计算机科学, 2017, 44(3):274-277,312.

DOI

[ Zhao Y Y, Liu J P, Yang Y, et al. Mixed geographically and temporally weighted regression and two-step estimation[J]. Computer Science, 2017, 44(3):274-277,312.]. DOI:10.11896/j.issn.1002-137X.2017.03.056

[24]
汤庆园, 徐伟, 艾福利. 基于地理加权回归的上海市房价空间分异及其影响因子研究[J]. 经济地理, 2012, 32(2):52-58.

[ Tang Q Y, Xu W, Ai F L. A GWR-based study on spatial pattern and structural determinants of Shanghai' s housing price[J]. Economic Geography, 2012, 32(2):52-58.]. DOI:10.3969/j.issn.1000-7636.2010.06.005

[25]
温海珍, 李旭宁, 张凌. 城市景观对住宅价格的影响——以杭州市为例[J]. 地理研究, 2012, 31(10):1806-1814.

[ Wen H Z, Li X N, Zhang L. Impacts of the urban landscape on the housing price: A case study in Hangzhou[J]. Geographical Research, 2012, 31(10):1806-1814. ] DOI: 10.3321/j.issn:1000-0933.2009.08.057

[26]
石忆邵, 张蕊. 大型公园绿地对住宅价格的时空影响效应——以上海市黄兴公园绿地为例[J]. 地理研究, 2010, 29(3):510-520.

[ Shi Y S, Zhang R. Temporal-spatial impact effects of large-scale parks on residential prices: Exemplified by the Huangxing Park in Shanghai[J]. Geographical Research, 2010, 29(3):510-520. ] DOI:10.3321/j.issn:0375-5444.2009.02.004

[27]
顾杰, 贾生华. 公共交通改善期望对住房价格及其价格空间结构的影响——基于杭州地铁规划的实证研究[J]. 经济地理, 2008, 28(6):1020-1024,1034.

[ Gu J, Jia S H. The effects of expected transport improvements on housing prices and price spatial distribution—Hangzhou-based research evidence of planning the mass transit railway[J]. Economic Geography, 2008, 28(6):1020-1024,1034. ] DOI:10.3321/j.issn: 1000-131X.2007.04.017

[28]
Liu J P, Yang Y, Xu S H, et al. A geographically temporal weighted regression approach with travel distance for house price estimation[J]. Entropy, 2016, 18(8):303. DOI: 10.3390/e18080303

Options
文章导航

/

〈 〉