A reinforcement learning framework includes an agent and an environment, as illustrated in Figure 1, which aims at maximizing a long-term reward through abundant interactions between the agent and the environment. At each step t, the agent observes states s t and executes action a t ; based on its observation and policy, the environment receives action a t , then emits states s t + 1 , and issues a reward r t + 1 to the agent. Compared with supervised learning, the actions of RL are not labeled, that is, the agent does not know what the correct action is during training and can only be trained through the trial and error approach to explore the environment and maximize its reward.

The interaction between the agent and environment can be modeled by a Markov decision process, which is a standard mathematical formalism of sequential decision problems. A typical Markov decision is denoted by a tuple < S , A , P , R , γ > , where S is the state space, and it is the complete description of the environment which is represented by a real-valued vector, matrix, or higher-order tensor. A is often called the action space that is also represented by a real-valued vector matrix or higher-order tensor, whereas different environments allow different kinds of actions, that is, discrete and continuous action spaces. P is the transition probability function, and P ( s ′ | s , a ) is the probability of transitioning into state s ′ by taking action a on state s . R is the reward function, and R ( s , a , r ) is the probability of receiving a reward r from action a and state s . γ ∈ [ 0,1 ] is the reward discount factor. The agent learns to find a policy to maximize the total discounted reward as presented in (11), and T is the number of time steps in each episode. G t = R t + 1 + γ R t + 2 + ⋯ + γ T − t − 1 R t + T . (10)

The policy is a rule used by an agent to decide what actions to take, which maps the action from a given state. A stochastic policy is usually expressed as π θ ( a t | s t ) , in which parameter θ denotes the weights and biases of the neural network in deep reinforcement learning algorithms.

The state value functions V π ( s ) is the expected return starting from state s following a certain policy as defined in (12), which is used to evaluate the state: V π ( s ) = Ε ( G t | s t = s ) . (11)

The action-value function Q π ( s , a ) is the expected return starting from state s , taking action a , and then following policy π , denoted as (13), which is utilized to evaluate the action: Q π ( s , a ) = Ε ( G t | s t = s , a t = a ) . (12)

The advantage function A π ( s , a ) corresponding to policy π measures the importance of each action in this state, which is mathematically defined as shown in (14): A π ( s , a ) = Q π ( s , a ) − V π ( s ) . (13)

In general, the DRL algorithms can be divided into the value-based, the policy-based, and the actor-to-critic (A2C) framework. The proximal policy optimization (PPO) algorithm follows the A2C framework with an actor network and a critic network.

The main advantage of applying PPO algorithm to the AGC optimization problem is that the new control action decision updates from the policy network does not change too much from the previous policy and can be restrained within the feasible region by the clipping mechanism. During the off-line training process, the PPO also converges faster than other DRL algorithms. Also, during the on-line operations, the PPO generates smoother, less variance, and more predictable sequential decisions, which is desired for the AGC optimization.

The overall structure of the PPO algorithm is presented in Figure 2, including an actor network and a critic network. The AGC training environment sends the experience tuples < s t , a t , r t + 1 , s t + 1 > to the trajectory memory pool to form finite mini-batches of samples and returns to the PPO algorithm.

The actor network and the critic network are realized by deep neural networks (DNNs) with the following equations: O i = f ( W i I i + b i ) ,         i = 2 , … ,   n l a y e r + 1 , (14) O D = O n l a y e r ( … O 2 ( O 1 ( s t ) ) ) , (15) where I i and O i represent the input array and output array of the i t h layer of the DNN. The layers are connected as (16), I i + 1 = O i ,   a n d   I 1 = s t . n l a y e r is the number of total layers and W i and b i are the weights and bias matrices of the i t h layer. The ReLU functions are used as the activation function f ( ⋅ ) .

The actor network contains the policy model π θ with network parameters < W θ ¯ , b θ ¯ > . It is responsible for the sequential decisions of AGC optimizations. The rewards are taken in by the critic network V μ with parameters < W μ ¯ , b μ ¯ > , which is a value function and maps the state s t to the expected future cumulative rewards.

The conventional policy gradient-based DRL optimizes the following objective function (Wang et al., 2020): L P ( θ ) = E ^ [ l o g π θ ( a t | s t ) A ^ t ] , (16) where E ^ [ ⋅ ] is the empirical average over a finite mini-batch of samples, π θ is a stochastic policy, and A ^ t is an estimator of the advantage function at time t . In this work, a generalized advantage estimator (GAE) is used to compute the advantage function, which is the discounted sum of temporal difference errors (Schulman et al., 2017). A ^ t = δ t + ( γ λ ) δ t + 1 V + ( γ λ ) 2 δ t + 2 V + … + ( γ λ ) U − t + 1 δ U − 1 V ,     (17) δ t = r t + γ V μ ( s t + 1 ) − V μ ( s t ) , (18) where γ ∈ [ 0,1 ] is the discount factor, λ ∈ [ 0,1 ] is the GAE parameter, U is the length of the sampled batch, and r t is the reward at time t . The objective function L V ( ⋅ ) can be formulated as: L V ( μ ) = E ^ [ L t V ( μ ) ] = E ^ [ | V ^ μ t a r g e t ( s t ) − V μ ( s t ) | ] , (19) V ^ μ t a r g e t ( s t ) = r t + 1 + γ V μ ( s t + 1 ) , (20) where V ^ μ t a r g e t ( ⋅ ) is the target value of time-difference (TD) error. The parameters of the critic network V μ can be updated by the stochastic gradient descent algorithm in Duan et al. (2020) according to the gradient ∇ L V ( μ ) with a learning rate η .

The input of the actor network is the observation state s t , and the outputs are the normal distribution mean value and standard deviation of the actions, that is, the strategy distribution π θ ( a t | s t ) . The importance sampling is used to obtain the expectation of samples gathered from an old policy π θ o l d ( a t | s t ) under the new policy π θ ( a t | s t ) . This process converts the PPO algorithm from an on-policy method to an off-policy method, which means that the actor network is updated asynchronously to further stabilize the performance of AGC actions. The following surrogate object function is being maximized:   L C P I ( θ ) = E ^ [ ( π θ ( a t | s t ) / π θ o l d ( a t | s t ) ) A ^ t ] = E ^ [ r t ( θ ) A ^ t ] , (21) s . t .   E ^ [ K L [ π θ o l d ( ⋅ | s t ) , π θ ( ⋅ | s t ) ] ] ≤ ξ , (22) where K L [ ⋅ ] is the Kullback–Leibler (KL) divergence, r t ( θ ) = π θ ( a t | s t ) / π θ o l d ( a t | s t ) denotes the ratio of the probability of action a t under the new and old policies, and ξ is a small number. In order to simplify the penalty by the KL divergence to a first-order algorithm and attain the data efficiency and robustness, a clipping mechanism, c l i p ( r t ( θ ) , 1 − ϵ , 1 + ϵ ) A ^ t , is introduced to modify the surrogate objective by clipping r t ( θ ) . It removes the incentive of moving r t outside of the interval ( 1 − ϵ , 1 + ϵ ) . The objective function with the c l i p ( ⋅ ) function is defined as, L C L I P ( θ ) = E ^ [ L t C L I P ( μ ) ] = E ^ [ min ( r t ( θ ) A ^ t ,   c l i p ( r t ( θ ) , 1 − ϵ , 1 + ϵ ) A ^ t ) ] . (23)

The importance sampling and clipping function help the PPO DRL algorithm achieve better stability and reliability for AGC online operations, better data efficiency and computation efficiency, and better overall performance.

State space S: the setting of the state space should consider the factors that may affect the decision as much as possible. In this work, the state space is determined as a vector of system information representing the current system condition at time t and prediction system information at time t + 1 . Specifically, the former includes the real power output of all units (AGC and non-AGC unit) P G , t r , the frequency deviation Δ f t r , the power deviation of the tie-line Δ P T , t r , and the area control error A C E t r . The latter includes the prediction system information of load P l , t + 1 f , wind power P w , t + 1 f , frequency deviation Δ f t + 1 f , power deviation of tie-line Δ P T , t + 1 f , and area control error A C E t + 1 f . It is set as follows: S : { P G , t r , Δ f t r , Δ P T , t r , A C E t r , P l , t + 1 f , P w , t + 1 f , Δ f t + 1 f , Δ P T , t + 1 f , A C E t + 1 f } . (24)

Action space A: action space is the decision variable in the optimization model, including the ramp direction and ramp power. In this article, to avoid the lack of generality, the action is defined as power increments of AGC units Δ P A G , t a at each optimization time, which are subjected to the ramp power limits of corresponding AGC units. A is set as A : { Δ P A G , t a } . (25)

The design of reward function is crucial in DRL. It generates reward r t at time t in each decision cycle, which evaluates the agent’s actions based on the AGC control performance under the impacts of uncertainties in the system variables. In this work, the values of load, wind generations, frequency deviations, the tie-line power deviations, and ACE are used to formulate the reward function, which consists of cost terms, punishment terms, and performance terms. The reward r t is calculated by the formula: r t = F c o s t + r p e n e l + f c p s , (26)

where the cost term F c o s t represents the total cost of the system. It includes AGC adjustment ancillary service cost and the load shedding cost, which are calculated as follows: F c o s t = − c 1 f A G C − c 1 P c , t 2 , (27) where c 1 and c 2 are the corresponding cost coefficients. The load shedding P c ,   t must be set reasonably as follows: P c , t = { 0 , | Δ f | ≤ 0.2 , ( | Δ f | − 0.1 ) ⋅ Δ P t , | Δ f | > 0.2. (28)

Here, the real power deviations Δ P t are utilized to reflect the stochastic process caused by the load and wind power fluctuations. At time t , the power deviations in the system are calculated as follows: Δ P t = ∑ i = 1 N P G , i , t + P w , t r − P L , t r − P T , t r − P l o s s , t r , (29) where N is the total number of thermal power units in the system, including AGC and non-AGC units, P G , i , t is the power output of thermal power unit i at t period, and P w , t r , P L , t r , and P T , t r are the power output of wind power, load, and tie-line power. Note that the power flowing out of the system is taken to be positive. P l o s s , t r is the system loss at time t .

The power output of any unit (AGC or non-AGC unit) relates to the frequency deviation, tie-line power deviation, and ACE. Taking an interconnection power grid of two areas as an example, the system contains region A and region B. The control strategy of the two areas is the tie-line bias frequency control. It is assumed that Δ P L A , t and Δ P L B , t , and Δ P G A , t and Δ P G B , t are the change of load and the change of power output of units in region A and region B at time t , respectively, and K A and K A are the frequency regulation constants of region A and region B. We define Δ P A , t and Δ P B , t as power imbalance of the two areas which can be calculated as follows: Δ P A , t = Δ P L A , t − Δ P G A , t , (30) Δ P B , t = Δ P L B , t − Δ P G B , t . (31)

Frequency deviation, tie-line power deviation, and area control error can be calculated as follows: Δ f t = − Δ P A , t + Δ P B , t K A + K B , (32) Δ P T , t = K A ⋅ Δ P B , t − K B ⋅ Δ P A , t K A + K B , (33) e A C E , t = Δ P T , t − 10 B ⋅ Δ f t , (34) where B is the equivalent frequency regulation constant for the control area in MW/0.1Hz and the value is negative.

The punishment term r p e n e l formulates the operation and control limits in AGC dynamic optimization, including generation unit power output limits, CPS1, frequency deviation limits, and tie-line power deviation limits and given as: r p e n e l = r 1 + r 2 + r 3 + r 4 . (35)

The AGC units participate in both primary and secondary frequency control; thus, the outputs of AGC units at time t + 1 are calculated as P A G , i , t + 1 ′ = P A G , i , t + 1 + Δ P A G , t a − K G i ( Δ f t + 1 r − Δ f t r ) , (36) where Δ P A G , t a is the regulated power of AGC unit i at time t , that is, the power increment of secondary frequency control; and K G i ( Δ f t + 1 r − Δ f t r ) is the primary frequency control power of AGC unit i , where K G i is the frequency regulation constant of unit i , and Δ f t + 1 r and Δ f t r are system frequency deviations at time t + 1 and t , respectively.

Accordingly, the power outputs of non-AGC units at time t + 1 are calculated as P N G , i , t + 1 ′ = P N G , i , t − K G i ( Δ f t + 1 r − Δ f t r ) . (37)

The outputs of AGC and non-AGC units are subjected to the corresponding maximum and minimum power limits: P A G , i , t + 1 = { P G , i , m i n , P G , i , t + 1 ′ < P G , i , m i n , P G , i , m a x , P G , i , t + 1 ′ > P G , i , max , P G , i , t + 1 ′ ,     else, (38) r 1 = { 0 , P A G , i , m i n < P G , i , t < P A G , i , m a x , k 1 , else, (39) where k 1 is the punishment coefficient. The CPS1-related punishment term is formulated as, r 2 = { 0 , K c p s 1 ≥ 200 % , − k 2 | e A C E − e A C E ∗ | , 100 % ≤ K c p s 1 ≤ 200 % , − k 3 | K c p s 1 − K c p s 1 ∗ | , K c p s 1 < 100 % , (40) where k 2 and k 3 are the punishment coefficients of ACE and CPS1 and e A C E ∗ and K c p s 1 ∗ are, respectively, ideal values of ACE and CPS1. In this article, the ideal values of ACE and CPS1 are 0 and 200%.

The following two functions denote the frequency deviation and tie-line power transfer deviation punishments, respectively: r 3 = { 0 , Δ f min ≤ Δ f ≤ Δ f max , k 4 , else, (41) r 4 = { 0 , Δ P T min ≤ Δ P T , t ≤ Δ P T max , k 5 , else, (42) where k 4 and k 5 are the corresponding punishment coefficients.

In this article, we added an additional performance evaluation term f c p s with a coefficent c 3 to the reward function which makes the PPO-based DRL algorithm has the capability of further improving the long-term AGC performance: f c p s = − c 3 ( 2 − K c p s 1 ) 2 . (43)

State transition probability P: in this work, the reinforcement learning algorithm based on the model-free method is utilized, so the state of the agent at the next time and rewards can be obtained by the interaction with the environment, and they make up state transition probability P including environmental stochasticity.

Discount factor γ ( γ ϵ [ 0,1 ] ) determines the importance of rewards in future to current reward. When γ = 0 , it means that the impact of current decisions on the future system operating status is not considered, and only the operating cost of the current control period is optimized; when γ = 1 , it means that the impact of current decisions on the operating status of the system at every moment in the future is equally considered. For AGC dynamic optimization control, the decision at the current moment will have an important impact on the future operating state of the system, and the closer the distance to the current decision period, the greater the impact.

Based on the aforementioned analysis, this article transforms the AGC dynamic optimization problem into a sequential decision issue and utilizes the PPO deep reinforcement learning algorithm to solve the proposed problem. The AGC dynamic optimization problem based on the PPO algorithm is shown in Figure 4. The specific process is described as follows:

1) Initialize the weight and bias of the neural network, actor and critic neural network learning rate, reward discount factor γ , and hyperparameters ε and other parameters. Set the number of episode M and decision cycle N .

In this article, the PPO agent for AGC dynamic optimization control is tested on the modified IEEE 39 bus system model which includes three AGC units and seven non-AGC units. The tie-line is connected to bus 29, and a wind farm with 130 MW installations is connected to bus 39. A single-line diagram of the system is shown in Figure 5.

The forecasting and actual data of load and wind power come from New England power grid¹. The basic parameters of three AGC units and test system are shown in Tables 1, 2. The control period is set at 15 min, and it is assumed that the deviation of frequency and tie-line transmission power at the initial time are 0.

The action space refers to the regulation power of AGC units at each optimization moment, which is determined by ramp power limits of each AGC unit. The per unit action space of the three AGC units is set as follows: A 1 = { − 0.3 , 0.3 } , A 2 = { − 0.45 , 0.45 } , A 3 = { − 0.6 , 0.6 } . (45)

In addition, the state space dimension is 35 according to the preceding description, which includes information on 19 forecasting loads at time t+1, actual output power of 10 units at time t, and actual and forecasting values of system frequency deviation, transmission power deviation, and ACE separately at time t and t+1. The dimensions of state space and action space, respectively, correspond to the neural numbers of input and output layers. Therefore, this work sets up three hidden layers both in actor and critic neural networks, and the number of neurons in each layer is 64, 128, and 32, respectively. The activation function in each hidden layer is the ReLU function. A larger learning rate α accelerates the convergence of the algorithm, while a smaller α tends to enhance the stability. In this article, learning rate α both in actor and critic networks is set to 0.0001.

Based on the preceding model and significant parameters, the PPO agent is coded using the TensorFlow framework with Python 3.7. The results of CPS1 index are shown in Figure 6.

The x-axis represents the number of episodes being trained, while the y-axis represents the value of CPS1 index in each episode. It can be observed that the CPS1 values of the first few hundreds of episodes are relatively low and unstable. As training episodes increase, CPS1 values are kept within a stable range around 191.3%, which fits CPS. In comparison, the deep Q learning (DQL) algorithm and the duel deep Q learning (DDQL) algorithm are also implemented. The average CPS1 values are 187.4 and 184.5%. This shows that the PPO architecture for AGC unit dynamic optimization proposed in this article can effectively learn the growing uncertainties in the power system. Once the agent is trained, it can make proper decisions based on its trained strategy combined with environmental observation data feedback. Specifically, the agent trained in this work receives data from the power system, including actual information of unit output power, frequency, tie-line transmission power, ACE, and forecasting information of load, wind power, frequency, tie-line transmission power, and ACE as its observation, and then makes decisions for the regulation power of AGC units at time t, that is, advanced control of AGC units, in order to reduce the frequency deviation at time t+1.

In addition, taking a typical control period of the system as an example, Figure 7 shows the actual load, the wind power generations, and the net load by subtracting the wind generations from the actual load. The load at each bus is allotted in proportion to the load of the original IEEE-39 node system.

Using PI hysteresis control and PPO algorithm for frequency control in this period, the results of system frequency deviation, transmission power deviation of tie-line, and ACE are represented in Figures 8A–C, respectively.

It is observed in Figure 8A that the outputs of both the PPO agent and optimization method can meet the requirements of frequency deviation (i.e., ±0.05 Hz). Moreover, maximum frequency deviation of the system controlled by the PPO agent is 0.0175 Hz, which is superior to -0.044 Hz that is controlled by the optimization method. This demonstrated that the dynamic optimization strategy of AGC units based on PPO algorithm is able to mitigate the frequency fluctuation of the system efficiently by advanced control of AGC units.

Figure 8B shows the transmission power deviation of tie-line. The power deviation under the optimization method is relatively large, which contains three times the off-limit conditions owing to the AGC resources in the system that are insufficient. While under PPO agent controller, the power deviation fluctuation is smaller without over-limit time. It proves that the system operation is more stable when using the PPO agent controller optimization method. In addition, as shown in Figure 8C, the agent performs much better than the optimization method when calculating the values of ACE. Figure 9 shows the AGC power regulation curve of all AGC units.

In the training process, the cumulative reward of each episode is recorded. Then, the results and the filtered curve are shown in Figure 10.

Due to the loads and wind power fluctuations being different in each episode, the needs of frequency regulation in each episode are also different. Therefore, it is normal that there exists slight oscillation of cumulative rewards for each episode. As the training process continues, the cumulative rewards tend to converge as shown in Fig.