Nonlinear Analysis on Transmission Line Parameters Estimation From Noisy Phasorial Measurements

This paper proposes a complete and novel approach for estimation of impedance and admittance parameters of transmission lines from phasorial voltage and current measurements. The work uses the rigorous modeling of the noise on the phasorial measurements, and the solution is obtained based on a nonlinear approach of the estimation problem. The distortion introduced into the classical nonlinear model from the noisy data is analyzed by using statistical techniques. Finally, the performance of the estimation method is studied considering various load conditions, time-varying loads, number of samples, and noise on the phasorial measurements. The proposed method shows high accuracy in the line parameters estimation, even with distortions due to high noise levels. In addition, a good performance is also observed with a small quantity of samples, which is not observed in the linear solutions.


I. INTRODUCTION
The previous knowledge of the electrical parameters of transmission lines is essential for the parameter adjustment of equipment in power systems. Such parameters are required to set digital relays, calculate load flow, control and stability analysis, determining the propagation characteristics of the transmission line, and many other applications. In this sense, the reliability of an entire power system depends on the previous and accurate knowledge of the electrical parameters of the transmission lines, i.e. the series impedance and shunt admittance.
Conventionally, the electrical parameters of overhead transmission lines are obtained based on the physical and geometric characteristics of the conductors in bundles, steel towers, height of each conductor, conductivity of the soil The associate editor coordinating the review of this manuscript and approving it for publication was Meng Huang . under standard conditions, and other approximations. [1]. Furthermore, the calculation of the parameters takes into consideration some phenomena present in alternating current transmission, such as the skin effect in the phase conductors and the return current through the soil [2]. However, common-place methods for obtaining such parameters neglect the dynamic behavior of some very important constants during the calculation, e.g. soil conductance and electric permittivity of the air. A major example is in Carson's corrections, the soil conductivity is constant throughout the line length. Such assumption is definitely not true, since it depends on the local characteristics of the soil, such as stratification of the soil layers, humidity, and different geological formation in each region in which the line is through [3], [4].
One can notice that poor estimations of transmission line parameters affect several aspects of a power system. For this reason, to increase the accuracy of the estimations, the problem must be considered as a dynamical mathematical problem. This becomes especially relevant in the new context of smart grids and modern transmission systems, which require a more efficient parameter estimation methodology, in which their dynamic behavior and time-varying nature are taken into account.
Some parameter-estimation methods have taken place in the last decade in which dynamical equations are used to model measurements of transmission lines. In some cases, the parameter estimation is carried out in the time domain by employing current and voltage measurements during transient state at both sending and receiving terminals [5], [6]. Despite these methods show a good performance, the parameters estimation is only possible after transient occurrences (fault, abrupt varying load, lightning etc). On the other hand, there are also estimation methods that are developed in the frequency domain by using phasor measurements at both ends of the line [7]- [11].
Several methods propose the estimation by using equations based on the well-known π model of the transmission line, and the Gauss-Newton method for non-linear equations [12]- [14]. Although these methods can be applied in a steady state at any time, i.e. it is not restricted only to eventual transient occurrences, the phasorial methods do not consider the noise in the phasor measurements. A relevant improvement was introduced where systematic errors in the synchrophasor measurements were considered in the estimation method [15]. In this work, the noise is properly modeled as a Gaussian term added to the phase and magnitude of the phasor measurement. Each noise term related to each phasorial measurement presents different modeling according to the system load. Although this work had considered a more complete model for the phasorial measures, the level used was quite unrealistic, since the current works [16], [17] have shown that the phasorial measures present a noise with at least 1 % in the magnitude and 0.5 rad in the phase. Considering the estimation methods that use a nonlinear approach and some noise modeling, the solutions are obtained only by the Gauss-Newton method [13], [14] and there is no statistical validation for the model.
A few estimation methods are based on multiple measurements selected in an iterative and adaptive form. The noise was considered and modeled by using a Gaussian distribution [8], though the correct noise modeling was not considered as presented in [17]. It was considered a simplified model defined just as a percentage of the amplitude. These measurements compose a system of nonlinear exponential equations, where each parameter has a standard deviation, and the solution is found using a numerical procedure based on Newton's method. In addition, a method based on moving-window Total Least Squares (TLS) and Kernel density estimation was proposed [18]. In this last research, a load flow approach is used to describe a system with 8 equations with non linearity in all terms. The solution is obtained by fitting a Kernel density for each parameter.
Another TLS variant is employed as well [19]. In this case, the authors focus the modeling on differences between the measurements obtained from current transformers and capacitor voltage transformers. This approach leads to a solution using the Weight Total least Squares method (WTLS). Both methods do not consider a complete noise modeling, using a simplified model in which the noise is added directly in the equations. Thus, this last approach does not consider the error propagation of the phase and the magnitude. In a summarized description, some terms are considered Gaussian and unbiased, which is definitely not true.
Some accurate statistical models for noise analysis were also proposed in the technical literature [17], [20]- [22]. For example, in reference [17], a model is developed based on realistic specifications of ITs and PMUs. Since the specifications are given in polar coordinates, the conversion to rectangular coordinates is required to compute the noise covariance matrix, which is therefore used in different estimation methods. In this work, the solution is obtained and compared by using different methods, more specifically: OLS, TLS, and WTLS methods. Another important issue is that this work analyzes a transmittance matrix model for a short-length transmission line. A similar method was previously proposed in the technical literature, with an exact linear reformulation of the problem, which can be solved in closed form [20]. In such case, the distributed-parameter model for long transmission lines is considered, and its parameters are estimated in a non-iterative manner.
A few estimation methods have been proposed by using Kalman-filter based solutions, where voltage and current measurements at both line terminals are required to track parameters of long transmission lines dynamically [23]. Furthermore, these methods are applicable to transmission lines with different series compensation configurations. However, the data model is simplistic. In [23], the noise projection into the real and imaginary axes is not considered. Moreover, the practical aspects were not developed, e.g, it was not mentioned if the measures are provided by PMUs or Supervisory Control and Data Acquisition (SCADA).
Besides the are methods for noise modeling, there are also works that model the inaccuracy in each equipment as a non-random term, and use the error propagation based on the error theory and Monte-Carlo statistical simulations [24], [25]. These approaches intend to deal with systematic errors in instrument transformers and how it affects the estimation process.
It is worth highlighting that nowadays the estimation of the line impedance of distribution systems corresponds to a relevant issue derived from the general transmission line parameter estimation problem. In such case, the equations are based on the classical state estimation problem, and the common measures are mixed with PMU data to compose the problem. From such approach, it is possible to mention several contribution in the technical literature [26]- [29].
Currently, most of the estimation methods, available in the technical literature, solves the line parameters estimation problem by means of a simplistic linear approach. Furthermore, some of these methods do not consider properly the presence of noise in the measurements, which represents a unrealistic situations [13], [30], [31]. There are also some methods which take into account the noise modeling by using very simplistic methods, which could also lead to estimation errors [14], [15]. In this context, we propose a complete statistical analysis for the nonlinear solution applied to the transmission line parameter estimation problem.
The main contributions of this research are listed as • This paper proposes a new set of phasorial equations to solve the π− model; • The method uses a complete and realistic approach for phasorial measurements and the solution is obtained as a nonlinear problem; • The proposed approach shows an improved accuracy when compared to linear methods with similar noise modeling and even with works that employ a simplified noise modeling; • The classical nonlinear model and its statistical hypothesis are verified after the application of the method; • The performance of the solution is analyzed under several conditions of load, number of samples and noise level; • The robust study allowed us to analyze the condition that leads to the worst accuracy and how to avoid this situation.

II. TRANSMISSION LINE MODEL
A symmetric and balanced three-phase transmission line, classified as medium-length in the literature, is modeled in this work [32]. The equivalent circuit (positive sequence) for this line is presented in Fig. 1, whereV s ,V r ,İ s ,İ r are respectively the voltage and the current at the sending and receiving ends. From the equivalent circuit in Fig. 1, it is possible to describe the phasor equations of the proposed estimation method for the transmission line parameters R, b, X . Being b and X defined as where ω is the angular frequency of the system. By applying such laws to the circuit, the following equations are obtained: whereV sr ,V rr ,Ī sr ,Ī rr are the real parts of the exact phasor measurements, i.e without noise, andV si ,V ri ,Ī si ,Ī ri are the imaginary parts of the exact phasor measurements. It is possible to write the set of equations as the following matrix equationȲ =f (θ).
In the Equation (7), the vector of parameters (θ) is given by The other elements of the matrix equation are presented in appendix X.

III. NOISE MODELING
Considering the works which deal with transmission line parameter estimation, there are some ways to incorporate the inaccuracy of the measures into the method. In general, the noise is totally neglected [12] or considered in a simple way as a additive term in each equation which composes the estimation method [8], [14] [18], [19]. In the latter cases, all measures are considered with normal distribution and zero mean, which is not verified in field case [17].
The complete model for the data obtained from a PMU was presented by [33]. Following the approach given by this work, a generic phasorial measurew =X eφ is corrupted by a zero-mean Gaussian noise in the magnitude and phase. The magnitude term is denoted by ρ,x and the phase is φ,x . Therefore, a generic field PMU measure can be approximated by The noise on the measures represents the inaccuracy of each component that composes the measurement system. The main resources of inaccuracies are the instrument transformers (ITs), i.e, current transformers and voltage transformers, and the PMUs. The total noise could be considered as the sum of the noise from ITs and the noise from PMUs [17]. Individually, the PMU has a precision better than 1% in the magnitude and 1µs for measures of time [34], for 60 Hz of rated frequency, this value is equal to 3.77 × 10 −4 rad. Concerning the ITs, the inaccuracy depends on the class of operation. The main class of operations is shown in Tab. 1.
As the noise has a Gaussian distribution, the maximum values presented in Tab. 1 must be converted to the standard deviation. Fixing the standard deviation as the maximum value divided by three, the measurements lie to the interval bounded for the maximum value with a probability of 0.9973. Thus, such relation was used to model the noise.
Considering noise in all terms, the equation (7) is not more where Y is the vectorȲ considering noise in the measures and f (θ) is the vector functionf (θ ) when the noise is incorporated. Now, the residual vector is defined as The optimization problem P1, which corresponds to the solution by using the nonlinear approach is represented by whereθ is the solution that minimizes the 2-norm of the residual vector.

IV. LEVENBERG-MARQUARDT METHOD (LMM)
The Gauss-Newton method is the most common method to solve nonlinear least-squares problems. This method is simple to implement, however, it presents a poor convergence when the Jacobin matrix is near to the singular matrix [35]. Furthermore, the solution in the Gauss-Newton method is strongly dependent on the initial guess. Thus, the method might be present a poor accuracy when the initial guess is far from the optimal solution [36].
An alternative to overcome these issues is the Levenberg-Marquardt method, in which the objective function is given by where, θ (k) is the vector of parameters in the iteration k, r(θ, θ (k) ) is the first order Taylor approximation of r(θ) centered in θ (k) , and λ (k) is a positive parameter updated in every iteration. It is worth highlighting that the Levenberg-Marquardt method is a multi-objective formulation, where the aims are to minimize the norm of the affine approximation of r( . ) and to choose new approximations that are near to previous iterations.
For this algorithm, the solution is calculated as The operator J(θ (k) ) is the Jacobian matrix of the function r( . ) calculated in θ (k) .

V. SIMULATION ASPECTS
In this section, we will present the simulation structure used to obtain the phasorial measures and the different scenarios considered in the work. The transmission line analyzed in this paper has the parameters presented in appendix X, obtained from [38]. Each lumped parameter R, L, and X is calculated by multiplying the distributed parameter R , L , and X to the line length l, i.e: In Fig. 2, R shunt is necessary to ensure numerical convergence of the simulation, and its numerical value should be small compared to the other lumped parameters of the line. Therefore, the shunt resistance was set to R shunt = 0.05 . We have analyzed the load as a percentage of the Surge Impedance Loading (SIL), defined in [38]. In our case, the SIL of the line is 140 MW.

A. ALGORITHM DESCRIPTION
The algorithm implemented in this work can be divided into three different parts. The first part corresponds to the numerical modeling for the real data from field applications of three phase transmission lines. In such part, the structure VOLUME 10, 2022 presented in Fig. 2 is used to create the phasorial measures without noise. Then, the noise is artificially added upon the absolute value and the phase of each measure, according to that it was described in Section III. It is important to emphasize that the noise has two components, the first originated by the ITs (α, β) and the second related to the PMUs (ρ PMU .φ PMU ). The algorithm 1 represents the pseudo-code of the real data generation process. The variable N represents the number of samples used on the estimation method. The vector function used on the estimation method, denoted by f conc (θ), is the concatenation of the function f (θ) for each set of phasorial equations, i.e the function has dimension 6N. It guarantees more stability for the method since the method has more linearly independent equations. Algorithm 1 Generation Data 0: Run MATLAB Simulink using Table 6 0: ObtainV s ,V r ,İ s ,İ r without noise The second part of the algorithm corresponds to the implementation of the Levenberg-Marquardt method to solve the nonlinear least squares problem. For the proposed study, the Levenerg-Marquardt implementation chosen was the common-place design presented in [39]- [41]. Such implementation can be accessed through the built-in function presented in MATLAB called lsqnonlin. A simplified flowchart for the proposed solution is shown in Fig. 3. The Levenberg-Marquardt part is based on the algorithm presented in [42].
Finally, the last part is dedicated to analyzing the results obtained from the estimation process. The quality of the calculated parameters is verified by using the absolute value of the relative error (9) wherep is the estimated parameter and p is the exact parameter. Furthermore, in order to complete the statistical description for the method, we proposed a residual analysis for the model based on the ordinary definition for the residuals in nonlinear least squares problems [43] where r can be described by six entrances, defined as r = [r 1 r 2 r 3 r 4 r 5 r 6 ] T .

VI. RESULTS
In this section, the application and performance of the estimation method will be discussed through computational simulations. First, the estimation method will be analyzed considering that the line carrying an active power corresponding to 0.9 SIL. For that case, the class of errors in ITs chosen was class 0.1, which leads to ρ = 1.1% and φ = 1.9 · 10 −3 rad. The number of samples (N ) was considered equal to 200. The other parameters, necessary to set the Levenberg-Marquardt algorithm, were considered equal to the default value present in the optimization toolbox from MATLAB [44]. In order to establish a fair comparison of the proposed method and others given by the literature, it was used the methodology presented in [45], [46]. The result of this analysis is shown in Tab. 4. It is worth emphasizing that the methods, chosen to compose the comparison, have a similar approach to estimate the transmission line parameters, i.e, these methods are using phasorial measures (from PMU or SCADA), consider some noise modeling upon the data, and evaluate the result by using the relative error.
The result in Tab. 4 shows that the proposed method has better accuracy compared with others in the literature, even with those with a very small noise level [15], [17]. Another relevant aspect is the number of samples used in the estimation. We reached better results using a reduced number of samples when compared with [8]. However, this last method does not consider the error propagation and neither the relationship between correlated measures, which becomes the solution more stable, although less realistic.

A. ANALYSIS OF THE RESIDUALS
The standard nonlinear regression model can be represented as [43] y where x i represents represents a vector of known variables associated with the ith observable response y i , θ is a p × 1 vector of unknown parameters, the response function f is assumed to be known, continuous and twice differential in θ, and the residual vector The operator E[ . ] represents the expected value of a random variable [47].
However, the noise modeling for the phasorial measures and the set of equations lead to a distortion in the standard nonlinear model. Basically, it occurs due to two reasons. First, the noise is independent and identically distributed (i.i.d) for the absolute value and the phase of each measure, but it is not true for the projections in the real and imaginary axis. Second, the set of equations mix the observations measures (x i ) and the output measures (y i ). Then, it is not possible to separate the residual term as an additive term in the equation. Therefore, we propose to analyze such modification introduced by the modeling in the residuals as a form to improve the statistical discussion for the method.
Using the estimated value for each parameter (p), it is possible to construct a histogram of the residual for the model. Remembering that the residual is a vector, then there are six different histograms for the residual. To illustrate the behavior of the residual in the method, we present the QQ-plot for r 1 , presented in Fig. 4.
The histogram and the QQ-plot show that the residuals has approximately a normal distribution, with a small perturbation for negative high values. The same behavior  is observed for the others entries. In order to establish if each entry has mean zero, it is possible to do a two tailed hypothesis test, using the following definition Using the definition of the p-value for each sampling distribution presents in the residual vector [48], it was obtained the results shown in Tab. 2.
Fixing the level of confidence in 95%, it is possible to affirm that the average value of the entries r 1 , r 3 , r 4 , r 5 , r 6 and the supposed mean (µ ri = 0) do not have significant statistical difference, while the average value of the entry r 2 have. Thus, the noise modeling and the estimation method introduced a bias in the second entry of the residual vector.
In order to verify the effect of the noise modeling and the estimation method on the correlation of each residual vector entry, we proposed to use the Pearson coefficient (ρ) of the random variables X and Y , defined as [49] where Cov(X , Y ) is the covariance between X and Y , Var(X ) represents the variance of the random variable X , and ρ ∈ [−1, 1] The result is presented in Tab. 3. As can be seen in Tab. 3, the residual vector presents a high distortion considering the correlation between each entry. This can be explained by analyzing the set of phasorial equations and noise modeling. If the model is ideal, the residual vector is directly related to the output vector (Ȳ ) and this vector has some identical entries, more precisely the second is equal to fourth and the second is equal to the sixth. Then, before applying the method, it was expected  that such entries might present a high correlation. Moreover, the noise modeling introduces on the equations terms are correlated, e.g, the terms V ei and V er are correlated because the phase for both has the same model. Therefore, the set of phasorial equations and the noise modeling give to the model an intrinsic correlation for the residual vector entries. As the proposed noise modeling is the most generic and rigorous form to include noise on the measures [16], [22], the behavior shown in Tab. 3 can be considered acceptable since the values of the estimated parameters do not present a significant deviation.

VII. INFLUENCE OF THE LOAD ON ESTIMATION METHOD
According to demand of the power system, the transmission line can carry different active powers during a day of operation. In order to study the influence of the load condition within the deployment time-window of the estimations, we have simulated different load conditions, as a percentage of the natural power. The result of the relative deviation for each parameter using a different load is presented in Fig. 5.
By observing Fig. 5, it can be seen that the relative error for the susceptance increase when the load increase. Such behavior can be explained by the problem formulation. Note that the first and second equations depend on the susceptance only. Basically, it is possible to write two simple equations to calculate this parameter where each projection of a phasorial measure on complex plan is a random variable as described in Section III. The deviation present in parameter b is related to the bias introduced in (11) by the noise modeling. Such bias (v) can where b is the actual value of the parameter b and E[ . ]is the expected value operator [48]. To prove the relation between the behavior present in Fig.5 and the bias, we performed a numerical analysis of the absolute bias for all load conditions. The result is present in Tab. 5.
Comparing the relative deviation for the parameter b present for each load and the corresponding bias, it is possible to establish a direct relationship between the bias and the relative deviation. The noise modeling leads to a bias on b calculation, such bias affects the method and it is reflected on the relative deviation. The remaining parameters present a relative deviation smaller than 1% for. Besides, with respect to these parameters, it is possible to affirm that they present an oscillating value, i.e, it is not increasing or decrease with the load. Such fact occurs because the equations used to calculate these parameters are nonlinear. Then, the behavior is more difficult to predict than the susceptance case.
Another explanation on the susceptance case is that the transmission line represents a capacitor circuit during operation with low load demand, i.e. the reactive power is mainly concentrated in the shunt capacitance in the equivalent circuit in Fig. 1. Otherwise, for a higher load profile, most part of the reactive power is concentrated in the series inductance because of the high current demand from the load. Thus, best estimation values for susceptance are obtained from lower power demands whereas more accurate series parameters are estimated from higher load values.
Besides the linear variation of the load represented in Fig. 5, it is possible to analyze the behavior of the method under a non-linear variation of the noise, i.e, considering an abrupt variation of the active power carried through the line. For such analysis, the load profile is shown in Fig. 6.
In order to verify the accuracy of the proposed method under the load profile presented in Fig. 6, it was performed a study of the relative deviation of the parameters. The result of such analysis is given by Fig. 7.
Observing the Fig. 5, it is possible to conclude that the method has a stable performance even when the load is abruptly varied. Moreover, it is worth highlighting that the susceptance and the resistance have an opposite behavior. The   resistance has poor accuracy for small loads, for example, when the load is 0.05 SIL the relative deviation for the resistance is approximately 1 %, whereas for the susceptance the relative deviation increases when the load increases. Nevertheless, for all load conditions the accuracy is better than 10%, the value presented as desirable in [50].

VIII. THE INFLUENCE OF THE NOISE STANDARD DEVIATION ON THE ESTIMATIONS
.
The standard deviation of the noise is related to the quality in which the phasorial measurements are obtained in a power system. The main parts of the errors present in phasorial measurements are associated with measurement instruments. Such errors are difficult to estimate. Then, the standard deviation can variate in a large range of values according to the situation [51].
Therefore, by varying the total standard deviation of the noise (sum of the noise from ITs and PMUs), it is possible to observe the performance of the proposed methods with respect to this factor. For this analysis, the number of samples was fixed in N = 10000 and the load condition fixed in 0.9 SIL. First, we fixed the noise related to the phase in 8 · 10 −3 rad and variate the standard deviation related to the magnitude, the result obtained in this analysis is shown in Fig. 8. After, the standard deviation related to the magnitude was considered fixed equal to 1% and the term related to the phase was modified as shown in Fig. 9.
Comparing with the load, the standard deviation of the noise has a higher influence on the relative deviation for each estimated parameter. Especially, when it is analyzed the results presented by the susceptance in the case in which the standard deviation for the magnitude was modified. All the parameters present a poor performance when the   standard deviation level assumes a high value, it can be explained by analyzing the noise modeling. As presented in [17], the bias present in the measures is affected by the standard deviation of the noise. Then, a high value of the standard deviation leads to a high bias on the measures, and, consequently, the residual vector presents a high norm. Therefore, the relative deviation for the estimated parameters presents a bad performance.
Considering that accuracy in line parameter estimation below 10% would be desirable [50], it is possible to affirm that the proposed method present a efficacy for ρ ≤ 3.0 % and φ ≤ 30 · 10 −3 rad.

IX. INFLUENCE OF THE NUMBER OF SAMPLES
We have simulated the line with a load of 90% of the SIL and varied the number of samples. The noise level was considered the same as presented in Section VI. The result of such analysis is shown in Fig. 10.
The values presented in Tab. 10 show that the relative deviation for each parameter is small affected for the number of samples. Analyzing the results for each parameter, the same order of magnitude is observed for the relative deviation. As all deviations are less than 10% [50], they can be considered acceptable. Another relevant aspect is that the relative deviation is oscillating with respect to the number of samples. It occurs for one main reason, the solution is calculated using a nonlinear least square approach. Then, the intrinsic correction given by the iterative method provides more stability for the final estimated parameter. It can be proved when such solution is compared with OLS, TLS, and WLS solution, where the relative deviation depends strongly on the number of samples [17]. For those estimators, the result present by the method using a large number of samples is compared with the present in this work when they are used less than 100 samples. Therefore, the proposed method presents an advantage when compared with the others present in the literature since it is possible to obtain the same relative deviation using small samples.
An important aspect related to interactive methods is the computational cost when compared with non-iterative methods. Indeed, if the number of samples is the same for the methods, a non-iterative one presents a lower computational cost. However, the non-iterative estimation methods present in the literature [17], [20], [22] use a large number of samples to obtain a satisfactory performance for the parameters, such methods employed more than 1000 samples in the estimation process. The presented method shows good performance even for less than 200 samples. Besides, the interactive method presented a fast convergence in the studied problem. For all conditions, the number of iterations was less than eight.

X. CONCLUSION
A novel nonlinear formulation has been presented in this paper to estimate the positive sequence parameters of a three phase transmission line. The complete and rigorous modeling for the noise was incorporated into the phasor measurements. The influence of such modeling on the classical nonlinear model was analyzed and it was possible to show numerical proof that some entries of the residual vectors present a distortion in the mean value and mainly the residual vector present entries strongly correlated. These facts lead the estimation to present a systematic error in the estimated parameters.
The simulation performed in this paper showed the importance of the accuracy in which the phasor measurements are obtained for the transmission line parameters estimation. Besides, we verify the efficacy of the method for many conditions of load, number of samples, and noise modeling. Comparing the proposed method with the others that consider the rigorous noise modeling as well, but using a linear approach, the method presents a better performance for many conditions. Also, the nonlinear approach presented a good performance even for a small number of samples.

APPENDICES A DESCRIPTION OF THE EXACT MATRIX EQUATION
The elements of the exact matrix equation are described bȳ