结合机载激光扫描与树种组成预测实现加拿大西部魁北克老龄北方混交林林龄（≤250年）精准制图

本研究创新性地融合机载激光扫描（ALS）与树种组成预测数据，成功构建了可精准预测北方老龄混交林林龄（最高达250年）的遥感模型。该模型在多种分辨率下均展现优异性能（R2达0.574–0.759），显著提升了对老龄林演替多样性的表征能力，为 boreal forest 生态系统保护与可持续管理提供了关键技术支撑。

　　

引言

准确绘制森林年龄图谱对实现可持续森林管理具有基础性意义。林龄（这里定义为上次林分替代干扰后的时间）是评估木材生产、碳固存、生态系统服务及高保护价值区域的关键属性。然而，野外获取林龄信息耗时费力，且传统方法（如树木年轮取样）存在选择性偏差和幼树抑制期遗漏等问题。遥感技术通过直接照片判读或与遥感变量结合建模，在林龄制图中扮演重要角色。机载激光扫描（ALS）等新技术的发展，更通过提供森林三维特征信息持续提升模型精度。

北方森林的林龄预测面临特殊挑战：树木相对矮小，达到峰值高度的时长远短于南方树种；从早期到晚演替树种的渐进变化虽是林龄的广义指标，但当林分进入耐荫树种主导阶段后，组成更多取决于次生干扰和局部环境条件而非火后时间。某些北方树种（如黑云杉、欧洲赤松）既可出现在早期也可存在于晚期演替阶段，其主导的森林随时间推移可能几乎无树种组成变化。次生干扰还可导致老龄林出现局部逆行演替（从晚期退回到早期物种）。耐荫树种对幼树抑制的高耐受性、树干/根际在土壤有机层的埋藏以及晚演替北方树种直径与年龄关系的弱相关性，严重限制了基于野外调查数据估算老龄林林龄的可靠性。因此，无论是通过遥感还是照片判读，北方森林的林龄制图往往存在较低的最大年龄阈值（100-150年），这通常对应这些生态系统老龄阶段的开始。超过此阈值，上述限制导致林龄预测出现显著噪声。

这种有限的时效深度给可持续森林管理带来实质性问题。工业化前，由于长火周期或地表火主导，北方森林相当大部分由老龄林构成。这些森林提供特定生境（通常与特化物种相关）和独特的生态系统服务（如高碳储存）。当前全球最后的大片老龄林主要位于北方生物群落，使其在全球保护政策中占据关键地位。其次，老龄林的结构、组成、生境和生态系统服务在超过老龄阈值后仍持续变化。同一地区两个不同年龄、干扰历史或环境特征的老龄林可提供截然不同的生境和服务。最后，采伐活动存在偏向性，最具生产力和可及性（位于独特环境条件）的老龄林更易被砍伐。现有年龄数据的有限精度意味着无法保证残留的老龄林在其功能和生境上能完全代表其工业化前的范围。因此，改进老龄林内部的林龄制图对于确保保护政策的成功至关重要。

加拿大东部北方混交林是准确估算林龄面临困难的典型例子，尤其是在可持续森林管理和老龄林保护的背景下。严重野火是主要原生干扰，但火周期通常长于300年，产生了以数百年历史森林为主导的景观。因此，森林演替是根据上次火后的时间定义的。此过程通常遵循从火后更新的不耐荫硬木树种到耐荫软木树种的过渡。然而，在干旱土壤或贫瘠有机土壤上，由于其防火和耐荫适应性，黑云杉可能在整个演替过程中保持优势种地位。类似地，东部云杉蚜虫爆发等次生干扰可导致逆行演替，在不耐荫硬木在较大林窗中定殖。云杉蚜虫爆发是研究区主要的次生干扰，复发周期约35年，可导致耐荫软木树种大量死亡，形成由未干扰森林连接的、不同大小的林窗镶嵌体。该地区标准的森林野外调查方法无法提供关于最老森林年龄的可靠信息，有时存在显著低估。欧洲人在此地的永久定居较晚（19世纪下半叶之后），限制了关于火历史和森林景观年龄结构的信息。与此同时，工业化的森林管理持续大幅减少该地区老龄林的面积和连通性。因此，亟需开发能够覆盖广泛演替梯度的、精确且广泛适用的工具来绘制这些森林的年龄图谱。

材料与方法

研究区域

本研究在杜帕奎湖研究与教学林（LDRTF）进行，该林地位于加拿大魁北克省西部东南部北方森林，属于魁北克和安大略省的北部克莱带的一部分，这个大的生理区域由威斯康星冰期后冰川湖巴洛和奥吉布韦最大扩展期留下的湖相沉积物形成。区域气候为寒冷大陆性气候，年均温约1.0°C，年均总降水量约985毫米，其中30%发生在生长季。降雪占年总降水量的25%。大部分液态降水发生在生长季，但6、7月的蒸发可能限制植物生长。平均无霜期为64天，但任何时间都可能发生冰冻。香脂冷杉是优势种，黑云杉、白云杉和白桦（纸皮桦）也是森林的重要组成部分。先锋树种（如北美短叶松、颤杨和纸皮桦）在林火干扰后占据大面积区域。基于生态系统的该区域森林管理规范将任何超过80年的森林视为“老龄林”。

该研究区域因其独特的生态意义而被选中，在此区域内，起源于1717年至1944年间连续13次火灾的森林共存于一个集中的地理区域，同时自1920年代欧洲人永久定居以来基本未受采伐影响。对该地区火历史的详细了解是连续研究项目的成果，为东部加拿大北方混交林动态提供了独特视角。火灾地图由Dansereau和Bergeron通过几个步骤确定：首先，在171个样地（取芯725棵树）对优势树基部取芯至髓心，并采样或分析包含火疤痕的横截面以补充早期分析；建立了北美短叶松的主年表以实现死亡个体的交叉定年；利用这些信息和20世纪该地区拍摄的航空照片，通过照片判读确定了火边界和日期。相比之下，魁北克北方森林可靠的火灾记录仅能追溯到1890年，且各地区可用性不同；而该省的判读森林清查在80或100年（取决于林分）后停止估算年龄。在此背景下，该研究区域代表了确定我们更准确评估老龄林年龄能力的重要资源，其火后时间梯度范围独特。

地面调查数据

本研究使用了Maleki等人先前提供的数据集。该数据集包含横跨7场野火（发生于1760年至1944年：1760, 1797, 1823, 1847, 1870, 1916, 1944年）起源林分的432个野外样地（16 m × 16 m, 256 m²)组成的样带网络。这些样带于1990-1991年依据Dansereau和Bergeron的火灾图设立，以覆盖不同年龄的火灾，并于2009年重新测量。样地间的目标距离为50 m。在样地内，对树木进行取芯以确保火灾图与样地位置树木年龄的一致性。该时间序列已成为众多北方混交林火后动态研究的基础。这些样地的数据随后用于校准SORTIE-ND（一个基于个体的、空间显式的北方混交林动态模型）。SORTIE-ND模拟数据与野外数据的优异一致性表明，这些样地可靠地代表了北方混交林演替。因此，这些地点火历史的时间准确性已得到充分证实，使其成为检验我们利用遥感数据估算北方混交林年龄能力的理想地面真值数据。

虚拟样地定义

研究使用野外样地中心点定义不同大小的虚拟样地用于提取遥感数据：20×20 m (400 m²), 40×40 m (1600 m²), 70×70 m (4900 m²), 和 100×100 m (10000 m²)。最小面积400 m²代表加拿大森林调查常用标准样地面积，而最大面积1公顷 (10000 m²)已被证明能有效捕捉老龄北方森林的结构异质性。两个中间面积（1600 m² 和 4900 m²)旨在提供介于其间的梯度，以评估样地大小如何影响预测。虚拟样地的火灾日期为其中心野外样地的日期。由于某些样地邻近以及某些样带局部交叉，部分虚拟样地重叠。为确保相同大小的连续虚拟样地不重叠，定义了虚拟样地中心点之间的最小距离：400 m²样地为25 m，1600 m²样地为50 m，4900 和 10000 m²样地为100 m。无法满足这些条件的虚拟样地从数据集中移除。虚拟样地边界也可能与中心野外样地关联火灾年龄不同的火灾区域或非林地重叠。对于最老的火灾，其确切边界也存在不确定性。因此，移除了所有中心点距离另一火灾边界小于150 m的虚拟样地。同时移除了部分覆盖非林地区域的虚拟样地，以及中心点距离杜帕奎湖岸小于50 m的虚拟样地。这些方法确保每个虚拟样地覆盖的森林位于一个火灾多边形内，其日期与中心野外样地相对应。此阶段最终虚拟样地数量为：400 m²面积263个，1600 m²面积220个，4900 m²和10000 m²面积各60个。

ALS预测变量

ALS勘测于2017年由魁北克政府进行，点密度为2–4点/m²，勘测条件经过优化以实现最小遮挡和研究区全覆盖。基于这些ALS勘测，政府已创建多个衍生指标并以开放获取的栅格格式共享，而非原始点云。本研究使用了从原始ALS点云衍生的冠层高度模型（1 m像素大小）、地形湿度指数（2 m像素大小）和坡度（2 m像素大小）产品。选择关注开放获取产品是为了便于本研究开发模型的应用，因为这些是省内森林利益相关者使用的产品。针对每个虚拟样地大小，提取了23个ALS预测变量，这些变量对应已确立的林龄和结构指标、森林演替阶段指标，尤其适用于具有不同冠层复杂性和物种组成程度的森林。

树种组成预测变量

基于这些ALS数据，以及卫星影像和政府森林调查数据，魁北克政府还制作了魁北克森林的树木测量图，提供了林分结构（如立木密度、断面积、商业材积）信息，以及以商业材积相对丰度表示的树种组成信息。与使用更粗分类并在更大多边形中定义的判读森林图相比，该产品能更精细地预测树种组成。在研究区，树种组成多边形面积为0.55 ± 0.59 ha (5^th百分位数: 0.02 ha, 95^th百分位数: 1.59 ha)。此处仅选择了与混交北方林相关的主要商业树种数值：香脂冷杉、北美乔柏、黑云杉、白云杉、北美短叶松、白桦以及杨树（颤杨和香脂杨归为同一类别）。由于树木测量图也以对应同质森林结构和组成斑块边界的多边形形式存在，将每个物种组成值（如香脂冷杉在商业材积中的百分比）转换为10×10 m分辨率的栅格格式。随后为每个虚拟样地提取了10个树种组成预测变量，并根据每个像素覆盖面积进行加权。

分析

对所有虚拟样地分辨率，林龄建模遵循相同步骤。首先，样带采样结构意味着相邻虚拟样地间存在潜在空间互相关。为确保虚拟样地空间独立性，对每个分辨率的预测变量进行了主成分分析（PCA），使用FactoMineR R包中的PCA()函数。PCA允许将所有33个预测变量汇总在有限数量的轴上。对所有分辨率，前两个PCA轴解释了超过50%的方差（400 m²: 56.11%；1600 m²: 57.15%；4900 m²: 59.12%；10000 m²: 61.67%）。随后使用spdep R包中的localmoran()函数计算了每个分辨率下所有虚拟样地第一和第二PCA轴值的局部莫兰指数（I），k=1。任何具有显著局部莫兰指数（p < 0.05）的样地从数据集中移除。因此移除了400 m²分辨率的11个虚拟样地，1600 m²分辨率的15个，4900 m²分辨率的5个，10000 m²分辨率的2个。移除的样地数量有限，这与先前排除过于接近样地的过程、限制连续样地数量的样带采样性质，以及由局部环境条件和次生干扰空间异质性塑造的北方混交林高度变异性是一致的。

第二步，使用增强回归树（BRT）模型通过遥感预测变量预测上次火后时间（观测林龄）。模型使用gbm R包中的gbm()函数创建，该函数应用增强回归树来捕捉预测变量间的复杂非线性关系和相互作用。BRT模型配置使用用户参数的默认值，包括5000棵树、0.001的学习率、3的树复杂度（交互深度），并在每次迭代中使用0.5的袋装比例进行随机抽样。模型测试采用“留一样带出”方法，即每次重复模型时排除一个样带（按Maleki等人数据集定义）的样地，此时将其用作独立测试样地。虚拟样地被分为31个样带（400 m²分辨率：每样带2至22个虚拟样地，平均10.3个；1600 m²分辨率：每样带1至19个样地，平均8.4个；4900 m²分辨率：26个样带，每样带1至6个虚拟样地，平均2.6个；10000 m²分辨率：26个样带，每样带1至6个虚拟样地，平均2.8个）。排除空间依赖虚拟样地的上一步确保测试虚拟样地与训练虚拟样地的邻近性不会影响预测。通过比较测试虚拟样地的观测和预测林龄评估模型性能。计算了每个模型的确定系数（R²）、均方根误差（RMSE）和偏差。使用线性回归的截距和斜率评估线性模型与理想1:1线的偏差，从而深入了解预测中的系统偏差或尺度错误。

对每个虚拟样地分辨率，拟合了三种不同模型：一个包含所有预测变量的完整模型，一个排除所有树种组成预测变量的模型，以及一个仅使用完整模型中相对重要性>1%的预测变量的简约模型。物种组成数据特定于魁北克省开发的方法学。无物种预测变量的模型允许评估，如果物种数据不可用，仅ALS预测变量在不同分辨率下提供准确年龄预测的程度。简约模型旨在确定更简单模型的可靠性以促进其大规模应用。总体上，用此方法定义了12个不同模型（4种虚拟样地分辨率 × 3种模型类型 [完整, 无物种, 简约]）。

最后，对LDRTF保护区内每个模型在其各自分辨率下进行了全域预测，该区域由起源于不同火灾的森林高度集中而定义。通过移除所有高度回返95^th百分位数低于7 m（该研究区成熟林阈值）的像素排除了非林地区域。将这些预测图与LDRTF火灾图和Maltman等人制作的加拿大森林年龄图进行视觉比较。此步骤不仅允许对我们的模型进行视觉比较，还能看到我们的方法与基于树木年代学数据判读火后时间制图以及加拿大尺度森林年龄模型相比的贡献。所有分析均使用R统计软件4.4.1版进行。

结果

模型性能与比较

完整模型的性能在不同样地大小间相对均匀，从400 m²到10000 m²略有增加，但4900 m²模型性能略有下降。R²值范围从0.566 (400 m²) 到 0.759 (10000 m²)，1600 m² (0.621) 和 4900 m² (0.579) 为中间值。RMSE呈现不规则模式，值范围从30.825 (10000 m²) 到 40.289 (4900 m²)，其他分辨率为中间值。偏差接近零，值范围从-1.093 (4900 m²) 到 3.933 (400 m²)，但在样带尺度上存在强烈变化。完整模型的斜率范围从1.03 (400 m²) 到 1.14 (10000 m²)，截距范围从-23.89 (10000 m²) 到 -0.29 (400 m²)，意味着观测与预测间关系在不同分辨率下具有稳定性。简约模型的性能略优于完整模型，意味着一些引起微小预测噪声的变量被移除。无物种预测变量的模型性能低于完整和简约模型，但随着分辨率降低而增加（R²从0.251到0.592，RMSE从52.58到39.300，对应400 m²和10000 m²分辨率），但4900 m²分辨率这些指标属于最低之列。偏差范围在-0.342 (400 m²) 和 -1.54 (4900 m²) 之间，但同样在样带间存在强烈变化。此处，与完整和简约模型相反，降低分辨率减少了空间变异性和噪声，导致更可靠的森林年龄预测。然而，最低分辨率（10000 m²)的模型性能仅略优于最高分辨率（400 m²)的完整和简约模型。

属性重要性

完整模型的预测变量重要性随样地大小而变化，反映了结构和物种组成预测变量相对影响的转变。总体而言，树种预测变量（按物种或按演替组）在所有模型中扮演关键角色。对于400 m²分辨率模型，树种属性占预测变量平均相对重要性的94.3%。这意味着在此尺度上，不同树种的组合而非三维预测变量是估算森林年龄最有效的。在10000 m²分辨率，物种属性仍占平均相对属性重要性的66.5%，而ALS属性占剩余的33.5%。这突显了在最低分辨率下，ALS属性对精确年龄预测的重要性增加，但树种预测变量（特别是区分早晚期演替物种）对精确预测仍起关键作用。在树种预测变量中，有六个在所有模型中的平均相对重要性超过1%：北美乔柏比例、香脂冷杉比例、纸皮桦比例、北美短叶松比例、早期演替物种比例和晚期演替物种比例。对于源自ALS产品的预测变量，粗糙度以及7至12 m间回返比例是唯一在所有模型中平均相对重要性超过1%的预测变量。

空间预测图

完整模型的预测在不同空间分辨率下显得相对一致，分辨率的逐步降低仅导致高分辨率预测的平均化。在所有分辨率下，预测与研究区的火灾图接近。然而，观察到空间上 distinct 的区域，其预测年龄高于火灾图，可能表明火灾严重度较低的区域；或预测年龄较低，此处可能表明次生干扰的影响。Maltman等人定义的森林年龄模型预测研究区大部分区域年龄小于100年，且野外样地的预测年龄也显著低估。Maltman等人使用的方法显著低估了我们研究区域森林的年龄。与完整模型相比，对于无物种属性的模型，分辨率的降低限制了在高分辨率下观察到的局部噪声，因为年龄可能从一个像素到相邻像素突然变化。然而，这些模型似乎仅能提供对规则、同龄结构森林与不规则、异龄结构森林的广泛区分，很少有像素的预测超过200年。它们在识别和区分不同年龄的老龄林方面仍比Maltman等人的模型更有效。最后，简约模型的预测与完整模型接近，表明移除重要性最低的变量对模型预测影响很小。

讨论

本研究开发的模型通过结合描述冠层树种组成和三维结构的遥感预测变量，提供了林分年龄的可靠估计。当包含树种预测变量时，模型精度随分辨率增加保持稳定，但仅使用ALS预测变量时显著下降。这些结果表明，仅凭ALS数据，需要捕获更大表面积以减少空间异质性带来的噪声，从而有效识别不同年龄混交北方老龄林的典型结构模式。相反，树种组成属性的纳入将ALS数据置于特定的森林演替背景中。

遥感数据允许精确估算混交北方林中上次严重干扰后时间（最高达250年）

我们所有的完整和简约模型，无论其分辨率如何，都提供了高达250年的混交北方林林龄的可靠预测。据我们所知，这些遥感模型可能是唯一能够覆盖老龄北方林最广年龄范围的模型。由于超过这些阈值后精度下降，林龄预测的遥感模型通常停止在100至150年的最大年龄。考虑到老龄北方林随时间在结构、组成、生境和生态系统服务方面的重大变化，能够区分早期和晚期老龄林对于可持续森林管理至关重要。关于最老森林年龄的信息常常缺失，部分原因是，从林业角度来看，关于已超过采伐成熟度森林年龄的详细知识被认为不相关。北方树木早期的树高生长平台期、演替过程中有时有限的树种变化以及随时间积累的次生自然干扰进一步复杂化了森林年龄估算。Maltman等人预测与我们模型预测在研究区的鲜明对比证明了这一挑战。Maltman等人的加拿大森林年龄预测模型使用年龄/优势高关系的渐近线来定义林龄。虽然这种方法可能适用于树高可超过50 m的温带森林，但在年龄/高度渐近线通常在100年届满前就达到的北方森林中，可能导致显著的年龄低估。使用遥感数据的大尺度制图方法日益被批评为掩盖了地方特异性，从而可能提供适得其反的信息。我们的结果支持这一论点，表明适应特定生态系统特征的方法可以准确识别具有高生态价值的森林，而这些森林会被更全局的模型遗漏。

所有分辨率的完整和简约模型在整个演替过程中表现出显著的稳定性。它们使用上次火后时间作为林龄指标进行训练，但在ALS数据捕获（2017年）与火灾（此处介于74至258年之间）之间可能发生了若干次生干扰，有时严重程度显著。间隙动态的开始在加拿大东部北方混交林演替中相对较早发生（上次火后70-150年）。云杉蚜虫爆发是研究区主要的次生干扰，复发周期约35年，20世纪初和末至少发生过两次严重爆发。严重爆发通常在完整森林矩阵中产生许多有限面积（<200 m²)的林窗，尽管偶尔也能观察到大规模连续死亡（>1 ha）。其他干扰，如风倒或冰暴，可能已在研究区发生，但其影响通常比昆虫爆发更轻且更随机。大林窗的局部开放可能有利于不耐荫物种（如白桦或颤杨）的更新，扰乱预期的演替动态，并形成小群丛的同龄群体。北方树木的低矮也限制了这些森林可以呈现的结构异质性梯度。北方混交林结构、组成和动态的这种固有复杂性可能解释了模型预测中的部分误差，一些样地显示出与其年龄通常预期特征局部不同的特征。这也解释了无物种属性模型的有限准确性，因为用ALS数据描述三维特征不足以在没有物种数据来情境化此结构的情况下区分不同年龄的森林。然而，完整和简约模型的总体结果显示出无论森林年龄如何，其预测总体稳健。模型的空间预测也与研究区的火灾图一致，尽管可以观察到一些局部不一致。该火灾图也有其局限性，例如在清晰识别火边界、火灾严重度变化以及后续次生干扰的影响方面。因此，我们的结果可用于细化长期干扰历史的制图，例如识别火灾严重度较低的区域或发生严重次生干扰导致逆行演替的区域。

本研究中定义的模型并未覆盖整个森林演替范围，仅适用于至少70年生的森林。此选择是由于研究区附近缺乏具有可靠地面真值数据的更近期火灾。火烧区域内严重度的强烈变化以及随之而来的更新和生长动态变化，限制了在没有可靠野外数据情况下定义虚拟样地的使用。火后严重度导致的结构差异在干扰后的几十年尤其显著，随后随森林演替逐渐减小。起源于皆伐的森林不能用作幼龄林的训练数据，因为其生长和物种组成与火后林存在实质性差异。然而，魁北克的森林调查在80年前是可靠的，只有超过此年龄后其可靠性才急剧下降。尽管如此，在此使用Maltman等人的森林年龄图来识别70年以下的森林似乎并不合适，因为Maleki等人样地网络中的大多数老龄林被识别为低于此年龄。因此，我们的方法学应被视为对当前使用判读方法的补充，使得将年龄估算扩展到最老的森林成为可能。此外，原始树种组成数据被聚合到多边形中（研究区平均面积0.55 ± 0.59 ha）。由于这些变量在我们模型中的决定性作用，预测因此将在这些多边形的区域内部分聚合，取决于其分辨率。然而，它们在加拿大东部北方景观和自然干扰机制尺度上的有限面积对预测的林龄结构空间格局仅有微小影响。

景观尺度的详细树种组成数据对于情境化混交北方林中描述森林结构的ALS属性是必要的

描述树种组成的属性在所有分辨率下都扮演关键角色，尽管其相对影响随分辨率降低而减小。火灾后，加拿大东部北方混交林可呈现广泛的树种组成，取决于环境条件和干扰前林分特征。即使到老龄林阶段，林分组成也可能因干扰历史和环境条件而高度可变。然而，这些森林遵循一个共同的全局演替模式，即不耐荫物种被耐荫物种逐步替代。树种组成属性在高分辨率下特别有用，因为样地面积太小无法捕获足够的结构变异性。在这种情况下，基于ALS属性且无物种属性的模型表现不佳。每个物种相对于其他物种的相对存在从而实现了高分辨率下森林年龄的有效预测。相反，降低分辨率通过更好地表征林分的三维结构增加了ALS属性的重要性，从而限制了物种和地形的局部效应。然而，树种组成属性在区分具有相似三维结构但处于不同演替阶段的样地方面仍扮演重要角色。有趣的是，完整和简约模型的精度随分辨率增加保持稳定，而无树种组成预测变量的模型精度显著下降。这些结果强调，对于ALS衍生产品（如冠层高度和数字地形模型），需要低分辨率来捕捉北方混交林的三维异质性，尽管存在一些限制。因此，描述森林演替的补充产品（本例中也部分基于ALS数据）对于情境化森林三维结构至关重要，尤其是在高分辨率下。

用于创建我们模型的地面真值数据也在预测质量中扮演重要角色。大多数使用遥感数据预测北方森林年龄的模型使用政府调查样地作为地面真值数据。虽然可在大面积区域大量获取，但这些数据通常不适合研究老龄林。年龄数据通常通过在高1.3米处对少数优势树取芯获得