High-Precision Urban Boundary Extraction and Urban Sprawl Spatial-Temporal Analysis in China's Prefectural Cities from 2000 to 2016
-
摘要: 城区边界和城区面积是城镇化的重要表征和扩展分析的基础。然而城区边界存在概念和提取标准不统一、精度较低、可比性较差等问题。为此，提出了基于高分辨率影像和地理信息资料辅助的城区半自动化提取方法，充分利用高分辨率影像上的城市景观特征、先验地理信息知识和一系列标准规则，以得到精度高、一致性强的数据。以中国337个地级以上城市为研究区，采用该方法得到了2000年、2005年、2010年、2016年4期城区边界成果，并开展了城区时空扩展及用地效率等相关分析。结果表明：①16年间城区扩展迅速，城区主要集中分布在东部和中部，东西部地区差异大；②城市用地效率与城镇化发展水平显著相关，城区扩展以外延型为主；③大多城市城区扩展超前于人口增长，少量城市城区扩展滞后于人口增长；④城区扩展以占用耕地为主。Abstract: Urban boundary and area is an important indicator for urban sprawl analysis. However, concept and extraction standard of urban boundary were inconsistent, and the resolution of images was coarse, resulting in lower precision and incomparable urban boundary. To solve above problems, a semi-automatic method of urban boundary extraction is proposed by using high-resolution image and geographic information data. Urban landscape and form characteristics, geographical knowledge, and a series of standardized rules are combined to generate high-precision and consistent urban boundary. Urban boundaries of China's 337 prefectural cities in year 2000, 2005, 2010 and 2016 are extracted with above-mentioned method. Together with other thematic data, urban sprawl is analyzed. Results show that: ① China's prefectural cities has undergone a rapid urbanization from year 2000 to 2015 with a high degree of concentration and regional differences. ②The efficiency of urban land use is significantly related to the development level of urbanization. The sprawl models of most cities are extensive. ③Most cities are faced with the problem that urban expansion is ahead of population growth, while a few cities expansion lag behind population growth. ④Cultivated land is the main source of urban expansion areas.
-
Keywords:
- urban boundary /
- urban sprawl /
- high-resolution images /
- China's prefectural cities
-
大地测量学科发展的总趋势是向地球科学纵深发展，深入到其他地学学科的交叉领域，其主要任务是监测和研究地球动力学现象，研究地球本体的各种物理场，认识与探索地球内部的各种物理过程并揭示其规律。以空间大地测量为标志的现代大地测量技术，不论在测量的空间尺度上还是已达到的精度水平方面，已经有能力监测地球动力学过程产生的运动状态和物理场的微变化[1]。利用GPS监测点复测资料所得的位移或速度研究区域地壳运动与应变是大地测量深入地球科学研究最常用的方法之一[2-7]，在数据处理过程中基本都采用最小二乘方法(least squares, LS)，没有考虑到反演模型中系数矩阵误差的影响。文献[8]采用总体最小二乘方法研究地壳应变参数的反演，考虑了系数矩阵的误差影响以及系数矩阵中含有非随机元素的情况，但是在公式推导时对公式进行了线性化处理，引入了模型误差，因此得到的参数结果存在偏差。目前关于地壳应变参数反演模型中系数矩阵误差对模型参数的具体影响情况，还没有明确的定论。
部分变量误差(partial-EIV，PEIV)模型作为变量误差(errors-in-variables，EIV)模型一般化的表达形式，既考虑了函数模型中观测向量的误差，又顾及了系数矩阵含有误差的情况，同时对系数矩阵的结构及组成无特殊限制，更适合用于EIV模型问题的求解。目前，国内外研究并使用PEIV模型的文献并不多。文献[9]在PEIV模型的基础上，推导出了加权总体最小二乘(weighted total least squares, WTLS)更一般化的解，该算法涵盖了一般WTLS算法需要特殊处理的各种情况，同时大大减少了模型的待估量，算法的简单形式使后续估计值的精度评定更加便利；文献[10]提出了基于PEIV模型的可靠性理论，并分析了模型中观测向量和系数矩阵可靠性的基本特点。本文以PEIV模型为基础，采用加权总体最小二乘算法进行地壳应变参数反演，进一步揭示GPS站点坐标误差对地壳应变模型参数估计的影响。同时本文算法可以很好地解决反演模型的系数矩阵中存在的随机元素、非随机元素和重复元素的情况，算法的推导公式中也避免了文献[8]中存在的问题，结合推导得到的最小二乘改正项公式，对地壳反演模型中坐标点误差对反演参数求解的影响进行了分析。通过模拟数据和川滇地区的实际数据，分析了系数矩阵误差对地壳应变参数反演的影响，展示了更多反演过程的细节，对深入理解地壳应变反演机理有一定的价值。
1 由坐标位移反演地壳应变参数的PEIV模型算法
在假定介质均匀且邻近测点间应变均匀的情况下，可以建立邻近点间相对形变量与地壳应变量的线性关系。设第j点两期纵坐标位移uj，横坐标位移vj，网中共有m点，由坐标位移反演应变参数的公式(二维)为[11]：
$$ \left\{ \begin{array}{l} {u_j} = u + {x_j}{\varepsilon _x} + {y_j}{\varepsilon _{xy}} - {y_j}w\\ {v_j} = v + {x_j}{\varepsilon _{xy}} + {y_j}{\varepsilon _y} - {x_j}w \end{array} \right. $$ (1) 式中，xj、yj(j=1, 2, …, m)为第j点近似坐标；u、v为网的平移量；w为其转动量；εx、εy、εxy为应变状态参数。
将式(1) 改写成矩阵形式：
$$ \left[{\begin{array}{*{20}{c}} {{u_1}}\\ {{v_1}}\\ \vdots \\ {{u_m}}\\ {{v_m}} \end{array}} \right] = \left[{\begin{array}{*{20}{c}} 1&0&{{x_1}}&0&{{y_1}}&{-{y_1}}\\ 0&1&0&{{y_1}}&{{x_1}}&{{x_1}}\\ {}&{}& \vdots &{}& \vdots &{}\\ 1&0&{{x_m}}&0&{{y_m}}&{-{y_m}}\\ 0&1&0&{{y_m}}&{{x_m}}&{{x_m}} \end{array}} \right]\left[{\begin{array}{*{20}{c}} u\\ v\\ {{\varepsilon _x}}\\ {{\varepsilon _y}}\\ {{\varepsilon _{xy}}}\\ w \end{array}} \right] $$ (2) 令y=[u1 v1 … um vm]T，y∈R2m×1；β= [u v εx εy εxy w]T, β∈R6×1; A∈R2m×6
$$ \mathit{\boldsymbol{A = }}\left[{\begin{array}{*{20}{c}} 1&0&{{x_1}}&0&{{y_1}}&{-{y_1}}\\ 0&1&0&{{y_1}}&{{x_1}}&{{x_1}}\\ {}&{}& \vdots &{}& \vdots &{}\\ 1&0&{{x_m}}&0&{{y_m}}&{-{y_m}}\\ 0&1&0&{{y_m}}&{{x_m}}&{{x_m}} \end{array}} \right] $$ 所以，式(2) 可以表示为：
$$ \mathit{\boldsymbol{y}} = \mathit{\boldsymbol{A\beta }} $$ (3) 基于PEIV模型的坐标位移反演应变参数的思想是：不仅考虑坐标位移的误差，即y的误差对应变参数的影响，以及系数矩阵A误差的影响，同时考虑系数矩阵的元素并非全都是随机的，而部分是由0或者1等常数组成的，所以有[9]:
$$ \left\{ \begin{array}{l} \mathit{\boldsymbol{y}} = \left( {{\mathit{\boldsymbol{\beta }}^T} \otimes {\mathit{\boldsymbol{I}}_{2m}}} \right)\left( {\mathit{\boldsymbol{h}} + \mathit{\boldsymbol{B\bar a}}} \right) + \mathit{\boldsymbol{\varepsilon }}\\ \mathit{\boldsymbol{a}} = \mathit{\boldsymbol{\bar a}} + {\mathit{\boldsymbol{\varepsilon }}_a}\\ {\rm{vec}}\left( {\mathit{\boldsymbol{\bar A}}} \right) = \mathit{\boldsymbol{h}} + \mathit{\boldsymbol{B\bar a}} \end{array} \right. $$ (4) 式中，y∈R2m×1表示观测向量；ε∈R2m×1表示观测向量y的随机误差；h∈R12m×1主要由系数矩阵A中非随机元素组成；B∈R12m×t是固定矩阵，其阶数取决于系数矩阵A中随机元素的数目(此处假设为t)；I2m ∈R2m×2m为单位阵；A∈Rt×1是系数矩阵A中随机元素组成的列向量(若随机元素重复出现，则只提取一次)，其真值用a表示；εa∈Rt×1是A中包含的随机误差；A·β=(βT⊗In)vec(A)，其中A表示系数矩阵A的真值，vec(·)表示矩阵拉直运算。
地壳应变参数反演的PEIV模型中各矩阵的具体形式表示为：
$$ \begin{array}{l} \mathit{\boldsymbol{h}} = {\left[{\mathit{\boldsymbol{h}}_1^{\rm{T}}\;\;\mathit{\boldsymbol{h}}_2^{\rm{T}}\;\;\mathit{\boldsymbol{h}}_3^{\rm{T}}\;\;\mathit{\boldsymbol{h}}_4^{\rm{T}}\;\;\mathit{\boldsymbol{h}}_5^{\rm{T}}\;\;\mathit{\boldsymbol{h}}_6^{\rm{T}}} \right]^{\rm{T}}}, \\ \mathit{\boldsymbol{B}} = {\left[{\mathit{\boldsymbol{B}}_1^{\rm{T}}\;\;\mathit{\boldsymbol{B}}_2^{\rm{T}}\;\;\mathit{\boldsymbol{B}}_3^{\rm{T}}\;\;\mathit{\boldsymbol{B}}_4^{\rm{T}}\;\;\mathit{\boldsymbol{B}}_5^{\rm{T}}\;\;\mathit{\boldsymbol{B}}_6^{\rm{T}}} \right]^{\rm{T}}}, \\ \;\;\;\;\;\;\;\;\;\;\;{\mathit{\boldsymbol{h}}_1} = {\left[{1\;\;0\;\; \ldots \;\;1\;\;0} \right]^{\rm{T}}}, \\ \;\;\;\;\;\;\;\;\;\;\;{\mathit{\boldsymbol{h}}_2} = {\left[{0\;\;1\;\; \ldots \;\;0\;\;1} \right]^{\rm{T}}}, \end{array} $$ $$ \begin{array}{l} {\mathit{\boldsymbol{h}}_3} = {\mathit{\boldsymbol{h}}_4} = {\mathit{\boldsymbol{h}}_5} = {\mathit{\boldsymbol{h}}_6} = {\left[{0\;\;0\;\; \ldots \;\;0\;\;0} \right]^{\rm{T}}}, \\ \;\;\;\;\;\;\;\;\;\;\;{\mathit{\boldsymbol{h}}_i} \in {R^{2m \times 1}}, \left( {i = 1, 2 \ldots 5, 6} \right), \end{array} $$ $$ \begin{array}{l} {\mathit{\boldsymbol{B}}_1} = {\mathit{\boldsymbol{I}}_\mathit{\boldsymbol{m}}} \otimes \left[{\begin{array}{*{20}{c}} 0&0\ 0&0 \end{array}} \right], {\mathit{\boldsymbol{B}}_2} = {\mathit{\boldsymbol{I}}_\mathit{\boldsymbol{m}}} \otimes \left[{\begin{array}{*{20}{c}} 0&1\ 0&0 \end{array}} \right], \\ {\mathit{\boldsymbol{B}}_3} = {\mathit{\boldsymbol{I}}_\mathit{\boldsymbol{m}}} \otimes \left[{\begin{array}{*{20}{c}} 1&0\ 0&0 \end{array}} \right], {\mathit{\boldsymbol{B}}_4} = {\mathit{\boldsymbol{I}}_\mathit{\boldsymbol{m}}} \otimes \left[{\begin{array}{*{20}{c}} 0&0\ 0&1 \end{array}} \right], \\ {\mathit{\boldsymbol{B}}_5} = {\mathit{\boldsymbol{I}}_\mathit{\boldsymbol{m}}} \otimes \left[{\begin{array}{*{20}{c}} 0&1\ 1&0 \end{array}} \right], {\mathit{\boldsymbol{B}}_6} = {\mathit{\boldsymbol{I}}_\mathit{\boldsymbol{m}}} \otimes \left[{\begin{array}{*{20}{c}} 0&{-1}\ 1&0 \end{array}} \right], \end{array} $$ $$ \begin{array}{l} \;\;\;\;\;\;\;\;\;{\mathit{\boldsymbol{B}}_i} \in {R^{2m \times 2m}}, \left( {i = 1, 2, \cdots, 5, 6} \right), \\ \mathit{\boldsymbol{a}} = {\left[{{x_1}\;\;{y_1}\;\; \cdots \;\;{x_m}\;\;{y_m}} \right]^{\rm{T}}}, \mathit{\boldsymbol{a}} \in {R^{2m \times 1}}; \end{array} $$ 其中⊗表示直积运算。
平差的随机模型为：
$$ \left\{ \begin{array}{l} \mathit{\boldsymbol{D}}\left( \varepsilon \right) = {\sigma ^2}{\mathit{\boldsymbol{Q}}_{\mathit{\boldsymbol{\varepsilon \varepsilon }}}}\\ \mathit{\boldsymbol{D}}\left( {{\varepsilon _a}} \right) = {\sigma ^2}{\mathit{\boldsymbol{Q}}_{{\varepsilon _a}{\varepsilon _a}}} \end{array} \right. $$ (5) 式中，D表示方差-协方差矩阵；σ2为单位权方差；Qεε=P1-1为观测值的协因数阵；Qεaεa=P2-1为系数矩阵中随机元素的协因数阵。
通过观测值y和系数矩阵随机元素A估计PEIV模型中a和β的加权总体最小二乘平差准则为[9]：
$$ \varphi = {\mathit{\boldsymbol{\varepsilon }}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{\varepsilon }} + \mathit{\boldsymbol{\varepsilon }}_\mathit{\boldsymbol{a}}^{\rm{T}}{\mathit{\boldsymbol{P}}_2}{\mathit{\boldsymbol{\varepsilon }}_\mathit{\boldsymbol{a}}} = \min $$ (6) 当初始权阵P1和P2已知时，可以采用下面的步骤进行未知参数求解[9]
第1步 初始值设定。给定观测向量y和系数矩阵A，提取矩阵A中随机元素组成列向量a，根据系数矩阵结构，构造常数向量h和固定矩阵B，计算参数初始值$ {{\mathit{\boldsymbol{\hat \beta }}}^{\left( 0 \right)}} = {\left( {{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}A} \right)^{ - 1}}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{y}} $；
第2步 迭代过程，利用公式$ \widehat {\mathit{\boldsymbol{\overline a}} }\mathit{\boldsymbol{ = }}{\left( {\mathit{\boldsymbol{S}}_\beta ^{\rm{T}}{\mathit{\boldsymbol{P}}_1}{\mathit{\boldsymbol{S}}_\beta } + {\mathit{\boldsymbol{P}}_2}} \right)^{ - 1}}\left( {{\mathit{\boldsymbol{P}}_2}\mathit{\boldsymbol{a}} - \mathit{\boldsymbol{S}}_\beta ^{\rm{T}}{\mathit{\boldsymbol{P}}_1}\left( {\sum\limits_{i = 1}^m {{\mathit{\boldsymbol{h}}_\mathit{\boldsymbol{i}}}{{\mathit{\boldsymbol{\hat \beta }}}_\mathit{\boldsymbol{i}}}} } \right) + \mathit{\boldsymbol{S}}_\beta ^{\rm{T}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{y}}} \right) $计算$ \widehat {\mathit{\boldsymbol{\overline a}} } $，式中$ {\mathit{\boldsymbol{S}}_\beta } = \sum\limits_{i = 1}^m {{\mathit{\boldsymbol{B}}_\mathit{i}}{{\mathit{\hat \beta }}_\mathit{i}}} $；利用公式$ \mathit{\boldsymbol{\hat \beta }} = {\left( {{{\mathit{\boldsymbol{\hat A}}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{\hat A}}} \right)^{ - 1}}{{\mathit{\boldsymbol{\hat A}}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{y}} $计算$ {{\mathit{\boldsymbol{\hat \beta }}}^{\left( {i + 1} \right)}} $，式中$ \mathit{\boldsymbol{\hat A}}\mathit{ = }{\rm{ve}}{{\rm{c}}^{ - 1}}\left( {\mathit{\boldsymbol{h}} + \mathit{\boldsymbol{B}}\widehat {\mathit{\boldsymbol{\overline a}} }} \right) $，其中vec-1(·)表示拉直运算的逆运算。
第3步 迭代终止。$ \left\| {{{\mathit{\boldsymbol{\hat \beta }}}^{\left( {i + 1} \right)}} - {{\mathit{\boldsymbol{\hat \beta }}}^{\left( i \right)}}} \right\| \le {\delta _0} $，其中δ0表示设定的阈值。
2 系数矩阵误差对WTLS估计值的影响
对基于PEIV模型的地壳应变参数加权总体最小二乘解的表达式进行相应调整可以得到[10]：
$$ \begin{array}{l} {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{WTLS}}}} = {\left( {{{\left( {\mathit{\boldsymbol{A}} + {{\mathit{\boldsymbol{\hat E}}}_A}} \right)}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\left( {\mathit{\boldsymbol{A}} + {{\mathit{\boldsymbol{\hat E}}}_A}} \right)} \right)^{ - 1}}\\ \;\;\;\;\;\;\;\;\;\;\;\;{\left( {\mathit{\boldsymbol{A}} + {{\mathit{\boldsymbol{\hat E}}}_A}} \right)^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\left( {\mathit{\boldsymbol{A}}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + {{\mathit{\boldsymbol{\hat e}}}_y}} \right) \end{array} $$ (7) 式中，$ {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{WTLS}}}} $表示加权总体最小二乘估值；$ {{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}}} $表示加权最小二乘估值；$ {{{\mathit{\boldsymbol{\hat E}}}_A}} $表示通过WTLS估计得到的系数矩阵误差的残差；$ {{{\mathit{\boldsymbol{\hat e}}}_y}} $表示通过WLS估计得到的观测向量的残差。
将式(7) 展开可以得到：
$$ \begin{array}{l} {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{WTLS}}}} = {\left( {{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}} + {\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat E}}}_A} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat E}}}_A}} \right)^{ - 1}}\left( {{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}} \right.{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left. {{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y}} \right) \end{array} $$ (8) 令$ {\mathit{\boldsymbol{N}}_1} = {\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}, {\mathit{\boldsymbol{N}}_2} = {\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat E}}}_A} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\hat E}}_A} $，根据分块矩阵求逆公式，式(8) 可以改写为：
$$ \begin{array}{l} {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{WTLS}}}} = \left( {\mathit{\boldsymbol{N}}_1^{ - 1} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{{\left( {{I_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)}^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}} \right)\\ \;\;\;\;\;\;\;\;\;\;\;\;\left( {{\mathit{\boldsymbol{N}}_1}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + {\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\hat e}}_y} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + {{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y}} \right) = \\ \;\;\;\;\;\;\;\;\;\;\;\;{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{P}_1}{{\mathit{\boldsymbol{\hat e}}}_y} + \mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + \mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\hat E}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y}\\ \;\;\;\;\;\;\;\;\;\;\;\; - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\\ \;\;\;\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} \end{array} $$ 对上式进行移项调整可以得到：
$$ \begin{array}{l} \mathit{\boldsymbol{b}} = {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{WTLS}}}} - {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} = \left( {\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)}^{ - 1}}} \right. - \\ \;\;\;\;\;\;\;\left. {\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)}^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}} \right){{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} + \\ \;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} - \\ \;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} \end{array} $$ (9) 令$ \begin{array}{l} {\mathit{\boldsymbol{H}}_1} = \mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}} - \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\\ \;\;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}\mathit{\boldsymbol{A}}, \;\;{\mathit{\boldsymbol{H}}_2} = \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} + \mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} - \\ \;\;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{A}}^{\rm{T}}}{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} - \\ \;\;\;\;\;\;\;\;\mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}{\left( {{\mathit{\boldsymbol{I}}_6} + \mathit{\boldsymbol{N}}_1^{ - 1}{\mathit{\boldsymbol{N}}_2}} \right)^{ - 1}}\mathit{\boldsymbol{N}}_1^{ - 1}{{\mathit{\boldsymbol{\hat E}}}^{\rm{T}}}_A{\mathit{\boldsymbol{P}}_1}{{\mathit{\boldsymbol{\hat e}}}_y} \end{array} $。
式(9) 可以改写为：
$$ \mathit{\boldsymbol{b}} = {\mathit{\boldsymbol{H}}_1}{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{LS}}}} + {\mathit{\boldsymbol{H}}_2} $$ (10) 式中，b表示WTLS估值和WLS估值之间的改正项，由于采用最小二乘求解地壳应变参数反演模型得到的是一个近似解，在求解过程中未考虑函数模型中系数矩阵误差的影响，且系数矩阵误差对WLS估计值的影响未知，因此解算结果的优劣无法判断。采用总体最小二乘进行地壳应变参数的求解，得到的是一个严密解，此时在WTLS估值和WLS估值之间需要添加一个改正项。在式(10) 的基础上进行模拟算例和实际算例的定量分析，可以对最小二乘准则和总体最小二乘准则之间的差异性有一个初步的了解，同时可以客观研究系数矩阵误差对地壳应变参数WTLS解的影响。
3 算例分析
3.1 模拟实验
对3个方案(见表 1)进行模拟实验，以方案1为例，在边长为100 m2的正方形区域，沿x和y轴每隔10 m进行取样，共得到121个坐标点，如图 1。给定该区域应变参数εx=200×10-8，εy=-80×10-8，εxy=1 000×10-8，ω=1 100×10-8(模拟算例中未考虑平移量)，将坐标数据代入式(2) 得到各点对应的位移量真值。分别给坐标点和位移量加入服从正态分布的随机误差，此时函数模型式(4) 中，观测向量y由加入误差后的位移量组成，其中D(ε)=σ2Q，σ2=1 cm2，Q=diag(0.1, 0.1, …, 0.1)；而系数矩阵中的随机元素A由含有误差的坐标组成，且满足D(εa)=σ2Qεaεa，Qεaεa=diag(0.3, 0.3, …，0.3)。采用加权最小二乘(WLS)、基于PEIV模型的加权总体最小二乘法(PEIV-WTLS)两种方法求解表 1中3种方案的应变参数，同时计算相应的改正项b用于衡量WLS与PEIV-WTLS法解算结果之间的差距。表 2为方案1、方案2、方案3的应变参数计算结果。图 2给出了利用基于PEIV模型的加权总体最小二乘法(PEIV-WTLS)求解得到的各方案观测向量和系数矩阵随机元素残差值。在算例模拟过程中，为了使模拟的位移量符合实际情况，在方案2和方案3中，当模拟区域扩大时，相应的应变参数需要进行调整。
表 1 各方案的情景设置Table 1. The Scene of the Three Schemes方案1 方案2 方案3 区域 100 m×100 m 200 m×200 m 1 000 m×1 000 m 采样间隔 沿x、y轴每隔10 m 沿x、y轴每隔20 m 沿x、y轴每隔100 m 坐标点总数(采样) 121个 121个 121个 模拟应变参数真值 εx=200×10-8
εy=-80×10-8
εxy=1 000×10-8
ω=1 100×10-8εx=200×10-9
εy=-80×10-9
εxy=1 000×10-9
ω=1 100×10-9εx=200×10-10
εy=-80×10-10
εxy=1 000×10-10
ω=1 100×10-10![]()
图 1 方案1中地壳应变参数反演的模拟站点Figure 1. Simulated Network Stations of the Inversion of Crustal Strain in Case 1表 2 各方案参数计算结果Table 2. Results of Parameters in Each Caseεx εy εxy w 方案1 真值 200×10-8 -80×10-8 1 000×10-8 1 100×10-8 WLS 189×10-8 -86×10-8 1 006×10-8 1 097×10-8 PEIV-WTLS 189×10-8 -86×10-8 1 006×10-8 1 097×10-8 b -1.5×10-18 2.6×10-17 -1.5×10-17 -1.7×10-17 方案2 真值 200×10-9 -80×10-9 1000×10-9 1100×10-9 WLS 195×10-9 -81×10-9 1005×10-9 1100×10-9 PEIV-WTLS 195×10-9 -81×10-9 1005×10-9 1100×10-9 b -3.2×10-21 8.9×10-21 -4.8×10-21 -7.1×10-21 方案3 真值 200×10-10 -80×10-10 1 000×10-10 1 100×10-10 WLS 202×10-10 -83×10-10 999×10-10 1 103×10-10 PEIV-WTLS 202×10-10 -83×10-10 999×10-10 1 103×10-10 b -2.9×10-23 9.8×10-23 -1.0×10-22 -1.2×10-22 ![]()
图 2 各方案观测向量和系数矩阵随机元素残差值(PEIV-WTLS)Figure 2. Residuals of Observation Vector and Coefficient Matrix in Each Case(PEIV-WTLS)从表 2中可以看出，采用WLS和PEIV-WTLS这两种算法求解得到的应变参数估值基本相同，无法看出两者之间的区别，因此通过分析表格中的改正项b对这种现象进行探索。从式(10) 中可以看出，影响改正项b大小的主要因素来自地壳应变参数的加权最小二乘估值$ {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{WLS}}}} $，以及矩阵H1和H2中的元素(主要为系数矩阵的残差$ {{\mathit{\boldsymbol{\hat E}}}_A} $和WLS估计得到的观测向量的残差$ {{\mathit{\boldsymbol{\hat e}}}_y} $)。在地壳应变参数反演模型中，系数矩阵的随机元素是由观测站点的坐标组成的，站点坐标所含误差的量级一般为10-3~10-2，观测向量元素由位移量组成，残差$ {{\mathit{\boldsymbol{\hat e}}}_y} $的量级一般在10-3~10-2，应变参数估值的量级一般与观测站点坐标的量级呈负相关。从表 1中可以看出，随着模拟区域扩大，站点坐标值增大，应变参数随之减小。比较方案1和方案2可以看出，当模拟区域扩大4倍，即坐标点观测值量级扩展到102时，相应的改正项减小10-4个量级；比较方案2和方案3，当模拟区域扩大25倍，即坐标点观测值量级扩展到103时，相应的改正项减小10-2个量级，产生这一现象的主要原因为随着坐标点观测值量级的增大，求解得到的应变参数变小，从式(3) 中可以看出，系数矩阵中随机元素的误差是通过与应变参数的乘积形式对函数模型产生影响的，因此当应变参数减小的同时也会削弱系数矩阵中随机元素对模型解算的影响。也可以解释为，由于反演模型系数矩阵基本是由监测点坐标组成的，数值一般在106量级，坐标观测值中含有的误差被坐标值本身淹没，导致计算结果与WLS结果基本相同。在许多研究坐标转换的文献中也存在相同的问题[12]。从图 2可以看出，采用严密的PEIV-WTLS算法可以求解得到观测向量的残差值，却无法恢复系数矩阵中随机元素的残差值。区域性形变研究是一个长期的过程，任何微小的变化经过长时间的积累都会达到一个可观的量级，因此在监测地壳形变和运动时，再微小的变化也需要考虑在内，地壳应变参数反演模型中系数矩阵的误差随着观测区域的扩大会逐渐减小其对整个模型参数估计的贡献，却不能忽略其微小的影响，需要采用严密的方法将其考虑在内。
3.2 川滇地区地壳应变参数反演
在实际情况中，GPS监测站点的分布一般比较分散，且数量较少，站点质量存在差异。这也给区域性地壳形变研究带来了更多的不确定性。本文采用GPS速度场数据对川滇地区8个块体的应变率参数进行反演，以此来探索在实际情况中是否也存在上述3个模拟方案中相同的情况。川滇地区主要包括北纬20°~34°、东经96°~106°的中国大陆西南区域的青藏高原东南部，分布于其中的大量活动断层将该区域分成了许多块体，主要的断层有红河断层、小江断层、腾冲-景洪断层、鲜水河断层、安宁河断层、丽江断层、金沙江断层和龙门山断层等，这些断层所分割成的块体有马尔康块体(Ⅰ)、川西北块体(Ⅱ1)、滇中块体(Ⅱ2)、保山块体(Ⅲ1)、景谷块体(Ⅲ2)、勐腊块体(Ⅲ3)、羌塘块体(Ⅳ)和华南块体(Ⅴ)[13]，各块体与断层之间的空间关系见图 3。
![]()
本文所用的GPS水平速度场的数据来自于文献[15]，即1998~2004年中国地壳运动观测网络的GPS水平速度场及其方差-协方差矩阵数据，主要来源于中国地壳运动观测网络。在应变率参数反演过程中，首先将每个GPS速度场数据对应的测站的经、纬度坐标经过高斯投影转换为高斯平面直角坐标(区域坐标)，在坐标转换过程中，选取每个块体中心单独投影[16]。由于GPS测站的经纬度坐标的测定必然含有误差，而且在坐标转换中也会引入其他误差，而本文采用的地壳应变参数反演的PEIV模型求解方法可以顾及这些坐标误差。利用式(10) 计算各个块体的改正项b，具体计算结果见表 3。此处需要明确，当采用的观测数据为位移时，计算得到的是地壳应变参数; 当采用的观测数据为速度时，计算得到的是应变率参数。图 4给出了利用基于PEIV模型的加权总体最小二乘法求解得到的各个块体反演过程中观测向量和系数矩阵随机元素残差值，图中实线表示观测向量残差，虚线表示系数矩阵残差。
表 3 各块体应变率参数计算改正项b结果Table 3. The Calculation Results of Each Block Strain Rate Parameter改正项 u/a v/a εx/a εy/a εxy/a ω/a bⅠ -5.6×10-16 -2.2×10-15 -4.4×10-22 8.7×10-24 1.9×10-23 6.9×10-23 bⅡ1 2.2×10-15 -1.2×10-15 3.1×10-23 -1.6×10-23 -3.4×10-23 -5.0×10-23 bⅡ2 2.2×10-15 7.5×10-18 1.2×10-24 -4.5×10-25 -2.7×10-24 -4.3×10-24 bⅢ1 -2.7×10-15 9.4×10-15 1.6×10-20 -1.1×10-21 -2.7×10-21 -2.5×10-21 bⅢ2 8.6×10-16 2.1×10-15 1.7×10-23 -2.4×10-23 -2.8×10-24 -1.2×10-24 bⅢ3 -3.4×10-13 5.4×10-13 5.3×10-22 2.2×10-22 3.0×10-22 1.8×10-22 bⅣ 1.8×10-13 1.8×10-13 -1.3×10-20 -4.8×10-22 -6.9×10-21 -6.4×10-21 bⅤ -5.6×10-16 -3.7×10-16 -1.9×10-24 7.1×10-25 5.4×10-26 1.2×10-24 ![]()
图 4 川滇地区地壳应变参数反演过程中观测向量和系数矩阵随机元素残差值Figure 4. Residuals of Observation Vector and Coefficient Matrix in the Region of Sichuan and Yunnan1) 在利用GPS监测点复测资料所得的位移、速度反演应变或应变率参数时，用最小二乘方法未考虑由监测点坐标组成的反演模型系数矩阵的误差影响，而PEIV-WTLS算法考虑了观测向量的误差，同时考虑了反演模型系数矩阵的扰动对反演应变参数的影响，从理论上，PEIV模型算法更加合理。
2) 结合图 3和表 3可以看出，保山块体(Ⅲ1) 区域内的GPS监测点分布较为集中，且覆盖区域较小，因此相比于其他块体，改正项b的量级更大些；滇中块体(Ⅱ2) 和羌塘块体(Ⅳ)区域内的GPS监测点覆盖范围较广，因此相比于其他块体，改正项b的量级更小些，这也与上述模拟方案中得到的结论相一致。
3)表 3中的改正项b主要用于衡量WLS与PEIV-WTLS算法解算结果之间的差距。通过比较不同块体地壳应变率参数反演改正项b的大小，可以更直观地分析影响地壳应变率参数反演过程的因素。地壳应变率参数反演模型中系数矩阵的随机元素由监测点坐标组成，含有随机误差，若不考虑其影响，则会导致解算结果存在偏差。
4) 从图 4中可以看出，利用PEIV-WTLS算法可以估算出各块体相应区域中GPS监测站点位移中包含的随机误差，而无法估算出GPS监测站点坐标含有的随机误差，即系数矩阵中的随机误差，这与模拟方案中的结果相同，主要是由于GPS监测站点坐标值量级过大，相应的应变参数偏小，较小的应变参数抑制了系数矩阵中的随机误差对模型估计的贡献度。
5) 区域性地壳应力应变研究需要极高的测量精度，因此在进行地壳应变率参数求解时必须采用理论严密的算法。虽然系数矩阵中随机元素的误差相对于随机元素本身量级过小，但在一定程度上也会影响模型估计过程，因此必须将其考虑在内。
4 结语
本文将以PEIV模型为基础的加权总体最小二乘算法应用于地壳应变参数反演，该算法考虑了以往一直被忽略的系数矩阵误差的影响，得到了合理的参数解。通过改正项b定量分析系数矩阵误差对参数估计的影响。通过对模拟数据和川滇地区的实际数据进行处理，可以得知系数矩阵误差对地壳应变参数反演的影响主要受到GPS站点坐标值量级以及应变参数量级的牵制。同时存在加权总体最小二乘算法无法恢复系数矩阵随机元素误差这一问题，对于这一问题的解决方法还有待进一步研究。
-
![]()
图 1 全国337个地级以上城市及地域分布
Figure 1. 337 Prefectural Cities and the Spatial Distribution of China
![]()
图 4 全国地级以上城市2000—2016年城区面积空间分布图
Figure 4. Spatial Distribution Map of Urban Area of China's Prefectural Cities 2000—2016
![]()
图 5 全国地级以上城市2000—2016年城区扩展面积空间分布图
Figure 5. Spatial Distribution Map of Urban Sprawl Area of China's Prefectural Cities 2000—2016
![]()
图 6 2000—2016年全国地级以上城市扩展类型空间分布图
Figure 6. Spatial Distribution Map of Urban Sprawl Types of China's Prefectural Cities 2000—2016
![]()
图 7 全国各省市2000—2016年地级以上城市用地增长弹性系数图
Figure 7. Map of Urban Sprawl Coordination Coefficient of China's Prefectural Cities 2000—2016
![]()
图 8 2000—2016年全国各省市空间扩展占用土地比例图
Figure 8. Proportion of Land Use Occupied by Urban Sprawl Areas of China's Prefectural Cities 2000—2016
表 1 本文用到的主要遥感影像数据
Table 1 Main Remote Sensing Images Used in This Paper
序号 时相 主要数据类型(分辨率: m) 1 2000年 航空影像(1)、IKONOS(1) 2 2005年 Quickbird(0.61)、航空影像(1)、IKONOS(1)、SPOT5(2.5) 3 2010年 Worldview-1(0.5)、Worldview-2(0.5)、Quickbird(0.61)、航空影像(1)、SPOT5(2.5)、ALOS(2.5)、CIRS-P5(2.2) 4 2016年 航空影像(0.5)、高景一号(0.5)、Pleiades(0.5)、北京二号(0.8)、高分二号(1)、SPOT-6(1.5)、SPOT-7(1.5)、天绘一号(2)、ZY-3(2.1) 覆盖面积占比 Worldview-1/2占23.14%，航空影像占20.06%，IKONOS占17.75%，其他占39.05% 
下载: 导出CSV
-
[1] 龙瀛, 韩昊英, 毛其智.利用约束性CA制定城市增长边界[J].地理学报, 2009, 64(8):999-1008 doi: 10.3321/j.issn:0375-5444.2009.08.010 Long Ying, Han Haoying, Mao Qizhi. Establishing Urban Growth Boundaries Using Constrained CA[J]. Acta Geographica Sinica, 2009, 64(8):999-1008 doi: 10.3321/j.issn:0375-5444.2009.08.010
[2] 邹德慈.发展中的城市规划[J].城市规划, 2010, 34(1):24-28 http://d.old.wanfangdata.com.cn/Periodical/csgh201001004 Zou Deci. Urban Planning in Progress[J]. City Planning Review, 2010, 34(1):24-28 http://d.old.wanfangdata.com.cn/Periodical/csgh201001004
[3] 李昕, 文婧, 林坚.土地城镇化及相关问题研究综述[J].地理科学进展, 2012, 31(8):1042-1049 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=QK201204380935 Li Xin, Wen Jing, Lin Jian. Review of Research on Land Urbanization and Related Studies[J]. Progress in Geography, 2012, 31(8):1042-1049 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=QK201204380935
[4] 车前进, 段学军, 郭垚, 等.长江三角洲地区城镇空间扩展特征及机制[J].地理学报, 2011, 66(4):446-456 http://d.old.wanfangdata.com.cn/Periodical/dlxb201104002 Che Qianjin, Duan Xuejun, Guo Yao, et al. Urban Spatial Expansion Process, Pattern and Mechanism in the Yangtze River Delta[J]. Acta Geographica Sinica, 2011, 66(4):446-456 http://d.old.wanfangdata.com.cn/Periodical/dlxb201104002
[5] 牟凤云, 张增祥, 迟耀斌, 等.基于多源遥感数据的北京市1973-2005年间城市建成区的动态监测与驱动力分析[J].遥感学报, 2007, 11(2):257-268 http://d.old.wanfangdata.com.cn/Periodical/ygxb200702018 Mou Fengyun, Zhang Zengxiang, Chi Yaobin, et al. Dynamic Monitoring of Built-up Area in Beijing Du-ring 1973-2005 Based on Multi-original Remote Sensed Images[J]. Journal of Remote Sensing, 2007, 11(2):257-268 http://d.old.wanfangdata.com.cn/Periodical/ygxb200702018
[6] 胡德勇, 李京, 陈云浩, 等.基于多时相Landsat数据的城市扩张及其驱动力分析[J].国土资源遥感, 2006, 18(4):46-49 doi: 10.3969/j.issn.1001-070X.2006.04.012 Hu Deyong, Li Jing, Chen Yunhao, et al. An Analysis of Urban Expansion and Its Dynamics Based on Multi-temporal Landsat Data[J]. Remote Sensing for Land & Resources, 2006, 18(4):46-49 doi: 10.3969/j.issn.1001-070X.2006.04.012
[7] 龚健雅, 彭飞飞, 王乐, 等.利用视差进行立体像对的建成区检测方法[J].武汉大学学报·信息科学版, 2016, 41(4):427-433 http://ch.whu.edu.cn/CN/abstract/abstract5410.shtml Gong Jianya, Peng Feifei, Wang Le, et al.Built-up Area Detection from a Stereo Pair Using Disparity Information[J]. Geomatics and Information Science of Wuhan University, 2016, 41(4):427-433 http://ch.whu.edu.cn/CN/abstract/abstract5410.shtml
[8] 邹进贵, 陈艳华, 丁鸽, 等.利用DMSP/OLS灯光影像提取城镇建成区的聚类阈值法[J].武汉大学学报·信息科学版, 2016, 41(2):196-201 http://ch.whu.edu.cn/CN/abstract/abstract3455.shtml Zou Jingui, Chen Yanhua, Ding Ge, et al.A Clustered Threshold Method for Extracting Urban Built-up Area Using the DMSP/OLS Nighttime Light Images[J]. Geomatics and Information Science of Wuhan University, 2016, 41(2):196-201 http://ch.whu.edu.cn/CN/abstract/abstract3455.shtml
[9] 陈基伟, 韩雪培.高分辨率遥感影像建筑容积率提取方法研究[J].武汉大学学报·信息科学版, 2005, 30(7):580-582 http://ch.whu.edu.cn/CN/abstract/abstract2224.shtml Chen Jiwei, Han Xuepei.Semiautomatic Extraction of Floor Area Ratio Based on Construction Shadow in High Resolution Remote Sensing Image[J]. Geomatics and Information Science of Wuhan University, 2005, 30(7):580-582 http://ch.whu.edu.cn/CN/abstract/abstract2224.shtml
[10] 段德罡, 黄博燕.中心城区概念辨析[J].现代城市研究, 2008, 23(10):20-24 doi: 10.3969/j.issn.1009-6000.2008.10.004 Duan Degang, Huang Boyan. Differentiation and Analysis of the Concept of Central City Region[J]. Modern Urban Research, 2008, 23(10):20-24 doi: 10.3969/j.issn.1009-6000.2008.10.004
[11] 许泽宁, 高晓路.基于电子地图兴趣点的城市建成区边界识别方法[J].地理学报, 2016, 71(6):928-939 http://d.old.wanfangdata.com.cn/Periodical/dlxb201606003 Xu Zening, Gao Xiaolu. A Novel Method for Identifying the Boundary of Urban Built-up Areas with POI Data[J]. Acta Geographica Sinica, 2016, 71(6):928-939 http://d.old.wanfangdata.com.cn/Periodical/dlxb201606003
[12] Schneider A, Friedl M A, Potere D. Mapping Glo-bal Urban Areas Using MODIS 500-m Data:New Methods and Datasets Based on 'Urban Ecoregions'[J]. Remote Sensing of Environment, 2010, 114(8):1733-1746 doi: 10.1016/j.rse.2010.03.003
[13] 王丰龙, 刘云刚.中国城市建设用地扩张与财政收入增长的面板格兰杰因果检验[J].地理学报, 2013, 68(12):1595-1606 doi: 10.11821/dlxb201312001 Wang Fenglong, Liu Yungang. Panel Granger Test on Urban Land Expansion and Fiscal Revenue Growth in China's Prefecture-Level Cities[J]. Acta Geographica Sinica, 2013, 68(12):1595-1606 doi: 10.11821/dlxb201312001
[14] Bagan H, Yamagata Y. Landsat Analysis of Urban Growth:How Tokyo Became the World's Largest Megacity During the Last 40 Years[J]. Remote Sensing of Environment, 2012, 127(1):210-222 http://www.sciencedirect.com/science/article/pii/S0034425712003653
[15] 徐惟诚.中国大百科全书-精粹本[M].北京:中国大百科全书出版社, 2002 Xu Weicheng. Encyolopedia of China[M]. Beijing:Encyclopedia of China Publishing House, 2002
[16] 中华人民共和国建设部.城市规划基本术语标准:Gb/T 50280-98[M].北京:中国建筑工业出版社, 1999 Ministry of Housing and Urban-Rural Development of P.R.C. Standard for a Basic Terminology of Urban Planning:Gb/T 50280-98[M]. Beijing:China Architecture & Building Press, 1999
[17] Li X, Gong P, Liang L. A 30-year (1984-2013) Record of Annual Urban Dynamics of Beijing City Derived from Landsat Data[J]. Remote Sensing of Environment, 2015, 166(1):78-90 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=fd9d7e89ed039a0ed285e9cade3bc474
[18] 李苗, 臧淑英, 吴长山, 等.哈尔滨市城乡结合部不透水面时空变化及驱动力分析[J].地理学报, 2017, 72(1):105-115 http://d.old.wanfangdata.com.cn/Periodical/dlxb201701009 Li Miao, Zang Shuying, Wu Changshan, et al. Spatial and Temporal Variation and Its Driving Forces of Urban Impervious Surface in Urban-Rural Continuum of Harbin[J]. Acta Geographica Sinica, 2017, 72(1):105-115 http://d.old.wanfangdata.com.cn/Periodical/dlxb201701009
[19] 张甜, 王仰麟, 刘焱序, 等. 1987-2015年深圳市主城区景观演变过程多时相识别[J].地理学报, 2016, 71(12):2170-2184 http://www.cqvip.com/QK/90059X/201612/671014318.html Zhang Tian, Wang Yanglin, Liu Yanxu, et al. Multi-temporal Detection of Landscape Evolution in Western Shenzhen City During 1987-2015[J]. Acta Geographica Sinica, 2016, 71(12):2170-2184 http://www.cqvip.com/QK/90059X/201612/671014318.html
[20] Schneider A, Woodcock C E. Compact, Dispersed, Fragmented, Extensive? A Comparison of Urban Growth in Twenty-five Global Cities Using Remotely Sensed Data, Pattern Metrics and Census Information[J]. Urban Studies, 2008, 45(3):659-692 doi: 10.1177/0042098007087340
[21] Griffiths P, Hostert P. Mapping Megacity Growth with Multi-sensor Data[J]. Remote Sensing of Environment, 2010, 114(2):426-439 doi: 10.1016/j.rse.2009.09.012
[22] Friedl M A, Mciver D K, Hodges J C F, et al. Global Land Cover Mapping from MODIS:Algorithms and Early Results[J]. Remote Sensing of Environment, 2002, 83(1-2):287-302 doi: 10.1016/S0034-4257(02)00078-0
[23] Wan B, Guo Q, Fang F, et al. Mapping US Urban Extents from MODIS Data Using One-Class Classification Method[J]. Remote Sensing, 2015, 7(8):10143-10163 doi: 10.3390/rs70810143
[24] Bartholome E M, Belward A S. GLC2000:A New Approach to Global Land Cover Mapping from Earth Observation Data[J]. International Journal of Remote Sensing, 2005, 26(9):1959-1977 doi: 10.1080/01431160412331291297
[25] Costa J N D, Bielecka E, Calka B, et al. Uncertainty Quantification of the Global Rural-Urban Mapping Project over Polish Census Data[C]. The 10th International Conference Environmental Engineering, Vilnius Gediminas Technical University, Lithuania, 2017
[26] Guindon B, Zhang Y, Dillabaugh C. Landsat Urban Mapping Based on a Combined Spectral-Spatial Methodology[J]. Remote Sensing of Environment, 2004, 92(2):218-232 doi: 10.1016/j.rse.2004.06.015
[27] Gao F, de Colstoun E B, Ma R, et al. Mapping Impervious Surface Expansion Using Medium-Resolution Satellite Image Time Series:A Case Study in the Yangtze River Delta, China[J]. International Journal of Remote Sensing, 2012, 33(24):7609-7628 doi: 10.1080/01431161.2012.700424
[28] Sun Z, Wang C, Guo H, et al. A Modified Norma-lized Difference Impervious Surface Index (MNDISI) for Automatic Urban Mapping from Landsat Imagery[J]. Remote Sensing, 2017, 9(9):942-959 doi: 10.3390/rs9090942
[29] Son N T, Chen C R. Urban Growth Mapping from Landsat Data Using Linear Mixture Model in Ho Chi Minh City, Vietnam[J]. Journal of Applied Remote Sensing, 2012, 6(1):100-106 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=169117ba7fb9e4921de718b69d6a1574
[30] Angiuli E, Trianni G. Urban Mapping in Landsat Images Based on Normalized Difference Spectral Vector[J].IEEE Geoscience & Remote Sensing Letters, 2013, 11(3):661-665 http://ieeexplore.ieee.org/document/6587128/
[31] Wang L, Li C C, Ying Q, et al. China's Urban Expansion from 1990 to 2010 Determined with Sate-llite Remote Sensing[J]. Science Bulletin, 2012, 57(22):2802-2812 doi: 10.1007/s11434-012-5235-7
[32] 陈洪, 陶超, 邹峥嵘, 等.一种新的高分辨率遥感影像城区提取方法[J].武汉大学学报·信息科学版, 2013, 38(9):1063-1067 http://ch.whu.edu.cn/CN/abstract/abstract2752.shtml Chen Hong, Tao Chao, Zou Zhengrong, et al. Automatic Urban Area Extraction Using a Gabor Filter and High-Resolution Remote Sensing Imagery[J]. Geomatics and Information Science of Wuhan University, 2013, 38(9):1063-1067 http://ch.whu.edu.cn/CN/abstract/abstract2752.shtml
[33] 陶超, 谭毅华, 蔡华杰, 等.面向对象的高分辨率遥感影像城区建筑物分级提取方法[J].测绘学报, 2010, 39(1):39-45 http://d.old.wanfangdata.com.cn/Periodical/chxb201001008 Tao Chao, Tan Yihua, Cai Huajie, et al. Object-Oriented Method of Hierarchical Urban Building Extraction from High-Resolution Remote-Sensing Imagery[J]. Acta Geodaetica et Cartographica Sinica, 2010, 39(1):39-45 http://d.old.wanfangdata.com.cn/Periodical/chxb201001008
[34] 倪欢, 林祥国, 宁晓刚.直角点特征引导的遥感影像居民地提取方法[J].武汉大学学报·信息科学版, 2017, 42(5):648-655 http://ch.whu.edu.cn/CN/abstract/abstract5732.shtml Ni Huan, Lin Xiangguo, Ning Xiaogang. A Method for Extracting Human Settlements from Remote Sensing Image Using Right Angle Corners Features[J]. Geomatics and Information Science of Wuhan University, 2017, 42(5):648-655 http://ch.whu.edu.cn/CN/abstract/abstract5732.shtml
[35] Ning X, Lin X. An Index Based on Joint Density of Corners and Line Segments for Built-up Area Detection from High Resolution Satellite Imagery[J].ISPRS International Journal of Geo-Information, 2017, 6(11):338-356 doi: 10.3390/ijgi6110338
[36] 林祥国, 宁晓刚.融合直角点和直角边特征的高分辨率遥感影像居民点提取方法[J].测绘学报, 2017, 46(1):83-89 http://d.old.wanfangdata.com.cn/Periodical/chxb201701010 Lin Xiangguo, Ning Xiaogang. Extraction of Human Settlements from High Resolution Remote Sensing Imagery by Fusing Features of Right Angle Corners and Right Angle Sides[J]. Acta Geodaetica et Cartographica Sinica, 2017, 46(1):83-89 http://d.old.wanfangdata.com.cn/Periodical/chxb201701010
[37] 焦利民, 吴苏.利用自组织网络分析1990-2010年中国主要城市扩展特征[J].武汉大学学报·信息科学版, 2014, 39(12):1435-1440, 1471 http://ch.whu.edu.cn/CN/abstract/abstract3137.shtml Jiao Limin, Wu Su.Analyzing the Characteristics of the Expansion of the Metropolises in China from 1990 to 2010 Using Self-organizing Neural Network[J]. Geomatics and Information Science of Wuhan University, 2014, 39(12):1435-1440, 1471 http://ch.whu.edu.cn/CN/abstract/abstract3137.shtml
[38] 吴晓青, 胡远满, 贺红士, 等. SLEUTH城市扩展模型的应用与准确性评估[J].武汉大学学报·信息科学版, 2008, 33(3):293-296 http://ch.whu.edu.cn/CN/abstract/abstract1565.shtml Wu Xiaoqing, Hu Yuanman, He Hongshi, et al.Accuracy Evaluation and Its Application of SLEUTH Urban Growth Model[J]. Geomatics and Information Science of Wuhan University, 2008, 33(3):293-296 http://ch.whu.edu.cn/CN/abstract/abstract1565.shtml
[39] 戴芹, 陈雪, 马建文, 等.基于遥感数据SOFM网络分类的五种城市增长方式鉴别方法应用研究[J].武汉大学学报·信息科学版, 2005, 30(6):525-528 http://ch.whu.edu.cn/CN/abstract/abstract2201.shtml Dai Qin, Chen Xue, Ma Jianwen, et al.Analysis of Five Forms of Urban Expansion Identification Mo-del Based on SOFM Classification Results of Remote Sensing Data[J]. Geomatics and Information Science of Wuhan University, 2005, 30(6):525-528 http://ch.whu.edu.cn/CN/abstract/abstract2201.shtml
[40] 刘耀林, 刘艳芳, 明冬萍.基于灰色局势决策规则的元胞自动机城市扩展模型[J].武汉大学学报·信息科学版, 2004, 29(1):7-13 http://ch.whu.edu.cn/CN/abstract/abstract4591.shtml Liu Yaolin, Liu Yanfang, Ming Dongping. Urban Expansion Model Based on Cellular Automata of Decision-Making with Grey Situation[J]. Geomatics and Information Science of Wuhan University, 2004, 29(1):7-13 http://ch.whu.edu.cn/CN/abstract/abstract4591.shtml
[41] 第一财经周刊. 《以后请不要再说"北上广深": 2018中国最新1-5线城市排名出炉!》[EB/OL].[2018-04-26]. https://www.cbnweek.com/articles/normal/20927 China Business News Weekly. 《Please Don't Say "Beijing, Shanghai, Guangzhou and Shenzhen" in the Future: The Latest Ranking of 1-5 Level Cities in China is Released》[EB/OL].[2018-04-26]. https://www.cbnweek.com/articles/normal/20927
[42] Gong P, Howarth P J. Land-use Classification of SPOT HRV Data Using a Cover-Frequency Method[J]. International Journal of Remote Sensing, 1992, 13(8):1459-1471 doi: 10.1080/01431169208904202
[43] Gong P, Howarth P J. Frequency-based Contextual Classification and Gray-level Vector Reduction for Land-use Identification[J]. Photogrammetric Engineering and Remote Sensing, 1992, 58(4):423-437
[44] 杨永春, 杨晓娟. 1949-2005年中国河谷盆地大型城市空间扩展与土地利用结构转型——以兰州市为例[J].自然资源学报, 2009, 24(1):37-49 doi: 10.3321/j.issn:1000-3037.2009.01.005 Yang Yongchun, Yang Xiaojuan. Research on Urban Spatial Expansion and Land Use Inner Structure Transformation of the Large Valley-basin Cities in China from 1949 to 2005-A Case Study of Lanzhou[J]. Journal of Natural Resources, 2009, 24(1):37-49 doi: 10.3321/j.issn:1000-3037.2009.01.005
[45] 阿斯耶姆·图尔迪, 阿里木江·卡斯木.基于多源遥感数据的奎屯市城市扩展分析[J].水土保持研究, 2013, 20(2):233-237 http://d.old.wanfangdata.com.cn/Periodical/stbcyj201302045 Asiyemu Tuerdi, Alimujiang Kasimu. Analysis of Urban Expansion in Kuitun City Based on Multi-source Remote Sensing Data[J]. Research of Soil and Water Conservation, 2013, 20(2):233-237 http://d.old.wanfangdata.com.cn/Periodical/stbcyj201302045
[46] Wang H, Ning X, Zhang H, et al. Urban Boundary Extraction and Urban Sprawl Measurement Using High-Resolution Remote Sensing Images:A Case Study of China's Provincial Capital[J]. International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, 2018, XLII-3:1713-1719 doi: 10.5194/isprs-archives-XLII-3-1713-2018
[47] 刘纪远, 王新生, 庄大方, 等.凸壳原理用于城市用地空间扩展类型识别[J].地理学报, 2003, 58(6):885-892 doi: 10.3321/j.issn:0375-5444.2003.06.012 Liu Jiyuan, Wang Xinsheng, Zhuang Dafang, et al. Application of Convex Hull in Identifying the Types of Urban Land Expansion[J]. Acta Geographica Sinica, 2003, 58(6):885-892 doi: 10.3321/j.issn:0375-5444.2003.06.012
-
期刊类型引用(4)
1. 王乐洋，邹传义. PEIV模型参数估计理论及其应用研究进展. 武汉大学学报(信息科学版). 2021(09): 1273-1283+1297 . 
百度学术
2. 钟光伟，符平贵，杨钢，张俊. 融合LSC和REHSM的地壳形变分析模型. 矿山测量. 2020(05): 27-30 . 百度学术
3. 韩杰，张松林. 附加一次和二次等式约束的Partial-EIV模型及相应算法. 测绘科学技术学报. 2019(01): 17-22+27 . 百度学术
4. 吕志鹏，隋立芬. 基于非线性高斯-赫尔默特模型的结构总体最小二乘法. 武汉大学学报(信息科学版). 2019(12): 1808-1815 . 百度学术
其他类型引用(3)
计量
- 文章访问数: 1444
- HTML全文浏览量: 183
- PDF下载量: 354
- 被引次数: 7
