冰雹是降落到地面的冰球，是由强对流天气引起的一种气象灾害，冰雹天气出现的范围较小，时间短，但来势猛、强度大，常伴有狂风骤雨等灾害性天气，可给局部地区的农林牧业、工矿企业、电力、通讯、交通以及人民生命财产造成较大危害。四川省地处青藏高原东侧，地形复杂多样，部分区域是冰雹频发地，如川西高原的理塘和石渠年平均降雹近30 d。2008年7月27日20：30—21：00(北京时，下同)，川南的泸州市叙永县观兴、石坝、水潦、麻城、摩尼等7个乡镇突然刮起了罕见大风，降下冰雹，其中冰雹直径最大约1.5 cm。这次冰雹灾害致使1.54万人受灾，损坏房屋916间，农作物受灾1527 hm²，绝收533.6 hm²，粮食减产0.2万吨，造成直接经济损失1410万元。随着社会经济的发展，冰雹的危害越来越严重，因此加强冰雹预报技术研究，建立预报方法和预警系统，提高四川省冰雹预报预警能力十分重要。

探空分析是利用探空资料分析中尺度对流系统发生发展的局地垂直环境场特征，在短时效内确定中尺度天气发生的潜势^[1]。国外气象工作者对冰雹的时空分布特征、天气模型和环境条件曾做过较详细的分析，并提出了天气模型和物理参数相结合的降雹预报方法^[2]；国内王笑芳等^[3]利用北京08：00单站探空资料预报强对流天气有无及强度的判断方法，为强对流天气的短时预报提供了有应用价值的思路和方法；李耀东等^[4-5]分析了多种能量和动力参数对强对流天气预报的应用状况；廖晓农等^[6]利用探空资料计算了对流有效位能、抬升指数、风切变等物理参量，分析了北京一次罕见的大雹事件；雷蕾等^[7]通过探空资料分析了多种物理参量在冰雹、雷暴大风和短时暴雨天气下的差异；高宾永等^[8]发现大冰雹出现在高空能量锋区内，落区在低层湿有效能量正中心附近；陈晓红等^[9]发现对流有效位能、归一化对流有效位能、风暴强度指数和密度加权平均垂直风切变等参数较好地反映出强风雹天气过程；刘玉玲^[10]用探空计算对流有效位能、粗理查森数、螺旋度等对流参数并应用到强对流天气的潜势预报。

目前的冰雹预报多为以相关分析为基础的回归、判别、聚类等预测方法，而作为小概率事件，降雹与气象因子之间是复杂的非线性关系，且这种关系也不是用一个高次多项式所能描述的。相关分析是寻找因子的主要方法，但不能完全解决因子寻找问题。灾害性天气小概率事件的出现实质上是前期某些气象要素有条件的因子组合作用的影响结果，且这种因子组合，是形式多样、繁杂的非线性组合，它在气象要素的变化中占有很重要的位置。针对因子之间的非线性特点，文献[11]提出了因子组合分析方法，寻找预报因子与预报对象的关系。本文基于因子组合方法，寻找四川盆地南部冰雹与探空资料计算物理参量的关系，提出一套应用于日常业务并具有明显物理意义的冰雹潜势预报的探空指标。

冰雹资料取1978—2007年四川省宜宾探空站200 km范围内地面气象观测站的地面气象记录月报资料。探空资料为1978—2007年MICAPS系统的T-lnp探空数据。

本研究建立了一个时间样本序列，共45例时间样本，降雹事件样本7例，非降雹事件样本38例。由于冰雹等强对流天气属于中尺度天气范畴，水平尺度一般为20~200 km，垂直尺度为10~20 km左右，属于β中尺度，β中尺度的生命史一般在几小时到十几小时^[12]。根据中尺度天气的特点，遵守邻近与临近的探空资料选取标准^[13]，选取宜宾探空站200 km范围内的降雹记录，根据冰雹直径大小与发生地点的海拔高度两个因素，剔除冰雹直径小于5 mm的样本7例，剔除海拔高度超过3000 m的峨眉山站点样本1例，最终确定7例降雹事件样本 (表 1)。

表 1

表 1 降雹事件样本 Table 1 Hail event samples 序号 探空资料时间 降雹时间 降雹地点 降雹地与探空站距离/km 降雹地海拔/m 冰雹直径/mm 1 1983-04-10T20：00 22:47 宜宾 0 341 5 2 1985-08-06T08：00 17:57 荣县 72 382 20 3 1994-07-18T08：00 16:45 隆昌 88 374 9 4 1994-07-27T08：00 15:12 井研 106 419 8 5 1998-04-28T08：00 18:42 珙县 48 368 14 6 2000-07-20T08：00 17:18 井研 106 419 20 7 2002-04-04T08：00 18:50 自贡 61 357 8

表 1 降雹事件样本 Table 1 Hail event samples

由于实际业务工作中，T-lnp探空数据为每日08：00和20：00两个时次，本文所做的工作是针对12 h冰雹预警，因此选用降雹事件发生前12 h内的数据，如降雹发生在17：00，则选取08：00的数据进行计算分析。

因子组合分析方法^[11]认为某些气象要素发生不稳定状态后，在特定的组合作用下，会引起另一气象要素发生突变，因此可以分步骤确定影响预报对象发生的主要因子与条件因子，并且找出临界值以及组合关系，建立判别式描述预报因子与预报对象的关系。

① 建立小概率事件序列：首先找出预报的气象要素在出现小概率事件或者灾害性天气的时间样本，通常以1表示异常样本，0表示非异常样本，将预报对象的气象要素排列为0，1序列，称Y序列。

② 寻找主要因子：将预报因子样本集中的每个预报因子通过迭代计算得到其临界值，变换为0, 1序列，由直接寻找因子变为寻找与Y序列中1相对应关系序列的因子。

③ 主要因子组合：主要因子序列组合方式包括交运算、并运算、混合运算3种。根据符合的条件，确定其中一种运算，使最后得到的主要因子组合V′序列与Y序列中的1全部对应。

④ 寻找条件因子：将V′序列中1与Y序列中相对应的1全部充零，而将无对应的1保留下来，则得到Y′序列，再将Y′序列中的1和0对调，得到一个与Y′序列相反的Y″序列。最后寻找一个与Y″序列对应关系的序列即条件因子序列。

⑤ 条件因子组合：条件因子序列组合方式同样也有交运算、并运算、混合运算3种，但运算上与主要因子序列组合方式不同。根据条件，确定其中一种运算，使最后得到的条件因子组合序列V″与Y″序列中的0对应。

⑥ 建立判别式：根据各因子的临界值进一步确定主要因子组合序列与条件因子组合序列的搭配关系，建立判别式。

由于条件因子的制约作用是随机的，不可能仅仅出现在小概率事件发生时，因此条件因子组合V″序列中的0个数比Y″序列中0的个数多。考虑到探空资料的冰雹12 h预报后将进一步通过雷达卫星进行短时临近监测预报，在不出现漏报的情况下允许一定的空报，在因子分析方法基础上，本研究对条件因子组合的满足条件进行适当调整，使理论化的方法更加适应实际业务需求。

对45例时间样本资料进行0，1序列构建，1表示降雹事件样本，0表示非降雹事件样本，即得到Y序列，见表 2。

表 2

表 2 Y序列与前10个主要因子序列 Table 2 Sequence Y and sequence of top 10 major factors 样本 Y X′1 X′2 X′3 X′4 X′5 X′6 X′7 X′8 X′9 X′10 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 3 1 1 1 1 1 1 1 1 1 1 1 4 1 1 1 1 1 1 1 1 1 1 1 5 1 1 1 1 1 1 1 1 1 1 1 6 1 1 1 1 1 1 1 1 1 1 1 7 1 1 1 1 1 1 1 1 1 1 1 8 0 1 0 1 1 1 1 1 1 1 1 9 0 1 1 1 1 1 1 1 1 1 1 10 0 0 0 0 0 0 0 0 0 0 0 11 0 0 0 0 0 0 0 0 0 0 0 12 0 0 0 0 0 0 0 0 0 0 0 13 0 0 0 0 0 0 0 0 0 0 0 14 0 0 0 0 0 0 0 0 0 0 0 15 0 0 0 0 0 0 0 0 0 0 0 16 0 1 1 1 1 1 1 1 1 1 1 17 0 0 0 1 1 0 0 1 1 1 1 18 0 0 0 0 0 0 0 0 0 0 0 19 0 0 0 0 0 0 0 0 0 0 0 20 0 0 0 0 0 0 0 0 0 0 0 21 0 0 0 0 0 0 0 0 0 0 0 22 0 0 0 0 0 0 0 0 0 0 0 23 0 0 0 0 0 0 0 0 0 0 0 24 0 0 0 0 0 0 0 0 0 0 0 25 0 0 0 0 1 0 0 1 1 1 1 26 0 1 1 1 1 1 1 1 1 1 1 27 0 1 1 0 0 1 1 0 0 0 0 28 0 0 0 0 0 0 0 0 0 0 0 29 0 0 0 0 0 0 0 0 0 0 0 30 0 0 0 0 0 0 0 0 0 0 0 31 0 0 0 0 0 0 0 0 0 0 0 32 0 0 0 0 0 0 0 0 0 0 0 33 0 0 0 0 0 0 0 0 0 0 0 34 0 0 0 0 0 0 0 0 0 0 0 35 0 0 1 1 1 1 1 1 1 1 1 36 0 0 1 1 1 1 1 1 1 1 1 37 0 0 0 0 0 0 0 1 1 1 1 38 0 0 0 0 0 0 0 0 0 0 0 39 0 0 0 0 0 0 0 0 0 0 0 40 0 0 0 1 0 1 1 0 0 0 0 41 0 1 0 1 1 1 1 1 1 1 1 42 0 1 1 1 1 1 1 1 1 1 1 43 0 0 1 0 1 1 1 1 1 1 1 44 0 0 1 1 0 0 0 0 0 0 0 45 0 0 1 0 1 1 1 1 1 1 1

表 2 Y序列与前10个主要因子序列 Table 2 Sequence Y and sequence of top 10 major factors

建立一个由3422个气象物理参数组成的预报因子集，包括温度、高度、水汽压、假相当位温、K指数、饱和温湿能等，其中表征热力、动力、水汽条件及特殊对流指数等参数参考文献[14-17]，利用T-lnp探空数据计算得出。根据本文介绍的方法计算3422个预报因子各自的临界值，将超过临界值的因子以1表示，未超过临界值的因子以0表示，得到预报因子的0，1序列X′₁~X′₃₄₂₂。

按σ′_j由大到小排序，挑选出10个排列靠前的因子，见表 2。

(1)

X′₁为400 hPa饱和湿静力温度与850 hPa湿静力温度之差的0，1序列，X′₂为400 hPa与地面垂直气压梯度之差的0，1序列，X′₃为400 hPa饱和比湿与850 hPa比湿之差的0，1序列，X′₄为850 hPa正压大气比能的0，1序列，X′₅为400 hPa与850 hPa温湿能之差的0，1序列，X′₆为400 hPa与850 hPa相当位温之差的0，1序列，X′₇为850 hPa位势高度的0，1序列，X′₈为850 hPa位势能的0，1序列，X′₉为850 hPa与地面位势高度之差的0，1序列，X′₁₀为850 hPa与地面位势高度能之差的0，1序列。

由于X′₁~X′₁₀均满足σ′₁=σ′₂=…=σ′₁₀=7，即7例降雹事件样本X′₁~X′₁₀均报出，因此X′₁~X′₁₀序列进行主要因子交运算。选取空报数较少的X′₁所对应的因子作为第1主要因子，将X′₁分别与X′₂~X′₁₀进行交运算V′_1∧j=X′₁∧X′_j=min (X′₁, X′_j)，j=2, …, 10, 分别对9个交运算组合的因子序列进行求和，分别计算n′_1∧j，

(2)

计算得出的n′_1∧j最小值所对应的因子即为第2主要因子，它们的交运算X′₁∧X′_j为最佳组合。经过计算发现n′_1∧2=12最小，得出X′₂所对应的因子为第2主要因子，X′₁∧X′₂为主要因子序列的最佳组合。

X′₁为400 hPa饱和湿静力温度与850 hPa湿静力温度之差的0, 1序列，X′₂为400 hPa与地面垂直气压梯度之差的0, 1序列，若T_σ400^*-T_σ850≤2.94℃, 则X′₁=1，否则X′₁=0；若G_z400-G_zsurface≥-57.18 hPa/km, 则X′₂=1，否则X′₂=0。

由于V′_1∧2=X′₁∧X′₂=min (X′₁, X′₂)>7，即X′₁∧X′₂的组合序列中1与Y序列中1全部对应外，还有一些无对应的1，如样本9、样本16、样本26、样本27和样本42，但这5例样本实际上并未发生降雹事件，因此推论在同时满足第1主要因子和第2主要因子的临界值条件外，还受到某种条件因子约束导致未出现降雹事件。

将V′_1∧2序列与Y序列中相对应的1全部充零，而将多余的1保留下来，则得到Y′序列，将Y′序列中的1和0对调，得到Y″序列。寻找除第1主要因子与第2主要因子以外的3420个预报因子的临界值，将超过临界值的因子以0表示，未超过临界值的因子以1表示，得到预报因子的0，1序列X″₁~X″₃₄₂₀。

按σ″_j由大到小排序，挑选出8个排列靠前的因子 (表 3)。

(3)

式 (3) 中，c满足a^c=1-a，其中0≤a≤1。X″₁为700 hPa水汽压与饱和水汽压之差的0，1序列，X″₂为700 hPa露点温度与饱和湿静力温度之差的0，1序列，X″₃为700 hPa湿球温度与饱和湿静力温度之差的0，1序列，X″₄为700 hPa饱和湿静力温度与地面干空气静力温度之差的0，1序列，X″₅为700 hPa虚温与饱和湿静力温度之差的0，1序列，X″₆为400 hPa与700 hPa饱和潜热能之差的0，1序列，X″₇为400 hPa垂直气压梯度的0，1序列，X″₈为400 hPa与700 hPa饱和温湿能之差的0，1序列。

表 3

表 3 主要因子组合序列及其前8个条件因子序列 Table 3 Sequence of major factor combination and its top 8 condition factors 样本 Y V′1∧2 Y′ Y″ X″1 X″2 X″3 X″4 X″5 X″6 X″7 X″8 1 1 1 0 1 1 1 1 1 1 1 1 0 2 1 1 0 1 1 0 0 0 0 0 0 0 3 1 1 0 1 1 1 1 1 1 1 1 1 4 1 1 0 1 1 1 1 1 1 1 0 0 5 1 1 0 1 1 1 1 1 1 1 1 1 6 1 1 0 1 1 1 1 0 0 0 0 1 7 1 1 0 1 0 1 1 1 1 1 1 1 8 0 0 0 1 0 0 1 0 1 0 1 0 9 0 1 1 0 0 0 0 0 0 0 0 0 10 0 0 0 1 1 1 1 1 1 1 1 0 11 0 0 0 1 1 1 1 1 1 1 1 1 12 0 0 0 1 1 1 1 1 1 1 1 1 13 0 0 0 1 1 1 1 1 1 1 1 1 14 0 0 0 1 1 1 1 1 1 1 1 1 15 0 0 0 1 1 1 1 1 1 1 1 1 16 0 1 1 0 0 0 0 0 0 0 0 0 17 0 0 0 1 1 1 1 1 1 1 0 1 18 0 0 0 1 1 1 1 0 1 1 1 0 19 0 0 0 1 1 1 1 1 1 1 1 0 20 0 0 0 1 1 1 1 1 1 1 1 1 21 0 0 0 1 1 1 1 1 1 1 1 0 22 0 0 0 1 1 1 1 1 1 1 1 0 23 0 0 0 1 1 1 1 1 1 1 1 1 24 0 0 0 1 1 1 1 1 1 1 1 1 25 0 0 0 1 0 0 0 0 1 0 1 1 26 0 1 1 0 0 0 0 0 0 0 0 0 27 0 1 1 0 0 0 0 0 0 0 0 0 28 0 0 0 1 1 1 1 1 1 1 1 1 29 0 0 0 1 1 1 1 0 1 1 0 1 30 0 0 0 1 1 1 1 1 1 1 1 1 31 0 0 0 1 1 1 1 1 1 1 1 1 32 0 0 0 1 0 1 1 1 1 1 1 1 33 0 0 0 1 0 0 0 1 1 1 1 0 34 0 0 0 1 0 0 0 0 0 0 0 1 35 0 0 0 1 1 1 1 1 1 1 0 1 36 0 0 0 1 1 1 1 0 0 1 0 1 37 0 0 0 1 1 1 1 1 1 1 0 1 38 0 0 0 1 1 1 1 1 1 1 1 1 39 0 0 0 1 0 0 0 0 0 1 0 1 40 0 0 0 1 0 0 0 1 1 1 1 1 41 0 0 0 1 0 0 0 0 0 0 0 0 42 0 1 1 0 0 0 0 0 0 0 0 0 43 0 0 0 1 1 1 1 1 1 1 0 1 44 0 0 0 1 1 1 1 1 1 1 1 0 45 0 0 0 1 0 1 1 1 0 1 0 0

表 3 主要因子组合序列及其前8个条件因子序列 Table 3 Sequence of major factor combination and its top 8 condition factors

由于X″₁~X″₈均满足式σ″₁=σ″₂=…=σ″₈=40，即均出现与Y″序列对应的0，但均出现个别0与Y序列中的1相对应，故对这些序列进行条件因子并运算。X″₁分别与X″₂~X″₈进行并运算V″_1∨j=X″₁∨X″_j=max (X″₁, X″_j)，j=2, …, 8, 且均满足Y_i∧(V″_i)^c=min (Y_i, (V″_i)^c)=0, i=1, …, 45，得出X″₁为第1条件因子。将V″_1∨2~V″_1∨8的7个并运算组合的因子序列分别与Y序列做比较，V″_1∨2无对应的1最少，则X″₂所对应的因子为第2条件因子，因此确定X″₁∨X″₂为条件因子的最佳组合。

X″₁为700 hPa水汽压与饱和水汽压之差的0, 1序列，X″₂为700 hPa露点温度与饱和湿静力温度之差的0, 1序列，若e₇₀₀-e_s700≤-5.74 hPa, 则X″₁=0，否则X″₁=1；若T_d700-T_σ700^*≤-69.75℃, 则X″₂=0，否则X″₂=1。

最终判别式为 (X′₁∧X′₂)∧(X″₁∨X″₂)=min (min (X′₁, X′₂), max (X″₁, X″₂)), 其值为1，则预测有冰雹；其值为0，则预测无冰雹。

第1主要因子是400 hPa饱和湿静力温度与850 hPa湿静力温度之差，即T_σ400^*-T_σ850。雷雨顺^[14]提出，T_σ^*-T_σ0表征了一般未饱和湿空气的静力稳定度，即潜在性稳定度，值越小表示越不稳定。这是一个极重要的判据，在许多情形下, 比条件性稳定度和对流性稳定度能更恰当刻画大气静力不稳定性质。在对暴雨与雷暴能量条件的差异进行对比时，运用T_σ500^*-T_σ850的不同阈值来区分，且称之为对流层下部潜在稳定度参数，雹暴时T_σ500^*-T_σ850≤2℃，而本文的结果T_σ400^*-T_σ850≤2.94℃就表征了这一判定模式，数值越小，越不稳定，越有利于强对流发生。

第2主要因子是400 hPa与地面垂直气压梯度之差，即G_z400-G_zsurface。显然低层大气气压随高度减少得快，干冷空气比暖湿空气气压随高度减少得快^[16]。本文结果G_z400-G_zsurface≥-57.18 hPa/km，说明高层气压垂直梯度越大，低层气压垂直梯度越小，越有利于冰雹天气的生成，即上层越具有干冷空气特质、下层越具有暖湿空气特质，冰雹云越容易发生。这种高层大气密度增加、低层大气密度减小的密度不连续性特征加速了动力不稳定性的生成，动力不稳定是对流能否发展和维持的关键因素，不同强对流天气现象的发展和移动多与动力稳定度有直接关系^[18]。

第1条件因子是700 hPa水汽压与饱和水汽压之差，即e₇₀₀-e_s700。水汽压代表大气中的绝对水汽含量，又被称为绝对湿度。700 hPa处于0℃层附近，该气压层的水汽丰沛与否关系冰雹是否生成，本文结果e₇₀₀-e_s700≤-5.74 hPa说明700 hPa水汽压越是小于饱和水汽压，足够接近饱和而产生冰晶的水汽越少，从而抑制冰雹产生。

第2条件因子是700 hPa露点温度与饱和湿静力温度之差，即T_d700-T_σ700^*。本文结果T_d700-T_σ700^*≤-69.75℃，说明700 hPa露点温度小于饱和湿静力温度，不利于湿空气中水滴的析出，进而抑制冰雹的生成。

本文得到的冰雹预报指标判别式中包含能量、温度、气压、水汽等多方面的物理参数，反映了强对流的不稳定机理，符合对冰雹云生成环境状态的物理解释。

为了检验基于因子组合分析法得出的冰雹预报指标判别式的预报准确率，选取2008年四川省宜宾探空站200 km附近降雹事件共10例，各地面作业点根据宜宾市气象局预警监测实施的防雹作业记录15例，以此检验预报效果。将2008年全年732例时次的探空资料代入判别式进行试报，共预报出65例降雹事件，其中报出实际发生的10例降雹事件和11例防雹作业记录。若只考虑降雹事件，则该预报结果对于降雹事件无一漏报。若将15例防雹作业视为潜在的降雹事件，通过两类天气预报评分方法^[19-24]，计算得探测率为84%，空报率为67.7%，成功指数为30.4%。由于防雹作业还要根据雷达监测等多种临近预报技术^[25]进一步判定与筛查，因此在利用探空资料所做预报不出现漏报的情况下，存在一定数量的空报是可以接受的。

1) 利用每日08:00和20:00 T-lnp探空数据计算了表征热力、动力、水汽条件、温度、高度等气象物理参数，采用因子组合分析方法，对四川省宜宾探空站200 km范围内的冰雹预报进行分析研究，获得冰雹预报的主要因子与条件因子，得到它们的阈值以及组合关系，建立了冰雹预报指标判别式。

2) 利用2008年T-lnp探空数据进行冰雹试报，报出65次降雹事件，实际发生的10例降雹事件全部报出，若将15例防雹作业视为潜在的降雹事件，则报出11例防雹作业记录，漏报4例。计算得到探测率为84%，空报率为67.7%，成功指数为30.4%。

3) 因子组合分析方法可以预报降雹这一类的小概率事件，通过一套数理统计方法寻找影响小概率事件的主要因子及抑制事件发生的条件因子，并定量计算其临界值。合理的物理解释说明这套探空数据计算出的预报指标有较好的预报能力和实际应用意义。

目前本预报指标还存在一些不足，如降雹资料的样本有限，降雹多不在气象台站内，未记录在册，影响了对当日是否发生降雹的判断，从而影响降雹事件样本的搜集。探空站空间分布与时间尺度不理想，对时空尺度均为中小尺度的冰雹灾害，200 km的覆盖面积和12 h的间隔探空监测显得分辨率不够，直接影响预测结果。由于降雹事件多发生在午后至半夜，因此14:00的加密探空对于冰雹的临近预报至关重要^[26]。本文结果主要给出了当时能量场的状态，未对整体的三维气流变化描述，需要结合雷达观测、卫星监测与数值预报等多种资料，在减少空报率方面做进一步努力。随着国内外新的对流指数的出现，需要考虑更多因子，考察其对冰雹的预报效果，对冰雹指标进行更新修正。同时，随着冰雹观测记录的规范及探空资料加密观测的实施，可逐步完善因子组合形式，由此提高冰雹预报的准确性和及时性，为人工防雹作业提供参考。