遥感和空间分析支持下的耕地、城镇动态变化空间采样方法研究，本文主要内容关键词为：遥感论文,空间论文,耕地论文,支持下论文,城镇论文，此文献不代表本站观点，内容供学术参考，文章仅供参考阅读下载。
一、引言
我国是一个地域广阔、资源类型多样、土地利用类型和土地利用结构复杂的国家，同时，又是一个发展中的国家，正处在工业化的过程中，土地资源的利用类型变化巨大，尤其是和国民经济发展密切相关的耕地、城镇等重要利用类型的动态变化程度极高，而国民经济的高速发展迫切要求及时掌握耕地、城镇的动态变化信息。如何获取耕地、城镇的动态变化尤其是年际变化数据，不仅涉及到相关技术手段，同时也涉及到各项社会经济指标，是一项复杂的综合工程。由于现阶段我国特有的社会发展状况以及复杂多变的自然条件，使得采用传统方法需要大量的人力、物力和财力，无法保证数据的客观和准确，更为重要的是需要花费大量的时间，从而降低耕地、城镇监测结果的时态精度。无法在满足相应调查数据准确性的同时，在时间周期上满足国家相关决策部门对耕地、城镇动态变化数据的需求。全覆盖的遥感方法从纯技术的角度来看，要满足上述要求是可行的，但由于遥感信息源的限制以及考虑到社会经济综合指标，要实现每年一次的监测也存在着较大的困难。统计采样监测方法在我国的资源监测中已得到广泛的应用，如国家统计局覆盖全国800个县的综合监测网以及农业部300个县的农情监测网等。此外经典统计采样方法在历次资源调查中均用于解决混杂在主体类型中的不同的细小对象问题，但经典统计采样方法用于空间对象的动态变化的监测研究，尤其是针对全部国土面积的监测，尚存在着采样样本太大和精度无法满足的难题。因此研究在现有本底数据基础上，利用遥感数据作为基本数据来源，以遥感数据获取并结合我国特有的资源分布和利用特征的空间采样方法，实现在耕地、城镇年际自然动态周期内的动态变化监测，以满足技术和社会经济综合指标，就具有十分现实的意义。一方面能提供针对国家耕地、城镇年际动态监测的实用技术方法，为国家决策部门提供相关的数据服务；同时也可为其他的自然社会领域的年际动态变化状况监测如土壤侵蚀、沙漠化、森林草地等，提供可供参考的技术方法。
二、耕地、城镇动态采样框架
耕地、城镇动态采样方法在于利用遥感作为基本的数据来源，利用遥感信息获取和空间分析技术，结合统计采样方法中的样本估计、最优样本选择、精度检验等。主要关键点在于构建地理区域规律控制下的耕地、城镇动态采样框架，包括基于全国耕地、城镇动态变化程度的土地利用动态度区划建立；基于耕地、城镇不同变化速度的多重动态空间信息采样框架布设和基本采样单元的耕地动态度划分。
1.基于全国耕地、城镇动态变化程度的土地利用动态度区划建立
基于全国耕地、城镇动态变化程度的土地利用动态度区划是整个采样样本布设的基础。一方面通过土地利用动态度的区划，确定我国在一定时期内的土地利用动态变化的基本发展趋势，为实施分层采样提供分层的基础；为采样先验精度前提下的样本数量的确定提供定量依据；另一方面，土地利用动态度的最小基本空间单元即为耕地、城镇动态采样的基本空间事件。
（1）最小基本空间单元的生成
最小基本空间单元作为采样的基本空间事件应具备：均一性、稳定性和空间区域的代表性。均一性要求空间内部的自然条件、土地利用类型以及其他自然社会条件的一致性；稳定性要求在一定时间段内，该空间单元的边界稳定少动；空间区域的代表性则要求该空间单元的面积大小适中，具备一定的空间代表性。考虑到最小基本空间单元的上述特性，最小基本空间单元的确定采用全国县级行政区划图和全国三级土地利用区划图，经空间叠加后形成的具有县级行政单元属性和土地利用结构区属性的最小级空间单元。具体的空间操作如图1所示：
图1 最小采样单元形成的空间操作示意图
（2）耕地、城镇动态度区划
在确定了最小采样单元之后，区划单元内的耕地、城镇动态度是保证分层准确性的基本依据。耕地、城镇的动态变化程度和动态度区划可根据1992年度1∶25万土地利用数据和1995年本底数据库1∶25万土地利用数据的直接空间分析生成，即：假设某一区划单元的耕地动态变化程度为I[,d]，城镇动态变化程度为I[,u]，则有：
I[,d]＝100%（│A[,fh]－A[,fn]│）／A[,fh]（1）
I[,u]＝100%（│A[,ch]－A[,cn]│）／A[,ch]（2）
其中A[,f]h为该区划内的1992年度耕地数据，A[,fn] 为本期耕地数据；A[,ch]为该区划内的1992年度城镇、居民点和建设用地数据；A[,cn]为该区划内的本期城镇、居民点和建设用地数据。
2.动态采样框架构成
耕地、城镇动态变化采样框架的作用在于通过对不同变化频率的区域进行分层处理，在此基础上，对不同区域、不同层次的采样区划单元中的样本数量和样本的空间分布进行控制，减少随机抽样所需的样本数，从整体上保证采样的精度。整个框架由主体框架和辅助框架构成，主体框架由耕地、城镇动态度区划构成，主要功能包括：（1 ）提供采样数据获取实施的基本空间事件单元；（2 ）提供采样前的精度先定和采样后精度检验与采样结果计算的空间事件单元。
辅助框架由土地利用区划、县级行政区划和样区TM分景分布图构成，其主要功能包括：（1）提供耕地、 城镇动态采样数据获取后的基本数据报告单元，并提供可供空间操作的空间关系控制；（2 ）提供基于不同遥感图像类型的采样样区的空间分布及其与基本采样空间事件单元和报告单元之间的空间关系；（3 ）为实现样区具有点源属性的采样数据向精度检验单元和采样数据报告单元的数据的扩散，提供逐级空间操作的基本空间关系控制框架。
采样框架由土地利用结构区划、县级行政区划和耕地、城镇动态度区划和样区TM数据分景分布图构成（图2）， 动态采样框架的时间有效性为5年，每5年根据遥感全覆盖更新数据进行一次调整。
图2 动态采样框架构成示意图
3.基本采样单元的耕地动态度类型的划分
根据1992和1996年两期耕地面积（未经特殊说明，均指毛面积）的变化率百分数β，划分为6类，分别为：
（1）持平当β≥－2.5并β≤2.5
（2）轻度增加当β＞2.5并β≤12.5
（3）轻度减少当β≥－12.5并β＜－2.5
（4）显著增加当β＞12.5
（5）显著减少当β＞－30并β＜－12.5
或β≤－30且所属行政单元不是城镇
（6）重点监测区 β≤－30且所属行政单元为城镇
三、耕地、城镇动态采样样本布设方案的确定。
耕地、城镇动态框架为耕地、城镇动态采样提供了最小采样空间单元、采样样本数量以及样本的空间分布。在此基础上，需要解决有效样本的确定、样本布设方案、采样精度先验以及耕地、城镇年际动态变化数据的汇总以及采样结果精度后验的方法。
1.有效样本的确定
在土地利用结构区划图与县级行政区划图做叠加生成基本采样单元过程中，由于自然环境条件的变化与行政界限划分的不一致性，会形成一些小图斑，对这些图斑进行采样会减少样本的代表性，从而影响采样效率。因此在确定有效样本时将其剔除，通过对每一级土地利用结构区内采样单元耕地变化率的统计分布进行分析，并根据样区的采样实验结果，确定满足如下条件的样本单元为有效采样单元：如属4 级（耕地增加），其耕地面积（1992 年或1996 年）必须超过该单元所属县（市）1996年耕地面积的1／5；其余级别的单元，要求为超过该单元所属县（市）1996年耕地面积的1／3。重点监测区实行年度全覆盖，因此在这里不做抽样误差分析。
2.空间采样方案的确定
（1）采样精度先验与采样样本布设方案的确定
本次采样采用分层随机抽样，具体实施方式为分层系统抽样。由于各有效采样单元均具有土地利用区划代码与县级行政区划代码所共同生成的索引代码，则可对其进行各层（即各级）的随机排列，然后分别按间隔为1、3和9进行采样，分别对应于50%、25%和10%的抽样。 各层样本单元为随机排列，因此采样效果与分层随机抽样的效果基本是相同的，对应以上三种采样频率，各采10次。运用采样误差的计算方法，得出各层在不同采样比情况下的精度状况，选取合适的采样频率。从这10次随机抽样中选用采样误差最小的样本单元集合作为对未来耕地动态进行监测的采样区。另外规定，最小采样样本数不得少于15个。
（2）采样信息源的使用和耕地、城镇动态数据的获取
采样样方内，采用监测当年的Landsat-TM数字图象作为基础数据源，经几何精纠正后形成分样区的图象资料，并和原有的本底土地利用数据中的耕地、城镇类型进行比较，提取耕地、城镇年际变化信息，形成单个样本的变化数。
3.采样误差的计算
在每一级土地利用结构区内， 以各级别耕地动态度基本采样单元1996年的实际耕地作为真值，应用1992年的耕地面积和根据采样样本得到的该层1992—1996年的耕地变化率估计值来反推1996年的耕地面积，所得耕地面积与实际耕地面积的相对差的绝对值即构成采样误差。利用采样单元求耕地变化率的估计值时，对各样本的变化率进行耕地面积加权平均。主要计算过程如下：
c[,p]＝（∑[,ip]S[,ip]－∑[,ip]S＃[,ip]）／∑[,ip]S＃[,ip]；
（3）
ip＝1，2，…，np；
e[,p]＝100%×│c[,p]－C[,p]│×S＃[,p]／S[,p]；（4）
p＝1，2，…，5；
c[,p]和C[,p]分别为某区划单元采样样本空间和总体空间的p级单元的耕地变化率；
S[,p]和S＃[,p]，分别为1996和1992年某区划单元内p 级单元的耕地面积；
np为p级单元采样样本数；
e[,p]为p级单元的采样误差。
采样误差的计算方法适用于采样样本布设方案的精度先验，同时也适用于动态变化数据采样获取后的精度后验。
四、耕地、城镇动态变化采样数据的汇总
1.土地利用二级区划单元内耕地、城镇年度变化的确定
在获取了各样区内耕地、城镇年度变化的绝对值后，根据各样区在动态采样框架中耕地、城镇动态度区划单元的耕地、城镇年度变化的绝对值，确定其变化率，即：
C[,fi]＝100%×（A[,t]／F[,α]）（5）
C[,Ci]＝100%×（A[,t]／C[,α]）（6）
其中C[,fi]为样区内的耕地年度变率；C[,Ci] 为样区内的城镇年度变率；A[,t]为样区内的耕地、城镇变化面积；F[,α]、C[,α] 分别为样区内耕地、城镇变化的绝对值。
得到样区耕地、城镇的年度变率后，根据耕地、城镇动态度区划的样区类型，将各样区内的变率根据分层赋予整个区划类型，并根据各样区所占整个区划类型的面积比重，得到整个耕地、城镇动态度区划类型的耕地、城镇年度变率，即：
F[,zi]＝100%（∑C[,fi]／A[,ti]） （7）
C[,zi]＝100%（∑C[,ci]／A[,ti]） （8）
其中F[,zi]、C[,zi]分别为耕地、 城镇动态度区划的样区类型内的耕地、城镇年度变率；C[,fi]为耕地、城镇动态度区划的样区类型内的任一样区的年变率；A[,ti]为耕地、城镇动态度区划的样区类型内的任一样区的耕地或城镇的面积。
2.分省耕地、城镇年度变化面积汇总
在确定了耕地、城镇动态度区划的样区类型内的耕地、城镇的年度变率后，利用全国行政区划边界图和全国耕地、城镇动态度区划的样区类型进行空间叠加，生成各行政区划内分耕地、城镇动态度区划的样区类型的年度变率分布图，在此基础上，根据年度变率分布图和分行政区资源层面叠加后得到带有年度变率信息的年度土地资源更新层面。利用本底耕地、城镇数据和年度变率得到分省耕地、城镇年度变化面积。
五、东南丘陵山地水田、林地区耕地动态采样方案布设示例
1.基本样本单元两期耕地变化的统计特性分析
该耕地利用区划一级单元覆盖了浙江、江苏、福建、江西、湖北、湖南、广东、广西和贵州等省区的部分区域。其各层基本样本单元两期耕地变化率百分数和各层两期耕地总面积的统计特性见下页表1。
根据有效样本单元的选取原则，在总体采样单元中抽取有效样本，其统计特性见表2。
2.有效样本单元与总体单元统计特性的分析
由于有效样本空间小于实际总体空间，造成c与C之间的不一致性，产生了固定的采样误差，这种误差与抽样比没有关系。但是可以通过对这种误差的先验分布作出概率统计，在进行动态采样结果分析时加入这种误差。另外，c或C与β的平均值相差很大，因为求β的平均值没有经过各采样单元面积的加权。这也说明，在用采样样本求耕地变化率时，面积加权平均的过程是必要的。对误差的估计也与一般的统计方法不同，不能单纯通过对β的分布预测采样误差，尽管能以它的分布情况大致估计误差范围和确定采样样本数或相对采样频率。
之所以对各层采用50%、25%和10%抽样比的抽样方案，是出于主体抽样框架为等比例抽样。在各区具体实施时，可进行适当的调整。根据经典采样误差理论，如果各层总体耕地变化率采用采样单元的耕地变化率的几何平均来估计，则该估计值的方差可以表示为：
V（C）＝（1－f）×S[2]／n（9）
f为抽样比，S[2]为β的方差，n为样本数。
可以看出，各层等比例抽样情况下，估计值方差之比为S[2]／n 之比，即与β的方差成正比，与采样样本量成反比。虽然实际的耕地变化率的估计值来自于样本的面积加权，其方差分布与上式并不完全相同，但是其基本趋势是相似的。为了对具体的采样方案进行误差估算，必须进行实际的抽样实验。利用前面所述的误差计算方法，得出了各层的十次不同采样比之下的采样误差分布情况。结果见表3，其中带*号者表示满足采样精度要求的可供选择的采样百分比。对于每一层，随着抽样比的减小，平均误差、95%置信区间的误差范围和25%—75%频率内的误差范围均有所增大，可以利用本表，根据事先规定的误差范围，确定各层合适的抽样样本量。由表3还可以看出，4类和5 类在相同采样比的情况下，其采样误差明显大于1、2和3类。这与表2中，4类和5类中样本单元较高具有方差是一致的。尽管样本量大于其他类，但是他们的方差更大，S[2]／n的值就比较高，抽样误差也就比较大。
六、问题与讨论
（1）遥感与空间分析技术支持下的耕地、 城镇动态采样方法的最大优势在于对遥感和采样方法进行了比较好的结合。从示例研究来看，和采用全覆盖的遥感监测相比，监测使用的遥感数据将减少3／4，监测的工作量相应减少3／4。同时，由于采用遥感和空间采样相结合，使得针对整个国土面积的大样本随机采样得以实现，克服了传统采样方法针对整个国土面积无法保证随机抽样的缺点，并保证了相应的采样精度。
表1 总体样本单元的统计特性
面积单位：亩
基本采样 S# S S占总体Cβ的平均值
层单元数N
比重W
1 25
9350935
9347873 0.049 -0.03 -0.2746
2 77 24413800 26177567 0.138 7.22
7.5776
3 67 27820362 25709616 0.136 -7.59 -7.6202
4 340 67053267 95805899 0.506 42.88 96.0161
5 158 45920076 31735862 0.168 -30.89 -33.3457
6
6857599527203 0.003 -38.53 -35.419
总体
673 175416040 189304020 1.000 7.52 40.4612
层β的均方差β的方差
1
1.39 1.95
2
2.91 8.25
3
3.20 9.17
4 256.60 65847.96
5 18.32335.74
6
5.92 35.05
总体
191.33 36608.47
表2有效样本单元的统计特性面积单位：亩
有效样 S#
S S占总体C β的平均值
层 本数n
比重W
116 81076828104739 0.05
-0.04-0.07
24321210474
22710302 0.147.07 7.30
33923538542
21752900 0.13
-7.95-7.62
4
17961688589
25707203 0.53 -41.3450.86
57236076942
25707203 0.16 -28.74
-28.85
总有效样本 349
150622229 165464861 1.009.8520.18
层β的均方差β的方差
1
1.41 1.98
2
2.93 8.60
3
2.93 8.60
4 50.71
2571.67
5 13.45180.94
总有效样本49.62
2462.23
表3各层不同采样比下抽析的误差分析
抽样
区划 -95%置信 95%置信 中值 最小值 最大值 标准差
百分比 分级 平均值 区间
区间 MEDIAN MINIMUM MAXIMUM STD-DEV
50
10.33 0.09
0.57 0.230.02
1.12 0.34
25
10.25 0.10
0.40 0.170.07
0.80 0.21
*10
11.09 0.57
1.60 1.220.11
1.91 0.72
50
20.48 0.26
0.70 0.520.02
0.86 0.30
25
20.50 0.29
0.72 0.550.00
1.05 0.30
*10
21.42 0.71
2.14 1.350.06
3.10 1.00
50
30.39 0.21
0.56 0.350.05
0.83 0.24
25
30.85 0.45
1.24 0.770.07
1.55 0.55
*10
31.15 0.94
2.08 1.580.46
2.57 0.80
50
41.71 1.08
2.33 1.600.53
3.15 0.87
*25
43.97 2.05
5.54 3.100.34
8.25 2.44
10
44.84 2.38
7.30 4.980.45 11.58 3.43
50
53.83 2.08
5.58 3.970.20
8.42 2.45
*25
53.57 2.30
4.84 3.111.35
7.29 1.77
10
55.82 2.51
9.12 5.770.09 13.34 4.62
（2）空间统计的采样方法的最初应用， 主要从地质矿产调查发展起来，如早期的克里金法等，现已发展成为一门完整的地质统计学，其主要依据即在于以空间插值为基础的空间分析方法的应用。此外，国外现已将空间统计学的方法应用于地表水文过程的研究，但对复杂地理区域上的空间对象的年际动态变化的采样方法尚待进一步的研究。本文在本底数据的基础上，采用框架控制方法，在技术路线上是可行的，但由于需要每5年更新一次本底数据，增加了整个采样方法实施的难度。
（3）本方法不仅适用于耕地、城镇的年际动态监测， 同时也适用于特定空间对象的年际动态监测。
（4）在不稳定对象的动态变化空间插值的空间分析方法上， 有待于进一步的研究，包括采样点的控制范围、插值区域的大小、误差传递模型以及精度分析等。