姜 元，茹晓雅，罗 琦，李美荣，王 钊，王景红，冯 浩，张 东，苏宝峰，于 强，何建强

·农业生物环境与能源工程·

基于Copula函数的陕西苹果晚霜冻特征分析

姜 元1,2，茹晓雅1,2，罗 琦1,2，李美荣3，王 钊3，王景红3，冯 浩2,4，张 东5，苏宝峰6，于 强4，何建强1,2※

（1. 西北农林科技大学旱区农业水土工程教育部重点实验室，杨凌 712100；2. 西北农林科技大学中国旱区节水农业研究院，杨凌 712100；3. 陕西省农业遥感与经济作物气象服务中心，西安 710015；4. 中国科学院水利部水土保持研究所黄土高原土壤侵蚀与旱地农业国家重点实验室，杨凌 712100；5. 西北农林科技大学园艺学院，杨凌 712100；6. 西北农林科技大学农业农村部农业物联网重点实验室，杨凌 712100）

为探讨利用Copula函数对苹果晚霜冻进行特征分析的适用性，该研究首先基于陕西省苹果主产区7个气象站点1971－2018年的逐日最低气温（min）数据集，提取出晚霜冻事件的历时和强度两个特征变量。然后，基于6种不同的Copula函数构建晚霜冻特征变量的联合分布，并进行拟合优度评价。最后，利用优选的Copula函数分析晚霜冻发生的概率及重现期。结果表明：陕西苹果产区各站点1971－2018年受晚霜冻的影响在空间分布上由东南向西北方向加重，各站点晚霜冻的历时和强度之间均具有显著的正相关关系。当晚霜冻强度和历时增大时，其联合累积概率也相应增大，且增大趋势变缓。各站点联合重现期代表的“或”事件比同现重现期所代表的“且”事件更容易发生。当单变量重现期取值较小时，可将联合重现期和同现重现期视为单变量重现期的两种极端情况，对其实际范围进行估计。总体而言，陕西苹果产区各站点发生长历时且高强度晚霜冻事件的概率较小，但是位于延安果区的站点相较于其他果区站点更容易发生高强度或长历时的晚霜冻事件，以及高强度长历时同时发生的晚霜冻事件，需要重点加以关注。该研究可为陕西苹果产区应对晚霜冻灾害提供理论依据。

气象；灾害；晚霜冻；Copula函数；频率分析；重现期；苹果

0 引 言

苹果在中国的栽培面积和产量均位居世界首位[1-2]。苹果的广泛栽培使得苹果生产更易受到各种气象灾害的影响，晚霜冻是中国北方苹果产区面临的首要极端气象灾害[3-5]，会对苹果产量和质量有严重影响[4,6]，给当地苹果产业造成重大经济损失[4]。目前，全球气候变暖的形势日趋严峻[7]，冬季增温幅度也更为明显[8]，暖冬导致果树春季萌芽及开花提前[9-10]，使其遭受晚霜冻的风险日益增加。陕西地区的大陆性季风气候特征十分明显，全年约有40%的强降温天气发生在3－4月份。然而，该时期正处于苹果开花期，频繁的冷空气过程导致晚霜冻灾害风险较高[11-13]。

目前，关于苹果晚霜冻的研究多集中在建立苹果晚霜冻指数等级并探究晚霜冻事件时空变化特征和成灾因子。例如，李健等[14]基于气象站最低气温构建了苹果果树花期冻害预警指标，并分析了主果区花期冻害的时空变化特征。屈振江等[13]构建了冻害风险灾损率对苹果花期冻害风险进行了评估。柏秦凤等[15]基于花期冻害风险指数对苹果花期冻害进行降尺度风险区划。王景红等[12]利用人工气候箱模拟低温并结合历史灾情调查对苹果花期霜冻指标进行了修订优化。张琪等[3]基于物候模型探究了苹果霜冻发生特征。邱美娟等[1]结合晚霜冻气象指标，对苹果花期晚霜冻气候风险进行了评估。王明昌等[11]采用线性倾向法和风险指数法分析了陕西省礼泉和旬邑两县苹果花期冻害风险情况。

上述研究大多基于晚霜冻指标进行苹果晚霜冻灾害风险的分析，但是评估霜冻等气象灾害事件风险的较佳方法是通过详细了解以持续时间和低温强度为特征的霜冻事件的发生频率[16]。李美荣等[17]建立了苹果花期严重冻害最低气温的概率分布模型及重现期预测，但是只使用了单变量频率分析方法且仅考虑了影响霜冻的一个变量，无法捕捉到极端气象灾害问题的复杂性。

另一种评估气象灾害事件风险更为有效的的方法是建立多元联合分布[18]。然而，许多传统多元联合分布要求各变量的边缘分布属于同一类型[19]。Copula函数作为一种构造灵活的模型不限制各变量边缘分布的选择，能够通过边缘分布和相关性结构两部分来构建多维联合分布，形式灵活多样，被广泛应用于降雨频率[20]、干旱频率[21]和洪水频率[22]等水文和灾害事件分析。Chatrabgoun等[16]首次采用Copula函数对葡萄园霜冻进行了概率评估，通过对霜冻持续时间和严重程度这2个特征变量建立联合分布，探究了位于伊朗马来尔地区葡萄园霜冻现象的风险和影响，表明Copula函数可以作为构建霜冻多维特征变量联合分布函数的有效工具。然而，将Copula函数应用于苹果晚霜冻特征分析的研究目前还鲜有报道。

本研究以陕西省苹果主产区7个地面气象观测站1971－2018年的气象数据集为基础数据，采用Copula函数对陕西苹果产区的晚霜冻灾害进行概率评估，研究目标包括：1）逐个气象站提取晚霜冻事件发生的连续天数和日最低温度最小值的绝对值分别作为表征晚霜冻历时和强度的2个特征变量；2）对2个特征变量的单变量边缘分布进行拟合，并分析变量间相关性；3）基于6种Copula函数，建立2个特征变量的多元联合分布模型并进行拟合优度检验；4）分析计算联合累积概率和重现期，从而对陕西苹果产区各站点的晚霜冻灾害风险进行评估。本研究的结果有望在站点尺度上为陕西苹果晚霜冻灾害防御和苹果产业生态管理提供科学依据。

1 材料和方法

1.1 研究区概况与数据来源

陕西苹果主产区位于黄土高原地区，是符合苹果生长气象指标的优质产区之一。根据气候特点及物候差异，将陕西苹果主产区划分为延安果区、渭北东部果区、渭北西部果区和关中西部果区[23]（图1）。本研究使用的1971－2018年陕西省气象站点的逐日最低气温（min）数据来源于中国气象科学数据共享服务网（http://cdc.cma.gov.cn/）。使用该气象数据前需对各个站点进行检查，删除数据序列长度过短或缺失的站点，并对数据异常值进行处理。最终选择研究区域内的延安市宝塔区、洛川、长武、旬邑、铜川、白水和礼泉等7个气象站作为研究站点（图1）。由于陕西苹果花期物候观测数据不够完整，所以本文使用的苹果花期物候数据主要来自王润红等[24]基于物候模型对陕西省苹果花期进行模拟所得的数据。

1.2 研究方法

1.2.1 苹果晚霜冻事件识别

刘映宁等[26]根据苹果晚霜冻气象指标，结合灾情调查资料，发现苹果花期内不同程度的低温会对苹果生长发育产生不同的影响：日最低温度低于−4 ℃时，出现严重晚霜冻，中心花受冻率高达80%，减产30%以上；日最低温度在−2至−4 ℃时，出现中度晚霜冻，中心花受冻率达60%至80%，对产量、品质、商品率产生严重影响；日最低温度在0至−2 ℃时，出现轻度晚霜冻，中心花受冻率达30%至60%，对产量、品质、商品率产生明显影响。此外，相关文献也给出了苹果花期冻害的临界温度。李美荣等[26]将日最低气温min为0、−2、−4 ℃作为苹果在花期内不同时期受冻的临界温度；刘映宁等[26]提出把日最低气温分别低于0、−2、−4 ℃作为陕西苹果花期冻害农业保险的3个等级；王景红等[27]则采用0、−2、−4 ℃作为陕西苹果不同等级花期冻害的临界指标。因此，本研究将日最低气温min为0 作为发生苹果晚霜冻事件的临界温度，每个晚霜冻事件中日最低温度最小值的绝对值被作为该期间的晚霜冻强度。苹果晚霜冻不仅与低温强度有关，也与低温的持续时间相关，强度越大或持续时间越长受冻率越高[12]，因此将晚霜冻发生的连续天数定义为晚霜冻历时。

图1 陕西苹果产区亚区划分及气象站点分布

1.2.2 构建晚霜冻历时和晚霜冻强度的边缘分布函数

在基于Copula函数构建晚霜冻历时和强度联合分布之前，首先要确定晚霜冻历时和强度这2个特征变量的边缘分布函数，同时要考虑它们之间的相依性。本文选用正态分布、指数分布、伽马分布、威尔布分布、柯西分布、逻辑斯谛分布、对数正态分布7种概率分布函数构建这2个特征变量边缘函数。使用R语言的“fitdistrplus”包，利用极大似然法估计相关参数，对晚霜冻历时和强度这2个特征变量的边缘分布函数进行拟合，并通过Kolmogorov-Smirnov（K-S）检验选取最优的单变量边缘分布函数。采取3种常用的相关系数来进行2个晚霜冻特征变量间的相关性度量，包括皮尔逊积矩相关系数（Pearson product-moment correlation coefficient）、斯皮尔曼秩相关系数（Spearman’s rank correlation coefficient）、肯德尔秩相关系数（Kendall rank correlation coefficient）。

以陕西省礼泉站为例，因为使用逐日最低气温来识别晚霜冻事件，所以晚霜冻历时数据均为整数天，导致晚霜冻历时存在较多的相同值，因此晚霜冻历时的经验分布呈条带状（图2）。在对晚霜冻历时边缘分布函数进行拟合过程中，这种条带状的离散分布会导致样本个数明显少于实际值，进而影响拟合效果。参考韩会明等[28-29]引入干旱历时离散变量连续化处理方法，本研究对整数的晚霜冻历时加上[−0.5, 0.5]的均匀随机分布变量，从而使离散化的数据序列连续化作为处理后的经验分布。此外，Michele等[30]证明，对样本添加均匀分布的扰动不会改变原始的统计量信息。

图2 礼泉站晚霜冻历时边缘分布拟合情况

1.2.3 构建晚霜冻历时和晚霜冻强度的Copula联合分布函数

Copula函数可用于构造边缘分布不同的多个变量的联合分布函数。根据Sklar定理[31]，设=(<)和=(<)分别为特征变量晚霜冻历时和强度的边缘分布函数，其联合分布函数为(,)（式（1））。

(,)=(<,<)=[(),()]（1）

式中为特征变量取值，为Copula函数，为边缘分布函数。

确定单个变量的边缘分布函数后，要优选一种Copula函数连接单变量分布函数。本文选用了6种常用的二维Copula联合分布函数（表1），包括Archimedean Copula函数中的Gumbel Copula，Frank Copula，Clayton Copula以及Joe Copula，和椭圆Copula函数中最常见的Normal Copula和T Copula[16,20-21]。其中，参数估计均选用极大似然法进行估计，拟合优度的检验选用赤池信息准则（Akaike Information Criterion, AIC；式（2））和贝叶斯信息准则（Bayesian Information Criterion, BIC；式（3）），以赤池信息量AIC值和贝叶斯信息量BIC值最小作为拟合最优的检验标准。

AIC=−2ln()+2（2）

BIC=−2ln()+ln()·（3）

式中为参数的数量；为似然函数；为数据数量。

表1 本研究中所采用的6种二维Copula函数

注：为特征变量取值，为边缘分布函数，为关联参数，为生成函数自变量。

Notes:are characteristic variable values;,are the marginal distribution functions;denotes to the association parameters;is the generator function argument.

1.2.4 基于晚霜冻历时和强度的晚霜冻事件重现期确定

重现期作为水文分析中评估风险的一个通用标准，被广泛用于确定重大灾害事件发生的概率，极具现实意义和应用价值。晚霜冻是一种对苹果生长极具威胁性的气象灾害，采用重现期对晚霜冻事件进行评估，有助于了解晚霜冻事件的发生规模与发生频率之间的关系，对晚霜冻事件的发生进行预测，有助于陕西苹果产区更好地防御晚霜冻灾害。一般情况下，重现期的分析仅限于单变量，所谓单变量重现期（式（4））是指变量超过某一特定值出现一次的间隔时间，然而由于极端气象灾害问题的复杂性，导致单变量重现期分析结果具有一定的不准确性[17]。晚霜冻事件是晚霜冻历时和强度两个因素相互作用的结果，所以在晚霜冻事件特征分析中应该重点考虑重现期的联合属性和条件属性。对于多特征变量的联合分布而言，联合重现期（Joint return period；式（5））代表多个特征变量中某一个特征变量大于给定阈值出现一次的“或”事件的间隔时间，而同现重现期（Co-occurrence return period；式（6））代表多个特征变量中每一个特征变量均大于给定阈值出现一次的“且”事件的间隔时间。

传统基于单变量重现期的计算公式为

式中F()为变量的边缘分布函数；单变量同理。

二维联合重现期和同现重现期的计算公式分别为[32]：

2 结果与分析

2.1 陕西省苹果主产区晚霜冻事件特征分析

本研究利用陕西省苹果主产区各个气象站点数据，提取出1971－2018年间该区域苹果花期发生晚霜冻事件的晚霜冻历时和晚霜冻强度（表2），以及每次晚霜冻发生的起始和终止日期。该区域不同果区的站点发生晚霜冻的频率存在明显差异，其中位于延安和渭北西部果区的4个站点共发生晚霜冻频率高，共计266次，且晚霜冻的历时和强度值也相应较大；渭北东部和关中西部的站点发生晚霜冻共89次，频次相对较少，且历时和强度值也较小。在陕西苹果产区各气象站中，延安站不仅严重晚霜冻和总晚霜冻发生次数最多，而且晚霜冻强度平均值也最大；旬邑站的晚霜冻历时和强度的最大值均为最大。

基于晚霜冻历时的起止时间进一步分析晚霜冻历时过程（图3），可见总体上陕西省苹果主产区晚霜冻事件历时多分布于3月31日至4月16日。关中西部果区的站点（图3a）明显发生晚霜冻的次数少、历时短、强度低；渭北西部果区（图3d和图3e）和延安果区（图3f和图3g）的站点发生的晚霜冻次数多、历时长、强度高。陕西苹果主产区各站点在1980年晚霜冻最多发生达17次；在2013年发生晚霜冻共12次，为近10年晚霜冻发生次数最多的年份。

表2 陕西省苹果主产区1971－2018年晚霜冻事件统计

图3 陕西苹果产区7个气象站的晚霜冻历时分布

2.2 晚霜冻历时和强度边缘分布函数的拟合

本研究采用正态分布、指数分布、伽马分布、威布尔分布、柯西分布、逻辑斯谛（Logistic）分布、对数正态分布等7种概率分布，对陕西苹果产区晚霜冻历时和强度的边缘分布函数进行拟合，得到晚霜冻历时和强度边缘分布的K-S检验结果以及变量之间的相关系数（表3）。选取拟合度最高的概率分布类型（显著性水平=0.05），最终确定各站点晚霜冻历时边缘分布函数为对数正态分布。对于晚霜冻强度这一特征变量的边缘分布拟合，各站点情况不一，其中礼泉站、白水站和旬邑站为Logistic分布，铜川站、洛川站和延安站为指数分布，长武站为伽马分布。相关性检验结果表明各站点Pearson、Spearman和Kendall这3个相关系数均大于0.4且都通过了显著性水平=0.05的显著性检验，说明各站点的晚霜冻历时和强度具有一定的相关性，可基于Copula函数建立联合分布函数进行进一步分析。

表3 单变量边缘分布的K-S检验和变量间的相关系数

注：*表示通过了显著性水平= 0.05的显著性检验，**表示通过了显著性水平= 0.01的显著性检验。

Notes: * is statistically significance level of α= 0.05; ** is statistically extreme significant level of= 0.01.

2.3 基于Copula函数晚霜冻事件联合分布函数的建立

本研究基于拟合的晚霜冻历时和强度的边缘分布函数，建立了晚霜冻历时和强度之间的6种Copula分布函数，利用极大似然法对其中未知参数进行了估计。此外，还计算了相应的AIC、BIC值（式（2）和式（3）），对Copula分布函数拟合优度进行评价（图4）。根据AIC和BIC准则判断可知，Normal Copula函数对礼泉站和旬邑站晚霜冻事件的拟合效果最好；Clayton Copula函数对白水站和铜川站晚霜冻事件的拟合效果最好；对洛川站晚霜冻事件拟合最优的是Frank Copula函数；对延安站和长武站晚霜冻事件拟合优度最高的是Joe Copula函数。

图4 陕西苹果产区晚霜冻历时和强度的6种Copula分布函数拟合优度检验

2.4 晚霜冻历时和强度的联合概率

进一步进行本研究中各站点晚霜冻特征变量二维联合概率分布分析（图5），可知当晚霜冻强度和历时值增大时，其联合累积概率也相应增大。总体而言，各站点所在区域发生短历时低强度、短历时高强度、长历时低强度的晚霜冻事件的概率较大，而同时满足长历时和高强度的晚霜冻事件的发生概率较小。具体地，在礼泉站（图5a），晚霜冻强度低于3 ℃、晚霜冻历时在0.5至2 d时，其联合累积概率随特征变量的增大而迅速增大；在白水站（图5b），晚霜冻强度高于3 ℃、晚霜冻历时高于3 d时，其联合累积概率增大趋势明显变缓；在铜川站（图5c），晚霜冻强度不超过4℃、晚霜冻历时在0.5至3 d时，其联合累积概率随特征变量的增大而迅速增大；在旬邑站（图5d），晚霜冻强度不超过2 ℃、晚霜冻历时不超过4 d时，其联合累积概率随特征变量的增大而迅速增大；在长武站（图5e），晚霜冻强度超过1 ℃、晚霜冻历时超过2 d时，其联合累积概率增大趋势明显变缓；在洛川站（图5f），晚霜冻强度不超过2 ℃、晚霜冻历时不超过3 d时，其联合累积概率随特征变量的增大而迅速增大；在延安站（图5g），晚霜冻强度超过1 ℃、晚霜冻历时超过2 d时，其联合累积概率增大趋势明显变缓。

图5 陕西苹果产区7个气象站晚霜冻历时和强度联合概率分布

2.5 晚霜冻历时和强度的重现期

以等值线图的形式对晚霜冻事件的联合重现期进行描绘（图6）。根据晚霜冻两个特征变量值垂直相交的交点所在等值线可查出联合重现期，可发现随着晚霜冻特征变量不断增大，两变量的联合重现期也随之增大。根据等线图的疏密情况可明显看到，不同站点间晚霜冻事件的联合重现期存在一定的差异，各站点中礼泉站最不容易发生高强度或长历时的晚霜冻事件（图6a），而白水站和铜川站较不易发生高强度或长历时的晚霜冻事件（图6b和6c）。对于其他站点，若霜冻强度达到8 ℃或者晚霜冻历时达到8 d时，在旬邑站其联合重现期约为150 a（图6d），在长武站其联合重现期约为280 a（图6e），在洛川站其联合重现期约为150 a（图6f），在延安站其联合重现期约为50 a（图6g）。这说明相较于旬邑站、长武站和洛川站，延安站更容易受到高强度或长历时晚霜冻事件的威胁，也进一步表明联合重现期等值线图可以表征不同区域间的晚霜冻事件发生概率，能够服务于减灾预防工作。

通过单变量重现期（式（4））以及构建的单变量边缘分布的逆函数，分别计算当单变量重现期为3、5、10、20、50、100 a时，各站点晚霜冻历时和强度值，再计算联合重现期（式（5））和同现重现期（式（6））对应值（表4）。结果发现联合重现期和同现重现期均随着晚霜冻特征变量的增大而增大，这与等值线图（图6）所示的联合重现期变化规律相一致。在单变量取值同等增幅条件下，同现重现期的增幅要明显高于联合重现期，这说明该区域各站点联合重现期代表的“或”事件比同现重现期所代表的“且”事件更容易发生。然而，不同站点的重现期变化存在一定的差异，例如，当单变量重现期增幅相同时，延安站对应的同现重现期增幅明显大于旬邑站，这说明延安站更易受到长历时且高强度晚霜冻事件的威胁。各站的单变量重现期、联合重现期和同现重现期具有相同的大小顺序，单变量重现期总介于联合重现期和同现重现期之间，因此可以通过计算联合重现期和同现重现期来估计单变量重现期的范围，即可得到不同晚霜冻特征变量所代表的晚霜冻事件发生的频率，但是仅当单变量重现期较小时估计的范围才比较准确。

图6 陕西苹果产区7个气象站晚霜冻历时和强度联合重现期等值线

表4 陕西苹果产区7个气象站晚霜冻历时和强度的单变量重现期和联合分布重现期

续表

3 讨 论

在对极端气象灾害事件进行特征分析时，首先要明确气象灾害指标以及事件相关的特征变量，要能够对极端气象灾害事件进行综合性评价。以往对晚霜冻灾害的研究重点是探究低温强度的影响，然而张晓煜等[33]提出低温和低温持续时间是对霜冻危害程度起决定作用的关键因子；王景红等[12]表示苹果晚霜冻不仅与低温强度有关，也与低温持续时间密切相关，温度越低或持续时间越长，苹果受冻率越高；袁佰顺等[6]发现0 以下温度持续时间和低温强度都需作为判断霜冻受灾和受灾程度的因素。因此，本研究选取晚霜冻历时和晚霜冻强度作为特征变量，基于Copula函数建立特征变量的多元联合分布，克服了单变量频率分析方法不能捕捉到极端气象灾害复杂性的不足，为陕西苹果产区晚霜冻特征研究提供了新的思路和方法。

本研究通过对陕西苹果产区晚霜冻事件的特征分析可知，在空间上陕西苹果产区受晚霜冻的影响总体上呈现东南向西北加重的趋势，这与前人对陕西苹果晚霜冻的探究结果基本一致[13-14,17,25-26,34]。该空间分布趋势是由于延安果区和渭北果区所处的经纬度、高海拔及地势起伏大所导致的。李美荣等[17]对旬邑站进行单变量重现期计算，得到10 a重现期对应的晚霜冻强度为5.2 ℃，而本研究为3.81 ℃，这可能是由于研究所选的时间尺度不同、识别晚霜冻所用的气象指标不同，以及拟合分布类型不同所导致的。此外，本研究采用Copula函数通过对晚霜冻历时和强度这2个特征变量建立多元联合分布评估陕西苹果晚霜冻风险，相比于李美荣等[17]仅考虑低温强度而使用单变量频率的分析方法，更能反映晚霜冻事件受晚霜冻历时和强度两个特征变量共同影响的复杂性，从而为精准评估陕西苹果晚霜冻风险提供更为全面和科学的信息。

准确的苹果物候数据是评估苹果花期冻害风险的重要依据。然而，目前国内关于苹果物候观测数据的积累十分有限，存在时间序列短、缺测严重等问题[11]，并且现有的驱动物候模型模拟的气象数据也是有限的，从中国气象科学数据共享服务网所获得的气象数据基本只能到县域尺度，缺乏空间密集的气象观测数据，难以对苹果晚霜冻事件进行区域精准识别，所以未来研究要加强对陕西苹果产区苹果花期的多站点和时序观测，采用质量较高的气象数据产品来驱动物候模型进行模拟。苹果属于多年生作物，不同地区、不同树龄、不同品种的苹果对气候响应有所不同，然而目前晚霜冻指标没有进行针对性的细化，使得目前对晚霜冻事件的识别工作具有一定的局限性。本文仅对晚霜冻历时和强度这2个特征变量进行联合，后续研究可以增加对晚霜冻特征量的选取，但是多元联合分布函数也会随着变量维数增多变得愈发复杂。所以，苹果花期模拟、晚霜冻指标细化、晚霜冻事件精准识别、基于Copula函数的更高维特征变量分析等都是今后进行苹果晚霜冻灾害研究需要关注的重点方向。

4 结 论

本研究根据日最低气温阈值提取1971–2018年间陕西苹果产区各气象站晚霜冻事件的历时和强度，基于Copula函数建立了晚霜冻历时和强度的联合分布，并分析了研究区内晚霜冻事件的重现期概况，得到以下主要结论：

1）位于延安和渭北西部果区的站点共发生晚霜冻266次，是主要的晚霜冻灾害频发区，并且历时长、强度高；渭北东部和关中西部的站点发生晚霜冻共89次，频次相对较少，而且历时短、强度低。陕西苹果产区各站点受晚霜冻的影响在空间分布上呈由东南向西北加重的趋势。

2）陕西苹果产区各站点晚霜冻历时的最优边缘分布均为对数正态分布，晚霜冻强度的最优边缘分布随站点而各异。各站点晚霜冻的历时和强度之间的相关系数均大于0.4，具有显著的正相关关系。

3）本研究采用各站点拟合优度最高的Copula 函数建立陕西苹果产区各站点晚霜冻历时和强度的联合分布函数，分析晚霜冻特征变量的联合累积概率可知，同时满足长历时和高强度的晚霜冻事件的发生概率较小。

4）不同站点间晚霜冻事件的重现期存在一定的差异，若霜冻强度达到8℃或者晚霜冻历时达到8 d时，延安站联合重现期约为50 a，相较于其他站点百年以上的联合重现期，延安站更容易受到高强度或长历时，以及高强度且长历时晚霜冻事件的威胁。陕西苹果产区各站点联合重现期代表的“或”事件比同现重现期所代表的“且”事件更容易发生。当单变量重现期的值较小时，可将联合重现期和同现重现期可看作单变量重现期的两种极端情况，对单变量重现期的实际范围进行估计。

总之，本文基于Copula函数建立的联合分布模型能够同时描述晚霜冻事件的历时和强度之间的关系，可以利用重现期评估陕西苹果产区站点尺度未来晚霜冻发生的风险，从而为研究该区域各站点晚霜冻事件的发生规律，提高该区域各站点晚霜冻的预报能力提供新的科学方法和依据。

[1] 邱美娟，刘布春，刘园，等. 中国北方主产地苹果始花期模拟及晚霜冻风险评估[J]. 农业工程学报，2020，36(21)：154-163.

Qiu Meijuan, Liu Buchun, Liu Yuan, et al. Simulation of first flowering date for apple and risk assessment of late frost in main producing areas of northern China[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2020, 36(21): 154-163. (in Chinese with English abstract)

[2] 万炜，师纪博，刘忠，等. 栖霞市苹果园氮磷养分平衡及环境风险评价[J]. 农业工程学报，2020，36(4)：211-219.

Wan Wei, Shi Jibo, Liu Zhong, et al. Nitrogen and phosphorus nutrient balance and environmental risk assessment of apple orchard in Qixia city[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2020, 36(4): 211-219. (in Chinese with English abstract)

[3] 张琪，李荣平，张晓月，等. 绥中苹果霜冻灾害特征及其对气候变暖的适应分析[J]. 北方园艺，2016(13)：30-33.

Zhang Qi, Li Rongping, Zhang Xiaoyue, et al. Analysis of apple frost disaster characteristics and adapt to the climate warming in Suizhong[J]. Northern Horticulture, 2016(13): 30-33. (in Chinese with English abstract)

[4] 马延庆，刘长民，朱海利，等. 陕西咸阳渭北旱塬地区优质苹果基地生态气候特征分析[J]. 干旱地区农业研究，2008(1)：146-153.

Ma Yanqing, Liu Changmin, Zhu Haili, et al. The analysis of ecosystem characteristic of high quality apple producing base on Weibei arid plain region in Xianyang of Shaanxi Province[J]. Agricultural Research in the Arid Areas, 2008(1): 146-153. (in Chinese with English abstract)

[5] 李化龙，赵西社. 陕西黄土高原果业气候生态条件研究及应用[M]. 北京：气象出版社，2010.

[6] 袁佰顺，许彦平，姚晓红，等. 基于春季晚霜冻害的天水苹果气候风险区划[J]. 干旱区资源与环境，2015，29(2)：185-189.

Yuan Baishun, Xu Yanping, Yao Xiaohong, et al. The compartment of later frost climate risk for apple in Tianshui[J]. Journal of Arid Land Resources and Environment, 2015, 29(2): 185-189. (in Chinese with English abstract)

[7] Organization W M. WMO Provisional Statement on the State of the Global Climate in 2019[EB/OL]. http://www.cma.gov.cn/en2014/climate/featutes/202003/P020200312816904145935.pdf, 2022-6-23.

[8] Cook B I, Wolkovich E M, Parmesan C, Divergent responses to spring and winter warming drive community level flowering trends[J]. 2012, 109(23): 9000-9005.

[9] Ge Q, Wang H, Rutishauser T, et al. Phenological response to climate change in China: A meta-analysis[J]. Global Change Biology, 2015, 21(1): 265-274.

[10] 郭佳，张宝林，高聚林，等. 气候变化对中国农业气候资源及农业生产影响的研究进展[J]. 北方农业学报，2019，47(1)：105-113.

Guo Jia, Zhang Baolin, Gao Julin, et al. Advances on the impacts of climate change on agro-climatic resources and agricultural production in China[J]. Journal of Northern Agriculture, 2019, 47(1): 105-113. (in Chinese with English abstract)

[11] 王明昌，刘布春，刘园，等. 陕西苹果主产县花期冻害风险评估[J]. 中国农业气象，2020，41(6)：381-392.

Wang Mingchang, Liu Buchun, Liu Yuan, et al. Assessment on the freezing injury risk during apple flowering in Liquan and Xunyi[J]. Chinese Journal of Agrometeorolog, 2020, 41(6): 381-392. (in Chinese with English abstract)

[12] 王景红，刘璐，高峰，等. 陕西富士系苹果花期霜冻灾害气象指标的修订[J]. 中国农业气象，2015，36(1)：50-56.

Wang Jinghong, Liu Lu, Gao Feng, et al. Revision on meteorological indices of florescence frost disaster for Fuji apple in Shaanxi Province[J]. Chinese Journal of Agrometeorology, 2015, 36(1): 50-56. (in Chinese with English abstract)

[13] 屈振江，刘瑞芳，郭兆夏，等. 陕西省苹果花期冻害风险评估及预测技术研究[J]. 自然灾害学报，2013，22(1)：219-225.

Qu Zhenjiang, Liu Ruifang, Guo Zhaoxia, et al. Study of risk assessment and prediction of apple blooming freezing injury in Shaanxi Province[J]. Journal of Natural Disasters, 2013, 22(1): 219-225. (in Chinese with English abstract)

[14] 李健，刘映宁，李美荣，等. 陕西果树花期低温冻害特征及防御对策[J]. 气象科技，2008(3)：318-322.

Li Jian, Liu Yingning, Li Meirong, et al. Low temperature and freezing injury to fruit trees at bloom stage in Shaanxi and countermeasures[J]. Meteorological Science and Technology, 2008(3): 318-322. (in Chinese with English abstract)

[15] 柏秦凤，王景红，梁轶，等. 基于县域单元的降尺度苹果花期冻害风险区划[J]. 中国农学通报，2013，29(16)：153-158.

Bai Qinfeng, Wang Jinghong, Liang Yi, et al. The apple flowering frost damage risk zoning down-scaling based on county level[J]. Chinese Agricultural Science Bulletin, 2013, 29(16): 153-158. (in Chinese with English abstract)

[16] Chatrabgoun O, Karimi R, Daneshkhah A, et al. Copula-based probabilistic assessment of intensity and duration of cold episodes: A case study of Malayer vineyard region[J]. Agricultural and Forest Meteorology, 2020, 295: 108150.

[17] 李美荣，李星敏，柏秦凤，等. 苹果极端气象灾害气温极值的分布及重现期预测[J]. 干旱地区农业研究，2012，30(3)：257-261.

Li Meirong, Li Xingmin, Bai Qinfeng, et al. Probability distribution and recurrence interval prediction of the extreme value of air temperature in extreme meteorological disasters of apples[J]. Agricultural Research in the Arid Areas, 2012, 30(3): 257-261. (in Chinese with English abstract)

[18] Pouteau R, Rambal S, Ratte J P, et al. Downscaling MODIS-derived maps using GIS and boosted regression trees: The case of frost occurrence over the arid Andean highlands of Bolivia[J]. Remote Sensing of Environment, 2011, 115(1): 117-129.

[19] Kurowicka D, Cooke R M. Uncertainty Analysis with High Dimensional Dependence Modelling[M]. Chichester: John Wiley & Sons, 2006.

[20] 郭爱军，畅建霞，王义民，等. 近50年泾河流域降雨-径流关系变化及驱动因素定量分析[J].农业工程学报，2015，31(14)：165-171.

Guo Aijun, Chang Jianxia, Wang Yimin, et al. Variation characteristics of rainfall-runoff relationship and driving factors analysis in Jinghe river basin in nearly 50 years [J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2015, 31(14): 165-171. (in Chinese with English abstract)

[21] 王晓峰，张园，冯晓明，等. 基于游程理论和Copula函数的干旱特征分析及应用[J]. 农业工程学报，2017，33(10)：206-214.

Wang Xiaofeng, Zhang Yuan, Feng Xiaoming, et al. Analysis and application of drought characteristics based on run theory and Copula function[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2017, 33(10): 206-214. (in Chinese with English abstract)

[22] 蔺文慧，宋松柏. 基于混合Copula函数的洪水联合概率计算方法[J]. 水力发电学报，2021，40(12)：40-51.

Lin Wenhui, Song Songbai. Study on calculation method of flood joint probability based on mixed Copula function[J]. Journal of Hydroelectric Engineering, 2021, 40(12): 40-51. (in Chinese with English abstract)

[23] 邬定荣，霍治国，王培娟，等. 陕西苹果花期机理性预报模型的适用性评价[J]. 应用气象学报，2019，30(5)：555-564.

Wu Dingrong, Huo Zhiguo, Wang Peijuan, et al. The applicability of mechanism phenology model to simulating apple flowering date in Shaanxi provice[J]. Journal of Applied Meteorological Science, 2019, 30(5): 555-564. (in Chinese with English abstract)

[24] 王润红，茹晓雅，蒋腾聪，等. 基于物候模型研究未来气候情景下陕西苹果花期的可能变化[J]. 中国农业气象，2021，42(9)：729-745.

Wang Runhong, Ru Xiaoya, Jiang Tengcong, et al. Based on the phenological model to study the possible changes of apple flowering dates under future climate scenarios in Shaanxi province[J]. Chinese Journal of Agrometeorology, 2021, 42(9): 729-745. (in Chinese with English abstract)

[25] 刘映宁，贺文丽，李艳莉，等. 陕西果区苹果花期冻害农业保险风险指数的设计[J]. 中国农业气象，2010，31(1)：125-129，136.

Liu Yingning, He Wenli, Li Yanli, et al. Chinese Journal of Agrometeorolog[J]. A study on the risk index design of agricultural insurance on apple florescence freezing injury in Shaanxi fruit zone, 2010, 31(1): 125-129,136. (in Chinese with English abstract)

[26] 李美荣，朱琳，杜继稳. 陕西苹果花期霜冻灾害分析[J]. 果树学报，2008(5)：666-670.

Li Meirong, Zhu Lin, Du Jiwen. Analysis of frost damage during apple bloom period in Shaanxi Province[J]. Journal of Natural Disasters, 2008(5): 666-670. (in Chinese with English abstract)

[27] 王景红. 陕西主要果树气候适宜性与气象灾害风险区划图集[M]. 西安：陕西科学技术出版社，2012：85-95.

[28] 韩会明，刘喆玥，刘成林，等. 基于Copula函数的赣江流域气象干旱特征分析[J]. 水电能源科学，2020，38(8)：9-13.

Han Huiming, Liu Zheyue, Liu Chenglin, et al. Analysis of meteorological drought characteristics in Ganjiang River basin based on Copula function[J]. Water Resources and Power, 2020, 38(8): 9-13. (in Chinese with English abstract)

[29] 张卓群，冯冬发，侯宇恒. 基于Copula函数的黄河流域干旱特征研究[J]. 干旱区资源与环境，2022，36(1)：66-72.

Zhang Zhuoqun, Feng Dongfa, Hou Yuheng. Study on drought characteristics in Yellow River basin based on Copula function[J]. Journal of Arid Land Resources and Environment, 2022, 36(1): 66-72. (in Chinese with English abstract)

[30] Michele C D, Salvadori G, Vezzoli R, et al., Multivariate assessment of droughts: Frequency analysis and dynamic return period[J]. Water Resources Research, 2013, 49(10): 6985-6994.

[31] Sklar A, Fonctions de repartition an dimensions et leurs marges[J]. Publ. inst. statist. univ. Paris, 1959, 8: 229-231.

[32] Salvadori G. Bivariate return periods via 2-Copulas[J]. Statistical Methodology, 2004, 1(1): 129-144.

[33] 张晓煜，马玉平，苏占胜，等. 宁夏主要作物霜冻试验研究[J]. 干旱区资源与环境，2001(2)：50-54.

Zhang Xiaoyu, Ma Yuping, Su Zhansheng, et al. Experiments on frost injuries of main crops in Ningxia Province[J]. Journal of Arid Land Resources and Environment, 2001(2): 50-54. (in Chinese with English abstract)

[34] 王景红，柏秦凤，梁轶，等. 2013年陕西苹果花期冻害气象条件分析及受冻指标研究[J]. 果树学报，2015，32(1)：100-107，174.

Wang Jinghong, Bai Qinfeng, Liang Yi, et al. Study on freezing index and weather conditions causing Shaanxi apple florescence freezing injury in 2013[J]. Journal of Natural Disasters, 2015, 32(1): 100-107, 174. (in Chinese with English abstract)

Analysis of the characteristics of apple later frost risks in Shaanxi Province based on Copula functions

Jiang Yuan1,2, Ru Xiaoya1,2, Luo Qi1,2, Li Meirong3, Wang Zhao3, Wang Jinghong3, Feng Hao2,4, Zhang Dong5, Su Baofeng6, Yu Qiang4, He Jianqiang1,2※

(1.712100,; 2.712100,; 3.710015,; 4.,712100,; 5.,712100,; 6.,,712100,)

Late frost is one of the most destructive meteorological disasters in the Loess Plateau of China. A great threat can be posed to the sustainable production of apples, leading to great economic losses in the apple industry. Thus, it is a high demand to explore the occurrence of late frost events for the prevention of apple late frost disasters in the decision-making of the local apple industry. In this study, the late frost return periods were investigated on the duration and severity of late frost events using the Coupla functions. The reliability of the model was then verified to analyze the characteristics of apple late frost. The study area was taken as the apple-producing areas of Shaanxi Province in western China. The meteorological datasets were collected from the seven weather stations from 1971-2018. The daily minimum temperature (min) of 0℃ was taken as the critical temperature for the occurrence of apple late frost events, in order to extract the two characteristic variables of duration and severity of late frost events. These characteristic variables of late frost events were first fitted by seven common distribution functions, respectively. Kolmogorov Smirnov (K-S) test was then carried out to verify the model. The joint distributions of late frost characteristic variables were constructed to evaluate the goodness-of-fit using six Copula functions. The occurrence probability and return period of late frost events were analyzed with the optimized Copula functions. The results showed that the severity of late frost risks generally increased from the southeast to northwest in the study area from 1971 to 2018. The optimal marginal distribution of late frost duration was in the log-normal distribution, while there was a great difference in the optimal marginal distributions of late frost severity. A significant positive correlation was found between the duration variables and the severity of late frost at each station. The joint cumulative probability increased significantly, as the severity and duration of late frost increased, but the increasing trend was much slower than before. A much more significant increase was observed in the co-occurrence return period under the same increase of univariate value, compared with the joint. The univariate return period was always between the joint and co-occurrence return period. Once the univariate return period was small enough, the optimal range of the actual univariate return period can be estimated, according to the joint and the co-occurrence return period. In general, a low probability was found in the late frost events with long duration and high severity at the weather stations in apple-producing areas in Shaanxi Province of China. However, the stations in the Yan’an area were more susceptible to the late frost events with high severity or long duration, as well as the late frost events with both high severity and long duration. Thus, more attention can be paid to late frost risks in the Yan’an area. This finding can provide a theoretical base to deal with the late frost disaster in apple production.

meteorological; disasters; late frost; Copula function; frequency analysis; return period; apple

10.11975/j.issn.1002-6819.2022.23.018

S166

A

1002-6819(2022)-23-0170-11

姜元，茹晓雅，罗琦，等. 基于Copula函数的陕西苹果晚霜冻特征分析[J]. 农业工程学报，2022，38(23)：170-180.doi：10.11975/j.issn.1002-6819.2022.23.018 http://www.tcsae.org

Jiang Yuan, Ru Xiaoya, Luo Qi, et al. Analysis of the characteristics of apple later frost risks in Shaanxi Province based on Copula functions[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2022, 38(23): 170-180. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2022.23.018 http://www.tcsae.org

2022-08-21

2022-11-26

国家重点研发计划（2021YFD1900700）；陕西省重点研发计划重点产业创新链（群）-农业领域项目（2019ZDLNY07-03）；陕西省气象局秦岭和黄土高原生态环境气象重点实验室2019年开放研究基金课题（2019Z-5）；西北农林科技大学人才专项资金（千人计划项目）；高等学校学科创新引智计划（111计划）(B12007)

姜元，博士生，研究方向为农业生态系统模拟。Email：yuanjiang59@nwafu.edu.cn

何建强，教授，博士生导师，研究方向为农业生态系统模拟。Email：jianqiang_he@nwafu.edu.cn