协同多源遥感数据的北亚热带森林蓄积量贝叶斯分层估测

协同多源遥感数据的北亚热带森林蓄积量贝叶斯分层估测 Modeling forest growing stock volume in a north subtropical region using the hierarchical Bayesian approach based on multi-source data 1 2 4 first-author 林文科，研究方向为林业遥感。E-mail： 18120825445@163.com 林文科，研究方向为林业遥感。E-mail： 18120825445@163.com 18120825445@163.com 3 1 2 1 2 1 2 1 2 corresp 陆灯盛，研究方向为林业遥感。E-mail： ludengsheng@fjnu.edu.cn 陆灯盛，研究方向为林业遥感。E-mail： ludengsheng@fjnu.edu.cn ludengsheng@fjnu.edu.cn 福建师范大学 湿润亚热带山地生态国家重点实验室培育基地, 福州 350007 福建师范大学 湿润亚热带山地生态国家重点实验室培育基地, 福州 350007 State Key Laboratory for Subtropical Mountain Ecology of the Ministry of Science and Technology and Fujian Province, Fujian Normal University, Fuzhou 350007, China State Key Laboratory for Subtropical Mountain Ecology of the Ministry of Science and Technology and Fujian Province, Fujian Normal University, Fuzhou 350007, China 福建师范大学 地理研究所, 福州 350007 福建师范大学 地理研究所, 福州 350007 Institute of Geography, Fujian Normal University, Fuzhou 350007, China Institute of Geography, Fujian Normal University, Fuzhou 350007, China 国家林业和草原局华东调查规划设计院, 杭州 310000 国家林业和草原局华东调查规划设计院, 杭州 310000 Institute of East China Inventory and Planning, National Forestry and Grassland Administration, Hangzhou 310000, China Institute of East China Inventory and Planning, National Forestry and Grassland Administration, Hangzhou 310000, China 福建省减灾中心, 福州 350001 福建省减灾中心, 福州 350001 Disaster Reduction Center of Fujian Province, Fuzhou 350001, China Disaster Reduction Center of Fujian Province, Fuzhou 350001, China 精确估算森林蓄积量是国家实现2060年前碳中和目标的迫切需求，而基于遥感的森林蓄积量定量反演是当前遥感应用领域面临的重要挑战和研究热点。光学遥感数据由于无法获取森林高度信息并存在信号饱和问题，反演森林蓄积量的精度较低，而机载Lidar数据能获取高度信息，但成本高、观测范围有限。本研究利用Sentinel-2多光谱、资源三号立体成像和机载Lidar数据，定量评价多源遥感数据在北亚热带估测森林蓄积量的能力，并探索小样本情况下贝叶斯模型在森林蓄积量估算中的优势，旨在解决单一数据集估测森林蓄积量精度低和样本数量不足的问题。研究结果表明，协同Sentinel-2和林分冠层高度模型（CHM）数据的贝叶斯分层建模能够显著提高森林蓄积量估测精度，其相对均方根误差（rRMSE）为27.6%，其估测误差比Sentinel-2数据减少了13.6 m 3 /ha，比单独CHM数据减少了7.4 m 3 /ha。多源数据可以有效改善低值高估和高值低估问题，相比单一数据集，多源数据在低值区估测误差减少了约1/4，在高值区减少了约1/3。分层贝叶斯模型有效缓解了因精细分层对样本数量的需求，减少了森林类型和地形差异对森林蓄积量估测精度的影响。为亚热带森林类型和地形复杂的地区，开展森林蓄积量估测制图提供关键技术手段。 Accurate estimation of Forest Growing Stock Volume (FGSV) is needed to achieve the goal of carbon neutral. Quantitative inversion of FGSV using remote sensing technologies is still a research challenge. Optical remote sensing technology is one of the most important means for FGSV estimation, but cannot provide sufficiently accurate estimates due to lack of canopy structure features and data saturation problem. Although airborne Lidar can overcome the shortcoming of optical sensor data, its high cost in data collection and limited observation area constrain its extensive application. This research employs integration of Sentinel-2, ZY-3 stereo, and airborne Lidar data to explore the performance of FGSV estimation in north subtropical regions, and examines the advantages of using the hierarchical Bayesian approach to develop FGSV estimation models under the condition of small population of sample plots. The objective is to solve low modeling accuracy caused by the single sensor data and insufficient number of sample plots. The results indicate that the hierarchical Bayesian approach based on combination of Sentinel-2 and Canopy Height Model (CHM) data (subtraction of Lidar-derived digital elevation model data from ZY-3 stereo-derived digital surface model data) provides the best estimation results with relative Root Mean Square Error (rRMSE) of 27.6%. The Root Mean Square Error (RMSE) using this approach reduced by 13.6 m 3 /ha comparing with the RMSE based on Sentinel-2 data alone, and reduced by 7.4 m 3 /ha based on CHM data alone. The research shows that use of multi-source data can effectively improve the problems of overestimation when FGSV is small and of underestimation when FGSV is relatively high, that is, use of multi-source data can reduce the overestimation by one forth and the underestimation by one third comparing with use of single data source alone. Comparing with traditional modeling approaches such as linear regression and random forest, the hierarchical Bayesian approach can effectively reduce the requirement of number of samples due to use of stratification strategy and reduce the impacts of forest types and terrain differences on FGSV estimation accuracy. This research provides new insights of using integration of different data sources to develop FGSV estimation models to achieve accurate estimates, and provides key technology for FGSV mapping in subtropical regions. 遥感 森林蓄积量 贝叶斯分层模型 Sentinel-2 资源三号 机载Lidar 多源数据 remote sensing forest growing stock volume hierarchical Bayesian approach Sentinel-2 ZY-3 airborne Lidar multi-source data 1　引言 1 森林是减缓气候变化和提供多种生态服务功能的重要陆地生态系统，是国家重要的战略资源。森林蓄积量是反映森林资源质量状况和实现2060碳中和的重要基础数据，精准估测森林蓄积量也是国家制定科学合理的碳中和方案的迫切需求。亚热带地区具有独特的自然特征和森林分布特征，其复杂的地表增加了样地调查的难度和蓄积量估测的不确定性。而中国具有广袤的亚热带森林分布，如何提高亚热带地区森林蓄积量估测精度，实现快速、准确的森林资源变化监测，对实现森林资源可持续发展具有重要意义。 遥感手段估测森林蓄积量主要通过光谱信息或垂直结构信息与样地实测数据建立估测模型（ 李德仁 等，2012 李德仁 等，2012 ；Dos Reis等，2019）。光学遥感数据，特别是Landsat系列，是森林蓄积量估测最常用的遥感数据源（ Barrett等，2016 Barrett等，2016 ； Babcock等，2018 Babcock等，2018 ）。Sentinel-2数据比Landsat数据具有更高的空间、光谱和时间分辨率，在森林蓄积量估测中更具优势（ Chrysafis等，2017 Chrysafis等，2017 ； Mura等，2018 Mura等，2018 ）。 然而，包括Sentinel-2在内的各种光学遥感都存在着穿透性差，只能记录植被水平结构特征，难以反映林分垂直结构信息等问题；另外，光学遥感数据的饱和问题，特别是在森林茂密的情况下，导致森林蓄积量的严重低估，限制了实际中的广泛应用（ Zhao等，2016 Zhao等，2016 ； Nuthammachot等，2022 Nuthammachot等，2022 ）。 激光雷达数据能获得森林垂直结构信息，可以有效的提高蓄积量估测精度，已成为森林资源估测研究的重要数据源（ 李旺 等，2015 李旺 等，2015 ； Jayathunga等，2018 Jayathunga等，2018 ），但数据获取成本昂贵，覆盖面积小，难以进行大范围的森林蓄积量反演研究。卫星立体成像数据覆盖范围广，重返周期短，数据成本低，可获取大范围的地面高程信息。 Li等（2019） Li等（2019） 利用生长季和落叶季的资源三号（ZY-3）立体成像数据提取北方落叶松林的冠层高度模型（CHM），进而估测森林地上生物量，显著地提高了森林地上生物量的估测精度。但这种方法仅仅适用于地形平缓的落叶林区，并不适合地形复杂的山区及亚热带常绿阔叶林区。 森林蓄积量估测通常采用传统的线性回归或机器学习算法（ 罗环敏 等，2011 罗环敏 等，2011 ； Lu等，2016 Lu等，2016 ）。机器学习算法能有效地拟合复杂的非线性关系，在森林蓄积量估测研究中得到广泛的应用（ 蒙诗栎 等，2017 蒙诗栎 等，2017 ； Feng等，2017 Feng等，2017 ； Gao等，2018 Gao等，2018 ）。但机器学习方法中的参数优化和模型过拟合等问题，常导致预测结果精度低。贝叶斯分层模型以低样本需求和良好的防过拟合能力成为当前研究的热点（ Junttila和Laine，2017 Junttila和Laine，2017 ；Ver Planck等，2018）。贝叶斯方法能够实现复杂的、具有层次结构的参数关联的有效估计，可以在较小训练样本集和高度相关的多维数据上取得良好的效果（ Junttila等，2015 Junttila等，2015 ）。 为解决基于单一数据集、有限的建模样本、亚热带地区复杂的地形和多样的森林类型导致森林蓄积量模型估测精度低的问题，本研究将联合Sentinel-2多光谱数据，ZY-3立体像对和机载Lidar数据，结合森林类型和地形等因素，采用贝叶斯分层方法构建森林蓄积量估测模型，解决建模对样本量的需求以及亚热带地区样本采集困难的矛盾，提高森林蓄积量估测制图精度，为大范围的森林资源遥感调查提供关键技术支撑。 2　材料与方法 1 2.1　研究区 2 金寨县位于安徽省西部（ 图1 图1 ），属于北亚热带湿润季风区，年均降水量1500 mm左右，梅雨显著、夏雨集中。夏季平均气温在22℃以上，冬季平均气温在10℃以下。该县总面积为3814 km 2 ，居安徽省各县之首。全县平均海拔500 m，大别山山脉由西南向东北贯穿全境，形成相对高差1600余米的垂直地势分布特征，使得金寨县南北的物候相差半月左右。全县森林资源丰富，森林覆盖率高达72.75% （ http：//www.ahjinzhai.gov.cn/ http://www.ahjinzhai.gov.cn/ ［2021-08-09］）。典型植被属亚热带常绿阔叶林，但由于温度、降水量和光照等环境因子受海拔高度的影响，植被类型具有明显的垂直地带性，形成了杉木、马尾松、板栗、竹类、油桐等为主的森林植被分布特征。 1,2;3; 2.2　研究数据及技术路线 2 本研究利用了Sentinel-2多光谱、ZY-3立体成像、机载Lidar和样地调查等多种数据（ 表1 表1 ）。为了评估不同遥感数据在森林蓄积量估测中的作用，制定了 图2 图2 的技术路线。主要步骤包括（1）数据的预处理：对机载Lidar、ZY-3、Sentinel-2和外业样地调查数据的收集和处理，分别得到数字地形模型（DTM）、数字地面模型（DSM）、多光谱波段和样地蓄积量；（2）变量的提取：基于Lidar-DTM和ZY-3 DSM生成的林分冠层高度模型数据（CHM），Sentinel-2多光谱波段以及各种植被指数；（3）变量的提取和筛选：利用统计方法提取基于CHM和Sentinel-2多光谱的变量，并将二者合并作为多源数据变量，形成3个数据变量集，然后分别使用逐步回归方法筛选变量，筛选的变量分为整体变量（对所有样本建模选中的变量）和分层变量（分层样本单独建模选中的变量）；（4）贝叶斯分层建模：基于3种变量数据集筛选出来的变量分别利用双因素贝叶斯分层方法进行蓄积量估测建模，这里的双因素是指坡向和森林类型两个分层因素；（5）结果分析：比较分析3种数据方案的建模效果，选出最优的数据方案结合森林分布图反演整个研究区的蓄积量空间分布。 数据源 获取时间 描述 Sentinel-2 2020-04-09 Sentinel-2多光谱数据有13个波段，包括4个（蓝、绿、红、近红）10 m波段，6个（4个红边和2个短波红外）20 m波段和3个（气溶胶、水蒸气、短波卷云）60 m波段；本研究选用10 m和20 m的波段 ZY-3立体成像 2020-04-09 ZY-3 02星立体成像数据包括一个2.1 m分辨率的正视波段和2个2.5 m的前视、后视波段 Lidar 2019-06 使用RIEGL-VQ-1560i激光雷达航摄仪获取的点云密度为2点/m 2 的Lidar 数据 样地调查数据 2019-09—2019-10 样地为长宽均25.82 m的正方形 森林资源调查数据 2019 来自安徽省林地“一张图”数据集，标注了地类及树种类型 2.3　数据收集与整理 2 以研究区的森林类型和年龄为分层依据，运用分层抽样法布设样地。样地大小均为25.82 m × 25.82 m（1亩）。样地调查于2019-09—2019-10进行，共调查样地71个，包括阔叶林（麻栎、板栗、杨树、枫香等）29个、马尾松26个、和杉木16个。野外作业时，对样地进行每木检尺，起测胸径为5 cm，详细记录单木的树种、树高、胸径等因子，并记录样地的优势树种、平均树高、龄组、郁闭度、坡度、坡向等因子。依据树种、胸径查找安徽省立木材积表计算单木材积（安徽省地方标准：杉木一元立木材积表（DB 34/T 1724—2012），马尾松立木材积表（DB34/T 3345-2019）， 硬阔立木材积表（DB34/T 3907-2021），杨树立木材积表（DB34/T 3620-2020）），对样地内的单木材积累加得到样地蓄积量（m 3 /亩），进而换算为立方米/公顷（m 3 /ha）。 表2 表2 为样地数据的统计结果。 森林类型 样地数 蓄积量/（m 3 /ha） 平均值/（m 3 /ha） 标准差/（m 3 /ha） 变异系数/% 阔叶 29 24.8—213.7 99.6 52.0 52.2 马尾松 26 21.4—237.0 140.5 57.4 40.9 杉木 16 3.2—362.9 219.6 104.8 47.8 总计 71 3.2—362.9 141.6 82.1 57.9 两景Sentinel-2 L1C级多光谱数据由欧洲空间局网站下载，获取时间为2020-04-09。使用Sen2cor插件进行辐射定标和大气校正，将L1C级产品处理成L2A级（地表反射率数据）。根据研究目的，去除3个60 m波段（气溶胶、水蒸气、短波卷云），使用SNAP融合工具将4个10 m和6个20 m分辨率的波段融合成10 m分辨率（ Brodu，2017 Brodu，2017 ）。为了消除地形对光谱反射率的影响，使用SCS+C校正模型和机载Lidar生成的2 m空间分辨率的DTM数据对融合后的Sentinel-2影像进行地形校正（ 高永年和张万昌，2008 高永年和张万昌，2008 ）， 然后计算植被指数，如归一化植被指数 （NDVI），比值植被指数（RVI）、土壤调整植被指数（SAVI）、增强型植被指数（EVI））。矢量格式的金寨县森林地图提取自2019年安徽省林地“一张图”。“一张图”中的每个多边形都标有详细的土地覆盖或森林类型。通过叠加Sentinel-2多光谱影像、CHM和土地覆盖边界，对森林区域进行了修改，将土地覆盖分为4类：阔叶、马尾松、杉木和其他地类，用于森林蓄积量空间分布的反演，经过调整后的金寨县森林类型分布如 图1（c） 图1（c） 。 2.4　联合机载Lidar和ZY-3立体成像数据提取林分冠层高度模型（CHM）数据 2 机载Lidar点云数据使用RIEGL-VQ-1560i激光雷达航摄仪获取，是具有双通道的机载激光扫描仪系统，系统的内部无缝集成高性能IMU/GNSS系统以及一亿像素的数码相机。点云获取时间是2019-06，密度为2点/ m 2 。对点云数据进行去噪处理后，将离散点云回波点分成地面点和非地面点两部分，利用不规则三角网插值算法对地面点进行插值运算生成数字地形模型（DTM），以2 m的空间分辨率导出研究区的DTM数据。 ZY-3立体成像数据采集时间为2020-04-09。从立体成像数据提取DSM的主要步骤包括输入数据、选择地面控制点和连接点、计算模型、生成核线和导出DSM（ Xie等，2019 Xie等，2019 ）。控制点所用的高程信息以Lidar DTM高程数据为参考。根据立体成像原理提取DSM的过程分为4步：（1）相对定向：确定两个不同影像的相对位置，获得成像瞬时两影像之间的位置关系和姿态参数；（2）绝对定向：通过地面选取的控制点（GCP）确定立体模型在地面测量坐标系中的位置；（3）连接点：通过连接左右影像建立对应关系生成核线影像；（4）输出DSM：根据核线影像提取DSM，以2 m空间分辨率输出。经分析发现正视和后视影像生成的DSM比正视和前视获取的DSM噪声少，精度高。因此，将前者生成的DSM和Lidar数据提取的DTM进行差值运算得到空间分辨率为2 m的CHM。对该CHM进行后处理，检查高于50 m和低于0 m数值的像元，对照Lidar影像确定是否对该像元进行滤波处理或替代处理。 2.5　蓄积量建模变量的提取和筛选 2 鉴于外业调查样地的大小为25.82 m× 25.82 m，针对10 m空间分辨率的Sentinel-2光谱波段及植被指数数据，以样地中心点提取3 × 3个像素的统计值作为新的变量赋予对应的样地，统计值包括最大值、最小值、平均值、标准差和方差。对于CHM数据，则以13 × 13个像素为统计单元，提取样地范围的最大值、平均值、标准差、方差和高度百分位（10 th ，20 th ，…，90 th ）作为变量。 对于这些以样地范围提取的统计变量，需要通过相关性分析进一步筛选出用于森林蓄积量建模的变量。通过将训练样本的蓄积量值与样本所在范围提取的变量建立逐步回归方程，将被选入方程中的变量作为建模变量。具体方法为，对于3种数据（CHM，Sentinel-2和多源数据），利用逐步回归方法对3种数据集分别重复以下操作：（1）利用所有训练样本建立逐步回归方程提取整体变量；（2）将训练样本按森林类型进行分组，每个组的样本各自建立逐步回归方程提取分层变量；（3）将整体变量和分层变量合并作为贝叶斯分层建模变量。 2.6　基于贝叶斯分层方法的森林蓄积量估测模型 2 贝叶斯分层模型是考虑了先验分布的多层模型，允许同时对不同层次的数据进行建模，从而有效地考虑数据的嵌套和依赖关系。多层模型的核心是通过对线性组合预测因子 <math id="M1"><mi>η</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474130&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474127&type= 2.20133328 3.38666677 进行变换来预测响应变量 <math id="M2"><mi>y</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474143&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474140&type= 1.69333339 3.47133350 ，并假设 <math id="M3"><mi>y</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474143&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474140&type= 1.69333339 3.47133350 服从于一个通过反函数 f 假设确定的分布 D， <math id="M4"><msub><mrow><mi>y</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>~</mo><mi>D</mi><mo stretchy="false">(</mo><mi>f</mi><mfenced separators="|"><mrow><msub><mrow><mi>η</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced><mo>,</mo><mi>θ</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474150&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474146&type= 23.36800003 5.16466665 ， <math id="M5"><mi>D</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474154&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474152&type= 2.53999996 2.62466669 是对模型中误差分布和反函数的描述，参数 <math id="M6"><mi>θ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474158&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474156&type= 1.77800000 2.70933342 表示不随数据点变化的特定参数， i 表示数据点的个数。线性组合的预测因子一般表达为 <math id="M7"><mi>η</mi><mo>=</mo><mi>X</mi><mi>β</mi><mo>+</mo><mi>Z</mi><mi>μ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474162&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474160&type= 19.55800056 3.38666677 ，式中， <math id="M8"><mi>β</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474166&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474164&type= 1.86266661 3.38666677 表示整体水平上的系数， <math id="M9"><mi>μ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474170&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474168&type= 2.03200006 3.38666677 表示分层水平上的系数， <math id="M10"><mi>β</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474177&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474174&type= 1.86266661 3.38666677 和 <math id="M11"><mi>μ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474191&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474189&type= 2.03200006 3.38666677 组成模型的参数集，分别代表固定影响和随机影响。 X 和 Z 为对应于系数 <math id="M12"><mi>β</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474199&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474197&type= 1.86266661 3.38666677 和 <math id="M13"><mi>μ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474191&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474189&type= 2.03200006 3.38666677 组成的变量矩阵。因此 <math id="M14"><mi mathvariant="bold-italic">y</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474195&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474193&type= 1.77800000 3.47133350 、 X 和 Z 组成了模型的变量数据集。获得整体水平回归参数 <math id="M15"><mi>β</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474199&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474197&type= 1.86266661 3.38666677 需要总体样本的先验分布，但先验分布需要前期大量数据的支撑，因此假定整体水平的参数不受正态先验的限制。分层水平的回归参数 <math id="M16"><mi>μ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474203&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474201&type= 2.03200006 3.38666677 假设来自均值为零且协方差矩阵 Σ 未知的多元正态分布： <math id="M17"><mi>μ</mi><mo>~</mo><mi>N</mi><mo stretchy="false">(</mo><mn mathvariant="normal">0</mn><mo>,</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474207&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474205&type= 12.10733318 3.55599999 Σ ），当假设不同分组参数之间的协方差为零时，矩阵可以分解为： <math id="M18"><msub><mrow><mi>μ</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>~</mo><mi>N</mi><mo stretchy="false">(</mo><mn mathvariant="normal">0</mn><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474210&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474209&type= 16.84866714 4.74133301 ）， <math id="M19"><mi>k</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474215&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474212&type= 1.69333339 2.62466669 指示分组的因素。并假设同组别内不同级别（由 j 索引）的参数是独立的，那么可以得到： <math id="M20"><msub><mrow><mi>μ</mi></mrow><mrow><mi>k</mi><mi>j</mi></mrow></msub><mo>~</mo><mi>N</mi><mo stretchy="false">(</mo><mn mathvariant="normal">0</mn><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">V</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474220&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474217&type= 17.10266685 4.65666676 ），将协方差矩阵 <math id="M21"><msub><mrow><mi mathvariant="bold-italic">V</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474225&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474222&type= 3.38666677 3.72533321 作为参数进行建模，令： <math id="M22"><msub><mrow><mi mathvariant="bold-italic">V</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>=</mo><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mo stretchy="false">(</mo><msub><mrow><mi>σ</mi></mrow><mrow><mi>k</mi></mrow></msub><mo stretchy="false">)</mo></mrow></msub><msub><mrow><mi mathvariant="bold-italic">Ω</mi></mrow><mrow><mi>k</mi></mrow></msub><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mo stretchy="false">(</mo><msub><mrow><mi>σ</mi></mrow><mrow><mi>k</mi></mrow></msub><mo stretchy="false">)</mo></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474229&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474227&type= 25.14599991 5.58799982 ，这里 <math id="M23"><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mo stretchy="false">(</mo><msub><mrow><mi>σ</mi></mrow><mrow><mi>k</mi></mrow></msub><mo stretchy="false">)</mo></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474233&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474231&type= 6.34999990 4.31799984 表示具有对角元素的对角矩阵 <math id="M24"><msub><mrow><mi mathvariant="bold-italic">σ</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474237&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474235&type= 3.64066648 3.72533321 ， <math id="M25"><msub><mrow><mi mathvariant="bold-italic">Ω</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474241&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474239&type= 3.80999994 3.89466691 代表 <math id="M26"><msub><mrow><mi mathvariant="bold-italic">V</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474245&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474243&type= 3.38666677 3.72533321 参数化后的相关矩阵。为 <math id="M27"><msub><mrow><mi mathvariant="bold-italic">Ω</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474250&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474247&type= 3.80999994 3.89466691 指定参数是通过 LKJ -相关先验参数 <math id="M28"><mtext> </mtext><mi>ξ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474253&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474251&type= 2.45533323 3.30200005 > 0来定义先验（ Lewandowski等，2009 Lewandowski等，2009 ）， <math id="M29"><msub><mrow><mi mathvariant="bold-italic">Ω</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>~</mo><mi>L</mi><mi>K</mi><mi>J</mi><mo stretchy="false">(</mo><mi>ξ</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474259&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474256&type= 17.27199936 4.65666676 ，这种先验通常比二分之一的柯西先验更能使模型收敛，但仍然是相对薄弱的信息。通过 LKJ 模型来得到 <math id="M30"><msub><mrow><mi mathvariant="bold-italic">Ω</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474263&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474261&type= 3.80999994 3.89466691 ，并最终可以求出未知协方差矩阵 <math id="M31"><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474269&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474266&type= 3.55599999 4.23333359 ，从而得到分层水平的回归参数 <math id="M32"><mi>μ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474275&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474271&type= 2.03200006 3.38666677 。 上述过程由R程序中的 brms 包完成，而 brms 通过调用 Stan 程序完成计算。 brms 通过扩展的 lme 4语法提供了高效和可读的 Stan 代码，通过R调用 Stan 程序在后台运行，运行结果以可视化的形式再返回R中，关于 brms 包中模型的一般性概述可参考（ Bates等，2015 Bates等，2015 ； Bürkner，2017 Bürkner，2017 ）。 brms 中应用的公式语法建立在 lme 4的语法之上，语法的基本模式与 lme 4一致，但在分层和参数设置上具有更大的灵活性。 lme 4基本语法的形式为 <math id="M33"><mi>r</mi><mi>e</mi><mi>s</mi><mi>p</mi><mi>o</mi><mi>n</mi><mi>s</mi><mi>e</mi><mo>~</mo><mi>p</mi><mi>t</mi><mi>e</mi><mi>r</mi><mi>m</mi><mi>s</mi><mo>+</mo><mfenced separators="|"><mrow><mi>g</mi><mi>t</mi><mi>e</mi><mi>r</mi><mi>m</mi><mi>s</mi><mo stretchy="false">|</mo><mi>g</mi><mi>r</mi><mi>o</mi><mi>u</mi><mi>p</mi></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474280&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=25474277&type= 48.68333054 4.57200003 ， response 代表响应变量， pterms 部分表示通过观测这部分的变量对整体水平的影响是相同的， gterms 部分包含了所谓的组级效应，这些效应假定随着组中指定的组变量而变化， group 表示分组的策略。而 brms 在此基础上进一步扩展使其能使用多组级效应模型和多成员模型。 本研究同时考虑了两种分组因子—森林类型（阔叶林、马尾松和杉木）和坡向（阳坡、半阳坡、半阴坡和阴坡）。然后联合两种分组因子将样本分成12组建立模型，命名为双因素贝叶斯分层模型。在贝叶斯分层模型中，模型变量由两部分组成：整体变量和分层变量。根据先验知识将一些具有整体效应的变量作为整体变量。为了获得整体和分层变量，利用所有的样本和变量进行逐步回归建模，筛选出整体变量，将样本依据分组进行逐步回归建模，筛选出每个分层的变量，最后的变量组成用于贝叶斯分层建模的变量数据集。 2.7　森林蓄积量估测结果评价 2 基于样本蓄积量和遥感变量建立的蓄积量估测模型需要进行精度评价。本文采用决定系数 R 2 评价模型回归系数拟合优度，采用均方根误差（RMSE）和相对均方根误差（rRMSE）评价模型预测能力。 R 2 越高，模型拟合效果越好；而RMSE和rRMSE越低，模型的预测效果越好。由于样本数量的限制，本研究采用留一交叉验证法，也就是说，对于71个样本，保留其中一个单独的样本作为验证，其余70个样本作为训练数据构建模型，重复操作直至每个样本被验证。每一次验证，训练数据几乎使用了所有的样本，可以保证产生可靠的预测，避免了随机因素的影响。通过比较基于不同数据集建立的贝叶斯分层模型的性能，选择最佳的数据方案，并将其应用于整个研究区的森林蓄积量预测。 3　结果分析 1 3.1　森林蓄积量估测模型分析 2 基于CHM、Sentinel-2以及多源遥感数据3种方案的建模变量筛选和模型精度验证结果（ 表3 表3 ）表明，多源遥感数据的结合具有较好的建模效果。多源数据集在整体变量和分层变量都选中了高度和光谱变量，两者共同在模型中发挥作用，模型的决定系数 R 2 从单一数据Sentinel-2的0.81和CHM 的0.83提高到了多源数据的0.93。RMSE和rRMSE数据中也体现了多源数据在改善森林蓄积量估测精度的作用。基于CHM数据的蓄积量估测误差比Sentinel-2数据的估测误差降低了4.3%，二者的结合再降低了5.2%。从森林蓄积量反演结果（ 图3 图3 ）上看，CHM反演的结果主要集中在50—150 m 3 /ha，而Sentinel-2数据反演的结果集中在50—200 m 3 /ha，更多的出现在150—200 m 3 /ha，很少出现低于50 m 3 /ha的预测。基于多源数据的反演结果则出现较多低于50 m 3 /ha的预测，以及高于250 m 3 /ha的预测结果，蓄积量较高的地方主要集中在西南部。 注： H 表示CHM变量，下标min、m表示统计最小值和平均值，数字10、20、…、90表示高度百分位；RE、RED、SWIR、EVI分别表示多光谱变量中的红边、红、短波红外和增强型植被指数，数字1、2、3表示波段编号。 2 数据 2 模型变量 2 模型 R 2 3 验证精度 整体变量 分层变量 R 2 RMSE/（m 3 /ha） rRMSE/% CHM H m H m ， H 90 ， H 30 0.83 0.68 46.4 32.8 Sentinel-2 EVI， RE3 min ， SWIR1 min SWIR1 min ， RED m ， SWIR2 min 0.81 0.58 52.6 37.1 多源数据 EVI， RE3 min ， SWIR1 min ， H 30 SWIR1 min ， RED m ， SWIR2 min ， H m 0.93 0.77 39.0 27.6 1,2;3; 3.2　贝叶斯分层模型估算结果分析 2 基于3种数据集对不同森林类型的蓄积量估测的精度（ 表4 表4 ）分析表明，3种数据集对所有森林类型估测效果呈现相同的趋势，估测精度从高到低都依次是多源数据、CHM和Sentinel-2数据。其中，多源数据对阔叶林蓄积量估测精度的作用最明显，比CHM数据的估测精度提高了6.7%，比Sentinel-2数据精度提高了12.7%。3种数据集的蓄积量预测值与实际值之间的关系图和残差（ 图4 图4 ）可以看出，对于多源数据的森林蓄积量估测，验证样本点多靠近1∶1线上，残差也从单一数据时（Sentinel-2 或CHM）的0—100 m 3 /ha变为0—75 m 3 /ha，且残差并没有随蓄积量的增大而出现明显波动；然而，使用单一数据都不可避免的出现残差值随实测蓄积量的增大而增大。此外，比较两种单一数据发现，总体上Sentinel-2数据的预测结果不如CHM数据，特别是对于阔叶林的预测和低蓄积量值的预测。 3 类型 6 数据 2 CHM 2 Sentinel-2 2 多源数据 RMSE/（m 3 /ha） rRMSE/% RMSE/（m 3 /ha） rRMSE/% RMSE/（m 3 /ha） rRMSE/% 阔叶 42.5 42.7 48.5 48.7 35.9 36.0 马尾松 47.8 34.0 51.8 36.9 42.3 30.1 杉木 50.8 23.1 60.4 27.5 39.0 17.8 1;2;3;4; 不同坡向的森林蓄积量估测结果（ 表5 表5 ）表明，CHM数据集和多源数据集对蓄积量估测的精度受坡向的影响。在朝阳面（阳坡、半阳坡）两个数据集的估测误差都小于背阴面（半阴坡、阴坡），其中CHM数据阳坡比阴坡精度高出5.5%，多源数据则高出6.9%。这与ZY-3立体像对提取DSM时不同坡向的精度差异有关（ Rahlf等，2014 Rahlf等，2014 ）。而在Sentinel-2数据集中则没有这种规律，阳坡比阴坡仅仅高出0.3%，这是因为前期对Sentinel-2影像做过地形校正，很大程度上消除了不同坡向引起的光谱值差异。3种数据的估测结果表明对于包含CHM的数据集，进行坡向分层是必要的；而在Sentinel-2数据中，除了半阳坡，其他各个坡向之间的估测精度没有明显差异。此外，将两种数据联合后，各个坡向的估测误差都显著下降，且在坡向上体现出差异，说明多源数据进行坡向分层建模有利于提高森林蓄积量的估测效果。 3 坡向 6 数据 2 CHM 2 Sentinel-2 2 多源数据 RMSE/（m 3 /ha） rRMSE/% RMSE/（m 3 /ha） rRMSE/% RMSE/（m 3 /ha） rRMSE/% 阳坡 39.7 28.2 55.8 39.7 31.0 22.1 半阳坡 42.4 33.0 38.5 30.0 33.8 26.3 半阴坡 51.8 34.8 57.6 38.8 46.2 31.1 阴坡 52.6 33.7 62.4 40.0 45.7 29.2 分析3种数据集对不同森林蓄积量等级的预测效果（ 表6 表6 ）可以看出：除区间200—250 m 3 /ha外，CHM的估测误差都小于Sentinel-2数据，但两种数据都不可避免的在高值区域出现了较高的估测误差；多源数据则展现出良好的估测结果，除100—200 m 3 /ha的区间外，多源数据都有最低的估测误差。此外CHM数据和Sentinel-2数据的估测误差都呈现出两头大中间小的特点，但这种现象在多源数据中则不明显，各个蓄积量区间的RMSE趋于相等，说明数据联合有利于改善高值低估和低值高估的现象。比如，在0—50 m 3 /ha和50—100 m 3 /ha的蓄积量区间中，多源数据的估测误差比Sentinel-2数据分别降低了18.0 m 3 /ha和16.4 m 3 /ha，在大于250 m 3 /ha的蓄积量区间也降低了近21.8 m 3 /ha，rRMSE只有14.1%，估测精度得到大幅度提升。 2 蓄积量/（m 3 /ha） 2 CHM 2 Sentinel-2 2 多源数据 RMSE/（m 3 /ha） rRMSE/% RMSE/（m 3 /ha） rRMSE/% RMSE/（m 3 /ha） rRMSE/% <50 43.8 144.2 58.1 191.1 40.1 132.1 50—100 44.4 57.0 51.0 65.4 34.6 44.4 100—150 25.6 20.2 36.9 29.1 32.5 25.7 150—200 34.4 20.3 52.0 30.7 38.5 22.7 200—250 67.3 31.2 62.0 28.8 51.4 23.9 >250 58.7 20.2 62.9 21.6 41.1 14.1 4　讨论 1 4.1　多源数据联合能提高森林蓄积量估测精度 2 Sentinel-2等多光谱数据一般只能提供林分冠层水平信息，无法获得垂直高度信息。因树高与森林蓄积量密切相关，基于光谱数据的蓄积量估测普遍不如基于高度数据的估测精度（ Feng等，2017 Feng等，2017 ）。如利用Landsat TM数据估测雪松蓄积量时，当雪松蓄积量高于200 m 3 /ha时，饱和问题使得拟合效果大幅度下降（ Lee和Phua，2010 Lee和Phua，2010 ），在亚马逊热带雨林的森林生物量和中国亚热带森林生物量估测时遇到同样的问题（ Lu，2005 Lu，2005 ； Zhao等，2016 Zhao等，2016 ； Gao等，2018 Gao等，2018 ），利用其他单一数据源，如RapidEye（ Feng等，2017 Feng等，2017 ）、SAR数据（Dos Reis等，2019）也不可避免的因饱和问题影响估测结果。本研究基于每个树种的Sentinel-2数据中出现明显的光谱饱和现象（ 图5 图5 ），但在蓄积量估测结果的散点图中（ 图4（b） 图4（b） ），没有出现明显的饱和现象。原因在于使用的贝叶斯模型对于蓄积量的估测考虑了森林类型，同时本研究的样本又具有特殊性，即高森林蓄积量样本（ > 250 m 3 /ha）全部为杉木（ 图5 图5 ），且杉木的样本也基本分布在高蓄积量区域，而杉木在250 m 3 /ha之前的光谱还未出现光谱饱和，因此高值区域的估测很容易分布在大于250 m 3 /ha的区域，使得原来在高值区域的蓄积量饱和现象不明显。另外，马尾松和阔叶在蓄积量大于200 m 3 /ha已经出现光谱饱和，但由于应用了分层建模，这两种树种的光谱饱和不会对杉木的估测造成影响，从侧面也说明了分层建模在减少饱和问题上的作用。 基于树高变量的森林蓄积量估测，不存在数据饱和问题，但树高变量与蓄积量之间的关系受其他因素，如树种、年龄、立地条件等的影响。如何从不同维度获取森林的信息，充分利用多源数据各自的优势，通过数据互补充分发挥数据潜力来弥补单一数据的不足，是提高森林蓄积量估测精度的重要手段（ Simard等，2011 Simard等，2011 ； Laurin等，2018 Laurin等，2018 ）。如在多光谱数据的基础上加入年龄、地形因子和SAR数据，桉树蓄积量的估测误差减少了近1/3（Dos Reis等，2019）。利用 Lidar的冠层穿透力结合Landsat OLI数据估测亚热带森林生物量，比单一数据明显提升了估测效果（ 曹林 等，2016 曹林 等，2016 ）。本研究进一步证实了以上结果，Sentinel-2数据联合CHM后，估测精度提高近10%。 4.2　多源数据对于高值和低值区域估测精度的提升 2 在蓄积量估测研究中，蓄积量的低值高估和高值低估常常是困扰研究者的难题（ Lu等，2016 Lu等，2016 ； Gao等，2018 Gao等，2018 ）。本研究发现Sentinel-2光谱和CHM数据的结合可以有效降低这种现象对估测精度的影响。本研究中总样本的平均蓄积量值为141.6 m 3 /ha，而不同蓄积量范围的估测精度验证发现，当蓄积量位于100—200 m 3 /ha时，多源数据的估测精度不如CHM，特别是100—150 m 3 /ha时，反而降低了5.5%，但当蓄积量偏离平均值附近时，即在0—100 m 3 /ha和大于200 m 3 /ha时，多源数据的估测误差都小于其他两种数据，说明数据联合对于模型的改进主要集中在高值和低值区域，而对平均值附近的改进不明显。因为当响应变量和预测变量相关性不高时，通过预测变量得到的响应变量的估测值更容易趋近于平均值，这在3种数据方案结果中都明显存在，因此平均值附近的预测误差往往最低。虽然在平均值附近多源数据的精度有所下降，但不影响整体精度的提高，因为蓄积量估测的问题主要出现高值区和低值区，而数据联合方法在这两个区域很好地实现了估测误差的下降，说明通过数据联合来提高蓄积量估测精度是切实可行的方法（ Li等，2019 Li等，2019 ）。 4.3　分层贝叶斯建模的依据和优势 2 研究中基于森林类型分层是基于不同森林类型的蓄积量和树高的关系不同。另外在森林类型分层的基础上，发现CHM和蓄积量的关系与坡向相关，这种相关性可能来自3方面的原因，（1）树木生长受坡向的影响；（2）资源三号DSM与机载Lidar DTM之间位置有错位，使得二者相减沿错位移动方向产生坡向相关特征；（3）光照条件随坡向发生变化，导致点云匹配出现差异。本文的前期实验中，曾对只使用Sentinel-2数据的不同建模方法做比较，发现分坡向的建模结果并不明显优于不分坡向的建模结果，说明这种相关性更多的不在于坡向间树木生长的差异，可能在于CHM数据在不同坡向上生产的差异。 Rahlf等（2014） Rahlf等（2014） 发现基于ZY-3立体像对提取的DSM高程会受到坡向的影响，背光面往往出现DSM数值偏大的现象，而向阳面与实际地表的状况较为吻合。这种坡向引起的差异，说明分层建模的必要性。然而分层建模意味着更多的样本需求，本研究只有71个野外调查样地，因此利用贝叶斯分层建模既可以减少不同森林类型和坡向差异对整体样本的影响，又能有效利用整体样本，从而提高估测精度，减少低值高估和高值低估问题。从分树种以及分坡向精度验证结果也可以发现，不同树种和坡向的估测精度之间存在明显的差异，为进一步提高森林蓄积量估测精度提供了方向。 5　结论 1 森林蓄积量遥感估测中林分高度信息比光谱信息更有优势，CHM数据的估测结果明显优于Sentinel-2光谱数据，而通过联合光谱信息和CHM数据的方法是提高森林蓄积量估测精度的有效手段。通过本研究的实验结果，得到以下结论： （1）3种数据集中，多源数据在蓄积量估测上具有明显的优势，其总体估测结果和分层结果评价都优于单独使用Sentinel-2数据集或CHM数据集。基于多源数据的森林蓄积量估测误差比Sentinel-2数据降低了13.6 m 3 /ha，比CHM降低了7.4 m 3 /ha； （2）多源数据可以有效改进高值低估和低值高估问题。多源数据在0—50 m 3 /ha区域的估测误差比单独使用Sentinel-2数据估测误差降低约1/4；在大于250 m 3 /ha区域的估测误差为比单独使用Sentinel-2或CHM数据估测误差也降低了近1/3。因此多源数据对于提高森林蓄积量估测精度特别是高值和低值区域具有显著作用； （3）在亚热带地区使用多源遥感数据进行森林蓄积量估测的结果表明，当样本数不足或数据随某个属性存在明显差异时，使用贝叶斯分层建模方法能发挥不同数据的优势，从而提高遥感手段估测森林蓄积量的精度。 参考文献（References） Babcock C , Finley A O , Andersen H E , Pattison R , Cook B D , Morton D C , Alonzo M , Nelson R , Gregoire T , Ene L , Gobakken T and Næsset E . 2018 . Geostatistical estimation of forest biomass in interior Alaska combining Landsat-derived tree cover, sampled airborne lidar and field observations . Remote Sensing of Environment , 212 : 212 - 230 [ DOI: 10.1016/j.rse.2018.04.044 http://dx.doi.org/10.1016/j.rse.2018.04.044 ] Geostatistical estimation of forest biomass in interior Alaska combining Landsat-derived tree cover, sampled airborne lidar and field observations Barrett F , McRoberts R E , Tomppo E , Cienciala E and Waser L T . 2016 . A questionnaire-based review of the operational use of remotely sensed data by national forest inventories . Remote Sensing of Environment , 174 : 279 - 289 [ DOI: 10.1016/j.rse.2015.08.029 http://dx.doi.org/10.1016/j.rse.2015.08.029 ] A questionnaire-based review of the operational use of remotely sensed data by national forest inventories Bates D , Mächler M , Bolker B and Walker S . 2015 . Fitting linear mixed-effects models using lme4 . Journal of Statistical Software , 67 ( 1 ): 1 - 48 [ DOI: 10.18637/jss.v067.i01 http://dx.doi.org/10.18637/jss.v067.i01 ] Fitting linear mixed-effects models using lme4 Brodu N . 2017 . Super-resolving multiresolution images with band-independent geometry of multispectral pixels . IEEE Transactions on Geoscience and Remote Sensing , 55 ( 8 ): 4610 - 4617 [ DOI: 10.1109/TGRS.2017.2694881 http://dx.doi.org/10.1109/TGRS.2017.2694881 ] Super-resolving multiresolution images with band-independent geometry of multispectral pixels Bürkner P C . 2017 . brms: Bayesian Regression Models using ‘Stan’ . R package version 1. 9 .0, URL [EB/OL]. https://CRAN.R-project.org/package=brms https://CRAN.R-project.org/package=brms URL Cao L , Xu T , Shen X and She G H . 2016 . Mapping biomass by integrating Landsat OLI and airborne Lidar transect data in subtropical forests . Journal of Remote Sensing , 20 ( 4 ): 665 - 678 Mapping biomass by integrating Landsat OLI and airborne Lidar transect data in subtropical forests 曹林 , 徐婷 , 申鑫 , 佘光辉 . 2016 . 集成Landsat OLI和机载Lidar条带数据的亚热带森林生物量制图 . 遥感学报 , 20 ( 4 ): 665- 678 [ DOI: 10.11834/jrs.20165210 http://dx.doi.org/10.11834/jrs.20165210 ] 集成Landsat OLI和机载Lidar条带数据的亚热带森林生物量制图 Chrysafis I , Mallinis G , Siachalou S and Patias P . 2017 . Assessing the relationships between growing stock volume and Sentinel-2 imagery in a Mediterranean forest ecosystem . Remote Sensing Letters , 8 ( 6 ): 508 - 517 [ DOI: 10.1080/2150704X.2017.1295479 http://dx.doi.org/10.1080/2150704X.2017.1295479 ] Assessing the relationships between growing stock volume and Sentinel-2 imagery in a Mediterranean forest ecosystem Dos Reis A A , Franklin S E , De Mello J M and Junior F W A . 2019 . Volume estimation in a Eucalyptus plantation using multi-source remote sensing and digital terrain data: a case study in Minas Gerais State, Brazil . International Journal of Remote Sensing , 40 ( 7 ): 2683 - 2702 [ DOI: 10.1080/01431161.2018.1530808 http://dx.doi.org/10.1080/01431161.2018.1530808 ] Volume estimation in a Eucalyptus plantation using multi-source remote sensing and digital terrain data: a case study in Minas Gerais State, Brazil Feng Y Y , Lu D S , Chen Q , Keller M , Moran E , Dos-Santos M N , Bolfe E L and Batistella M . 2017 . Examining effective use of data sources and modeling algorithms for improving biomass estimation in a moist tropical forest of the Brazilian Amazon . International Journal of Digital Earth , 10 ( 10 ): 996 - 1016 [ DOI: 10.1080/17538947.2017.1301581 http://dx.doi.org/10.1080/17538947.2017.1301581 ] Examining effective use of data sources and modeling algorithms for improving biomass estimation in a moist tropical forest of the Brazilian Amazon Gao Y K , Lu D S , Li G Y , Wang G X , Chen Q , Liu L J and Li D Q . 2018 . Comparative analysis of modeling algorithms for forest aboveground biomass estimation in a subtropical region . Remote Sensing , 10 ( 4 ): 627 [ DOI: 10.3390/rs10040627 http://dx.doi.org/10.3390/rs10040627 ] Comparative analysis of modeling algorithms for forest aboveground biomass estimation in a subtropical region Gao Y N and Zhang W C . 2008 . Comparison test and research progress of topographic correction on remotely sensed data . Geographical Research , 27 ( 2 ): 467 - 477 Comparison test and research progress of topographic correction on remotely sensed data 高永年 , 张万昌 . 2008 . 遥感影像地形校正研究进展及其比较实验 . 地理研究 , 27 ( 2 ): 467- 477 [ DOI: 10.3321/j.issn:1000-0585.2008.02.024 http://dx.doi.org/10.3321/j.issn:1000-0585.2008.02.024 ] 遥感影像地形校正研究进展及其比较实验 Jayathunga S , Owari T and Tsuyuki S . 2018 . The use of fixed–wing UAV photogrammetry with LiDAR DTM to estimate merchantable volume and carbon stock in living biomass over a mixed conifer–broadleaf forest . International Journal of Applied Earth Observation and Geoinformation , 73 : 767 - 777 [ DOI: 10.1016/j.jag.2018.08.017 http://dx.doi.org/10.1016/j.jag.2018.08.017 ] The use of fixed–wing UAV photogrammetry with LiDAR DTM to estimate merchantable volume and carbon stock in living biomass over a mixed conifer–broadleaf forest Junttila V , Kauranne T , Finley A O and Bradford J B . 2015 . Linear models for airborne-laser-scanning-based operational forest inventory with small field sample size and highly correlated LiDAR data . IEEE Transactions on Geoscience and Remote Sensing , 53 ( 10 ): 5600 - 5612 [ DOI: 10.1109/TGRS.2015.2425916 http://dx.doi.org/10.1109/TGRS.2015.2425916 ] Linear models for airborne-laser-scanning-based operational forest inventory with small field sample size and highly correlated LiDAR data Junttila V and Laine M . 2017 . Bayesian principal component regression model with spatial effects for forest inventory variables under small field sample size . Remote Sensing of Environment , 192 : 45 - 57 [ DOI: 10.1016/j.rse.2017.01.035 http://dx.doi.org/10.1016/j.rse.2017.01.035 ] Bayesian principal component regression model with spatial effects for forest inventory variables under small field sample size Laurin G V , Balling J , Corona P , Mattioli W , Papale D , Puletti N , Rizzo M , Truckenbrodt J and Urban M . 2018 . Above-ground biomass prediction by Sentinel-1 multitemporal data in central Italy with integration of ALOS2 and Sentinel-2 data . Journal of Applied Remote Sensing , 12 ( 1 ): 016008 [ DOI: 10.1117/1.JRS.12.016008 http://dx.doi.org/10.1117/1.JRS.12.016008 ] Above-ground biomass prediction by Sentinel-1 multitemporal data in central Italy with integration of ALOS2 and Sentinel-2 data Lee J S and Phua M H . 2010 . Estimation of stand volume of conifer forest: A Bayesian approach based on satellite-based estimate and forest register data . Forest Science and Technology , 6 ( 1 ): 7 - 17 [ DOI: 10.1080/21580103.2010.9656352 http://dx.doi.org/10.1080/21580103.2010.9656352 ] Estimation of stand volume of conifer forest: A Bayesian approach based on satellite-based estimate and forest register data Lewandowski D , Kurowicka D and Joe H . 2009 . Generating random correlation matrices based on vines and extended onion method . Journal of Multivariate Analysis , 100 ( 9 ): 1989 - 2001 [ DOI: 10.1016/j.jmva.2009.04.008 http://dx.doi.org/10.1016/j.jmva.2009.04.008 ] Generating random correlation matrices based on vines and extended onion method Li D R , Wang C W , Hu Y M and Liu S G . 2012 . General review on remote sensing-based biomass estimation . Geomatics and Information Science of Wuhan University , 37 ( 6 ): 631 - 635 General review on remote sensing-based biomass estimation 李德仁 , 王长委 , 胡月明 , 刘曙光 . 2012 . 遥感技术估算森林生物量的研究进展 . 武汉大学学报(信息科学版) , 37 ( 6 ): 631- 635 [ DOI: 10.13203/j.whugis2012.06.001 http://dx.doi.org/10.13203/j.whugis2012.06.001 ] 遥感技术估算森林生物量的研究进展 Li G Y , Xie Z L , Jiang X D , Lu D S and Chen E X . 2019 . Integration of ZiYuan-3 multispectral and stereo data for modeling aboveground biomass of Larch plantations in North China . Remote Sensing , 11 ( 19 ): 2328 [ DOI: 10.3390/rs11192328 http://dx.doi.org/10.3390/rs11192328 ] Integration of ZiYuan-3 multispectral and stereo data for modeling aboveground biomass of Larch plantations in North China Li W , Niu Z , Wang C , Gao S , Feng Q and Chen H Y . 2015 . Forest above-ground biomass estimation at plot and tree levels using airborne LiDAR data . Journal of Remote Sensing , 19 ( 4 ): 669 - 679 Forest above-ground biomass estimation at plot and tree levels using airborne LiDAR data 李旺 , 牛铮 , 王成 , 高帅 , 冯琦 , 陈瀚阅 . 2015 . 机载LiDAR数据估算样地和单木尺度森林地上生物量 . 遥感学报 , 19 ( 4 ): 669- 679 [ DOI: 10.11834/jrs.20154116 http://dx.doi.org/10.11834/jrs.20154116 ] 机载LiDAR数据估算样地和单木尺度森林地上生物量 Lu D . 2005 . Aboveground biomass estimation using Landsat TM data in the Brazilian Amazon . International Journal of Remote Sensing , 26 ( 12 ): 2509 - 2525 [ DOI: 10.1080/01431160500142145 http://dx.doi.org/10.1080/01431160500142145 ] Aboveground biomass estimation using Landsat TM data in the Brazilian Amazon Lu D S , Chen Q , Wang G X , Liu L J , Li G Y and Moran E . 2016 . A survey of remote sensing-based aboveground biomass estimation methods in forest ecosystems . International Journal of Digital Earth , 9 ( 1 ): 63 - 105 [ DOI: 10.1080/17538947.2014.990526 http://dx.doi.org/10.1080/17538947.2014.990526 ] A survey of remote sensing-based aboveground biomass estimation methods in forest ecosystems Luo H M , Chen E X , Li Z Y and Cao C X . 2011 . Forest above ground biomass estimation methodology based on polarization coherence tomography . Journal of Remote Sensing , 15 ( 6 ): 1138 - 1155 Forest above ground biomass estimation methodology based on polarization coherence tomography 罗环敏 , 陈尔学 , 李增元 , 曹春香 . 2011 . 森林地上生物量的极化相干层析估计方法 . 遥感学报 , 15 ( 6 ): 1138- 1155 [ DOI: 10.11834/jrs.20110391 http://dx.doi.org/10.11834/jrs.20110391 ] 森林地上生物量的极化相干层析估计方法 Meng S L , Pang Y , Zhang Z J , Li Z Y , Wang X Q and Li S M . 2017 . Estimation of aboveground biomass in a temperate forest using texture information from WorldView-2 . Journal of Remote Sensing , 21 ( 5 ): 812 - 824 Estimation of aboveground biomass in a temperate forest using texture information from WorldView-2 蒙诗栎 , 庞勇 , 张钟军 , 李增元 , 王雪琼 , 李世明 . 2017 . WorldView-2纹理的森林地上生物量反演 . 遥感学报 , 21 ( 5 ): 812- 824 [ DOI: 10.11834/jrs.20176083 http://dx.doi.org/10.11834/jrs.20176083 ] WorldView-2纹理的森林地上生物量反演 Mura M , Bottalico F , Giannetti F , Bertani R , Giannini R , Mancini M , Orlandini S , Travaglini D and Chirici G . 2018 . Exploiting the capabilities of the Sentinel-2 multi spectral instrument for predicting growing stock volume in forest ecosystems . International Journal of Applied Earth Observation and Geoinformation , 66 : 126 - 134 [ DOI: 10.1016/j.jag.2017.11.013 http://dx.doi.org/10.1016/j.jag.2017.11.013 ] Exploiting the capabilities of the Sentinel-2 multi spectral instrument for predicting growing stock volume in forest ecosystems Nuthammachot N , Askar A , Stratoulias D and Wicaksono P . 2022 . Combined use of Sentinel-1 and Sentinel-2 data for improving above-ground biomass estimation . Geocarto International , 37 ( 2 ): 366 - 376 [ DOI: 10.1080/10106049.2020.1726507 http://dx.doi.org/10.1080/10106049.2020.1726507 ] Combined use of Sentinel-1 and Sentinel-2 data for improving above-ground biomass estimation Rahlf J , Breidenbach J , Solberg S , Næsset E and Astrup R . 2014 . Comparison of four types of 3D data for timber volume estimation . Remote Sensing of Environment , 155 : 325 - 333 [ DOI: 10.1016/j.rse.2014.08.036 http://dx.doi.org/10.1016/j.rse.2014.08.036 ] Comparison of four types of 3D data for timber volume estimation Simard M , Pinto N , Fisher J B and Baccini A . 2011 . Mapping forest canopy height globally with spaceborne Lidar . Journal of Geophysical Research : Biogeosciences , 116 ( G4 ): G 04021 [ DOI: 10.1029/2011JG001708 http://dx.doi.org/10.1029/2011JG001708 ] Mapping forest canopy height globally with spaceborne Lidar Ver Planck N R , Finley A O , Kershaw Jr J A , Weiskittel A R and Kress M C . 2018 . Hierarchical Bayesian models for small area estimation of forest variables using LiDAR . Remote Sensing of Environment , 204 : 287 - 295 [ DOI: 10.1016/j.rse.2017.10.024 http://dx.doi.org/10.1016/j.rse.2017.10.024 ] Hierarchical Bayesian models for small area estimation of forest variables using LiDAR Xie Z L , Chen Y L , Lu D S , Li G Y and Chen E X . 2019 . Classification of land cover, forest, and tree species classes with ZiYuan-3 multispectral and stereo data . Remote Sensing , 11 ( 2 ): 164 [ DOI: 10.3390/rs11020164 http://dx.doi.org/10.3390/rs11020164 ] Classification of land cover, forest, and tree species classes with ZiYuan-3 multispectral and stereo data Zhao P P , Lu D S , Wang G X , Wu C P , Huang Y J and Yu S Q . 2016 . Examining spectral reflectance saturation in Landsat imagery and corresponding solutions to improve forest aboveground biomass estimation . Remote Sensing , 8 ( 6 ): 469 [ DOI: 10.3390/rs8060469 http://dx.doi.org/10.3390/rs8060469 ] Examining spectral reflectance saturation in Landsat imagery and corresponding solutions to improve forest aboveground biomass estimation 10.11834/jrs.20221545 1007-4619（2022）03-0468-12 林文科，陆亚刚，蒋先蝶，李桂英，李登秋，陆灯盛.2022.协同多源遥感数据的北亚热带森林蓄积量贝叶斯分层估测.遥感学报，26（3）： 468-479 林文科，陆亚刚，蒋先蝶，李桂英，李登秋，陆灯盛.2022.协同多源遥感数据的北亚热带森林蓄积量贝叶斯分层估测.遥感学报，26（3）： 468-479 Lin W K，Lu Y G，Jiang X D，Li G Y，Li D Q and Lu D S. 2022. Modeling forest growing stock volume in a north subtropical region using the hierarchical Bayesian approach based on multi-source data. National Remote Sensing Bulletin， 26（3）：468-479 Lin W K，Lu Y G，Jiang X D，Li G Y，Li D Q and Lu D S. 2022. Modeling forest growing stock volume in a north subtropical region using the hierarchical Bayesian approach based on multi-source data. National Remote Sensing Bulletin， 26（3）：468-479 国家自然科学基金（编号：32171787）；国家重点研发计划（编号：2017YFD0600900） Supported by National Natural Science Foundation of China (No.32171787) ;National Key Research and Development Program of China (No.2017YFD0600900) 2021-08-09 2022-03-07 2022-03-28 2021-09-25 遥感学报 3 26