Online Monitoring and Early Warning of Subsynchronous Oscillation Using Levenberg–Marquardt and Backpropagation Algorithm Combined with Sensitivity Analysis and Principal Component Analysis

Over the past few years, with the access of large-scale new energy sources, the problem of subsynchronous oscillation (SSO) in power systems has presented a novel multisource and multitransformation form, which may be significantly threatening. Conventional control and protection methods primarily give rise to device protection actions in the presence of severe oscillation. On the whole, online monitoring only identifies the frequency and amplitude, whereas it cannot identify the attenuation factor. Moreover, the determination of the warning threshold is more dependent on human experience, so the reliability and rapidity of the early warning cannot be ensured. This study conducts an in-depth investigation of the wind-thermal power bundling and extreme high-voltage alternating current- (AC-) direct current (DC) hybrid transmission system. The major factors of SSO using this system are unclear, which brings difficulties to effective monitoring. Given the mentioned problems, a method combining Levenberg–Marquardt- (LM-) Backpropagation (BP) machine learning and Sensitivity Analysis (SA) and principal component analysis (PCA) is developed. First, the sensitivity analysis of each factor in the system is conducted to identify the major factors of SSO. Subsequently, the historical sample data are reduced with the principal component analysis to reduce the redundancy, which is adopted to train the regression model to determine the attenuation factor and frequency and then send them to the classifier for classification to complete the task of the assessment model. When a novel data signal is uploaded, the assessment model identifies the attenuation factor and frequency and subsequently determines the presence of SSO. Accordingly, an early warning is conducted. The system's refined simulation model and machine learning model verify the effectiveness of the method.


Introduction
Over the past few years, new energy power generation systems represented by wind power have been leaping forward worldwide. As impacted by the high volatility of wind farms, the combination of wind farms and thermal power plants has been commonly adopted, which can reduce the effect of wind farms fluctuations on the power grid by regulating the output of thermal power units [1]. e majority of new energy power plants are far away from the load center and require long-distance power transmission. Series capacitance compensation technology and high-voltage direct current transmission (HVDC) technology [2] have been extensively employed in new energy base power delivery occasions. With the continuous expansion of the scale of the power grid and grid-connected delivery of HVDC control systems for new energy power generation, dynamic reactive power compensation devices and other types of equipment have been installed in the power system to enhance the transmission capacity and the stability of the system. With the advancement of the mentioned devices and technologies, a novel method is developed for the power grid to enhance its transmission capacity under the existing conditions and underpin large-scale new energy grid-connected transmission, whereas it also poses certain risks to the stable operation of the power grid. Improper design and operation may trigger SSO in the power system [3].
A novel way of classifying SSO is proposed in the literature [4], helping to gain insights into the emerging SSO in practical power systems. Literature [5] presented an overview of potential subsynchronous control interaction (SSCI) mitigation techniques at various stages of the power system. Literature studies [6,7] proposed a simple practical grid-side subsynchronous damping controller (GSDC) and a rotorside subsynchronous damping controller (RSDC), respectively. e proposed scheme has verified its effectiveness and robustness by performing controller-hardware-in-the-loop tests to mitigate SSCI over various operating conditions. Literature [8] proposed a fractional-order sliding mode control (FOSMC) based on feedback linearization technique to mitigate SSCI in DFIG. A robust subsynchronous damping controller is proposed with the FOSMC method to mitigate SSCI induced by DFIG-based wind farm in existing studies [9,10]. e timely detection and accurate alarms of the power system SSO are critical to ensure the safety of the unit itself and the stable operation of the entire power system [11,12]. Literature [13] developed an SSO detection method based on the PRONY algorithm; the changing law of instantaneous current in continuous time was concluded to determine the presence of SSO. Literature [14] adopted the PMU to monitor the subsynchronous harmonic frequency characteristics of the active power oscillation waveform of the tie line to determine the presence of SSO. In addition, literature [15] proposed an SSO online monitoring method with both the synchronous extraction transform (SET) and the Naive Bayes method. An artificial intelligence method was presented to classify the parameter identification results.
In some areas of North China, a wind-thermal bundling power supply and extreme high-voltage alternating currentdirect current (EHV AC-DC) hybrid transmission system have been installed, whereas rare studies in the literature have solved SSO of the system. In the meantime, the accurate monitoring and early warning of SSO of this novel complex power supply system did not fully consider the role of the main factors of SSO, which has caused insufficient comprehensive early-warning criteria and the decrease in earlywarning accuracy. Specific to online parameter identification, the existing methods can only infer the occurrence of the oscillation phenomenon by identifying SSO frequency and amplitude. ey are required to exploit human experience to determine the threshold value and cannot assess the attenuation factor through online identification directly. e generation of warning information for consecutive occurrences will hinder the process of making judgment, and the harm caused by oscillation is difficult to avoid. To solve the mentioned problems, a machine learning method, capable of extracting patterns from the original data, should be introduced.
rough the extraction and training of critical information in considerable data, an assessment model is built, the autonomous identification of SSO frequency and attenuation factor is conducted, and a classifier is developed. Information is transmitted in time given different classifications, and online monitoring and early warning are conducted. In the meantime, by combining sensitivity analysis (SA) and principal component analysis (PCA), referred to as SA-PCA, the main factors of the SSO in the studied system are extracted, and their dimensions are reduced; besides, the dimensionality-reduced data are exploited in the training of the assessment model, which simplifies the model and increase its efficiency. us, the assessment results can be more precise.

Background
A study on the actual system in an area of China is conducted here [16]. is area has abundant wind resources, acting as a critical wind power base in China. On the one hand, wind-thermal power supply in this area is delivered by EHV AC line; on the other hand, it employs an EHV converter station to make it converted into DC and then transmitted to the load center via the transmission line. To enhance the transmission capacity of the EHV AC channel, fixed series capacitance compensation is also configured in the line. e system diagram and simulation diagram are illustrated in Figures 1 and 2.

The Determination of Main Factors by SA-PCA
Specific to a specific machine-network transmission system with series compensation configured, its internal components and their parameters may affect energy exchange in the subsynchronous frequency range and subsequently affect SSO characteristics. However, some components and their parameters have a greater impact. In actual engineering, its parameters can be regulated by predesign to underpin the reduction of the risk of SSO [17]. us, it is critical to know the effect of the mentioned components and their parameters on SSO.
In this study, the SA-PCA method is adopted to determine the main factors. Eigenvalue sensitivity refers to the corresponding variation of the characteristic root [18] when a system parameter is changing. It can indirectly reflect the effect of various parameter changes on the stability of the system, and it is critical to parameter setting and system stability analysis [19].
A linearization model of the entire system is built, consisting of a wind power plant model, a thermal power plant model, an EHV DC model, an EHV AC system model, a transmission model with series compensation, and an interface conversion model between the dq coordinate system and the xy coordinate system. A DFIG is installed in wind power plants, and its linearization model covers an asynchronous motor, a shafting system, a back-to-back converter, and its control strategy model. And the location of damping control in the DFIG system is presented in Figure 3. ermal power plant models consist of synchronous motors, excitation systems, and shafting models. EHVDC transmission systems include converter transformers, 2 Mathematical Problems in Engineering converters, reactive power compensation, and filter branches, as well as DC transmission. Set the A matrix in the linearization model of the system as a function of the parameter α; then it yields where μ i denotes the right eigenvector corresponding to the characteristic root λ i . Since μ i and λ i are implicit functions of α, in accordance with the principle of the derivative of the implicit function, the partial derivative of α in both sides is calculated.
ere is a premultiplication by transpose ] T i of the left eigenvector of the λ i matrix. Subsequently, since ] T i A � λ i ] T i , the terms on both sides are deleted. It yields where zλ i /zα denotes a complex number, which reflects the direction and magnitude of λ i 's movements with small changes in α. Its absolute value suggests the impact of parameter changes on system stability, and positive and negative indicate the change directions. Since zA(α)/zα is hard to obtain directly, this study uses the method of numerical differentiation, where Δα denotes a small microincrement. In this study, Δα � 0.02α. By analyzing the sensitivity of the system's eigenvalues, the primary factors of SSO of the system and the extent of the influence are listed in the following table. Table 1 lists the primary factors and influencing modes of SSO in the system. To facilitate the subsequent model training, the mentioned major factors are transformed into relevant input parameters. ere are 13 main input parameters, that is, the thermal power output active and reactive power P th and Q th , the wind power output active and reactive power P DFIG and Q DFIG , the HVDC output active and reactive power P HVDC and Q HVDC , the system current I, the rotate speed of DFIG ω r , the degree of series compensation K C , and two wind power control parameters (the rotor converter to speed stator voltage directional control integral coefficient K i_DFIG and the rotor converter to speed stator voltage directional control proportional coefficient K P_DFIG ), as well as two HVDC control parameters (i.e., the rectifier terminal constant current control integral coefficient K i_HVDC and rectifier terminal constant current control proportional coefficient K P_HVDC ). eir changes will to different extents affect the recognition of attenuation factor and the frequency in the later model.
In addition, Table 1 shows that the magnitude of system series compensation is the critical factor affecting SSO, and other factors show direct or indirect relationships to series compensation. If all relevant factors are directly inputted into the model training, it will cause information redundancy, reduce the running speed, and adversely affect the training and running of the model. PCA transforms multiple relevant variables into a small number of less related comprehensive indicators, thereby simplifying data and eliminating redundancy; however, the original information remains accurately reflected. is method aims to exploit the idea of dimensionality reduction to transform spatial coordinates, process the original interrelated information volume by data processing, and represent valid information in a novel form independent of each other. Assume that x 1 , x 2 , ..., x p are variable indicators, with a total of n samples, through which an n × p-order matrix can be formed. It is therefore suggested that there are 13 variable indicators, so p � 13, and the number of samples n is 10000. Of the mentioned related variables, 9000 sets of data are employed for the training model, and the other 1000 are adopted for the assessment model.  In principal component analysis, arithmetic processing is performed on the original variable data, and the specific process is as follows: (1) Standardize the raw data: the standardization process aims to draw a comparison in the characteristics' value between different dimensions, as an attempt to enhance the convergence speed and accuracy of the model. First, the mean μ and standard deviation σ of each information variable are calculated, and then the following formula is used to calculate the value of each information variable after normalization: (2) Calculate the sample correlation matrix: first, the correlation coefficient between the original information variables is calculated as follows: where r i,j (i, j � 1, 2, . . . , p) denotes the correlation coefficient between x i and x j as an attempt to obtain the correlation matrix Σ: (3) Solve the eigenvalues: a characteristic equation |λI − R| � 0 is set, the eigenvalue λ i (i � 1, 2, . . . , p) is solved, and then the first k eigenvalues and the unit vector a i (a i � [a i1 , a i2 , . . . , a ip ] T ) of the corresponding eigenvector are selected. Since the original information variables are to be compressed into the new characteristic subspace, a subset consisting of feature vectors that contain the most information should be taken. e size of the eigenvalue determines the significance of the eigenvector. By comparing the magnitude of the variance contribution rate, a suitable k value is taken. e variance contribution rate refers to the ratio of the eigenvalue λ i to the sum of all eigenvalues, i.e., λ i / p i�1 λ i . For the model adopted, the variance contribution rate of each input factor is presented in Figure 4. Figure 4 shows that the variance contribution rate of the 7th to the 13th input complies with the order of 10 −30 , so it can be ignored. Accordingly, K value here is obtained as 6. (4) Map the p-dimensional data to the k-dimensional subspace to build a novel sample matrix Z, where z 1 , z 2 , . . . , z k express the novel variable indicators.
us, the 10000 × 13 matrix is reduced into a 10000 × 6 matrix, which is convenient for training the following model. e correlation between the novel variables is presented in Figure 5. According to Figure 5, the novel variables are independent.

LM-BP Machine Learning
Online identification and classification of attenuation factor and frequency require the use of machine learning technology based on a neural network framework, exhibiting high nonlinear performance, good generalization ability, and the unique advantages and potential in automatically extracting characteristics and pattern recognition [20,21]. Accordingly, a novel research method is developed to study the fault diagnosis of complex systems [22]. is study also adopts Backpropagation neural network optimized by Levenberg-Marquardt algorithm to train and learn historical data to enhance the SSO online monitoring capabilities.

Principle of LM Algorithm.
According to the analysis of Section 3, 6 factors play a critical role in SSO, which can be represented by sets Z � z 1 , z 2 , . . . , z k ; the nonlinear relationship between the assessed value Y and the main factors Z is expressed as follows: where W 0 denotes the correlation matrix of the space on which the main factors depend; Z is the input vector of main factors. Substituting it into the LM algorithm, the neural network LM-BP prediction mathematical model is built. Set w (k) as the vector consisting of weights and thresholds of kth iteration. e vector of w (k+1) of (k + 1)th iteration is written as follows: where Δw calculated by Newton's method is defined as

Mathematical Problems in Engineering
where ∇ 2 E(w) denotes Hessian matrix of the output error and ∇E(w) is gradient. e output error is defined as follows: where e 2 k (w) is the error, d k is the predicted output, and y k is the actual output.
For the Hessian matrix, it satisfies the following relation: , where S(w) � l k−1 e k (w)∇ 2 e k (w) and J(w) is Jacobian matrix which is expressed as follows: In accordance with the rules of Gauss-Newton's law, is study adopts the LM algorithm, i.e., an optimization algorithm based on Newton's method; thus, equation (16) is converted as where the proportionality coefficient is a constant (μ > 0) and I is an identity matrix. According to the LM algorithm, if μ � 0, it represents a Gauss-Newton method; if μ > 0, it refers to a gradient descent method. When the LM algorithm is close to the practical value and μ approaches 0, the algorithm can converge efficiently and locally. When it is far away from practical value, μ increases, and the algorithm can converge globally. Furthermore, since [J T (w)J(w) + μI] − 1 is the positive definite matrix, so it is always reversible; in this regard, the LM algorithm is better than the Gauss-Newton method. e flowchart of the LM algorithm is illustrated in Figure 6.
On the whole, the mathematical model in intelligent online monitoring can be built.

LM-BP Training Process of LM-BP Neural Network
Model. To build intelligent online monitoring and predicting model for subsynchronous oscillation faults, the first step is to form training samples by collecting and screening historical load data. Subsequently, after normalization, sample data can act as input variables for the assessment model. e training steps are elucidated as follows: (1) Collect and select historical data information (2) Normalize and reduce the dimensionality of sample data to reduce the dimensionality (3) Determine network structure, including the number of hidden layers and nodes in each layer (4) Set network training parameters, which cover weights and thresholds of initial network, learning rate, expected error, and so on (5) Model assessment and actual output results output

Experiments and Results
e simulation in this study adopts the simulation model shown in Section 2. is model has been approved by the State Grid Electricity Regulatory Institute of China to simulate wind-thermal bundling and EHV AC/-DC hybrid transmission systems in a certain area of North China. It exhibits high accuracy and reliability.
First, the main factors extracted by the SA are verified by the simulation to impact SSO. Taking the rotor-side converter current inner loop proportional coefficient K i_DFIG as an example, by keeping other parameters unchanged, K i_DFIG is regulated to study its effect on the system's SSO characteristics.
As indicated in Table 2 and Figure 7, with the increase in the proportionality coefficient, the attenuation factor increases progressively, and the system tends to be unstable. As suggested from the electromagnetic transient simulation diagram, when the proportionality coefficient is 0.1, the system converges. However, when it increases to 0.3, the system diverges. When the proportionality coefficient continues to increase to 0.9, the system diverges faster than the situation when the proportionality factor is 0.3, and the amplitude is larger. It is therefore proved that this factor significantly impacts SSO.
Using control variables to verify other factors, it is indicated that the greater the system compensation, the greater the risk of subsynchronization and the more unstable the system will be; the greater the output power of DC transmission, the smaller the risk of SSO in the system and the stronger the stability will be; the greater the output power of thermal power plants and wind farms is, the stronger the system's SSO stability will be; DC transmission control parameters also have a certain effect on SSO stability. It is therefore demonstrated that the selected input parameters are of training significance. 10000 sets of sample data of this experiment are acquired with the simulation model. Here, 9 representative input parameters are taken, one of which is specific to all the units to run in the rated state, and the system parameter settings are appropriate. e other 8 groups, respectively, alter different input parameters. Moreover, the value of ω r and system current I from 2.3 s to 2.5 s are taken as the input parameters. e six major factors after dimensionality reduction in the sample data are determined by PCA. In Table 3, power is expressed in per-unit value.
When the neural network is being trained with dimensionality reduction parameters, the parameters of the neural network are optimized; it is also found that the training effect of the neural network with two hidden layers is enhanced when the number of hidden nodes in the hidden layer is 10 and 4, respectively. e neural network structure is illustrated in Figure 8. Mathematical Problems in Engineering e hidden layer transfer function is the Sigmoid function whose initial value is selected in the interval [0, 1], the learning rate generally ranges from 0.01 to 0.8, and the expected error is in the interval on the whole [0.001, 0.01]. In this experiment, the input layer node is taken as 6, the hidden layer node is taken as 10 and 4, and the output layer node is taken as 2. On that basis, a 3-layer BP network is formed. Its parameters are set to learning rate 0.1 and target error 0.001. Table 4 lists the assessment results and error numerical analysis results.
As revealed from Table 4, the assessed result based on LM-BP is significantly close to the practical value, and the average absolute errors of attenuation factor and frequency are 0.0473 and 0.026, which are not more than 5%, and the model classification is all correct, thereby proving its accuracy and reliability.  Figure 6: LM algorithm flowchart.   Original input parameters Z 4 (2) Figure 8: Neural network structure diagram.
Mathematical Problems in Engineering e fifth set of data in the table is selected to set the parameters of the simulation model, and the correctness of the assessed values was verified with the simulation model. e system current change chart is given in Figure 9. In this figure, the system current diverges from 2.0 sec. By 2.7 s, the current unit value reaches 4.0, and the divergence is significant. e Fourier analysis of this waveform is performed, and the analysis result is presented in Figure 10. e frequency of the oscillating current is 38 Hz, identical to the assessed value of 38.44, thereby verifying the reliability of the model again.
Likewise, the sixth set of data is selected for verification. e current change graph and spectrum analysis graph are illustrated in Figures 11 and 12. e system current diverges, and the current divergence frequency is 43 Hz, consistent with the data in the table as well. As indicated from the eighth set of data, the 42 Hz frequency current in the system takes up a large proportion, whereas the attenuation factor is negative, so the system converges.
e simulation results are illustrated in Figures 13 and 14, complying with the predicted results. e simulation results of attenuation factor and frequency are consistent with those predicted, thereby demonstrating the accuracy of the prediction.
All 1000 sets of data undergo model assessment and classification, and the classified results are presented in Figure 15. Tag 1 represents SSO, tag 0 expresses no SSO, and the accuracy rate takes up 99%.

Conclusions
In the present study, a method combining LM-BP machine learning and PCA-SA is proposed to solve the online monitoring and early-warning problem of the subsynchronous oscillation of wind-thermal power bundling and EHV AC-DC hybrid transmission system. By conducting the sensitivity analysis and the principal component analysis, the input signals of SSO in the complex system are screened and then dimensionalized, while their original characteristics are retained. By training the prediction model with the processed data, the attenuation factor and the frequency can be predicted more accurately and reasonably, and the accuracy of monitoring and early warning can be significantly enhanced. e attenuation factor and frequency are identified and then predicted with the LM-BP machine learning method. Accordingly, the problem that the warning threshold is largely determined by human experience is solved, and the real-time online monitoring of SSO is achieved. e accuracy of the proposed method is verified by simulation research. As revealed from the results, the proposed method has the accuracy for attenuation factor and frequency prediction reaching 99% and exhibits high performance when monitoring SSO. e subsequent work will highlight the use of artificial intelligence in the suppression measures of SSO. e occurrence of SSO will be mitigated by modifying the damping control strategy automatically in accordance with the precise positioning of the oscillation source.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare no conflicts of interest.