分位数回归模型的构建及评估

1、利用地磁ap指数建立电网GIC 评估模型 基于分位数回归理论华北电力大学 马学俊、金炜、王博目录摘 要3一、问题的提出3二、研究现状及存在的问题4三、数据的描述4四、经典回归模型的构建及评价54.1经典回归模型的构建54.2经典回归模型的评估6五、分位数回归模型的构建及评估75.1分位数回归模型的构建75.2分位数回归模型的评估9六、GIC的量表制作106.1不同ap指数下的GIC风险量表106.2高频ap指数下的GIC分布量表11七、总结12参考文献13摘 要空间天气的变化引起的地磁剧烈扰动称为磁暴，磁暴引起电场变化产生的地磁感应电流（简称GIC）对电力系统的安全运行产生不利的影响。随着我国

2、高压及特高压的长距离输电线路相继建成运行，GIC对电网的影响越来越受到关注，GIC的评估也变得越来越重要。对于GIC的评估主要分为物理模型和统计模型。建立物理模型考虑到很多因素的影响，这些影响因素的数据往往很难得到，因此在实际模型的建立中假设因素比较多，如平面波方法假设大地为电导率均匀的无限大半空间，这些都影响到模型的评估效果。统计模型不需考虑其复杂的物理机制，利用已掌握的数据，通过分析数据与GIC关联特征可以建立评估模型，因此受到GIC研究学者的重视。但是国内外对于GIC统计模型研究的进展依然缓慢，目前的模型主要是经典回归模型。但是经典回归模型只能描述ap对GIC均值的影响，不能体现对电网产

3、生严重影响上的GIC值，模型反映实测GIC的信息非常有限，特别是对于GIC分布尾部信息几乎没有表现出来。本文利用中国科学院空间环境研究预报中心数据共享平台的地ap数据和2004-2005年几次强磁暴发生时广州岭澳核电站监测GIC数据，建立GIC关于阿婆的分位数评估模型。克服了经典回归只能描述自变量对于因变量均值的关系的缺点，分位数回归模型能够全面描述ap和GIC之间的关系，且具有良好的稳健性，不易受到离群点的影响。本文还根据电网安全运行关注的问题，制作了建立GIC的风险值量表和高频ap指数的GIC分布量表。量表直观易懂，实际操作性强，能够方便的查询GIC风险值，对于电网的安全运行提供重要的参考

4、。关键词：磁暴；地磁感应电力（GIC）；ap指数；分位数回归；GIC风险值一、问题的提出太阳活动形成的太阳风经过地球时，有着巨大能量的太阳粒子流与地球磁层发生相互作用从而对地球磁场产生的不规则剧烈扰动称为磁暴。磁暴感应出的地面电场在输电线路、中性点接地的变压器和大地形成的回路中产生地磁感应电流(Geomagnetically Induced Currents,简称GIC )。GIC的产生对于石油和天然气管道、通信线路以及电力系统的安全稳定运行带来很大影响，例如造成电力变压器严重半波饱和、 发电机过热与振动、 继电保护装置误动、系统无功消耗增加和静止无功补偿电容器过载等1。1989 年 3 月

5、13 日，加拿大魁北克的大停电事故是历史上遭受GIC影响最为严重的一次停电事故，使魁北克电力系统瘫痪近9个小时，近600 万人无法正常用电，造成了巨大的社会影响和严重的经济损失。2001年以来，江苏上河、广州岭澳等地的变压器多次出现不明原因的强烈振动和噪声增大事件，后经证明是磁暴在电网产生的GIC所为2。随着我国500kv、750kv的长距离输电线路相继建成运行，GIC对电网的影响越来越受到关注，其中2004年11月7日广州岭澳核电站监测GIC最大峰值达到55.8A。由于500KV、750KV不能满足需要，我国需要建设1000kv特高压电网向东部沿海地区送电，而特高压电网网络结构和电气参数的特

6、点将使电网受GIC危害的风险大大增加3，因此评估GIC具有十分重要的意义。二、研究现状及存在的问题对于GIC的评估主要分为两种：物理模型、统计模型。对于物理模型，国外研究的相对比较早，加拿大、芬兰等从20世纪七十年代末开始对电网中的GIC进行检测。国内中国科学院地质与地球物理研究所马晓冰4等从地球物理的角度出发研究GIC的作用及评估；华北电力大学刘连光、刘春明2,3,5结合电网研究磁暴对中国的影响，并利用平面波和复镜像法、电网模型等评估电网中的GIC。由于物理模型受到很多因素的影响，这些影响因素的数据往往很难得到，因此在实际模型的建立中假设因素比较多，如平面波方法假设大地为电导率均匀的无限大半

7、空间，这些都影响到模型的评估效果。统计模型不需考虑其复杂的物理机制，利用已掌握的数据，通过分析数据与GIC关联特征可以建立GIC评估模型，因此受到GIC研究学者的重视。慕学玲6，刘连光 2分析了GIC与地磁数据、地磁Dst指数的相关性，加拿大Trichtchenko7 等用统计的方法研究了GIC与地磁ap指数之间的关系，并建立GIC与ap的一元线性回归模型。但是国内外对于GIC统计模型研究的进展依然缓慢，这一方面由于经典一元回归对残差要求是的要求非常高（要求独立同分布），实际数据很难满足这一要求，另一方面经典一元回归只能描述ap对GIC均值的影响，不能体现对电网产生严重影响上的GIC值，当GI

8、C的离散程度相当大时，模型反映实测GIC的信息非常有限，特别是对于GIC分布尾部信息几乎没有表现出来。本文主要基于分位数回归理论建立GIC与地磁指数ap的分位数回归模型，该模型能够评估各不同分位概率下的GIC，反应GIC的全部信息，能刻画GIC分布的尾部特征。而且分位数回归模型对于随机扰动项要求很低（仅要求分位数），数据很容易满足，模型不易受异常点的影响，比较稳健。三、数据的描述地磁指数是描述每个时间段内地磁扰动的总体强度或某类磁扰强度的指数，常用的指数有K，Kp，ap指数等。其中K指数是描述单个地磁台3小时时段内地磁磁情的指标，Kp是描述全球3小时段内地磁活动的指标，它从全球地磁台网中选择1

9、2个台站，先求每个台站的标准化指数，然后求这些标准指数的均值。按Kp指数的大小，磁暴分为三类：中常磁暴（Kp=5，6）、中烈磁暴（Kp=7，8）、强烈磁暴（Kp=9）。由于Kp指数与地磁扰动幅度的非线性关系使得它的运算很不方便，为此有必要把它们转化为线性幅度，这样得到的指数叫作ap ，也称为“等效的行星性3小时幅度”。ap指数也是描述全球3小时段内地磁活动的指数8。 ap与Kp的对应关系见表1：Kp000+1-101+2-202+3-303+4-404+ap02345679121518222732Kp5-505+6-606+7-707+8-808+9-90ap39485667809411113

10、2154179207236300400表1 Kp与ap的对应表本文ap数据来源于中国科学院空间环境研究预报中心数据共享平台。GIC数据来自2004-2005年几次强磁暴发生时广州岭澳核电站监测数据（表2）。事件Kp指数磁暴类型起止时间(UT)18急始04/11/07 10:30-11/08 12:3028急始04/11/09 18:30-11/10 22:3038急始05/01/21 17:00-01/23 00:3049急始05/05/15 02:30-05/16 19:3057缓始05/05/29 20:30-05/31 00:3069急始05/08/24 06:30-08/25 21:3

11、077缓始05/08/31 06:30-09/01 18:30表2 2004-2005年强磁暴事件本文中所有的数据的统计分析均是用R软件实现。四、经典回归模型的构建及评价4.1经典回归模型的构建经典的一元回归分析模型： 它的回归方程为模型假设： 其中。模型参数估计：其中arg表示最小二乘函数的极值点。需要注意的是模型描述自变量x对于因变量y均值的影响，它的假设条件：（1）残差之间不相关且等方差；（2）残差来自均值为0的正态分布。对于ap与GIC，它们的回归方程为： 表3为参数估计值：Estimatet valuePIntercept3.089073.4370.00094*ap0.058037.

12、4030.00000*、*分别表示在 0.001、0.01下显著表3模型参数估计得到的回归方程为： 4.2经典回归模型的评估我们首先检验模型的显著性。检验统计量F =54.81, 其对应的P值：p-value=0.00000<0.01，因此在0.01模型通过检验。再做系数显著性检验。表3显示，对于回归系数的检验， ap系数的检验统计量t=7.0403， P=0.000<0.01，说明ap影响显著。值得注意的是，此时F 及t检验基于残差满足假设的得到的。如果残差存在异方差及自相关，可能导致回归系数的t检验值高估，可能造成本来不显著回归系数变成显著，如果残差不满足模型条件，那么F、t检

13、验将失效9。所以下面进行拟合和效果观察及检验残差。通过拟合图（图1），发现存在较多数据点离图中拟合直线偏差较大。 图1 经典回归模型拟合图比如ap=179时，拟合GIC为13.49，实际的GIC均值是13.8240，实际数据中GIC对应值为4.02、4.44、11.11、16.59、32.96值。显然13.49体现的信息非常有限。这条直线只能反映GIC的平均信息。对于有可能对电网系统造成严重影响的GIC高尾信息反映不出来。我们进一步对于模型的残差进行分析： 自相关性验检用DW检验，它的原假设残差不存在一阶自相性。经计算得到DW的P值2.192e-07<0.001，可见显著拒接原假设，即残

14、差存在一阶自相关。异方差检验用等级相关系数法，其原假设是残差不存在异方差。通过计算得到其P值为0.003273<0.01，可以得到残差存在异方差。正态性检验采用采用常见的Shapiro-Wilk检验得到P为1.22e-08，从而拒绝残差来自正态分布。通过对残差的检验，发现残差不满足模型的假设，所以经典的一元回归模型效果不理想，同时也不能满足我们对于电网中产生严重影响的GIC的评估需求。为此我们引入分位数回归模型。五、分位数回归模型的构建及评估5.1分位数回归模型的构建分位数回归是Roger Koenker10提出，它依据因变量的条件分位数对自变量进行回归。其比传统的经典回归优点：首先它对

15、模型的随机扰动项不需要做任何分布假设，这样回归模型有很强的稳健性；其次它是对所有概率的分位数进行回归，对于数据出现的异常点具有耐抗性，最后它能够精确的描述自变量对于因变量的条件变化范围以及对条件分布形状的影响，捕捉分布的尾部特征，更加全面的刻画分布的特征11。设为实随机变量，其分布函数为，对于任意的，(有时也用tau表示)则称为的分位数。给定，分位数回归模型为： 分位数回归函数：模型参数估计： （1）其中（1）式可以转化为求进而转化为线性规划问题求解。最常用的求模型参数方法有单纯性法、内点法。单纯性法是由 Koenker and d'Orey12将最优化问题的单纯性算法应用到

16、分位数回归中，该算法适合样本量不大和不多的变量，当数据中存在大量离群点时，单纯算法估计出来的参数具稳定性比较好，但是在处理大量数据时运算的速度会显著的降低11。内点法是由Portnoy 、Koenker13提出，适合样本量比较大，并且自变量个数少的样本数据。实际中除了上面两种方法，还有预处理内点法及罚函数的方法。针对GIC与ap的数据特点，本文选择一元线性分位数回归。其回归方程为本文选择单纯性方法估计参数，得到不同tau下分位数模型的系数： tauValuet valuePr(>|t|)0.05(Intercept)1.332722.609740.01081*ap0.015012.387

17、440.01933*0.1(Intercept)1.474555.302960.00000*ap0.016572.332230.02220*0.2(Intercept)1.285363.415520.00100*ap0.029973.286670.00151*0.3(Intercept)1.554054.182690.00007*ap0.037493.653790.00046*0.4(Intercept)1.603625.251760.00000*ap0.053115.843170.00000*0.5(Intercept)1.947614.828490.00001*ap0.055577.371

18、430.00000*0.6(Intercept)2.401043.590720.00057*ap0.060836.927590.00000*0.7(Intercept)3.589343.14880.00230*ap0.056874.378190.00004*0.8(Intercept)5.754553.967760.00016*ap0.066934.289940.00005*0.9(Intercept)7.798516.784380.00000*ap0.085392.495920.01462*0.95(Intercept)7.943667.566840.00000*ap0.139762.553

19、870.01255*、*、*分别表示在 0.001、0.01、0.05下显著表4 分位数模型的参数估计5.2分位数回归模型的评估首先对模型参数进行检验。参数检验常用方法有三种：残差为独立同分布、残差为非独立同分布、Bootstrap法。Bootstrap方法是通过重新抽样来估计参数估计，常用的方法有X-Y Bootstrap、残差Bootstrap及马尔可夫链边际Bootstrap。其中X-Y Bootstrp是Bootstrp最一般的形式，不要求残差和自变量X独立，它是对X和Y进行重新抽样（样本数m不必等于原序列的样本数n），然后计算分位数回归的系数估计值，重复进行B次抽样得到B个系数估计值

20、。渐近分布协方差可以由此估计出来： （1）其中是B个参数估计值的均值。 表4中的第四、五列分别是参数检验的t值和相应的P值，可以看到所有的系数在0.05水平下都是显著的。再对模型进行检验。模型检验采用Koenker提出的拟似然比（quasi-likelood-ration） 14。构造的统计量为：其中，为稀疏函数。在原假设成立时，它服从自由度为q（自变量的个数）的卡方分布。利用R软件计算得到每个模型检验的拟似然比检验量值和相应的P值（表5），说明每个分位数回归模型均可以通过。tauQuasi-LR statisticProb(Quasi-LR stat)0.059.456360.00210*0

21、.117.032580.00003*0.224.361180.00000*0.343.279430.00000*0.447.433130.00000*0.555.537030.00000*0.642.609870.00000*0.736.595220.00000*0.831.893320.00000*0.912.745790.00036*0.958.780820.00304*、*分别表示在 0.001、0.01下显著表5 模型的显著性检验通过模型的拟合图观察：图2 分位数回归模型拟合图图2中自下而上的直线分别是分位数为0.05，0.5，0.7，0.8，0.95的回归直线。可以看出，GIC随着a

22、p的增长而呈增长趋势,经典一元回归线与tau=0.6的分位数回归线比较接近。在ap<110时，GIC值密集在较低的位置，少数GIC值偏高，此时我们可用不同的分位数回归拟合GIC。分位数回归可得到不同tau的GIC拟合值，分位数回归可以反映GIC的全部信息。图2中tau不同，得到的分位数回归方程斜率也不同，这也说明经典回归中确实存着异方差性，这种情形下建立分位数回归模型是合理的和必要的。六、GIC的量表制作6.1不同ap指数下的GIC风险量表表4给出了不同tau的分位数模型系数，根据系数表，既可以计算确定时不同ap下的GIC,也可以计算ap确定时不同下的GIC。根据表4： =0.95时，分

23、位数模型为：ap=179时分位数回归的预测值为44.567，也就是：或 说明，ap=179时，GIC的值超过44.567的概率是0.05，如果这个值对电网是安全的，那么就可以正常运行；如果这个值对电网安全是危险的，那么就应该引起注意，采取措施。所以我们把44.567称为0.05水平下的风险值，GIC风险值是在电网运行中特别重要的数字。表明了对于安全性要求的高低，这个数越小，安全性越高。但是也不能无限降低，因为那样会导致成本的无限提高，一般选择。由此，本文制作水平下的风险值量表：apKp风险值39515.1396348517.0313356518.7128467621.0249180623.75

24、73694626.7111730.2732132734.68716154739.31131179844.56602207850.4513236856.54677300969.99883表6 GIC的0.0.5风险值量表该表在实际中方便好用，根据磁暴发生前的ap，可以方便的查询到GIC风险值。6.2高频ap指数下的GIC分布量表我们对磁暴发生时ap的频率进行统计发现：ap值为48、67、94、179、236出现最为频繁，因此对于这几个指数下GIC的特征需要更详尽的了解。本文制作了GIC的条件分布表： aptau= 0.05tau= 0.50tau= 0.70tau= 0.90tau= 0.954

25、81.25013.2569765.16136712.54694117.03133671.5288854.2602026.72610515.24414821.02491941.9250535.6858398.9496819.07702226.71793.1722510.17395515.94982331.14347844.566022364.00860513.18363220.64403739.23510156.54677表7 常见ap对应的GIC的分布（低、中、高尾）表7中，tau=0.05的拟合值反映GIC低尾信息；tau=0.5的拟合值反映GIC的平均水平（不一定是均值），tau=0.7的

26、拟合值,体现的是GIC中高尾信息，tau=0.9、0.95反映GIC的高尾信息。对于给定的ap我们可以得到，GIC不超过某个范围的概率。如ap=179，我们可以也相信GIC不大于3.17概率为5%、不大于10.17概率为50%、不大于15.94概率为70%、不大于31.14概率的为90%、不大于44.57概率的为95%。七、总结本文首先对于ap与GIC数据，建立经典一元回归，但是该模型不满足模型假设条件，同时不能反应可能对电网产生严重影响的GIC信息，所以该模型效果并不理想。为此本文利用分位数回归理论建立了ap与GIC模型，该模型检验通过，并且对于不同概率下的GIC分位数能够很好的预测，并根据

27、电网运行安全关注的问题，给出常见ap（Kp>4）相对应的GIC的分布（低、中、高尾拟合值），还给出GIC的风险值表，该表的实际操作性极强，根据磁暴发生前的ap，该表可以方便的查询GIC风险值，对于电网的安全提供重要的参考作用。参考文献1 D. H. Botelera, R. J. Pirjolab, a and H. Nevanlinna,The effects of geomagnetic disturbances on electrical systems at the Earth's surface. Advances in Space Research

28、. 1998, 22(1), pp.17-272 刘连光、刘春明、张冰等.中国广州电网的几次强磁暴影响事件J.地球物理学报，2008,51(4),pp976-9813 刘连光、刘春明、张冰.磁暴对我国特高压电网的影响研究J.电网技术，2009,33(11)4马晓冰，Ian J.Ferguson.孔祥儒等地磁感应电流( GIC)的作用与评估，2005,48（66）5 ChunMing Liu, LianGuang Liu, Risto Pirjola, ZeZhong Wang. Calculation of geomagnetically induced currents in mid to lowlatitude power grids based on the plane wave method: A preliminary case study. Space Weather J, Vol. 7, S04005, 9 PP., 2009doi:10.1029/2008SW000439. 6慕学玲、何凤霞、刘连光.电网GIC数据与磁暴地磁数据相关性的分析J.中国地球物理学会第二