A Novel Detection Algorithm to Identify False Data Injection Attacks on Power System State Estimation

Abstract

This paper provides a novel bad data detection processor to identify false data injection attacks (FDIAs) on the power system state estimation. The attackers are able to alter the result of the state estimation virtually intending to change the result of the state estimation without being detected by the bad data processors. However, using a specific configuration of an artificial neural network (ANN), named nonlinear autoregressive exogenous (NARX), can help to identify the injected bad data in state estimation. Considering the high correlation between power system measurements as well as state variables, the proposed neural network-based approach is feasible to detect any potential FDIAs. Two different strategies of FDIAs have been simulated in power system state estimation using IEEE standard 14-bus test system for evaluating the performance of the proposed method. The results indicate that the proposed bad data detection processor is able to detect the false injected data launched into the system accurately.

1. Introduction

Integration of renewable energy resources into the power grid causes major deregulation in the electricity market. Keeping the generation cost down in the deregulated power system as well as providing a reliable load to customers necessitates the use of remote measurements and smart metering in the power system, which leads to smart grid development. A smart grid is able to establish a two-way communication by means of advanced sensors to improve the reliability and performance of a system. While these cyber resources enable market operators to enhance the performance of a system, they bring new challenges to the infrastructure. Cyber-physical systems (CPSs) render smart grids more vulnerable to cutting-edge undetectable cyber-physical attacks. Therefore, the impact assessment of cyber-attacks on smart grids is considered as one of the most serious issues in recent years [1].

Several types of cyber-attacks have been investigated, which are divided into three different categories, namely, data availability, confidentiality, and integrity. Cyber-attacks on availability or named denial-of-service (DoS) attacks effort to delay, deny, or alter transferring information at different layers of a communication network. Attacks on confidentiality aim to obtain unauthorized information from networks in smart grids. Attacks on integrity intend to deliberately and illegally modify or disrupt data exchange in a smart grid. Such attacks attempt to stealthily manipulate data to corrupt critical information exchange in smart grids [2]. Recently, the method of false data injection attack (FDIA) has been attracting the attention of researchers. The false data injection attack was designed to impact the state estimation [3] and real-time electricity market [4] by manipulating the data in the SCADA system. Another type of cyber-attack on a smart grid’s integrity, namely, alter-and-hide (AaH) is used to change the true values of digital variables such as substation switches and breakers. This kind of attack affects the topology processor and brings catastrophic consequences [5]. However, owing to the basic difference in the modeling of the various cyber-attacks, and the vast literature on this type of attack, the FDIA is modeled and formulated in this study to be investigated.

FDIAs aim to modify data at the control center or measurement units to reach a predetermined end. FDIAs target analog measured data from the power system, i.e., the nodal voltage magnitudes and angles, nodal power injections, line power flows, and digital data such as the status of breakers and switches [6]. Such data are essential to the power system operator (PSO) to perform state estimation (SE) in order to monitor and control the power grid operation. The result of SE is also exploited in the security constraint economic dispatch (SCED) to settle in the real-time power electricity market. Thus, the data integrity guarantees an accurate SE as well as effective monitoring and control of the power system operation.

The concept of FDIAs on SE in smart grids was first raised in [7]. The authors in [7] illustrated that by having adequate knowledge of the power system configuration, an adversary can launch FDIA without being detected. That investigation was an inspiration for other researchers to work on impacts of this type of attack on the power system and appropriate techniques to defend it. The authors in [8] presented analytical approaches to analyze the influence of FDIAs on AC state estimation. The authors in [9,10,11,12,13] discussed various methods of FDIA on SE in smart grids considering complete and incomplete knowledge about the topology of the considered power system.

As mentioned earlier, the estimated variables are utilized in SCED and real-time electricity market settlement. This issue can be considered as a motivation for the adversary to stealthily attack the SE by altering locational marginal price (LMP) values in the real-time electricity market. Wide areas of research have been carried out to investigate this matter. The economic consequences of FDIA on the real-time ex-ante market were discussed in [6]. FDIA on analog and digital measurements were studied in [14], and an optimized attack problem was presented. Given the Australian real-time electricity market mechanism, the authors in [15] proposed a compact FDIA-based approach to cyber topology. The researchers in [16,17] provided a thorough review of FDIAs on power system state estimation and their influence on the electricity market up to the year 2017.

Several methods in the literature have been developed to mitigate and detect FDIAs on SE in smart grids [18,19,20,21,22]. The authors in [18] developed a heuristic algorithm to find an undetectable vector of FDIA and correspondingly, an infinite norm bad data detection (BDD) processor is defined to detect the attack vector efficiently. Liu et al. [19] presented an approach utilizing the correlation and history of measurements as well as the sparse attribute of malicious attacks to detect FDIAs in smart grids. By applying the Kullback–Leibler distance (KLD) and assuming a probability distribution function for historical data, a statistical method for FDIA was introduced in [20]. Another statistical outlier idea was proposed in [21] to estimate the density ratio by using machine learning. By utilizing load forecasting, generation schedules and Phasor Measurement Units (PMUs), the author in [22] proposed an online bad data detection algorithm.

Recently, machine learning techniques have been widely applied for bad data detection and identification in smart grids. The authors in [23] proposed a hybrid machine learning method to detect stealthy attacks in smart grids. Data was labeled in the preprocessing stage by a supervised algorithm followed by training the labeled data with a support vector machine (SVM). Similarly, a deep-learning-based approach for detection of FDIA on smart meters’ data was investigated in [24]. The model used a state vector estimator (SVE) and a deep-learning-based identification algorithm to prevent electricity theft. Multiple hidden layers were considered in this method to enhance the detection accuracy. A sparse principal component analysis and approximation approach were developed in [25] to identify a stealthy FDIA on a smart grid. The accuracy of the method depended on the sparsity of data and it was less efficient for sparse vectors. Wang et al. [26] investigated a data-centric paradigm to identify FDIAs in smart grids and the developed model employed a margin setting algorithm to categorize a huge amount of data. It was illustrated in [26] that, by increasing the false data rate, the detection accuracy is decreased. In order to detect sparse cyber-attacks in smart grids, the authors in [26,27] developed a deep-learning-based interval state estimation algorithm. A multi-layer stacked autoencoder was used in that method, in which the last layer detected anomalies in electric load forecasting.

To the best knowledge of the authors, nonlinear autoregressive exogenous (NARX) (a special configuration of ANN) has not been utilized to detect stealthy optimized FDIA on contingency analysis. In spite of the distribution system, a transmission system is facilitated with adequate measurement devices; hence, analytical methods reveal accurate results of the state of the system. Still, measurements at the transmission level are exposed to FDIAs. There are BDD processors in the state estimator procedure even though, under specific conditions, FDIAs can pass through them. This paper aims to develop a novel BDD procedure based on a recurrent architecture of ANN, namely, nonlinear autoregressive exogenous (NARX) model, to identify FDIAs that could not be observed by conventional BDD processors.

The NARX neural network (NARXNN) is a robust engine for times series prediction by using historical data. Thus, it is practical to apply this learning engine to predict step-ahead values of the state by considering measurement values and historical data as input variables. However, owing to the nonlinear characterization of the forecast, the output would contain subtle errors in comparison to real values and SE results by analytical methods. Meanwhile, the results of analytical methods can be utilized without any restriction to validate and evaluate the estimated state.

This paper represents potential FDIAs that can be launched on smart grids without being detected by conventional BDDs. Moreover, by introducing NARXNN and its prominent features, which renders it an attractive predictor engine for estimation of the states, a novel FDIA detector is introduced and evaluated. The rest of this paper is organized as follows. Section 2 briefly presents the main concept of state estimation in power systems as well as bad data identification tests. Section 3 provides a review and formulation of undetectable FDIAs by conventional methods. A novel ANN-based bad data detection method is introduced and formulated in Section 4. Finally, numerical analyses are provided in Section 5, followed by conclusions in Section 6.

2. Preliminaries

In this section, the state estimation process in power systems and real-time power market operation are briefly introduced.

2.1. DC State Estimation

A state estimator receives measurements and other necessary information from supervisory control and data acquisitions (SCADA) system in order to estimate the power system state. Based on the nonlinear SE, the relation between measurements $z\in {ℝ}^m$ and state variables $x\in {ℝ}^n$ is as follows [28]:

$$z=h\left(x\right)+e$$

(1)

where $h\left(x\right)$ is the nonlinear measurement function of $x$ and $e\sim N\left(0,R\right)$ ($e\sim N\left(0,R\right)$ indicates the normal distribution with zero mean and variance $R$). Even though developing an estimation model based on AC nonlinear measurements, i.e., $h\left(x\right)$, is accurate, in many cases, it is a time consuming, costly procedure and might not result in an optimal solution for large scale power systems. On the other hand, using linearized DC state estimation may be less accurate, but simpler and more feasible, and practical for real-time LMP calculation [16]. In addition, the undetectable yet impressive FDIA constraints are too complicated to be formulated in AC state estimation. Hence, in order to formulate the optimization problem for the most consequential FDIA the DC load flow model is utilized in this paper. With an acceptable accuracy, a nonlinear measurement is linearized around an operation point. Hence, the linearized state estimation is formulated as follows:

$$z=Hx+e$$

(2)

where $H$ is a Jacobian matrix that depends on the impedance of the network topology.

Assuming that the control center employs the standard weighted least squares (WLS) state estimator, the estimated state of the system $\widehat{x}$ is computed by minimizing an objective function.

$$j\left(\widehat{x}\right)={\left(z-H\widehat{x}\right)}^T{R}^{-1}\left(z-H\widehat{x}\right)$$

(3)

Hence, the estimated state is as follows:

$$\widehat{x}={\left(H^T{R}^{-1}H\right)}^{-1}{H}^T{R}^{-1}z$$

(4)

where $G=\left(H^T{R}^{-1}H\right)$ is defined as the gain matrix. Thereby, the estimated value of measurements is as follows:

$$\widehat{z}=H\widehat{x}=Kz$$

(5)

and correspondingly, the measurement residual vector is defined as follows:

$$r=z-H\widehat{x}=\left(I-K\right)z$$

(6)

The state variable obtained from (4) is used for computing the power flow vector of transmission lines $f$ via the following relationship:

$$\widehat{f}=Y\widehat{x}$$

(7)

where $Y$ is the admittance matrix. The estimated state of power flow is employed for contingency analysis. The congested lines are defined as a set of $\ell =\left\{l:{\widehat{f}}_l>f_l^{max}\right\}$, where $f_l^{max}$ indicates the maximum constraint of power flow line l. The result will be used for SCED and LMP computation in the real-time market [15,16].

2.2. Bad Data Detection

Bad data detection is one of the most important processes of a state estimator, by which unacceptable measurement errors are identified.

In this paper, the ${\ell }_2-$norm-based residual test and the largest normalized residual (LNR) are utilized for BDD mechanism [28].

The ${\ell }_2-$norm-based residual test is a common test to detect possible erroneous or abnormal data in the measurement set. The following inequality is applied for this test:

$${\Vert r\Vert }_2^2\lessgtr {\tau }^2$$

(8)

Provided that ${\Vert r\Vert }_2^2>{\tau }^2$ the measurement vector contains bad data and if ${\Vert r\Vert }_2^2<{\tau }^2$ the result can pass the BDD processor without any problem. Random measurement errors follow the normal distribution with zero mean. Therefore, ${\Vert r\Vert }_2^2$ follows Chi-Square distribution with $m-n$ degree of freedom and probability of success $p \left({\chi }_{\left(m-n\right),p}^2\right)$. Correspondingly, $\tau$ is predetermined according to the false alarm probability, namely, $q=1-p$.

Therefore, the following inequality should be true in order for the BDD processor to confirm the result of SE:

$${\Vert r\Vert }_2\le \tau$$

(9)

Another test to detect bad data in measurements is the largest normalized residual (LNR), which can be more accurate in some cases and is formulated as follows [28]:

$$r_i^N=\frac{\Vert r_i\Vert }{\sqrt{{\Omega }_{ii}}}\le \epsilon$$

(10)

where $\epsilon$ is a determined identification threshold and is the same for all measurements and $\Omega$ is defined as the residues covariance matrix, which is defined as follows:

$$\Omega =R-H{G}^T{H}^T$$

(11)

Provided that Equation (10) is not satisfied for each measurement, it is suspected as bad data. Passing both aforementioned tests can assure the state estimator that the measurement vector and the estimated state vector are dependable and the result can be utilized for the next procedure. The first test checks that the sum of all residual values for estimated variables is less than the predetermined value; consequently, the second test checks the residual value for each measurement individually.

3. Undetected False Data Injection Attack on State Estimation

In order to attack the state estimation, the adversary would inject an attack vector $a$ $\left(a\in {ℝ}^m\right)$ to the measurement vector $z$. Thus, the residual vector would be altered as follows:

$$r_a=z_a-{\widehat{z}}_a=\left(z+a\right)-K\left(z+a\right)=r+\left(I-K\right)a$$

(12)

Correspondingly, the power flow vector is manipulated according to the following equation:

$${\widehat{f}}_a=Y\left(\widehat{x}+c\right)=\widehat{f}+Yc$$

(13)

where $c$ indicates the deviation of estimated state vector as a result of the attacked vector $a$. It is indicated in [14] that the congested line has a direct impact on dual variables in the LMP equation. Consequently, each congested line set is associated with a unique LMP vector. Considering this theorem, most publications focused on FDIAs that target congested line sets without being detected by the conventional BDD. Similar to [11,14,15,16,17], this paper assumes that FDIAs are launched without being detected by changing congested lines from their original situations. According to ${\ell }_2-$norm-based residual test, in order to launch an unidentified attack, the following constraint should be satisfied:

$${\Vert r_a\Vert }_2\le \tau$$

(14)

In practice, ${\Vert r_a\Vert }_2\le {\Vert r\Vert }_a$; therefore, provided that ${\Vert r\Vert }_a\le \tau$, it can be guaranteed that the BDD process cannot detect the injected data. Equation (14) is extended by using the previous equations as follows:

$$\Vert r_a\Vert =\Vert z_a-H{\widehat{x}}_a\Vert =\Vert z+a-H\left(\widehat{x}+c\right)\Vert =\Vert z-H\left(\widehat{x}\right)+a-Hc\Vert \le \Vert r\Vert +\Vert a-Hc\Vert \le \tau$$

(15)

where $c\in {ℝ}^n$ is an arbitrary non-zero vector. It is shown in [19] that, in cases $a=Hc$, a perfect and unidentifiable attack is launched to the state estimation. However, it is not feasible to assume that the system topology is completely known to an attacker. The attacker has to ensure $\Vert a-Hc\Vert \le \delta$ to launch an imperfect attack to avoid detection with high probability.

On the other hand, according to the LNR test, the following condition should be satisfied for each array of the residual vector to pass a BDD process:

$${\Vert r\Vert }_i\le \epsilon \sqrt{{\Omega }_{ii}}$$

(16)

hence

$$\left|r\right|\le \xi$$

(17)

where $\left|r\right|$ is the vector of absolute values of residual vector and $\xi$ is the threshold vector. Equation (17), at the presence of an attack, is rewritten as follows:

$$\left|r_a\right|=\left|r+a-Hc\right|\le \xi$$

(18)

$$\left|r_i+a_i-H_{i}c\right|\le \left|r_i\right|+\left|a_i-H_{i}c\right|\le \xi$$

(19)

therefore, in order to have a hidden attack, the attacker should consider the following constraint:

$$\left|a_i-H_{i}c\right|\le {\vartheta }_i$$

(20)

where $\vartheta$ is the threshold vector for the residual vector to have an undetected attack.

As shown in Figure 1, the attacker may either attack an individual measurement throughout the power system or gathered measurement vectors in the control center. Therefore, two different attack strategies can be considered. The attacker may try to attack the control center directly and alter the aggregated measurement vector in the first scenario. In this case, the number of available measurements that can be changed by the attacker is not restricted. However, in order to keep the attack unidentified, the size of the attack vector must be as small as possible. On the other hand, the attacker can launch an unidentified attack by changing the measurements values directly. In this case, limited measurements should be altered locally and simultaneously. Therefore, the lowest number of measurements that should be attacked to have an undetected and effective attack would be determined in this scenario. Both scenarios are introduced and explained in the following sections:

3.1. Attack Strategy I

In order to attack a SCADA system in a control center, before sending the gathered measurements to the state estimator, there is no limitation for the adversary to attack the number of measurements; whereas, limiting the ${\ell }_2-$norm of the attack vector can effectively impact on finding a feasible attack vector. In this case, the optimization problem is as follows:

$$min {\Vert a\Vert }_2^2=a^{T}Ia$$

(21)

s.t.

$$Y_{l}c+\beta <f_l^{max}-{\widehat{f}}_l \forall l\notin c{\ell }^{\prime }$$

(22)

$$Y_{l}c+\beta \ge f_l^{max}-{\widehat{f}}_l \forall l\in c{\ell }^{\prime }$$

(23)

$$\Vert a-Hc\Vert \le \delta$$

(24)

$$\left|a_i-H_{i}c\right|\le {\vartheta }_i$$

(25)

$$\beta >0$$

(26)

where ${\Vert a\Vert }_2^2$ is the sum of squared arrays of $a$, $a^T$ is the transposed vector of $a$, and $I$ is the identity symmetric matrix representing the quadratic function. $Y_l$ is lth row of the admittance matrix, $\beta$ is the uncertainty relaxation constant (because of the random distribution of measurements errors and the uncertainty of the real-time power generation and demand, in order to assure the effectiveness of the attack this parameter is needed to be adjusted), ${\widehat{f}}_l$ is the lth estimated power flow line in the original state, $c{\ell }^{\prime }$ is the attacker’s desired congested line, and $\delta$ is the imperfect attack threshold. Equation (21) is the objective function of the optimization problem that the attacker uses to manipulate the gathered measurements in the SCADA system. Equations (22) and (23) lead the state estimation to a situation that the attacker can change the congestion pattern. Equations (24) and (25) restrict the residual vector for each measurement to have an unidentified attack vector in the ${\ell }_2-$norm-based residual test and the LNR test, respectively. Equation (24) restricts and controls the sum of residual values for all measurements. Hence, the attacker cannot alter the value of measurements freely and there is a limitation in the aggregated value of the residual vector. Consequently, based on Equation (25) the attacker cannot alter the value of an individual measurement. In this case, there is an upper limit for each array in the attack vector. This constraint leads to change a large number of measurements to launch an undetectable attack. Finally, (26) guarantees that the uncertainty relaxation constant acquires a positive value.

3.2. Attack Strategy II

As illustrated in Figure 1, the attacker may want to attack individual measurements before being gathered by the SCADA system. In this case, the minimum sparsity objective function is considered by the attacker as the following optimization problem:

$$min {\Vert a\Vert }_0,$$

(27)

s.t. Equations (22)–(26).

Where ${\Vert a\Vert }_0$ indicates the number of non-zero arrays or sparsity of $a$. Equation (27) is the objective function of an optimization problem, by which the minimum number of manipulated measurements are found to achieve a successful attack. The constraints are the same as the constraints of the first attack strategy.

Finding the minimum sparsity of the attack vector is a nonlinear and nonconvex problem. Hence, in order to make the optimization problem solvable, the ${\ell }_0-$norm optimization problem can be estimated with ${\ell }_1-$norm or ${\ell }_2-$norm problems [27,29]. Here, in order to have a convex problem with a positive infinite symmetric matrix, a weighted ${\ell }_2-$norm optimization problem is substituted with the original problem. In addition, an iterative algorithm similar to the method in [30] is utilized.

The equivalent ${\ell }_2-$norm optimization problem is formulated as follows:

$$min a^T\omega a,$$

(28)

s.t. Equations (22)–(26).

Where $\omega$ is a weight matrix representing a quadratic function.

An iterative algorithm is performed to compute the optimum solution. The details are illustrated in Algorithm 1. The problem is solved several times to reach the stop criterion constraint.

Algorithm 1. Iteratively reweighted ${\ell }_2-$norm minimum sparsity.

Input: congestion pattern, imperfect attack relaxation, residual vector relaxation constraints

Output: attack vector $a^i$, $c^i$

$i=0$, initial weights ${\omega }^0=I$

while stop criterion false do

$a^i$, $c^i$$\leftarrow$$min a^T\omega a$ s.t. $\left(22\right) -\left(26\right)$

${\omega }^{i+1}\leftarrow {\left(\left|a^i+\zeta \right|\right)}^{-1}$ where $\zeta$ is a small positive value

$i\leftarrow i+1$

end while

4. Methodology

State variables of a power system have repetitive daily behavior in a midterm viewpoint. This periodic behavior is a remarkable indicator with valuable information in historical data. Moreover, an integrated structure of a power network creates identifiable dependencies between state variables and measurements across the grid.

NARX is an important class of discrete-time nonlinear system that employs historical values of output variable and exogenous input variables into the model. The equation for NARX is expressed as follows:

$$y\left(t+1\right)=f\left[y\left(t\right), \dots , y\left(t-n_y\right), u\left(t\right), \dots , u\left(t-n_u\right)\right]$$

(29)

where $u\left(t\right)$ and $y\left(t\right)$ indicate input and output variables for nonlinear time series modeling, respectively; $n_u$ and $n_y$ represent the input and output delay, respectively; and $f$ is a nonlinear function. The output variable is a function of the previous values of the output signal and past values of independent input variables. The function $f$ can be approximated by a neural network.

ANN has been widely employed in electrical engineering for forecasting undetermined parameters. ANN is a robust engine for prediction, classification, clustering, pattern recognition, and process modeling. It provides nonlinear mapping solutions between input features and output, due to its learning, generalization, error tolerance, and parallel processing abilities. The advantage of not requiring a clear relationship between inputs and outputs renders a neural network an attractive choice for forecasting [31].

A further developed neural network configuration is the recurrent neural network (RNN), namely, a supervised learning engine that contains feedback loops. The memorizing ability in RNN renders it interesting for time series prediction and modeling with historical data. The recurrent NARX neural network (RNARXNN) is a proper engine for chaotic nonlinear forecasting problems. RNARXNN uses feedback from the output layer, instead of hidden layers in the conventional RNN structure. This configuration improves the learning capability, convergence, and generalization performances of the learning engine [32].

A schematic diagram of the RNARXNN model with one time-series is depicted in Figure 2. The formulated RNARXNN is shown as follows.

$$y\left(t+1\right)=f_0[b_0+{\displaystyle \sum }_{h=1}^{Nh}{w}_{h0}{f}_h(b_h+{\displaystyle \sum }_{i=1}^{d_u}{w}_{ih}u\left(t-i\right)+{\displaystyle \sum }_{j=1}^{d_y}{w}_{jh}y\left(t-j\right))]$$

(30)

where $w_{ih}$, $w_{jh}$ and $w_{h0}$ are weights between independent input to the hidden layers, weights between input with past values of output to the hidden layer and weights between the hidden layer to the output layer, respectively. $b_0$ and $b_h$ are biases, and $f_0$ and $f_h$ are activation functions (linear) of the hidden layer and output layer, respectively.

A recent formulation is used to predict state variables of a power system. In this case, the input variable is the historical value of the measurement vector and the output is the state variables of the power system. Accordingly, in order to develop RNARXNN, Equation (30) is restated as follows:

$$\overline{x}\left(t+1\right)=f_0[b_0+{\displaystyle \sum }_{h=1}^{Nh}{w}_{h0}{f}_h\left(b_h+{\displaystyle \sum }_{i=1}^{d_u}{w}_{ih}\overline{z}\left(t-i\right)+{\displaystyle \sum }_{j=1}^{d_y}{w}_{jh}\overline{x}\left(t-j\right))\right]$$

(31)

where $\overline{z}\left(t\right)$ indicates the measurement vector at the time interval $t$ and $\overline{x}\left(t+1\right)$ is the predicted state variable by using the RNARXNN method.

By applying this method, state variables of a power system can be forecast with previous measurement values as input vectors. The result of forecasting may contain uncertainty and approximation, which is not as accurate as analytical methods like WLS for the estimation of states. However, it can be exploited to evaluate the accuracy of measurements and detect false data.

Therefore, by using the forecasted states, the residual vector is computed as follows:

$$\overline{r}\left(t\right)=\overline{z}\left(t\right)-H\overline{x}\left(t\right)$$

(32)

where $\overline{r}\left(t\right)$ is the estimated residual vector by using the predicted state variables. Instead of the conventional computation method for the residual vector (Equation (6)), the RNARXNN-based residual vector is utilized to test the measurements. Hence, the ${\ell }_2-$norm-based residual test equation will be updated as follows:

$$\Vert \overline{r}\left(t\right){\Vert }_2\le \tau$$

(33)

consequently, the LNR test is changed to the following:

$$\overline{r}{\left(t\right)}_i^N=\frac{\Vert \overline{r}{\left(t\right)}_i\Vert }{\sqrt{{\Omega }_{ii}}}\le \epsilon$$

(34)

Therefore, the FDIA will be identified using the updated value of the residual vector in the conventional tests; however, this type of attack could alter the power system characteristics without being detected by the conventional BDD procedure. The flowchart of the proposed method is illustrated in Figure 3.

5. Data Preparation and Application

5.1. Database Generation

According to Figure 4, the IEEE 14-bus test system is utilized to examine the proposed method [33]. In order to investigate the method more clearly, some modifications are undertaken in the system. Accordingly, the generation capacity is altered to 330 MW for unit 1, the line power flow limitation between buses 2–3, 4–5, and 6–11, are assumed to be 50 MW, 50 MW, and 30 MW, respectively. The power injection and power flow transmitters are considered for all buses and line ends. The measurement vector $z$ consists of 40 active power flow and 14 active power injection. The normal distribution is assumed for measurement errors. The standard deviation of injection and power flow measurements are determined as 0.01 and 0.008, respectively. The state variable vector $x$ consists of voltage angles of 13 nodes.

DC power flow and DC state estimation are performed by using the MATPOWER toolbox in MATLAB to examine the proposed model.

The load profiles illustrated in Figure 5 are presumed to perform DC power flow and state estimation in 5-min intervals during a day [34]. According to

Table 1. Load profile proportion for each node.

Node	Nominal Active Load (MW)	Residential (%)	Commercial (%)	Industrial (%)
1	0	0	0	0
2	21.7	0	50	50
3	94.2	0	50	50
4	47.8	100	0	0
5	7.6	0	50	50
6	11.2	0	50	50
7	0	0	0	0
8	0	0	0	0
9	29.5	100	0	0
10	9	0	50	50
11	3.5	0	50	50
12	6.1	0	100	0
13	13.5	0	100	0
14	14.9	100	0	0

, various proportions of residential, commercial, and industrial load are determined for each node in the power system.

In order to train and test the introduced ANN, a dataset of historical data is provided. As mentioned, the state estimation and state variables are computed by using load profiles and generation units by DC Optimal Power Flow (DCOPF) program at 5-min intervals. In order to consider the variation and uncertainty in the load and generation, load profiles of all nodes are generated by a Gaussian distribution with a 5 percent deviation ($\sigma =5\%$) in comparison to the original load profile. A subset of data containing 6048 measurement vectors and state variables is provided to train and test the network.

5.2. Structures of RNARXNN

There are two different architectures for RNARXNN, namely, series-parallel architecture, in which real output values are directly fed into the engine from the dataset; and parallel architecture, wherein estimated outputs are feedbacked into the neural network. Therefore, the series-parallel architecture is utilized for training and testing, whereas the parallel architecture is performed for real-time estimation. Sigmoid and linear activation functions are selected for the input neurons and hidden layer, respectively. A random selection of 70% of the historical data is used for training the network, 15% is exploited for testing the network, and 15% is used for validation. The maximum number of epochs is set to 100 and the performance goal is set as 0.0000001 after trial and error procedure in order to attain the best result. In order to examine the accuracy and precision of estimated states, mean square error (MSE) is adopted as an evaluation index, which is defined as follows:

$$MSE=\frac{{{\displaystyle \sum }}_1^N({\overline{x}}_i-{\widehat{\overline{x}}}_i)}{N}$$

(35)

where $N$ is the number of sampling observations and ${\overline{x}}_i$, ${\widehat{\overline{x}}}_i$ are measured and predicted values, respectively.

6. Numerical Results

This section presents the simulation results of the proposed method. In this case, FDIA strategies are discussed; consequently, the results of ANARXNN are illustrated, and, finally, the efficiency of the proposed method is indicated.

6.1. FDIA Strategies

Two different strategies of FDIA are formulated in Section 3. Both problems are quadratic optimization problems in the form of $a^{T}Qa$ with convex constraints. In this case, provided that the symmetric matrix $Q$ is a positive definite matrix and the constraints are convex, the problem has a unique answer. In strategy I, the identity matrix $I$ is positive and definite; hence, the problem has a unique answer. In strategy II, the $Q$ equivalent symmetric matrix is checked in an iteration to be positive and definite. MATLAB software is utilized to solve these problems. The nominal load profile is considered for this simulation. In order to consider the impact of the Gaussian distribution error of measurements, the Monte Carlo method is implemented, with the aim of attaining the expected results. Accordingly, the state estimation is performed 1000 times by generating and adding an accidental error vector to the measurement vector. The original congestion pattern is $c\ell =\left[1\right]$ in this case, which means line 1 is congested throughout all transmission lines. According to the assumption, lines 1, 4, and 5 have a power transmission limitation. Hence, in order to change the congestion’s pattern, the attacker has seven options. The simulation is performed for these seven different options and the results are illustrated in

Table 2. FDIA results for various congestion patterns.

$\mathit{\varsigma }{\mathit{\ell }}^{\prime }$	Attack Strategy I			Attack Strategy II
	IL (%)	$\mathit{E}\left({\Vert \mathit{a}\Vert }_0\right)$	$\mathit{E}\left(\raisebox{1ex}{$\Vert \mathit{a}\Vert $}\!\left/ \!\raisebox{-1ex}{$\Vert \mathit{z}\Vert $}\right.\right)$	IL (%)	$\mathit{E}\left({\Vert \mathit{a}\Vert }_0\right)$	$\mathit{E}\left(\raisebox{1ex}{$\Vert \mathit{a}\Vert $}\!\left/ \!\raisebox{-1ex}{$\Vert \mathit{z}\Vert $}\right.\right)$
[]	5.3	34	0.13	5.3	13	0.19
[4]	6.7	34	0.39	8.5	14	0.45
[5]	6.4	34	0.42	7.3	14	0.49
[1,4]	5.9	34	0.37	8	16	0.43
[1,5]	7.3	34	0.43	8.2	15	0.50
[4,5]	7.2	34	0.41	6.7	16	0.61
[1,4,5]	4.3	34	0.42	9.7	17	0.63

.

The first column illustrates various congestion patterns selected by the attacker. The percentage values of IL in the 2nd and 5th columns indicate the identification likelihood by the conventional BDD processor for both attack strategies. $E\left({\Vert a\Vert }_0\right)$ is the expectation value of sparsity of the attack vector. $E\left(\raisebox{1ex}{$\Vert a\Vert $}\!\left/ \!\raisebox{-1ex}{$\Vert z\Vert $}\right.\right)$ is the proportion of the size of the attack vector to the measurement vector. According to the results, in order to attack the state estimation and change the congestion pattern, all measurements are altered in strategy I; however, fewer measurements are manipulated in strategy II. On the other hand, the proportion of the attack vector in strategy I is lower in all attacks. Actually, the result describes the objective functions in the attack strategies. The optimization problem tries to minimize the size of the attack vector and there is no limitation in the number of measurements; hence, all the measurements are altered to have an attack with the minimum size of attack vector. On the contrary, the number of measurements to be altered is less in strategy II; meanwhile, the size of the attack vector is not optimized. The identification likelihood is not dependent on the attack strategy and congestion pattern. However, the conventional BDD processor only identifies less than 10% of the attacks.

6.2. RNARXNN Training and Testing

By using historical data and the ANN toolbox in MATLAB software, the proposed estimation of power system states is performed for different values of fundamental parameters and the results are illustrated in

Table 3. RNARXNN results for different delays and hidden layer neurons.

Input Delays	Feedback Delays	Hidden Layer Neurons	MSE	Epochs	Time Duration
5	5	25	$4.8779\times {10}^{-8}$	22	37 min
5	5	15	$3.1571\times {10}^{-8}$	42	28 min
4	4	20	$3.9819\times {10}^{-8}$	30	20 min
3	3	15	$8.66\times {10}^{-8}$	28	7 min
2	2	25	$3.1814\times {10}^{-8}$	50	15 min
2	2	15	$2.3324\times {10}^{-8}$	91	10 min
1	1	25	$3.8557\times {10}^{-8}$	95	8 min

.

According to the results illustrated in Table 3, the difference in the values of the input and feedback delay as well as hidden layer neurons do not impact remarkably on the performance of the network. In fact, owing to the high correlation between input and output variables, the network can predict the states reliably. However, the number of hidden layer neurons can desirably impress the number of epochs and the time duration for training the network. Accordingly, the higher the number of hidden layer neurons, the higher the time duration and the lower the number of epochs in the same input and feedback delays. In addition, it is quite obvious that, by increasing the value of delay, the time duration is increased to train the network. Figure 6 illustrates the training, validation and test performance of the state estimation for various delay times and different sizes of hidden layers. According to the results, the proposed method can predict the state vector effectively in all cases.

6.3. Proposed Strategy for Large Scale Power System

The real power system contains a large number of buses with several state variables and measurements. Hence, the number of input and output variables in NARX network would be much more than the proposed standard system. However, the large number of variables makes the results more accurate, because of including all correlations and links, the network would be more complicated and the training time duration may be intolerable and a more powerful computer system is needed [35]. In order to tackle this drawback, a simple still effective strategy is implemented here to make the method practical for real-world power systems. The power system is divided into several subsystems to reach to an appropriate subsystem with a smaller number of state variables and measurements. Accordingly, an individual network is considered for each zone where the number of inputs and outputs is considerably smaller than the whole power system. For instance, consider the IEEE 14-bus test system in Figure 4. The system has been divided into two different subsystems. Zone 1 contains five state variables and ten measurements. Consequently, zone 2 involves nine state variables and ten measurements. In order to examine the performance of the proposed method, two different networks have been considered for each zone.

The results in

Table 4. Performance and time duration of RNARXNN for each zone.

ZONE	Hidden Layer Size	MSE	Epochs	Time Duration
1	5	$3.9349\times {10}^{-8}$	40	6 s
1	10	$2.9264\times {10}^{-8}$	108	55 s
2	5	$9.5893\times {10}^{-8}$	55	2:18 min
2	10	$1.4916\times {10}^{-7}$	32	30 s

illustrate that the network training time has been decreased dramatically in comparison to the results for the whole system listed in Table 3. Meanwhile, the mean square error is near the same as the total system. Therefore, using two separate networks the state variables can be predicted with the same accuracy. This method can be implemented for large power systems. Hence, considering the results, it can be illustrated that the proposed method can be easily utilized for various power systems with a different number of measurements and state variables. In addition, Figure 7 depicts the performance of the network for the data from both zones with different hidden layers.

6.4. FDIA Detection Using RNARXNN

In order to examine the efficiency of the proposed method, a typical day is selected and a specific time-spot is opted for an FDIA on the state estimation without being identified by the conventional BDD processor. The results for the two introduced attack strategies are illustrated in Figure 7.

As shown in Figure 8, the value of state variables is changed in a specific time spot to alter the congestion pattern. Actually, in the time interval, the line between bus number 1 and 2 is congested; hence, LMP values are not the same in all buses. However, FDIA is considered here to alter the status of the congested line to a non-congested line. This attack can alter LMPs in the real-time power market virtually. The measurement vector is altered in such a way that the residual vector is still a small value and the attack cannot be realized by using the conventional BDD processor.

Instead of using the conventional residual vector and bad data detection processor that is not efficient in the case of optimal FDIA, the proposed method can detect similar false data at a higher rate. In order to examine the proposed method, a specific network is selected, in which the input and feedback delay is 4 and the size of the hidden layer is 20. Figure 9 illustrates the results of the state variables predicted by the trained network. In this case, the residual vector computed by Equation (32) exceeds the threshold and the presence of bad data is identified in measurement vectors. Table 4 illustrates the results of various attack scenarios and the performance of the proposed method to identify the attack considering the uncertainties of measurement error.

The second column in

Table 5. Results of FDIA identification likelihood by using the proposed method.

$\mathit{c}{\mathit{\ell }}^{\mathbf{\prime }}$	Attack Strategy I			Attack Strategy II
	${\mathit{\ell }}_{\mathbf{2}}-\mathbf{norm}\mathbf{IL}(\%)$	LNR IL (%)	LNR IMN	${\mathit{\ell }}_{\mathbf{2}}-\mathbf{norm}\mathbf{IL}(\%)$	LNR IL (%)	LNR IMN
[]	78.67	88	7	78.67	92	6
[4]	92.12	96	19	96.63	100	10
[5]	96.14	100	20	97.26	100	10
[1,4]	97.18	100	19	97.78	100	10
[1,5]	97.93	100	20	98.87	100	10
[4,5]	98.63	100	21	99.01	100	10
[1,4,5]	99.12	100	21	99.51	100	10

demonstrates the identification likelihood of the ${\ell }_2-$norm-based residual test by using RNARXNN-based BDD processor. Accordingly, the third and fourth columns depict the identification likelihood and the number of identified measurements of the LNR test for attack strategy I, respectively. Similarly, the results for attack strategy II are illustrated in the rest of the columns. According to the results, the RNARXNN-based BDD processor is able to detect FDIAs with remarkable accuracy. Except for the first scenario, all other scenarios are detected with more than 92% precision. The exception in the first scenario happens because the altered pattern is extremely similar to the original congestion pattern without considerable variations; thus, the estimation is performed with lower accuracy compared to other scenarios. The 100% percentage values in columns 3 and 6 indicate that the residual value of at least one measurement is more than the threshold preset value for all runs. Comparing two attack strategies, the number of manipulated arrays is lower in attack strategy II; nevertheless, the majority of the manipulated measurements are detected by using the proposed method.

7. Conclusions

This paper proposes a NARX-based BDD processor to identify FDIAs on power system state estimation. The desired state variables are predicted by using a NARX network and are compared to computed state variables. Hence, by using the residual vector, potentially bad data are identified. The computed state variables passing the BDD processor can be applied as accurate results in power system operation and market clearance.

Two different FDIA strategies are considered to formulate and simulate a cyber-attack on power system state estimation. It is illustrated that this kind of attack can be launched on state estimation and can alter the line congestion pattern without being detected by conventional bad data detection processors. The results depict that only less than 10% of the formulated attacks are detected by conventional processors.

The proposed method is implemented in the IEEE 14-bus standard test system. After having prepared adequate historical data for power system measurements and state variables, an RNARXNN-based model is developed to train, test, and validate the data in order to identify FDIAs on state estimation. According to the results, owing to the correlation between input and feedback variables, the network is trained and tested with great performance and can predict the state variable vector appropriately. The results illustrate the high efficiency of the proposed method for bad data detection of FDIA on state estimation.

The proposed method can be implemented on large scale power systems by dividing the power system into several subsystems and training individual network for each subsystem. This method can decrease the training time duration by keeping the performance of the network as nearly accurate as the results of the primary system.

In future work, the presented method can be extended for AC state estimation by adding reactive power measurements and PMUs to the measurement vector and voltage magnitude to state variables.