一种区域尺度上森林植被碳储量的估算方法与流程

1.本发明涉及一种森林生态、统计学和3s的技术领域，尤其涉及一种区域尺度上森林植被碳储量的估算方法。


背景技术：

2.目前森林植被碳储量与碳汇的精确估算一直是学界的一个难点，也限制了森林碳汇的交易。目前森林植被碳储量与碳汇的估算方法主要有生物量换算因子法、碳通量法、过程模型法、ippc法、遥感模型法等。而生物量换算因子连续函数法是国内应用较多的方法，但这种方法忽略了灌木林、疏林、未成林等碳汇资源，造成区域森林碳汇估算结果偏低。而基于涡度相关技术的碳通量法观测范围只有几百米，由于受站点数量限制，受环境因子影响较大等因素，该方法测定的数据还存在很大不确定性，仍处于探索阶段。模型法的缺点是大多模型都是国外开发的，所需参数较难获取，模型的结果也难以得到验证等，而遥感模型优点是获取速度快，运算简单，可以动态、可重复评估，其缺点林下部分碳汇量无法估算。本发明就是利用遥感技术、空间插值技术结合地面样地调查数据，估算区域森林植被地上、地表和地下部分碳储量与碳汇量。此方法的优点是能够快速、动态、精确估算区域尺度上森林植被碳储量与碳汇量，可适用于不同空间尺度。此方法的缺点是缺少植被生物学生态生理过程。
3.森林碳密度是评估森林碳汇功能的一个重要指标，由于遥感波段信息和森林碳密度具有一定的相关性，已经成为区域森林碳储量与碳汇估算的主要方法。蒋九华等利用遥感影像数据和标准样地调查数据，结合invest模型碳储量模块，评估分析了北京山区森林生态系统的碳储量，结果为北京山区森林生态系统的平均碳密度为99.95mg/hm
21.。张沁雨等基于基准样地法和国产高分数据对湖南省森林植被碳储量进行估测，湖南省碳储量总量为22.28mt
2.。张桂莲采用样地调查数据与landsat oli遥感影像，构建基于多元逐步回归模型和普通克里格残差矫正相结合的估算模型评估上海城市森林总碳储量为2.87mt
3.。liu等利用改进的biomebgc模型，结合样地调查和遥感数据，估算浙江省竹林地上部分平均碳密度，估算结果为浙江省竹林碳密度从2000年的6.75mg/ha增加到2014年的19.07mg/ha，平均年增加为0.88mg/ha
4.。虽然已有学者对区域森林碳储量与碳密度进行了研究，但大多数研究只估算地上部分森林碳储量，没有从地上、地表、地下部分分别估算，估算结果往往偏低。有些模型结构复杂，参数难于获取，造成估算结果偏差较大。因此，为了科学评估森林生态系统碳储量，需要研究新的、可复制、可推广的估算方法。


技术实现要素：

4.本发明旨在解决上述缺陷，提供一种区域尺度上森林植被碳储量的估算方法。
5.为了克服背景技术中存在的缺陷，本发明解决其技术问题所采用的技术方案是：这种区域尺度上森林植被碳储量的估算方法，该估算方法包括，
6.第一步、进行标准样地数据采集：
7.a、设置乔木林、灌木林研究区域标准样地，根据研究区域尺度大小设置标准样地个数，样地数量要尽可能均匀分布于研究区域森林植被覆盖区；
8.b、计算乔木林标准样地碳密度，包括：
9.b1、乔木林生物量计算：对乔木林标准地进行每木检尺，记录胸径和树高数据,把记录的胸径或胸径和树高代入一元生物量模型或二元生物量模型，计算出乔木林标准地生物量；对于缺少生物量模型的标准样地，采用平均标准木法，伐倒标准木，当场称干、枝、叶、根鲜重，分别采集干、枝、叶、根300-500g样品送回实验室置于100℃的烘箱中烘干至恒重，计算标准木生物量，并推算标准地生物量；
10.b2、林下灌木、草本、枯落物生物量计算：采用全部收获法，收获样方内所有灌木枝、灌木叶、灌木干、灌木根、草本、枯落物，称其鲜重，分别记录数据值，在每个样方内灌木枝、灌木叶、灌木干、灌木根、草本、枯落物各取100-200g鲜重，将其送回实验室置于100℃的烘箱中烘干至恒重，计算样方内灌木、草本、枯落物生物量；
11.b3、挖出该样方内所有根系，去掉泥土，称其鲜重，并采集100-200g样品送回实验室置于100℃的烘箱中烘干至恒重，计算样方内根系生物量；
12.b4、根据样方和标准地面积，分别计算出标准地内灌木、草本、枯落物和根系生物量；
13.b5、在实验室测定乔木、灌木、草本、枯落物和根的平均含碳率，根据乔木林标准地生物量和平均含碳率，分别计算标准样地乔木林、林下灌木、林下草本、枯落物和根系碳密度；
14.c、计算灌木林标准样地碳密度，采用收割法计算生物量，测定各部分含碳率，计算样地碳密度；
15.d、计算标准样地土壤碳密度，在标准地的样方内挖土壤剖面，土壤分0～10cm、10～20cm、20～30cm、30～40cm、40～60cm，共5层分别采样；
16.e、统计各标准样地数据，导入arcgis中，建立标准样地数据库geodatabase；
17.第二步、下载遥感数据与预处理：
18.a、下载研究区域遥感数据；
19.b、利用envi5.6的遥感数据处理软件，对下载的遥感数据进行预测处理；
20.第三步、乔木林、灌木林地上的部分碳储量估算：
21.a、提取生态遥感因子，在arcgis中，利用标准样地点坐标，提取样地点位单波段、植被指数、图像变换、纹理特征、地形、高程、坡度坡向等遥感生态因子；
22.b、相关性分析，利用spss软件，对样地乔木林、灌木林碳密度与提取的遥感生态因子做相关性分析，选取具有相关性，具有统计学意义的遥感生态因子参与构建碳密度遥感反演模型；
23.c、构建碳密度遥感反演模型，利用spss软件，以样地点乔木林、灌木林碳密度为因变量，样地点遥感生态因子为自变量，进行多元回归分析，建立多元回归模型，并对回归模型进行精度验证评价，直到精度满足要求；
24.d、根据构建满足精度的遥感反演模型和森林植被面积，反演与统计研究区域乔木、灌木林地上部分碳储量，并做相关分析与制图，或者做乔木、灌木林碳储量动态变化分析；
25.第四步、林下灌木层、草本层、枯落物、根系和土壤碳储量估算，通过利用统计学的空间插值技术估算，包括：
26.a、空间插值，包括三个重要理论为区域化变量、半变异函数和协方差函数：
27.a1、区域化变量是一个具有空间信息特征或现象的随机的变量；
28.a2、半变异函数，kriging插值首先要确定区域化变量的半变异函数，半变异函数是kriging插值的重要组成部分，其反映的是区域化变量z(x)与z(x+h)的相似程度的量，其函数决定了未采样点的值，半变异函数公式如下：
[0029][0030]
公式中z(xi)、z(xi+h)为区域化变量，二者相距为h的距离；
[0031]
a3、协方差函数，协方差函数反映是区域化变量z(x)与z(x+h)的差异程度的量，该函数同样也是决定未采样点的值，其与变异函数都是反映空间信息特征或现象相关系数的表现，采用变异函数：
[0032][0033]
公式中z(xi)、z(xi+h)为区域化变量，二者相距为h的距离，分别为z(xi)、z(xi+h)的样本平均数；
[0034]
b、插值数据预处理，kriging插值要求插值数据服从正态分布，如果插值数据不服从正态分布，采用kriging插值之前必须对数据进行预处理，数据预处理的主要方法是对数变换，可以将不具有正态分布的数据经过对数变换后，服从正态分布或基本上服从正态分布，满足插值要求即可；
[0035]
c、插值模型选择，kriging插值的类型主要有普通、简单、协同、贝叶斯、概率等插值，一旦kriging插值类型选择之后，就要选择相应的变异函数模型，变异函数模型的主要参数值有块金值、变程、基台值、偏台值、标准平均值；
[0036]
d、插值结果交叉验证和精度验证，根据arcgis中自动交叉验证结果，利用spss软件，对全体实测值与预测值进行曲线拟合验证，分析二者之间的相关性是否显著，同时利用检验样地提取预测数据，并对检验样地的实测数据和预测数据进行精度验证，分析二者之间的相关性是否显著，插值结果的精度。
[0037]
根据本发明的另一个实施例，进一步包括所述第一步的a步骤中设置乔木林、灌木林研究区域标准样地为在对角线上设置3个林下灌木样方为1m
×
1m，3个草本样方为1m
×
1m，3个枯落物样方为1m
×
1m，3个根系样方为1m
×
1m
×
1m，3个土壤剖面同根系样方向，并纪录每个标准样地中心点位经纬度坐标。
[0038]
根据本发明的另一个实施例，进一步包括所述第一步中d步骤的土壤采样用环刀、铝盒和自封袋分别取样，编号记录，送回实验室，测定土壤容重、含水率和土壤有机质，计算土壤碳密度。
[0039]
根据本发明的另一个实施例，进一步包括所述第二步中a步骤的研究区域遥感数据根据研究尺度大小下载landsat系列、sentinel2、gf1或者下载某一类系列遥感数据，做动态分析和研究区域dem。
[0040]
根据本发明的另一个实施例，进一步包括所述第二步中b步骤的预测处理包括几何精校正、辐射定标、大气校正、镶嵌、研究区域裁剪和监督分类的程序，并进行单波段、植被指数、图像变换、纹理特征的图像信息与增强处理。
[0041]
本发明的有益效果是：这种区域尺度上森林植被碳储量的估算方法操作简单，可以快速、动态估算区域森林植被生态系统碳储量，预测值精度较高。可以适用于不同尺度下森林植被生态系统碳储量的估算，并且该方法可推广、可负值。
附图说明
[0042]
下面结合附图和实施例对本发明进一步说明。
[0043]
图1是本发明标准样地的3个土壤剖面同根系样方向结构示意图；
[0044]
图2是本发明的实施流程图的结构示意图；
[0045]
图3是秦皇岛市森林植被碳密度分布图；
[0046]
图4是林下灌木、草本、枯枝落叶和土壤层碳密度空间分布图。
具体实施方式
[0047]
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0048]
如图2所示，为这种区域尺度上森林植被碳储量的估算方法的实施流程图，具体包括：
[0049]
第一步、进行标准样地数据采集：
[0050]
a、设置乔木林、灌木林研究区域标准样地，根据研究区域尺度大小设置标准样地个数，样地数量要尽可能均匀分布于研究区域森林植被覆盖区；
[0051]
b、计算乔木林标准样地碳密度，包括：
[0052]
b1、乔木林生物量计算：对乔木林标准地进行每木检尺，记录胸径和树高数据,把记录的胸径或胸径和树高代入一元生物量模型或二元生物量模型，计算出乔木林标准地生物量；对于缺少生物量模型的标准样地，采用平均标准木法，伐倒标准木，当场称干、枝、叶、根鲜重，分别采集干、枝、叶、根300-500g样品送回实验室置于100℃的烘箱中烘干至恒重，计算标准木生物量，并推算标准地生物量；
[0053]
b2、林下灌木、草本、枯落物生物量计算：采用全部收获法，收获样方内所有灌木枝、灌木叶、灌木干、灌木根、草本、枯落物，称其鲜重，分别记录数据值，在每个样方内灌木枝、灌木叶、灌木干、灌木根、草本、枯落物各取100-200g鲜重，将其送回实验室置于100℃的烘箱中烘干至恒重，计算样方内灌木、草本、枯落物生物量；
[0054]
b3、挖出该样方内所有根系，去掉泥土，称其鲜重，并采集100-200g样品送回实验室置于100℃的烘箱中烘干至恒重，计算样方内根系生物量；
[0055]
b4、根据样方和标准地面积，分别计算出标准地内灌木、草本、枯落物和根系生物量；
[0056]
b5、在实验室测定乔木、灌木、草本、枯落物和根的平均含碳率，根据乔木林标准地
生物量和平均含碳率，分别计算标准样地乔木林、林下灌木、林下草本、枯落物和根系碳密度；
[0057]
c、计算灌木林标准样地碳密度，采用收割法计算生物量，测定各部分含碳率，计算样地碳密度；
[0058]
d、计算标准样地土壤碳密度，在标准地的样方内挖土壤剖面，土壤分0～10cm、10～20cm、20～30cm、30～40cm、40～60cm，共5层分别采样；
[0059]
e、统计各标准样地数据，导入arcgis中，建立标准样地数据库geodatabase；
[0060]
第二步、下载遥感数据与预处理：
[0061]
a、下载研究区域遥感数据；
[0062]
b、利用envi5.6的遥感数据处理软件，对下载的遥感数据进行预测处理；
[0063]
第三步、乔木林、灌木林地上的部分碳储量估算：
[0064]
a、提取生态遥感因子，在arcgis中，利用标准样地点坐标，提取样地点位单波段、植被指数、图像变换、纹理特征、地形等遥感生态因子；
[0065]
b、相关性分析，利用spss软件，对样地乔木林、灌木林碳密度与提取的遥感生态因子做相关性分析，选取具有相关性，具有统计学意义的遥感生态因子参与构建碳密度遥感反演模型；
[0066]
c、构建碳密度遥感反演模型，利用spss软件，以样地点乔木林、灌木林碳密度为因变量，样地点遥感生态因子为自变量，进行多元回归分析，建立多元回归模型，并对回归模型进行精度验证评价，直到精度满足要求；
[0067]
d、根据构建满足精度的遥感反演模型和森林植被面积，反演与统计研究区域乔木、灌木林地上部分碳储量，并做相关分析与制图，或者做乔木、灌木林碳储量动态变化分析；
[0068]
第四步、林下灌木层、草本层、枯落物、根系和土壤碳储量估算，通过利用统计学的空间插值技术估算，这部分碳储量对维持土壤肥力、森林生态系统能量流动和物质循环、土壤c库的平衡等具有重要意义。由于landsat8等光学影像波长原因，其穿透能力有限，无法估算林下灌木层、草本层、枯枝落叶层、土壤等碳储量，这部分碳储量估算主要利用统计学的空间插值技术，包括：
[0069]
a、空间插值，包括三个重要理论为区域化变量、半变异函数和协方差函数：
[0070]
a1、区域化变量是一个具有空间信息特征或现象的随机的变量，区域化变量是一个具有空间信息特征或现象的随机的变量。假设研究区域为a，z(x)是其中的一个变量，x表示空间位置，z(x)可以用一个空间点函数来表示：z(x)＝z(x
x
，xy，xz),x
x
，xy，xz为三维空间信息，也可以是二维空间信息。z(x)在已知采样点xi(i＝1、2
……
，n)处的信息特征值为z(xi)。区域化变量还要具有某种程度的空间相关性，即变量在点x处与相距距离为h的点x+h处的值z(x)与z(x+h)具有空间相关性，同时还要满足正态分布和二阶平稳性；
[0071]
a2、半变异函数，kriging插值首先要确定区域化变量的半变异函数，半变异函数是kriging插值的重要组成部分，其反映的是区域化变量z(x)与z(x+h)的相似程度的量，其函数决定了未采样点的值，半变异函数公式如下：
[0072]
[0073]
公式中z(xi)、z(xi+h)为区域化变量，二者相距为h的距离；
[0074]
a3、协方差函数，协方差函数反映是区域化变量z(x)与z(x+h)的差异程度的量，该函数同样也是决定未采样点的值，其与变异函数都是反映空间信息特征或现象相关系数的表现，采用变异函数：
[0075][0076]
公式中z(xi)、z(xi+h)为区域化变量，二者相距为h的距离，分别为z(xi)、z(xi+h)的样本平均数；
[0077]
b、插值数据预处理，kriging插值要求插值数据服从正态分布，如果插值数据不服从正态分布，采用kriging插值之前必须对数据进行预处理，数据预处理的主要方法是对数变换，可以将不具有正态分布的数据经过对数变换后，服从正态分布或基本上服从正态分布，满足插值要求即可；
[0078]
c、插值模型选择，kriging插值的类型主要有普通、简单、协同、贝叶斯、概率等插值，一旦kriging插值类型选择之后，就要选择相应的变异函数模型，变异函数模型的主要参数值有块金值、变程、基台值、偏台值、标准平均值，插值类型和模型的选择对于预测结果至关重要。目前arcgis中常用的变异函数模型主要有指数模型、球面模型、高斯模型、圆形模型等。最优模型的选择主要是通过辨别模型的参数值大小来决定。变异函数模型的主要参数值有块金值、变程、基台值、偏台值、标准平均值等。根据最优模型的评价标准：块金值最小，结构比最大，标准平均值最接近于0，平均标准误差最接近于1，均方根最小；
[0079]
d、插值结果交叉验证和精度验证，根据arcgis中自动交叉验证结果，利用spss软件，对全体实测值与预测值进行曲线拟合验证，分析二者之间的相关性是否显著，同时利用检验样地提取预测数据，并对检验样地的实测数据和预测数据进行精度验证，分析二者之间的相关性是否显著，插值结果的精度。
[0080]
优选例，第一步的a步骤中设置乔木林、灌木林研究区域标准样地为在对角线上设置3个林下灌木样方为1m
×
1m，3个草本样方为1m
×
1m，3个枯落物样方为1m
×
1m，3个根系样方为1m
×
1m
×
1m，3个土壤剖面同根系样方向，并纪录每个标准样地中心点位经纬度坐标，见图1。
[0081]
优选例，第一步中d步骤的土壤采样用环刀、铝盒和自封袋分别取样，编号记录，送回实验室，测定土壤容重、含水率和土壤有机质，计算土壤碳密度。
[0082]
优选例，第二步中a步骤的研究区域遥感数据根据研究尺度大小下载landsat系列、sentinel2、gf1或者下载某一类系列遥感数据，做动态分析和研究区域dem。
[0083]
优选例，第二步中b步骤的预测处理包括几何精校正、辐射定标、大气校正、镶嵌、研究区域裁剪和监督分类的程序，并进行单波段、植被指数、图像变换、纹理特征的图像信息与增强处理。
[0084]
实施例：
[0085]
以秦皇岛市森林植被为研究对象，共建立了181个标准地，覆盖了整个研究区域，能够代表研究区域主要森林植被类型。按照本发明的技术方案：
[0086]
第一步、开展标准样地调查，核算标准样地碳密度，建立标准样地geodatabase数据库；
[0087]
第二步、采集研究区域landsat8遥感数据和dem地形数据，进行数据几何精校正、辐射定标、大气校正、镶嵌、研究区域裁剪等预处理，对预测后的研究数据进行监督分类，提取研究区域森林植被信息。对预测后的研究数据进行单波段、植被指数、图像变换、纹理特征等图像信息与增强处理，共78个自变量因子，分别为单波段因子b1、b2、b3、b4、b5、b6、b7；植被指数因子ndvi、rvi、div、evi、savi、msavi；k-t图像变换因子b、g、w；k-l图像变换因子p1、p2、p3；8种纹理特征值(56个特征值)，分别是平均值m1、m2、m3、m4、m5、m6、m7、方差v1、v2、v3、v4、v5、v6、v7、均一性h1、h2、h3、h4、h5、h6、h7、对比度c1、c2、c3、c4、c5、c6、c7、相异性d1、d2、d3、d4、d5、d6、d7、熵e1、e2、e3、e4、e5、e6、e7、角二阶矩a1、a2、a3、a4、a5、a6、a7、相关性c1、c2、c3、c4、c5、c6、c7；地形因子e、s、p；
[0088]
第三步、构建地上部分碳密度遥感反演模型与模型精度验证。随机选取80％左右的标准样地，即选择146个样地碳密度值与自变量因子进行相关性分析，选择相关系数达到显著水平或极显著水平因子参与遥感反演模型构建，并对模型精度进行验证。利用spss软件，选择相关性高的自变量因子参与逐步回归建模。构建的回归模型为：
[0089]
y＝-35.319+0.023
×
elevation-0.068
×
pac3+20.119
×
e6+20.921
×
a3+10.277
×
c2+28.223
×
a2，剩余20％左右的标准样地进行验证。模型的预测精度为80.53％，模型精度满足要求。
[0090]
构建森林植被地上部分碳密度遥感反演模型。利用构建的遥感反演模型，反演得出研究区域总的森林植被总碳储量为11.88
×
106t，平均碳密度为30.49t/hm2，图3所示。乔木碳储量为11.46
×
106t，平均碳密度为33.49t/hm2；灌木林为0.42
×
106t，平均碳密度为8.87t/hm2。
[0091]
第四步、林下灌木层、草本层、枯落物、根系和土壤碳储量估算，根据插值的原理、方法和步骤，分别估算森林植被下灌木层、草本层、枯枝落叶层和土壤层的碳储量和碳密度。首先利用co-kriging插值法估算研究区域内的森林植被层下灌木层的碳密度，然后利用森林植被的边界提取出林下灌木碳密度的分布图，最后统计出林下灌木层的碳储量为：0.83
×
106t，平均碳密度为2.13t/hm2,图4中a所示。同理，利用该方法估算草本层、枯枝落叶层和土壤层的碳储量和碳密度，估算结果分别为林下草本层碳储量为0.23
×
106t，平均碳密度为0.58t/hm2，图4中b所示，林下枯枝落叶层碳储量为1.03
×
106t，平均碳密度为2.63t/hm2，图4中c所示，林下土壤层有机碳储量为23.08
×
106t，平均碳密度为58.72t/hm2，图4中d所示。
[0092]
以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，根据本发明的技术方案及其发明构思加以等同替换或改变，都应涵盖在本发明的保护范围之内。