Robust Power System Dynamic State Estimator with Non-Gaussian Measurement Noise: Part II--Implementation and Results

IEEE TRANSACTIONS ON PO WER SYSTEMS, VOL. , NO. , 2017 1 Rob ust Po wer System Dynamic State Es timator with Non-Gauss ian Measurement No ise: P art II–Implementa tion and Results Junbo Zhao, Studen t Member , IEEE , Lamine Mili, F ell ow , IEEE Abstract —This p aper is the second of a two-part series t h at discusses the implementation issu es and test results of a rob ust Unscented Kalman Filter (UKF) for power system dynamic state estimation with non-Gaussian synchrophasor measurement noise. The tuning of the parameters of our Generalized Max imum- Likelihood-type robust UKF (GM-UKF) is presented and dis- cussed in a systematic way . Using simulations carried out on the IEEE 39-bus system, its performance is evaluated under different scenarios, includ ing i) the occurrence of two different types of noises following thick-tailed distributions, namely the Laplace or Cauchy probability d istributions for real and reacti ve power measureme nt s; ii) the occurrence of observ ation and innovation ou t liers; iii) the occurrence of P MU measurement losses due to communication failures; iv) cyber attacks; and v) strong system n onlinearities. It is also compared to the UKF and th e Generalized M aximum-Likelihood-type robust iterated EKF (GM-IEKF). S imulation results revea l that the GM-UKF outperfor ms the GM-IEKF and the UKF in all scenarios consid- ered. In particular , wh en the system is operating u nder stressed conditions, indu cing system nonlin earities, th e GM-IEKF and the UKF diver ge while our GM-UKF does conv erge. In addi t ion, when th e power measurement noises obey a Cauchy d istribution, our GM-UKF con ver ges to a state estimate vector th at exhibits a much higher statistical efficiency than th at of the GM-IEKF; by contrast, the UKF fails to conv erge. Finally , potential applications and future work of th e proposed GM-UKF ar e d iscussed i n concluding remarks section. Index T erms —Rob ust dynamic state esti mation, unscen t ed Kalman filter , phasor measurem ent un it, state tracking, Laplace noise, Cauchy noise, outliers, cyber attacks, strong nonlinearity . I . I N T R O D U C T I O N R ELIABLE and fast d ynamic state e stima to r (DSE) plays a vital role in power system mon itoring an d control. In the literature, both the p rocess and the observation noises of th e system nonlinear dy namic models are assumed to be Gaussian when developing a DSE. Furthermo re, the dyna mical system mo del is suppo sed to be accu rate and the PMU measuremen ts are secure. Howe ver , these assum ptions d o not hold true for practical power system s as elaborated in the first p art of this two-p art series. T o address th ese pro blems, se veral robust E xtended Kalm an Filter (EKF) an d Unscen ted Kalman Filter (UKF)-based D SE s hav e been pro posed [1]–[5]. In [ 1]–[3], while filters ba sed on the Huber M- estimator are propo sed to suppress observation o utliers, they are vulner able to m easurement noise obeying thick- tailed distribution and innovation outliers that ar e in duced by model param eter errors. Junbo Zhao and Lamine Mili are with the Bradl ey Department of E lec- trical Computer E nginee ring, V irgi nia Polytech nic Institute and State Uni- versi ty , Northern V irginia Center , Falls Church, V A 22043, US A (e-mail: zjunbo@v t.edu, lmili@vt.edu ). These lim itatio ns are m itigated in [4], [5] via a Gener a lized Maximum- Likelihood-ty pe robust iterated EKF (GM-IEKF). While the latter is resistant to ob servation and innovation outliers and se veral ty pes of cyb e r attacks, it has poor statisti- cal efficiency u n der thick-tailed n o n-Gaussian measuremen t noises. I n ad dition, du e to the inherent limitation of the EKF , th e GM- IEKF m ay fail to c o n verge when the system is operating under stressed conditions, tha t is, when exhibiting strong model nonlinear ities. In Part I [6], we develop a Genera lized Maximu m- Likelihood-ty pe ro bust UKF-b ased DSE, termed th e GM-UKF for short. W e show that our GM-UKF is able to han d le thick - tailed non-Gaussian measuremen t noises with g ood statistical efficiency and is robust to o bservation a nd innovation ou tliers. In this second part, we will fo c us on its implementation . Specifically , we will discuss how to set and tune the pa r ameters of the GM-UKF alon g with th e choice of state initializatio n of the algo rithm that solve for the filter . W e will th en evaluate the perfor mance o f our GM-UKF th rough the following scenarios: i) the occ u rrence of two different types of n oises f ollowing Laplace or Cauchy pro bability distributions for the real and reactive power measurem ents; ii) the occu rrence of ob ser- vation a n d innovation ou tliers; iii) the o c currence of PMU measuremen t losses due to occasional com munication failures; iv) cyber attacks; and v) strong system n onlinearities. WE perfor m compar iso n s between our GM-UKF and the UKF and the GM-IEKF p roposed in [4]. W e show that our GM-UKF outperf orms the GM-IEKF and the UKF for all the scenarios being con sidered. When the system is operating under stressed condition s, our GM-UKF converges while the GM- IEKF a nd the UKF diver ge. Furtherm ore, if the power measureme n t noises follow a Cauchy distribution, the UKF fails to co n verge; although the GM -IEKF can ha n dle that case, it has a muc h lower statistical relative ef ficiency with respect to the GM- UKF . In this paper, the statistical ( asymptotic) efficiency of an estimator at a g i ven proba bility distribution (e.g ., Gaussian, or L a placian, or Cauchy distribution) is defined as the ratio between the inverse o f the Fisher infor m ation evaluated at that distribution and the (asympto tic) variance of the no rmalized estimator when all th e assumptions un derlying th e model are exact. As for the robustness of an estimator to outliers, Hampel [7] pro poses to in vestigate how th e asymp totic bias an d the asymptotic variance of the estimator in crease with the f raction of contaminatio n, ǫ , termed bias- and v ar ian ce-robustness analysis. T o quantif y bias-r o bustness, he intro duces the (b ias- IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 2 Dynamical System Model • Nonlinear system equation Robust prewhitening • T empora l correlation analysis • PS value ca lculations • Statistical test to PS values • Downweight outliers Robust regression • Huber -type GM-estimator • Iteratively reweighted so lution Robust error covariance matrix updating • T otal influence function • Asymptotic variance matrix • Asymptotic normal distribution of the GM-estimator Next snapshot Initialization Batch-mode regression form • State/measurement predictions • Statistical linearization of observation equation Fig. 1. Flo wchart representing the four main steps of the proposed GM-UKF . )breakd own point, the influence fun c tion and the asympto tic maximum bias curve while to quantify variance-robustness, he introdu c es the chan ge-of- variance f unction, wh ich measur e s the sensitivity of the asym ptotic variance to an infin itesimal change in ǫ abo ut zer o. A r o bust estimator h as a finite bias and a finite variance when subjec t to contam ination up to the breakd own point. Using simulations car r ied ou t o n the IEEE 39-bus system, we will show th a t o ur G M -UKF satisfies all these robustness and efficiency pro p erties. The paper is organized as follows. Section II deals with the modeling of power system d y namics, imple m entation issues and a systematic way to tune parameters. Section III shows the test results u nder se veral scenarios using the deta iled two- axis generator m o dels. Finally , Section IV co ncludes the paper and presents some interesting future research directions. I I . A L G O R I T H M I M P L E M E N TA T I O N The developed ro bust GM-UKF is a g eneric technique fo r online mo nitoring of many dynam ical cybe r-phy sical systems, including smart grids, au tonomo us vehicles, airc raft tracking, GPS trackin g and navigation, rad ar systems, to nam e a few . In this paper, we take power system dyna mic state estimation as an illustrativ e example to demo nstrate the capabilities of our GM-UKF for supp ressing observation and innovation o utliers while being ab le to filter o ut th ick-tailed non -Gaussian noise. A. Non linear Discr ete-T ime P ower System Dyna mical Mod e l For an electr ic power system, its discrete- time state space representatio n can be expressed as x k = f ( x k − 1 , u k ) + w k , (1) z k = h ( x k , u k ) + v k , (2) where the state vector x k contains the rotor an gle, the r o tor speed, the d - and q- axis state variables of the synchro nous generato r, th e exciter , the voltage regulator , and the governor . Here, u k represents the input vecto r; h ( · ) is the vector-v alued measuremen t function while f ( · ) is the vector-v alued fun ction that relates x k to x k − 1 ; z k is the measurement vector that contains a c o llection of voltage phasors, curren t phasors, real and reactiv e power flows an d power injections so that the system dynamical mo del is observable. The noises w k and v k , which may be non-Gau ssian n oise, ar e assumed to b e white and indepen dent of each other . In this pap er , the detailed 9th order two-axis generator mode l with IEEE - DC1A exciter and TGO V1 turbine-g overnor is assumed and tested [8]. B. I mplementation o f the GM-UKF The flowchart of the pr o posed GM-UKF is sh own in Fig. 1. It consists of fou r majo r steps, n amely a batch -mode regression form step , a robust pre-whiten ing step, a r obust regression and ro bust error covariance matr ix up d ating steps. Specifically , after state initialization and the application of statistical linearization to the n onlinear system pro c ess model, we calculate th e pred icted state and its associated covariance matrix. Next, applying statistical linearization to the nonlinear observation fu nction around the predic ted state, we d eriv e the expression o f the predicted measurement and its cov arianc e matrix. Th e n, by proc e ssing the obser vations an d pr e d ictions simultaneou sly , we obtain the batch-m ode regression form. Next, we apply the PS to a matrix that consists of two-time sequence of the pred icted state an d innovation vectors to detect the presence of any observation and innovation outliers. This in turn allows u s to carry out a robust prewhitening of the regression model. T o sup p ress outliers an d filter ou t thick- tailed non-Gau ssian measur ement noise, the GM-estimator is used and so lved by m eans of the Iteratively Reweighted Least Square s (IRLS) alg orithm. Finally , th e total influ ence function of our GM-U KF is d eriv ed and utilized to derive the asymptotic state estimation err or cov arian ce matrix. Note that during the iterative solution of the GM -UKF , we advocate to use the W eighted Least Squares (WLS) estimation for th e first iteration an d then switch to the IRLS algorithm . By d o ing so we improve the con vergence speed of that algorith m . Remark. Cor ollary 3.1 in P art I of the two-part series states that the matrix Z r oughly follows a bivariate Gaussian distribution. Using th at pr operty , we ca rry ou t extensive Monte Carlo simulations an d QQ-p lots to determine the pr o bability distribution of the PS. F r om Fi g. 3 in P a rt I, we infer that they IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 3 appr oximately o bey a chi-squa r es distribution with 2-degr ees of freedom. Consequ ently , we set the outlier detection thr esh- old of th e sta tistical test applied to the PS to η = χ 2 2 , 0 . 975 at a signifi cance level of 97.5 %. The deta iled implementatio n pr ocedures of th e PS algorithm ar e pr esented in App endix A . C. T uning the P arameters of the GM-UKF The tun in g o f the GM-UKF inv olves th e settings of the breakp o int λ o f the Huber ρ - function , the par ameter d o f the weighting fu nction, an d the c o n vergence toler ance o f th e IRLS algorithm. λ determines the trade-off tha t we wish to achieve be twe en a least-squar es and a least-absolu te-value fit. Indeed , when λ → 0 , the Huber ρ -function tends to the least-absolute-value ρ -fun ction and when λ → ∞ , it tend s to the least-sq u ares ρ -fun ction. Regarding th e parameter d , it determines the statistical efficiency of the PS at the assumed probab ility distribution alo ng with th e ro bustness of the GM- estimator [9]. Decre asing this param e ter too much shrin k s the dimen sions o f the 97.5% confid ence ellipse. As a result, good measuremen ts may be unduly downweighted, which yields a d ecrease in th e statistical efficiency . On the othe r hand, increasing d will incr ease th e bias of the GM-estimato r . Extensive simulations have shown th at th e paramete r s λ and d can be set to 1.5 to achiev e a good statistical efficiency at the Gaussian, the Laplacian, and the Cauchy distributions while achieving a go o d ro bustness to outliers. Regarding the conv ergenc e toleranc e threshold of the IRLS algorith m, a typical value is 0 . 0 1 ; decre asing this value results in sm all incrementa l ch a nges of the state estimates wh ile incr easing the compu ting time of the algorithm . Proposition 1. The coefficient of the estimation err o r covari- ance ma trix of o ur GM-UKF expr essed as E Φ [ ψ 2 ( r S ) ] { E Φ [ ψ ′ ( r S )] } 2 is equal to 1.0369 for the Huber cost function with λ =1.5. Pr oof. See the p roof in Appendix B. I I I . N U M E R I C A L R E S U LT S The performa n ces of our GM-UKF to handle thick -tailed non-Gau ssian noise a long with innovation and observation out- liers ar e assessed on the IEEE 39-bus system. The UKF a n d the GM-IEKF [4] are in c luded for co mparisons. The time-domain simulation results are used to generate a collection of sample s of th e nodal voltage magn itu des a n d ph ase ang le s as well as of the real an d reactiv e power injections at the termina l buses of all the gen erators. A samp ling rate of 50 sample s/seco n d is assumed . A syn th etic noise is added to the tru e values following th e pr obability distributions displayed in Fig. 1 of Part I. Specifically , a ze r o mean Gau ssian n oise is assum e d for the v oltag e angles, a bimo d al Gaussian mix ture d istribution is assumed for the noise of the voltage m agnitudes, and either a Laplace or a Cauch y distribution is assumed f or th e no ise of the r eal and reactive power measurem ents. Note that the random variable that follows a Gaussian mixture distribution is ge n erated via Matlab functio ns; the Cauchy ran dom variable  is ob tained b y sampling the inverse cumulative distribution function of the distribution given by  = β + α · tan ( π ( U 1 − 0 . 5)) , (3) 0 2 4 6 8 10 0 10 20 30 40 50 δ 5−1 in degree 0 2 4 6 8 10 0.9995 1 1.0005 ω 5 in pu 0 2 4 6 8 10 1.7 1.8 1.9 E fd5 in pu 0 2 4 6 8 10 4.95 5 5.05 time/s TM 5 in pu true value GM−IEKF UKF GM−UKF (a) 0 2 4 6 8 10 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 time [seconds] Mean Absolute Error [per unit] GM−IEKF UKF GM−UKF (b) Fig. 2. Case 1: Tracking performanc e of the GM-IEKF , UKF , and GM-UKF without outlie rs; (a) the estimate d rotor angle and s peed, field vol tage and mechanic al power of Generator 5 are used for illustra tion purpose s; (b) mean absolute error of each of the three filters. where β is the lo c a tion and α is the scale parameter ; U 1 are values r a ndomly sampled from the uniform distribution on the interval (0,1); sample s obeying the Laplace ran dom variable ζ with m ean µ and scale b are genera ted using ζ = µ − b sgn ( U 2 ) ln (1 − 2 | U 2 | ) , (4) where U 2 is a rand om variable drawn from the unif o rm distribution in the interval (1 /2, 1/2 ]. The two-axis generato r model is assumed and tested, whose parameters are ta ken fr o m [8]. A d isturbance is applied at t =0.5s by o pening the transmission line between Buses 15 IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 4 and 16. The maximal number of iterations allowed for the IRLS a lg orithm is 20 . For th e state initialization, the steady - state values with 10% e r rors ar e used. Due to space limitation, not all th e 9 state variables of each gener ator are shown; instead estimated values of the rotor angle and spe e d, the field voltage and th e mecha n ical power of Generato r 5 are utilized for illustration pu rposes. The mean absolu te err or ( MAE) is utilized as the index to ev aluate the overall p erforman ce of each method. A. Case 1: Thick-tailed Non - Gaussian Mea sur ement Noise without Outliers In this section, we first e valuate the performan ce of th e GM- IEKF , the UKF , and the GM- UKF und er norm al conditio ns. Specifically , a zero mean Gaussian noise with a standard deviation o f 10 − 2 is add ed to the voltage angles; the noise of the voltage mag nitudes follows a bimodal Gau ssian mixtu re with zero mean, variances of 10 − 4 and 1 0 − 3 and we ights of 0.9 and 0 .1, respectively; Laplace n oise with zero mean and scale 0.2 is a d ded to the real and reactive power in jections measuremen ts. The test results are d isplayed in Fig . 2. It is observed th at the UKF is not ab le to cope with Laplacian noise ev en in the absence of outliers. By con trast, the GM-IE KF and the GM- UKF ca n filter out such a n o ise while achieving good tracking pe rforman ce. Ho wever , the GM-I EKF has much lower relativ e statistical efficiency with respect to our GM-UKF . In particular, the GM-IEKF poorly estimates the field voltage and the roto r speed . By o bserving Fig. 2, it is interesting to note that the turbine mechanical power is changing during the transient p rocess. Acco r ding to the CIGRE report [ 1 0], it can significantly vary when co n trol f eatures such as fast valving or special protection sch emes are u sed to lim it th e output of the steam driven g e n erator du ring tra n sients. Consequ ently , it is of vital impo rtance to not assume it to be fixed at a con stan t steady-state value as commonly don e in most of the literatu re, but to obtain accur ate dynamic state variables of the governor for controls and stability analysis. B. Case 2: Thick-tailed Non - Gaussian Mea sur ement Noise with Observation Outliers The settings are the same as those of Case 1 except fo r the pr esence of observation outliers f rom t =4s to t =6s. The latter are simulated by add ing 20% errors to th e real and the reac ti ve power measurements of Generator 5. The r e sults are presented in Fig. 3. Fro m th is fig ure, we observe that the UKF is no t robust to o bservation outliers since it y ie ld s significantly b iased results. Althoug h the GM-I EKF can handle them, it produces increased biases o n the e stimates at th e time when ob servation outliers occur ( see the estimated ro tor speed and the field voltage for example). By con trast, the GM-UKF suppresses th e o utliers and prod uces much less bia s than the GM-IEKF . Note th a t the Gau ssianity o f the GM-estimator u sed in the estimation step of the GM- UK F allows that method to filter out th ick-tailed no ise while its statistical r obustness enables it to sup p ress the outlier s, hence ach ieving very go od estimates. 0 1 2 3 4 5 6 7 8 9 10 15 20 25 30 35 δ 5−1 in degree 0 1 2 3 4 5 6 7 8 9 10 0.998 0.999 1 1.001 ω 5 in pu 0 1 2 3 4 5 6 7 8 9 10 1.75 1.8 1.85 1.9 1.95 E fd5 in pu 0 1 2 3 4 5 6 7 8 9 10 5 5.2 5.4 5.6 5.8 time/s TM 5 in pu 5 6 7 8 5.02 5.04 5.06 true value GM−IEKF UKF GM−UKF (a) 0 2 4 6 8 10 0 0.05 0.1 0.15 0.2 time [seconds] Mean Absolute Error [per unit] GM−IEKF UKF GM−UKF (b) Fig. 3. Case 2: Trackin g performance of the GM-IEKF , the UK F , and the GM-UKF in the presence of observ ation outliers from t =4s to t =6s, where the real and reacti ve po wer measurements of Generat or 5 are corrupted with 20% errors; (a) the estimated rotor angle and spee d, field voltage and mech anical po wer of Generato r 5 are used for illustrati on purposes; (b) m ean absolute error of each of the three fi lters. C. Case 3: Thick-tailed Non-Ga ussian Measu r ement Noise with Innovation Ou tliers The settings are the same as those of Case 1 except fo r the presence of in novation outliers from t =4s to t = 6s. They are simu lated b y ad ding 2 0% er r ors to the predicted rotor angle of Gen erator 5. Th is innovation outlier is induced by a gro ss par a m eter v alu e in th e model. The co mparison results are shown in Fig. 4. As expe c ted, du e to the non -robustness of the UKF , it is un able to handle innovation outliers. On the other han d, th e GM-I EKF can h andle them, b ut it p roduce s larger biases comp ared with Case 2. Th is can be e xp lained by the fact tha t the mo d el errors will not on ly affect the predicted state vector but also will produce a smearing effect throug hout the Jacobian matr ix. As a result, it downweights IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 5 0 1 2 3 4 5 6 7 8 9 10 20 30 40 50 δ 5−1 in degree 0 1 2 3 4 5 6 7 8 9 10 0.99 0.995 1 ω 5 in pu 0 1 2 3 4 5 6 7 8 9 10 1.2 1.4 1.6 1.8 2 E fd5 in pu 0 1 2 3 4 5 6 7 8 9 10 5 5.5 6 time/s TM 5 in pu true value GM−IEKF UKF GM−UKF 4 5 6 0.9995 1 1.0005 4 5 6 7 8 5.02 5.04 5.06 true value GM−IEKF UKF GM−UKF (a) 0 2 4 6 8 10 0 0.05 0.1 0.15 0.2 0.25 time [seconds] Mean Absolute Error [per unit] GM−IEKF UKF GM−UKF (b) Fig. 4. Case 3: Trackin g performance of the GM-IEKF , the UK F , and the GM-UKF in the presence of innov ation outliers from t =4s to t =6s, where the predic ted rotor angle of Generator 5 is corrupte d with 20% errors; (a) the estimated rotor angle and speed, field volt age and mechanical power of Generat or 5 are used for illustratio n purposes; (b) mean absolute error of each of the three filters. se veral g ood m easurements. By contrast, ou r GM-UKF is capable of handling both observation and innovation outliers, yielding compar able performan ces. Because only the sigma points associated with the mode l errors will be a ffected an d downweighted, this filter o btains better estimates than the GM- IEKF in presence of model errors. D. Case 4: Thick-tailed Non -Gaussian Measur emen t No ise with Measur ement Losses The settings are th e same as those of Case 1 except for the losses of PMU measur ement from t =5s to t =8s; specifically , all the PMU m e asurements at the terminal bus of Generato r 5 are lost d ue to co mmunica tio n failures or cyber attacks. The test results are presen ted in Fig. 5. Due to its lack o f robustness, the UKF is unable to hand le th e loss of PMU measuremen ts 0 1 2 3 4 5 6 7 8 9 10 20 30 40 50 δ 5−1 in degree 0 1 2 3 4 5 6 7 8 9 10 0.998 0.999 1 1.001 ω 5 in pu 0 1 2 3 4 5 6 7 8 9 10 1.75 1.8 1.85 1.9 1.95 E fd5 in pu 0 1 2 3 4 5 6 7 8 9 10 5 5.2 5.4 5.6 time/s TM 5 in pu true value GM−IEKF UKF GM−UKF 5 6 7 8 9 10 5 5.05 5.1 (a) 0 2 4 6 8 10 0 0.05 0.1 0.15 0.2 time [seconds] Mean Absolute Error [per unit] GM−IEKF UKF GM−UKF (b) Fig. 5. Case 4: Tracking performance of the GM-IEKF , the UKF , and the GM- UKF in the presence of PMU measurement losses from t =5s to t =8 s, where all the termina l m easurements of Generator 5 are lost; (a) the estimated rotor angle and s peed, field voltage and m echanical po wer of Generator 5 are used for illustrat ion purposes; (b) mean absolute erro r of each of the three filters. because in that ca se, only noise is recei ved an d taken as PMU measuremen ts. As for th e GM-I E KF an d the GM-UKF , thanks to their ro bustness an d to the batch-mo de regression fo rm that provid es en hanced data redu ndancy , they will rely on the majority of pred icted states and the good measuremen ts to filter o ut the n oise and strongly downweigh t the lost PMU measuremen ts, which are flag g ed as outliers by the PS. As a result, they b oth achieve reasonab le state estimates, with an clear advantage for the GM -UKF since it exhibits smaller mean absolute e r ror . E. Case 5 : Ha ndling Cauchy P ower Measur ement Noises From Part I, it is obser ved that th e real and the rea c tive power measurem ent noises may ob ey a Cauchy distribution, which is a very thick tailed distribution with n o moments IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 6 0 1 2 3 4 5 6 7 8 9 10 14 16 18 20 22 δ 5−1 in degree 0 1 2 3 4 5 6 7 8 9 10 1 1.0002 1.0004 ω 5 in pu 0 1 2 3 4 5 6 7 8 9 10 1.8 1.85 1.9 E fd5 in pu 0 1 2 3 4 5 6 7 8 9 10 5 5.05 5.1 time/s TM 5 in pu true value GM−IEKF GM−UKF (a) 0 2 4 6 8 10 0 0.005 0.01 0.015 time [seconds] Mean Absolute Error [per unit] GM−IEKF GM−UKF (b) Fig. 6. Case 5: Trackin g performance of the GM-IEKF , the UK F , and the GM-UKF in the presence of Cauchy po wer measurement noise. The estimat ed rotor angle and speed, field voltage and mechanica l po wer of Generator 5 are used for illustra tion purposes; (b) mean absolute error of each of the three filters. Since the UKF div erges, its resul ts are not sho wn in the fi gure. being defined . T o test the capability of our GM-UKF to handle that case, we assum e that the simulation settings ar e the same as those of Case 1 except that the Cauchy n oises with zero media n and a scale of 0.0 05 is added to the real and the reactive power injection measur e m ents. The obta in ed test results are displaye d in Fig . 6. Note th at in pr esence of Cauchy measureme nt noise, the UKF has a no n-positive definite covariance matrix , resu lting in its divergence. By contrast, thank s to the Gaussian n ormality and robustness of our GM - UKF , the total influence functio n-based covariance matrix up dating appro ach ca n always g u arantee its p ositiv e- definiteness. On the other ha n d, compar ed w ith the r e sults obtained when u sing Laplacian power m easurement noises, the GM-IEKF prod uces larger biases of the state estimates an d 0 1 2 3 4 5 6 7 8 9 10 −50 0 50 δ 5−1 in degree 0 1 2 3 4 5 6 7 8 9 10 0.98 0.985 0.99 0.995 1 ω 5 in pu 0 1 2 3 4 5 6 7 8 9 10 −4 −2 0 2 4 E fd5 in pu 0 1 2 3 4 5 6 7 8 9 10 −5 0 5 10 time/s TM 5 in pu true value GM−IEKF UKF GM−UKF (a) 0 2 4 6 8 10 0 1 2 3 4 5 6 time [seconds] Mean Absolute Error [per unit] GM−IEKF UKF GM−UKF (b) Fig. 7. Trackin g perfo rmance of the GM-IEKF , the UKF , and the GM-UKF in the presence of strong system nonli nearity . (a) The estimated rotor angle and speed, field volt age and mechanical power of Generato r 5 are used for illustra tion purposes; (b) mean absolute error of ea ch of the three filters. takes much longer time to ap p roach the true system states. This is not the case for our GM-UKF since it achieves a compa rable perfor mance and tracks the sy stem state at the very beginnin g of the transient process. F . Ro bustness to Str o ng System Nonlinearity Practical systems may be he avily loaded, resulting in strong nonlinear dynamics. T o illustrate the c apability of o u r GM- UKF to handle that case, we assume tha t the load at Bus 7 is in c reased from 2 33.8 MW to 1 500 MW to stress th e system bef ore switching Lin e 15- 16 while the other simulation settings are the same as those of Case 1. Note that the steady- state m a x imum loadability at Bus 7 is around 2 000 MW . After the line switching, the system op erates u n der even g reater stressed conditions. The test results are disp lay ed in Fig. 7. It is observed fr om these two figur es that th e GM-IEKF fails IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 7 T ABLE I A V E R AG E C O M P U T I N G T I M E S O F T H E T H R E E D S E M E T H O D S F O R E V E RY P M U S A M P L E , W H E R E N A R E P R E S E N T S N O T A P P L I C A B L E . Cases UKF GM-IEKF GM-UKF Case 1 6.28m s 9.64m s 9.52m s Case 2 6.31m s 9.68m s 9.55m s Case 3 6.38m s 9.72m s 9.63m s Case 4 6.36m s 9.70m s 9.59m s Case 5 N A 9.80m s 9.70m s to conv erge at the very b eginning of the transien t pr ocess wh ile the UKF div erges around t =8s. By contrast, ou r GM-UKF is able to handle this scenario while achie vin g excellent tracking perfor mance. The und erlying reasons are as follows: • Under stressed system o perating conditions, the first- order T aylor series expansion used in the EKF and the GM-IEKF is too app roximate and is un able to account for strong system nonlinearities. As a resu lt, these two filters produ ce very la rge approx imation error s and eventually div erge; • Thanks to the sigma-po ints-based unscen ted tran sfor- mation a n d its approx imation accuracy up to at lea st third-or der T ay lo r series expan sion, both the UKF and the GM-UKF ar e able to handle stron g system nonlin - earities. However , due to the accumulative estima tio n error induced by th e non- Gaussian measuremen t error, the estimation error covariance matrix o f the UKF is close to no n-positive defin iti veness. T herefor e , it pro - duces large estimation error s an d finally d i verges. By contrast, our GM - UKF le verages the stren gth of th e unscented transfo r mation to handle system no nlinearities while the ro bustness of the GM-estimator allows it to filter out thick-tailed non- Gaussian measurement noise, yielding good estimation results. G. Breakdown P oint of the GM-UKF to Cyber Attacks W ith the strong reliance o f smar t grid fu nctions on com- munication networks, cyber attacks have becom e a m ajor concern . T ypically , they are classified as b ias injection attack, denial of service attack, and rep lay attack [11], [12]. Bias injection attack occurs wh en an adversar y attempts to cor rupt the content of either the m easuremen t o r the c o ntrol signals; for examp le, the man-in -the-midd le inter c e pts the PMU mea- surement sign a ls and cor rupts them with large b iases. Denial of ser v ice attack occur s when the actuator and sensor data are prevented from r e a ching their respective d estinations, resulting in the absence of data fo r th e DSE; f or instance , this will be the case if the PMU metered values do not r each the phasor data concen trator . Replay attac k occurs when a hacker first perfor ms a disclosure attack from a certain time period , gathering seq uences of data, and then begins replay ing the data during a certain p eriod; for in stance, the cu rrent PMU measuremen ts pr o cessed b y a dynamic state e stima to r are replaced b y past values. In other words, those attacks ind u ce observation or innovation outliers. T o in vestigate the break down point o f the GM- U KF to cyber attacks, which is defined as the maximum nu mber of outliers that the filter can handle without yield ing unreliable estimates, we carry o u t extensive simulations on the IEEE 39- bus test system using th e concept of finite sample break down in n onlinear r egre ssion introduce d by Stromberg and Rup pert [13]. By replacing a varying percentag e of o bservations b y outliers in th e vector y k , it is observed th a t the GM -UKF can handle at least 25% of corrupted observations. It is worth noting that the breakdown point o f the GM-estimator in n onlinear regression is still u n known. This pro blem will be investi g a te d as a f u ture work. An other intere stin g problem is the deter mination of the maximu m breakdown po int that any regre ssion estimator m ay have in structur e d non linear regression such as power system state estimatio n problem s; this will b e an interestin g extension of the r esults proved in Mili and Coakley [14] in the linear case. H. Comp utationa l Efficienc y T o v alidate th e applicability o f the p roposed GM-UKF to online estimation with a PMU samp lin g rate of 30 o r 60 samples per second, its comp utational efficiency is analyzed and compar ed to that of the UKF and th e GM-IEKF in Cases 1 -5. The test is perfo rmed on a PC with Intel Core i5, 2.50 GHz, 8GB of RAM. The average computin g tim e of eac h metho d for every PMU samp le is displayed in T able I. W e observe from th is ta b le that the UKF has the best computatio nal efficiency , exhib iting compu tin g times m uch lower than the PMU sampling period , which ar e 33 .3ms and 16.7ms for 30 sample/s and 60 samples/s, respectively . Although the execution times o f the GM-IEKF and the GM- UKF is lon ger , they are still smaller than the PMU sampling period, demonstratin g their ab ility to track system real-time dynamic states. I V . C O N C L U S I O N A N D F U T U R E W O R K In this second part of the two-paper series, the prop osed GM-UKF is im plemented , tested, an d validated. V ario us sce- narios have been con sidered and explo red to evaluate its perfor mance with Laplacian an d Cauchy measure m ent n o ises, observation and innovation outliers, and strong system non- linearities. Its brea kdown poin t to cyber attacks has also b een in vestigated. Compar iso n results with existing meth o ds show that the GM-UKF outperfor ms the GM-IEKF and the UKF in all the simulated scenar ios. I t is interesting to note that when the system is operating u nder stressed condition s, the GM-IEKF and the UKF fail to conv erge while ou r GM-UKF conv erges. Furthermor e, if the power mea surement no ises fo l- low a Cauch y distribution, the UKF fails to converge while th e GM-IEKF achieves much lower relativ e statistical ef ficiency with r espect to our GM-UKF . There are different possible a venues to further inv estigate the stud y consider ed in this pap er . The propo sed centralized GM-UKF can h andle observation and the innovation outliers, but pr ovides p o or results in presence of stru ctural outliers. The latter may b e induced by gr oss errors in circuit b reaker statuses or in the pa rameters o f the turbine-ge n erators and the transmission lines. T o address this problem , we will inv estigate IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 8 a dece ntralized GM-UKF , wh ich will b e imple m ented at the generating unit level usin g local voltage and curren t phaso r measuremen ts. Here, ad ditional av ailable measurements on rotor speed , ter minal-bus real and r eactiv e power , and field voltage an d curren t can be utilized for improving the measu r e- ment redunda ncy . Further m ore, we will d evelop a ge neralized GM-UKF for simu ltaneously e stima tin g the system states and mod el par ameters wh ose values are either inaccurate o r incorrect. Furthermo re, we will inves tigate the case where the GM-UKF fails to produ ce g ood results due to very strong system nonlinear ities. The dev elo p ment of the GM-particle filter can b e a g ood cand idate to hand le that case. Finally , we will extend the prop osed GM-UKF to the generator model calibration and validation, the estimation o f dyna mic lo ad model parameter s and power system oscillatory mo des. A P P E N D I X A P R O J E C T I O N S TA T I S T I C S A L G O R I T H M The main steps of impleme nting the projec tio n statistics algorithm are shown as follows: • Step 1 : For a point l i in an n -dimension al spa ce, calculate the coordinate- wise median giv en by M =  med j =1 ,...,m ( l j 1 ) , ..., med j =1 ,...,m ( l j n )  , (5 ) where m is the numbe r o f points; • Step 2 : Calculate the directions f o r projec tio ns u j = l j − M , j = 1 , ..., m ; • Step 3 : No rmalize u j to g et ℓ j = u j k u j k = u j q u 2 j 1 + ...u 2 j n ; j = 1 , ..., m ; (6) • Step 4 : Calcu late the standard ized projections of the vectors { l 1 , ..., l m } on ℓ j , which are giv en by ζ 1 j = l T 1 ℓ j ; ζ 2 j = l T 2 ℓ j ; ..., ζ mj = l T m ℓ j ; (7) • Step 5 : Calcu late the m edian o f { ζ 1 j , ..., ζ mj } = ζ med,j ; • Step 6 : Calculate th e med ian absolute d eviation (MAD) M AD j = 1 . 4826 · b · me d i | ζ ij − ζ med,j | , where the correction factor is b = 1 + 15 / ( m − n ) ; • Step 7 : Calculate the stan dardized p rojections P ij = | ζ ij − ζ med,j | M AD j f or i = 1 , ..., m ; (8 ) • Step 8 : Repea t step s 4 –7 f or a ll vectors { ℓ 1 , ..., ℓ m } to get the standar dized pr ojections { P i 1 , ..., P im } f or i = 1 , ..., m ; • Step 9 : Calcu late the p rojection statistics P S i = max { P i 1 , ..., P im } f or i = 1 , ..., m. (9) A P P E N D I X B P R O O F O F T H E P RO P O S I T I O N 1 Pr oof. In the previous companio n pap e r, the estimation erro r covariance matrix P xx k | k has b een shown to follow asympto tic a Gaussian distribution and since the covariance o f ξ k is an identity matrix, the standard ized residual r S therefor e follows asymptotic a norm al distribution. As a consequen ce, the pr ob- ability d istribution fu nction of the stand ardized re sid u al can be expressed as φ ( r S ) = 1 √ 2 π e − r S 2 2 . On the o ther han d, from the Huber function with λ = 1 . 5 , we can calculate ψ ( r S ) =  r S for | r S | ≤ λ λsig n ( r S ) for | r S | > λ , (10) ψ ′ ( r S ) =  1 for | r S | ≤ λ 0 for | r S | > λ . (11) Then, we can f urther obtain E h ψ ′ ( r S ) i = Z ∞ −∞ ψ ′ ( r S ) φ ( r S ) dr S = 1 √ 2 π Z ∞ −∞ e − r 2 S 2 dr S = 2Φ ( λ ) − 1 = 0 . 8 664 , (12) E  ψ 2 ( r S )  = Z ∞ −∞ ψ 2 ( r S ) φ ( r S ) dr S = λ 2 √ 2 π Z − λ −∞ e − r 2 S 2 dr S + 1 √ 2 π Z λ − λ r 2 S e − r 2 S 2 dr S + λ 2 √ 2 π Z ∞ b e − r 2 S 2 dr S = λ 2 Φ ( − λ ) − 2 λ √ 2 π e − λ 2 2 + 2Φ ( λ ) − 1 + λ 2 (1 − Φ ( λ ) ) = 0 . 778 4 . (13) Finally , we can calculate E  ψ 2 ( r S )  ( E [ ψ ′ ( r S )]) 2 = 0 . 7784 (0 . 8664) 2 = 1 . 0369 . (14) A C K N OW L E D G M E N T The autho rs would like to thank Dr . Zhenyu Huan g from PNNL for p roviding us with real PMU d ata for an alyzing the statistical prob ability distributions o f the PMU measur e ment errors. R E F E R E N C E S [1] A. Rouhani, A. Abur , “Linear phasor estimator assisted dynamic state estimati on, ” IEEE T rans. Smart Grid , 2016. [2] X. W ang, N. Cui, J. Guo, “Huber- based unscente d filtering and its applic ation to vision-ba sed relati ve navigati on, ” IET radar , sonar and navigat ion , vol. 4, no. 1, pp. 134–141, 2010. [3] L . Y . Chang, B. Hu, G. K. Chang, et al, “Huber -based nove l robust unscente d Kalman filter , ” IET Scien ce, Measur ement & T echno logy , vol. 6, no. 6, pp. 502–509, 2012. [4] J. B. Zhao, M. Netto, L. Mili, “ A robu st iter ated extend ed Kalman filter for po wer system dyna mic state estimat ion, ” IEE E T rans. P ower Syst. , DOI: 10.1109/TPWRS.2016.2628344, 2016. [5] J. B. Zhao, L. Mili, A. Abdelhadi, “Robu st dynamic state estimator to outlie rs and cyber attacks, ” in Pr oc. of IEEE PES General Meeting , July 16-20, Chicago 2017. [6] J. B. Zhao, L. Mili, Z. Huang, “Po wer system dynamic state estima- tion with non-Gaussian synchrophasor measurement noi se using robust unscente d Kalman filter –Part I: Theory , ” IEEE T rans. P ower Syst. , 2017. [7] F . R. Hampel, E. M. Ronchetti , P . J. Rousseeuw , and W . A. Stahel , Robust Stati stics: The A ppr oach Based on Influence F unction s . New Y ork: John W ile y & Sons, Inc., 1986. [8] IE EE PES T F on Ben chmark System for Stability Controls, “Benc hmark systems for small-signal stabilit y analysis and con trol, ” Aug. 2015. [9] M. Gandhi, L. Mili, “Robust Kalman filter based on a generalized maximum-lik elihood-typ e estimator , IEE E T rans. Signal Pr ocessing , vol. 58, no. 5, pp. 2509–2520, 2010. [10] “Modeling of gas turbines and steam turbine s in combined cycle power plants, ” CIGRE repor t, T ask Force C4.02.25, 2003. IEEE TRANSACTIONS ON POWER SYSTEMS, V OL. , NO. , 2017 9 [11] A. T eixeira, D. Perez , H. Sandberg, K. H. Johansson, “ Attack models and scenari os for network ed control systems, ” Pr oceedings of the 1st ACM interna tional conf erence on H igh Confidenc e Network ed Systems , pp: 55-64, 2012. [12] J. B. Zhao, G. X. Zhang, M. La Scala, Z. Dong, C. Chen, J. W ang, “Short-te rm state forecasting-a ided m ethod for dete ction of smart grid general false data injecti on attacks, ” IEEE T rans. Smart Grid , 2016. [13] A. J. Stromber g, D. Ruppert, “Break down in nonli near regressio n, ” J ournal of the A merican Statistic al Association , V ol. 87, no. 420, pp. 991–997, 1992. [14] L. Mili and C. W . Coakle y , “Rob ust esti mation in s tructu red linear regre ssion, ” Annals of Statist ics , vol. 24, no. 6, pp. 2593–2607, 1996.