北京地区春季木本植物花粉起始期预报模型研究

王春玲，叶彩华，姜 江

（北京市气象服务中心，北京 100089）

0 引言

花粉是有花植物的雄性生殖细胞，按照传播方式可以分为虫媒花粉和气传花粉。气传花粉相比虫媒花粉具有数量大、体积小、质量轻、外壁光滑易于被风吹散的特点，易进入人类呼吸器官，造成过敏等一系列疾病[1-3]。每年春季，北京地区各类木本植物先后进入开花期，并且大部分属于气传花粉[4-6]。在全球气候变暖的背景下，由于每年气候的差异性增大，导致各类植物花粉期年际变化较大，部分年份能相差20天以上，这给开展花粉观测研究以及公众气象服务带来了困难。

以积温为基础的各类物候模型是开展植物花期预报的主要手段，国内外很多学者开展了相关研究。Chuin等[7]基于法国南部城市蒙彼利埃1974—1992年的花粉观测数据，对比了4种主要物候模型对不同花粉植物花粉期的模拟效果，认为需要根据不同植物的特点选择适合的物候模型。Chung等[8]基于改进的积温模型，利用美国华盛顿地区1991—2010年的樱花始花期历史数据建立积温物候模型，模拟了2020—2080年IPCC 2种排放模式下(A1B、A2)樱花始花期的变化。国内方面，仲舒颖等[9]利用中国42个站点4种观赏植物1962—2013年的始花期数据建立积温物候模型，所建模型能够准确模拟多个站点和年份的始花期、末花期变化，均方根误差为4～6天。陆佩玲[10]利用北京地区1950—2001年4种木本植物（山桃、杏树、紫丁香、刺槐）的物候观测数据，建立了基于温度的非线性植物开花期预报模型，模拟结果较好，均方根误差为3～6天。

北京市全年气传花粉含量呈现出春季、秋季2个高峰期，其中春季主要是木本植物花粉。笔者以北京市春季4类主要花粉植物（柏科、杨柳科、松科、榆科）为研究对象，利用5个站点的花粉浓度观测数据确定每种植物的开花期，基于4种主流物候模型（Spring Warming模型、Alternating模型、Sequential模型、Parallel模型），对不同植物种类的花粉起始期进行模拟研究，确定每种植物的最适物候模型以及最佳模型参数组合，以期实现对主要植物花粉起始期的精确预报，为花粉观测研究、花粉公众气象服务等提供依据。

1 资料与方法

1.1 研究区概况和数据来源

1.1.1 研究区概况 北京市位于中国华北地区北部，地形复杂，三面环山，东南部为平原，属于典型的北温带半湿润大陆性季风气候，夏季高温多雨，冬季寒冷干燥。得益于良好的气候条件以及植树造林等城市绿化工作，北京植物资源丰富，截至2018年底，北京市城市绿化覆盖率达到48.44%，人均公共绿地面积达到16.3 m2。丰富的植被资源也带来了种类多样的花粉，欧阳志云等[4]研究认为，北京城区（五环内）共有19科32属99种植物会产生致敏性花粉。其中春季花粉主要来源于柏科、松属、杨柳科、榆科等木本植物[5-6]。

1.1.2 数据来源 本研究的花粉浓度数据来源于北京市气象局与北京同仁医院2012年联合开展的花粉监测项目，在全市总共设立了12个花粉观测点。根据每年的气候情况（通常通过2月平均气温确定）确定当年的花粉观测开始日期，一般比当年的花粉季起始日早7～10天左右。开始观测后，在每天的20:00在观测装置内放置花粉采样片，通过重力沉降法，即利用空气中花粉颗粒本身的重力自然沉降到载玻片上进行观测取样，载玻片每24 h更换一次，经碱性复红液染色后在显微镜下由专业人员根据花粉颗粒的形态特征进行人工分类和计数，计数结果换算成单位截面积(1000 mm2)下的花粉数量。在开展花粉读数时，由于需要通过人工开展花粉分类和读数，受限于观测人员的经验和精力等原因，无法实现精确分类，只能将花粉种类精确到科。

笔者选取北京地区5个花粉观测站点作为研究站点，分别属于昌平、海淀、丰台、石景山、朝阳5个行政区（图1），这5个站点均位于北京平原人口稠密区域，具备较好的地理、气候一致性。根据文献[11-12]，当某种花粉累计量达到其全年总量的2.5%时，该植物花粉期开始，到达全年花粉含量的97.5%时，该植物花粉期结束。笔者依据此方法统计了2012—2019年花粉浓度观测数据，确定了每年各科植物的花粉期起始日期，并得到一共140个有效花粉期数据（榆科32个、柏科33个、杨柳科35个、松科40个），笔者选取距离花粉观测站点最近的5个气象观测站2011年秋季—2019年的日平均气温资料，用于计算相应的积温数据（表1）。以上数据经过北京市气象局信息中心初步质量控制，部分缺测数据采用插值方法补齐。

表1 北京花粉监测站点资料说明（2012—2019）

图1 北京地区花粉监测站分布图

1.2 研究方法

采用积温物候模型的方法研究模拟植物花粉期。1735年，法国学者德列奥米尔最早发现积温对于植物生命周期的影响，1928年，前苏联科学家李森科发现冬季低温刺激可以促进部分作物开花。此后积温逐渐成为植物物候预测的主要方式，并逐渐发展为2类模型。一类模型只考虑春季积温，例如Spring Warming模型[13-14]，另一类则同时考虑积温和前一年冬天低温刺激对植物开花期的促进作用（春化作用）[15-17]，例如Alternating 模型[18]、Sequential模 型[19]、Parallel模型[20-21]。笔者用4种模型分别对北京地区春季主要花粉的起始期进行模拟。

针对每一科植物，通过计算其从前一年9月1日直至花粉期开始期间经历的低温刺激日数和积温总量，分别用4种模型建模，并通过机器学习的方法对积温阈值(Tb)、低温刺激阈值(To)等参数进行调参优化，使计算值达到最优。

1.2.1 Spring Warming模型 Spring Warming模型认为植物的生长发育与能量累积（积温）呈线性关系[22]，当积温达到一定的量，植物即进入某一发育期，模型计算如式(1)～(2)。

式中，y为开花期，t0为积温计算开始日期，GDDforc为到达开花期需要的积温阈值，xt为日平均气温，Tb为植物发育临界温度，F(xt)为当日有效积温。

1.2.2 Alternating模型 Alternating模型通过累积日数计算植物经历的低温刺激量[22]，当累积低温刺激日数达到一定阈值后植物才能进入下一个发育周期，模型计算如式(3)～(7)。

式中，y为开花期，t1为低温刺激起始日期，t2为低温刺激完成日期，GDDchill为低温刺激完成的阈值，a、b为无量纲参数，Sc为累积低温刺激日数，R(xt)表示当日低温刺激，xt为日平均气温，F(xt)为当日有效积温，Tb为植物发育临界温度。

1.2.3 Sequential模型 Sequential模型认为植物只有经过一定量的低温刺激后才能进入下一个发育周期[22]，模型计算如式(8)～(11)。

式中，y为开花期，t1为低温刺激起始日期，t2为低温刺激完成日期，GDDchill为低温刺激完成的阈值，GDDforc为到达开花期需要的积温阈值，R(xt)表示当日低温刺激，To为最优低温刺激阈值(-3.4～10.4℃)，xt为日平均气温，F(xt)为当日有效积温，Tb为植物发育临界温度。

1.2.4 Parallel模型 Parallel模型假定植物的低温刺激和积温促进可以独立运行[22]，当低温刺激阈值没有达到时，积温对植物发育的促进效率受到负面的影响，模型计算如式(12)～(15)。式中，y为开花期，t1为低温刺激起始日期，GDDforc为到达开花期需要的积温阈值，R(xt)表示当日低温刺激，To为最优低温刺激阈值，GDDchill为低温刺激完成的阈值，Km为模型调节参数(0.1～0.001)，xt为日平均气温，Tb为植物发育临界温度，F(xt)为当日有效积温。

1.2.5 参数率定及模型检验方法 使用偶数年资料进行模型参数率定，通过最小二乘法寻找最佳匹配参数，即当观测值和模拟值误差的平方和最小时，判断参数为最优。由于使用的4组模型中，参数最多达到4个，在限定区间后，参数组合仍然超过4.5亿组，如果通过遍历试参法，运算耗时巨大，因此使用同类研究常用的模拟退火算法进行调参[23-25]。其原理是模拟物理学中固体物质退火过程的随机寻优算法，可以以一定概率接受非最优解，从而避免陷入局部最优解，并最终达到全局最优解。当获得最优参数组合后，采用偶数年资料进行模型内部检验，采用奇数年的资料对模型进行外部检验，模型检验选择均方根误差(RMSE)和决定系数(R2)，均方根误差表示模拟值与实际值的偏离程度，值越小表示模拟结果越好。R2表示模拟值与实际值的拟合优度，介于0～1之间，值越大表示模拟结果越优，计算如式(16)～(17)。式中，n为样本量，xi为模拟值，yi为实测值，为实测值的平均值。

2 结果与分析

2.1 模型参数率定

为了提高运算效率，在参照前人研究的基础上[22,26]，对部分参数进行了简化，假定Spring Warming模型中积温起始计算日期为每年的1月1日，假定Alternating、Sequential、Parallel模型低温刺激起始计算日期为前一年的9月1日。同时对参数范围进行了限定，参数范围及率定结果见表2。

表2 参数率定结果

2.2 模型模拟结果评价

采用偶数年资料率定模型参数并进行内部验证，采用奇数年资料对模型进行外部验证。从模拟结果（表3）可以看到，榆科植物在内部检验阶段，4个模型的模拟结果相差不大，RMSE均在1.7～1.9天之间，R2介于89%～93%之间，没有明显的差别；但是在外部检验阶段，Alternating模型表现出明显的优异性，其RMSE只有1.6天，R2为83.1%，而其他模型RMSE均大于3天，R2小于60%。对于柏科植物，在内部验证阶段，Spring Warming模型结果较差，RMSE为2.8天，R2为69.7%。其他3个模型RMSE相差不大，介于2.2～2.3天，但是Alternating模型R2明显较优，为84.2%，在外部检验中，Alternating、Sequential模型明显优于其他2个模型，RMSE为2.3～2.4天，R2均为89.4%。对于杨柳科植物，在内部验证阶段，4个模型表现相差不多，RMSE介于1.8～1.9天之间，R2介于85%～87%；在外部验证阶段，Sequential模型明显优于其他3个模型，RMSE只有1.7天，R2为83.7%。对于松科植物，在内部检验阶段，Alternating模型表现较好，RMSE为1.3天，R2为90.6%；在外部验证阶段，Alternating、Sequential、Parallel模型表现相差不大，RMSE介于1.3～1.4天，R2介于95%～96.1%。

表3 模型结果检验

图2～5分别是榆科、柏科、杨柳科和松科植物所有年份的花粉起始期观测值与模拟值散点图。从拟合结果和表3的检验结果看，认为榆科、柏科、松科植物拟合最优模型为Alternating模型，R2分别为92.1%、91.2%、95.0%；杨柳科拟合最优模型为Sequential模型，R2为92.8%。

图2 榆科植物花粉起始期模拟值与实测值散点图(2012—2019)

图3 柏科植物花粉起始期模拟值与实测值散点图(2012—2019)

图4 杨柳科植物花粉起始期模拟值与实测值散点图(2012—2019)

图5 松科植物花粉起始期模拟值与实测值散点图(2012—2019)

3 结论与讨论

笔者以北京城区榆科、柏科、杨柳科以及松科4类春季主要木本花粉植物为研究对象，利用2012—2019年的花粉浓度观测数据计算出每种植物的花粉起始期。进一步建立了每种植物的花粉期起始日预报模型，并对模型进行了横向对比，得到以下结论：（1）总体来看，本研究4类木本植物花粉期起始日均受到前一年冬季低温刺激的影响，即前一年冬季温度越低，当年开花期所需要的积温越少，同等气候条件下花粉期起始日会提前。（2）经过对Spring Warming、Alternating、Sequential、Parallel 4种模型进行横向对比发现，榆科、柏科、松科植物的最适模型为Alternating模型，RMSE介于1～3天，R2超过91%。杨柳科植物的最适模型为Sequential模型，RMSE小于2天，R2为92.8%，说明基于积温原理的物候模型，对北京地区春季木本植物的花粉期起始日有较好的模拟预测能力，具备应用于业务保障的潜力，同时对同类研究具备参考价值。（3）由于北京市气象局的花粉浓度观测项目起始于2012年，目前积累的数据样本偏少，会对模型结果产生不可避免的影响，未来待资料积累较多后，可对模型进行再次率定和检验。另外，由于花粉期结束日不仅受积温影响，同时受到风、降水等要素的影响，影响过程较为复杂，本研究中暂未涉及花粉结束期，未来可以考虑将其他气象要素如降水、日照等引入模型，从而实现对花粉期结束日期的模拟预报。