See

discussions, stats, and author profiles for this publication at: [Link]

Identification of multiple harmonic sources in

power systems using independent component
analysis and mutual information

Article March 2010

CITATIONS READS

6 328

3 authors:

Masoud Farhoodnea Azah Mohamed

National University of Malaysia National University of Malaysia
39 PUBLICATIONS 193 CITATIONS 573 PUBLICATIONS 3,638 CITATIONS

SEE PROFILE SEE PROFILE

Hussain Shareef
United Arab Emirates University
324 PUBLICATIONS 1,677 CITATIONS

SEE PROFILE

Some of the authors of this publication are also working on these related projects:

Impact study of PHEV in distribution network with considering renewable energy View project

Islanding detection and load shedding scheme for radial distribution systems View project

All content following this page was uploaded by Hussain Shareef on 17 March 2017.

The user has requested enhancement of the downloaded file. All in-text references underlined in blue are added to the original document
and are linked to publications on ResearchGate, letting you access and read them immediately.
Identication of multiple harmonic
sources in power systems using
independent component analysis
and mutual information
Masoud Farhoodnea, Azah Mohamed, Hussain Shareef

Department of Electrical, Electronic and Systems Engineering, University Kebangsaan, Malaysia

E-mail: azah@[Link], hussain_ln@[Link]

This paper presents a novel technique for accurate determination of the probable locations of multiple harmonic sources in power
systems. In the proposed methodology, independent component analysis (ICA), which is one of the blind source separation techniques,
is applied for estimating the proles of injected harmonic currents produced by nonlinear loads. By using the ICA algorithm, the
reconstruction of the injected harmonic currents can be implemented easily without prior knowledge about the system components.
From the reconstructed harmonic current proles, the mutual information theory is then applied for determining the exact location of
harmonic sources by computing the pair-wise mutual information between the extracted current trends and bus voltages. Numerical
simulations were made to prove the accuracy of the proposed method in locating multiple harmonic sources in the IEEE 30 and 34
bus test distribution systems. Results proved that the proposed method can accurately determine the probable location of harmonic
sources in radial and non radial power systems with mean square error less than 0.00035.

Keywords: Harmonic distortion; harmonic source location; power quality; multiple harmonic sources.

or lessen the problems related to harmonic distortion. Meth-

1. INTRODUCTION ods for identifying the location of harmonic sources in power
systems are generally categorized as single point and multiple
Due to the penetration of various distributed nonlinear and har-
point strategies. The real power direction method [1] is one of
monic producing loads, the diffusion of harmonics in power
the earliest single point method proposed in locating harmonic
systems is escalating and has become an important power
sources. However, the real power direction method has been
quality problem. In addition harmonics are the main rea-
found to be only 50 percent reliable and therefore the accuracy
son for causing resonance problems, overheating in capacitor
of the method is questionable [2]. Other single point meth-
banks, transformers and conductors, and telecommunication
ods for harmonic source localization include incentive-based
conflictions. These problems may decrease the reliability and
method [34], the critical impedance method [5] and voltage
increase the maintenance costs of power systems. From the
magnitude comparison method [6]. These methods require
above-mentioned issues, it is important for power providers
complete knowledge of system parameters in harmonic fre-
to identify the source of harmonic distortions so as to solve
quency or implementation of switching tests for obtaining the
E-mail: farhoodnea_masoud@[Link] harmonic impedances of the system. However, such data are
 IDENTIFICATION OF MULTIPLE HARMONIC SOURCES IN POWER SYSTEMS

not available in most of the time, and implementing switching In addition, unlike the above mentioned techniques, the pro-
tests is not possible in practical power systems. More recent posed method only requires voltage measurements of all buses
single point based methods such as total harmonic distortion and eliminates other kinds of measurements. In addition, it is
(THD) method [7], harmonic vector method (HVM) [8], and assumed that there are no prior information about the system
IEEE14592000 standard based method [9] tried to solve the parameters, harmonic impedance matrix and characteristics
above mentioned problems in the previous methods. Nonethe- of harmonic sources. Hence estimation will be done com-
less, these current methods still suffer from some problems pletely in blind situation. This assumption is significant in
such as the use of unsuitable index harmonic modeling and the deregulated power systems where complete power system
inability to determine the location of harmonic sources for parameters, load behavior and utility response is not known
individual harmonic orders. In general, implementing single because of economic and security restrictions in electricity
point based methods are very economical in comparison with market.
other methods, but such methods are unable to work in practi-
cal systems because harmonic sources exist at different points
2. INDEPENDENT COMPONENT
in a power network. To tackle this problem, several methods
have been developed for locating multiple harmonic sources ANALYSIS
in power systems. Harmonic state estimation (HSE) [1011]
is one of the most popular multiple point based methods. But, Independent component analysis (ICA) is one of the blind
the HSE based methods require complete knowledge about the source separation techniques (BSS) based on the statistical
system parameters for each individual harmonic order, which independence between the source signals. Due to the gen-
are usually unknown or difficult to calculate in practice. In ad- erality of the ICA model, it has found many applications in
dition, the measurement placement technique [1213] which different areas such as audio separation, econometrics, brain
is widely applied in HSE for optimizing the number of mea- imaging, and telecommunication data separation. Basically,
surements in the system needs to have prior knowledge about ICA is a technique that transforms the observed signals into
harmonic sources and their locations. Another concern with a linear transform of source signals that are statistically and
HSE is that this approach needs various types of harmonic mutually independent from each other. To reconstruct the
measurements such as voltage, real and reactive powers which source signals from the observed signals, certain assumptions
are costly for large systems. In [14], the weighted least square are considered such that the source signals are statistically in-
estimation technique and Euclidean norm have been applied dependent, all but one source signals must have non gaussian
to recognize the location of multiple harmonic sources. In distributions and the number of observations should be greater
this method, the weighted least square estimation technique or equal to the number of sources [16].
and HSE are used to locate the possible buses with harmonic Assuming that there are N sources and M measurements,
sources while the Euclidean norm determines the exact buses the linear mixing model of ICA can be written as
with harmonic sources. However, the disadvantage of this
xj = aj 1 s1 + aj 2 s2 + ... + aj n sn for all
method is the use of all types of harmonic measurements at
all buses which is not economical. Another multiple point j = 1, 2, ..., M (1)
based method is the Independent component analysis (ICA)
The matrix representation of Eq. (1) is
method [15]. In this method, measured voltages at the selected
buses, determined by the measurement placement technique X = AS (2)
are used for estimating the harmonic impedance matrix of
the system, and harmonic current profiles of the harmonic where,
sources at each harmonic frequency. The minimum electri-
cal distances between the estimated impedance matrix and the S(t) = [s1 (t), s2 (t), ..., sn (t)]T : N-dimensional vec-
actual impedance matrix are then obtained by using the ex- tor of unknown source signals
haustive search technique. The buses with minimal distances
X(t) = [x1 (t), x 2 (t), ..., xm (t)]T : Mdimensional
are considered as the location of harmonic sources. However,
vector of observed signals
this method is not quite practical because it needs to determine
the actual impedance matrix of a system at each harmonic fre- A : M N full column coefficient matrix aij called
quency. In addition, this method requires prior information as the mixing matrix
about system parameters and historical records of the loads
for solving indeterminacies related to ICA algorithm. In ICA, the objective is to find estimates of the S and A
This paper presents a new method based on ICA and mutual from the available observation vector X. The unknown source
information theory (MI) to determine the accurate location of matrix S can be estimated as
multiple harmonic sources in a power system. In this method,
ICA is used to rebuild or estimate the harmonic currents pro- Sest = W X (3)
duced by the harmonic sources while the MI theory is used
where,
for determining the locations of harmonic sources, respec-
tively. The proposed method for identifying the location of Sest : estimate of the sources S with N T dimension
harmonic sources is different from the previous ICA method
[15], in the fact that it uses the mutual information theory for W : N Mseparating matrix which is the pseudo inverse
final decision making stage in locating the harmonic sources. of the mixing matrix A
M. FARHOODNEA ET AL

For estimating the model of ICA given by Eq. (3), sev- where D and E are the diagonal matrix of eigenvalues and
eral optimization techniques can be used such as maximiza- the orthogonal matrix of eigenvectors, respectively. They are
tion of nongaussianity, minimization of mutual information, obtained from the covariance matrix of the observed vector
and maximum likelihood estimation [17]. Here, the ICA al- X.
gorithm with maximization of nongaussianity of the source For maximizing the contrast functions in Eq. (7), the Fas-
signal has been applied. Based on this approach, the entropy tICA algorithm [18] is applied and it is described in terms of
of a random variable which is related to the information that a flowchart as shown in Figure 1.
the observation of the variable gives is considered. The more Note that for convergence of the algorithm, the old and new
random the variable is, the larger is its entropy [17]. The en- values of w should be in the same direction.
tropy, H of a random vector, y with density function, py ()
is defined as:


3. MUTUAL INFORMATION THEORY
H (y) = py () log py ()d (4)
Basically, mutual information (MI) is a measure of the infor-
where y is a random variable with mean and standard devi- mation that members of a set of random variables have on
ation . the other random variables in the set. Using the concept of
It should be noted that a gaussian random variable has the entropy defined for continuous random variables in Eq. (4),
largest entropy among all random variables with equal vari- the discrete form of entropy is derived as,
ance. Therefore, entropy can be used as a measure of non-
gaussianity of random variables, which is always nonnega-

N
tive. An important measure of nongaussianity is given by the H (y) = p(y) log p(y) (9)
normalized version of differential entropy which is called as y=1
negentropy J and is defined as:
J (y) = H (ygauss ) H (y) (5) The entropy in Eq. (9) can be interpreted as a measure of
uncertainty of the events with probability function p . In other
where ygauss is a gaussian random vector with the same co- words, entropy is small for events with probability close to
variance matrix as y. 0 or 1, and large for probabilities between 0 and 1. In fact,
Negentropy is always nonnegative, and it is zero if y has the entropy of a random variable can be defined as the degree
a gaussian distribution. Estimating negentropy of Eq. (5) of information that the observation of the variable gives. To
is very difficult because it needs to estimate the probability determine uncertainty about variable A in subsequence of tri-
functions of the source signals. Hence, by using only one als in which B occurs, the conditional entropy is used [19].
quadratic function, the approximation of J becomes [18]: Assuming random variables A and B consisting of the N ele-
J (y) [E {G(y)} E {G(v)}]2 (6) ments of ai and bi for i = 1, 2, ..., N , the conditional entropy
can be defined as:
where G is practically any nonquadratic function with slow
variation such as tanh(u), v is a gaussian variable with zero 
N
mean and unit variance and y is a random variable with zero H (A |bi ) = p(ai |bi ) log p(ai |bi ) (10)
mean and unit variance. i=1
To estimate the independent components, maximization of
the estimated negentropy (6) is necessary. By maximizing the Assuming that B is the average of H (A |bi ), the conditional
sum of N one unit contrast functions and taking into account entropy of variables A and B can be written as:
the constraint of decorrelation, an optimization problem is
obtained as follows: 
N
H (A |B ) = p(bi )H (A |bi ) (11)

n
i=1
Maximize J (wiT x)
i=1
  By using Eq. (10) and (11), the mutual information I of
under constraint E (wkT x)(wjT x) = j k (7) variable A and B can be obtained as:

where wi for i = 1, 2, ..., N are the rows of the separating

I (A,B) = H (A) H (A |B ), I (A,B) 0 (12)
matrix, W which is inverse of coefficient matrix A.
Before implementing ICA, the observed vector X is prepro-
Equation (12) is symmetric and can be interpreted as the
cessed by centering and whitening. Preprocessing in ICA is
information about A contained in B and vice versa. If A
known to be a suitable technique for simplifying the ICA esti-
and B are independent, their mutual information should be
mation .[16]. By centering, the mean of vector X is subtracted
zero. In addition, mutual information is a better function for
and transformed to zero-mean variable, while by whitening;
measuring dependency between two random variables rather
the observed vector X is transformed linearly to a new vector
than the correlation function [20]. The reason is that mutual
X where its components are uncorrelated and their variances
information is able to measure general dependency between
are unity. In whitening, Eq. (8) is considered,
two variables, while the correlation function can only measure
X = ED 1/2 E T X (8) linear dependency between them.
 IDENTIFICATION OF MULTIPLE HARMONIC SOURCES IN POWER SYSTEMS

tion under non-sinusoidal condition which can be written in

Start matrix form as:
V h = ZhI h (13)
where V , Z and I are the bus voltage, impedance and bus
injected current vectors, respectively and h is the harmonic
Choose the number of
order of the frequency.
independent components
If voltage V and impedance Z are identified, then the har-
monic current I can be computed using Eq. (13). However,
it is difficult to form the impedance matrix because of com-
Centering of the observed plexity of obtaining the values of harmonic impedance of the
vector X system for individual harmonic orders. Assuming that in the
absence of any information of the system parameters and char-
acteristics of non linear loads, the only available information
Whitening of the observed is the bus voltages which are obtained by distributed harmonic
vector X using Eq. (8) voltage measurement devices for each harmonic order.
The linear system of equations for the harmonic load flow
(13) under non sinusoidal condition is in close similarity with
the model in Eq. (2) which makes ICA a useful instrument for
Initializing weight vector w multiple harmonic source estimation. In addition, fast fluctu-
randomly with unit norm ations of the loads consumption during specific time that are
related to second to second or minute to minute variations,
have non gaussian characteristics and act independently [21
23], which satisfy the preliminary conditions for applying ICA
Update separating matrix w using gradient algorithm. Applying ICA to Eq. (13), X represents the har-
optimization algorithm to maximize monic voltage measurement vectors V, the mixing matrix A
negentropy by taking the gradient of Eq. represents the harmonic impedance matrix Z in the harmonic
(7) domain, and the estimated signal S represents the bus injected
current vector. In Eq. (13), the current matrix I is zero for
linear loads and non zero for harmonic or non-linear loads.
From (13), the injected harmonic currents generated by the
Normalizing separating matrix w to a new non-linear loads at the respective buses have a relationship
matrix with unite norm with the bus voltages. By using the mutual information the-
ory, the produced harmonic current at each bus has maximum
mutual information with its own bus voltage and the mutual
information between this harmonic current and other bus volt-
No ages is reduced due to the current division between branches.
If converged The procedure in implementing the proposed method using
ICA and mutual information theory for estimating the location
of multiple harmonic sources is summarized as follows:
Yes Measure harmonic voltages at all buses.
I. Estimate current traces of the harmonic sources using the
Repeat the algorithm for each independent
fast ICA algorithm (Figure 1).
component
II. Calculate MI between the reconstructed harmonic cur-
rent traces and harmonic bus voltages.

End III. Estimate the location of harmonic sources by using com-

puted mutual information.

Figure 1 Implementation procedure of the FastICA.

5. SIMULATION RESULTS
To validate the performance and accuracy of the proposed
4. HARMONIC SOURCE method in radial and non radial distribution systems, two test
LOCALIZATION USING ICA AND MI systems include the non radial IEEE 30 bus and radial IEEE 34
bus test systems are used. In both test systems, three harmonic
For estimating the probable location of multiple harmonic sources are considered at the selected buses. In addition, for
sources in power systems firstly, consider the system equa- both cases it is assumed that all loads are with constant power
M. FARHOODNEA ET AL

Figure 2 IEEE 30 bus non radial test system.

Figure 3 Actual and reconstructed harmonic currents at bus 7, 16, and 30 for the 5th harmonic order for the 30 bus test system.

factor and the harmonic loads are modeled as harmonic cur- 5.1 Results of the Non Radial 30 Bus Test
rent injection sources with harmonic spectrums given in [24]. System
The procedure of estimating the location of harmonic sources
using ICA begins with harmonic voltage measurements. Har- In this subsection, to verify the accuracy of the proposed
monic power flow simulations are performed to compute bus method, a non radial IEEE 30 bus test system shown in Figure
harmonic voltages and these voltages are considered as mea- 2 is used [25]. In this test system, three harmonic sources
surement vector, X in the ICA algorithm. The injected har- containing the 5th , 7th , 11th , 13th , and 17th harmonic orders
monic currents are then reconstructed as illustrated in Figure are placed at bus 7, 16, and 30. The PCFLO software [26] is
1. As we use random number for initialization of the sepa- employed to perform harmonic power flow simulations so as
rating matrix w, an average of ten runs of the algorithm is to determine the harmonic voltages at all the buses. To create
considered. The locations of harmonic sources are then de- variations in the harmonic measurement vectors, harmonic
termined by using the mutual information theory. The test power flows are simulated for different loading conditions.
results for the two test systems are presented in the following The loads are then multiplied by the Laplace distributed ran-
subsections. dom variables to create different operating conditions for each
 IDENTIFICATION OF MULTIPLE HARMONIC SOURCES IN POWER SYSTEMS

Figure 4 Actual and reconstructed harmonic currents at bus 7, 16, and 30 for the 7th harmonic order for the 30 bus test system.

Figure 5 Actual and reconstructed harmonic currents at bus 7, 16, and 30 for the 11th harmonic order for the 30 bus test system.

load at each sampling time. Harmonic bus voltages are com- To evaluate the accuracy of the ICA algorithm, correlation
puted by solving the linear system equation (13) for each har- coefficients of the actual and the reconstructed signals are
monic frequency of interest. In this case, 350 harmonic volt- computed as shown in Table 1. From the table, it is clear that
age samples have been generated to represent the harmonic all the correlation coefficients are within the acceptable range,
measurement vector, X. Since the number of harmonic mea- that is, close to 1. This implies high accuracy between the
surements is greater than the number of harmonic sources, estimated and the actual harmonic currents. Table 2 shows the
preprocessing procedures are taken into account to reduce the mean square error of the signals so as to quantify the difference
dimension of the harmonic measurement vectors, X [17]. For between the estimated and the actual harmonic currents. The
reconstructing variations in the harmonic currents produced results of Table 1 and 2 show that the Fast ICA algorithm can
by the harmonic sources, the Fast ICA algorithm programmed rebuild signals with high accuracy in which the mean square
in Matlab is used. Figures 3 to 7 show the actual and the ex- error of the estimated and the actual harmonic currents is less
tracted current traces at bus 7, 16, and 30 for the 5th , 7th , 11th , than 0.0015.
13th and 17th harmonic orders, respectively. In the figures, To identify the location of harmonic sources, the mutual
the red lines represent the actual or original harmonic currents information theory has been applied using Eq. (15) to calcu-
while the blue lines represent the reconstructed harmonic cur- late the pair-wise mutual information between the extracted
rents over the 350 samples. It should be noted that for ease of current trends and bus voltage for each harmonic frequency.
comparison, both signals in each figure are normalized with The results of the mutual information are shown in Table 3,
their largest values in per-unit (p.u). in which each column and row represents the estimated cur-
M. FARHOODNEA ET AL

Figure 6 Actual and reconstructed harmonic currents at bus 7, 16, and 30 for the 13th harmonic order for the 30 bus test system.

Figure 7 Actual and reconstructed harmonic currents at bus 7, 16, and 30 for the 17th harmonic order for the 30 bus test system.

Table 1 Correlation coefficients between the actual and the reconstructed Table 2 Mean square error between the actual and reconstructed signals for
signals for the 30 bus test system. the 30 bus test system.

Harmonic order Bus 7 Bus 16 Bus 30 Harmonic order Bus 7 Bus 16 Bus 30
5 0.9811 0.7789 0.9878 5 4.61e-4 1.52e-3 5.69e-4
7 0.9829 0.9316 0.9930 7 9.76e-5 3.90e-4 4.00e-5
11 0.9847 0.8851 0.9887 11 8.76e-5 6.56e-4 6.45e-5
13 0.9846 0.8422 0.9885 13 8.81e-5 9.01e-4 6.54e-5
17 0.9696 0.9389 0.9842 17 1.73e-4 3.48e-4 9.05e-5

rent and the bus harmonic voltage for individual harmonic estimated harmonic currents and voltages and therefore this
order, respectively. The intersection between each column bus is identified as the location of harmonic source. For ex-
and row represents the computed mutual information between ample, in Table 3 the mutual information of 2.387 between the
the three estimated currents at each harmonic frequency, and first estimated current of the 5th harmonic order (first column)
the specified bus harmonic voltage. As shown in Table 3, and the harmonic voltage at bus 7 (7th row) has the greatest
the underlined-bold mutual information values imply that the value in the related column. Therefore, the location of one of
related bus has the greatest mutual information between the the possible harmonic sources (for 5th harmonic order) is at
 IDENTIFICATION OF MULTIPLE HARMONIC SOURCES IN POWER SYSTEMS

Figure 8 IEEE 34 bus radial test system.

Table 3 Mutual information between estimated currents and bus voltages for the IEEE30 bus test system.

Harmonic order
5 7 11 13 17
1 2 3 1 2 3 1 2 3 1 2 3 1 2 3
Bus 1 1.644 0.747 1.026 0.739 0.820 1.575 0.836 1.909 0.903 1.649 0.769 0.920 0.869 1.793 0.831
Bus 2 1.990 0.700 0.946 0.754 0.793 1.732 0.830 2.169 0.911 1.763 0.773 0.850 0.861 1.844 0.854
Bus 3 1.628 0.755 1.036 0.739 0.820 1.575 0.836 1.909 0.903 1.649 0.769 0.920 0.868 1.790 0.824
Bus 4 1.635 0.754 1.031 0.733 0.835 1.581 0.833 1.918 0.905 1.643 0.766 0.918 0.853 1.790 0.822
Bus 5 2.335 0.753 0.856 0.782 0.726 1.730 0.768 2.269 0.766 1.183 0.766 0.779 0.762 1.746 0.785
Bus 6 2.338 0.776 0.864 0.798 0.775 1.918 0.785 2.327 0.861 1.952 0.729 0.812 0.850 1.767 0.847
Bus 7 2.387 0.753 0.856 0.782 0.726 2.330 0.768 2.469 0.766 2.182 0.767 0.778 0.762 1.846 0.785
Bus 8 2.385 0.757 0.872 0.786 0.756 1.788 0.838 2.417 0.841 1.846 0.770 0.809 0.856 1.770 0.860
Bus 9 0.966 0.728 1.353 0.759 1.007 1.066 0.746 0.779 1.438 1.115 0.755 1.213 1.344 0.992 0.810
Bus 10 0.754 0.625 1.581 0.711 1.218 0.865 0.741 0.933 1.692 0.962 0.734 1.543 2.015 0.922 0.792
Bus 11 0.966 0.728 1.353 0.759 1.007 1.066 0.746 0.779 1.438 1.115 0.755 1.213 1.346 0.994 0.807
Bus 12 0.612 0.564 1.603 0.682 1.493 0.787 0.705 0.737 1.671 0.922 0.709 1.779 2.065 0.890 0.738
Bus 13 0.612 0.564 1.603 0.682 1.493 0.787 0.705 0.909 1.767 0.922 0.709 1.779 2.065 0.890 0.738
Bus 14 0.782 0.704 2.413 0.698 1.476 0.794 0.705 0.909 1.767 0.897 0.699 1.835 2.084 0.892 0.742
Bus 15 0.907 0.836 1.499 0.684 1.449 0.800 0.705 0.943 1.708 0.726 0.638 1.684 2.050 0.898 0.734
Bus 16 0.850 0.688 2.541 0.679 1.726 0.738 0.667 0.930 1.784 0.899 0.713 1.856 2.421 0.909 0.850
Bus 17 0.875 0.791 1.360 0.718 1.464 0.821 0.696 0.860 1.721 0.900 0.689 1.738 1.529 0.959 0.882
Bus 18 0.828 0.742 2.392 0.655 1.357 0.835 0.692 0.943 1.720 0.920 0.775 1.791 2.002 0.892 0.752
Bus 19 0.776 0.683 2.205 0.710 1.338 0.846 0.704 0.951 1.727 0.911 0.740 1.786 1.999 0.896 0.770
Bus 20 0.746 0.717 2.163 0.714 1.305 0.822 0.741 0.976 1.751 0.906 0.759 1.710 2.022 0.907 0.774
Bus 21 0.760 0.651 1.799 0.715 1.216 0.886 0.766 0.918 1.684 0.946 0.748 1.610 1.917 0.921 0.757
Bus 22 0.735 0.615 1.820 0.709 1.176 0.828 0.747 0.944 1.655 0.895 0.763 1.600 1.905 0.901 0.758
Bus 23 0.860 0.720 1.867 0.586 1.172 0.638 0.686 0.838 1.298 0.900 0.740 1.344 2.208 0.896 0.680
Bus 24 0.915 0.674 1.856 0.749 1.068 0.875 0.705 0.968 1.436 0.943 0.815 1.502 2.088 0.889 0.734
Bus 25 0.818 1.350 0.929 1.463 0.728 0.728 1.162 1.033 1.028 0.832 1.196 0.966 1.826 0.914 0.854
Bus 26 0.792 1.302 0.860 1.410 0.762 0.713 1.141 1.002 1.003 0.801 1.105 0.944 1.787 0.906 0.858
Bus 27 0.895 1.791 0.876 2.330 0.729 0.700 1.474 0.951 0.911 0.748 1.635 0.778 1.557 0.913 0.892
Bus 28 2.059 0.817 0.984 0.779 0.705 1.441 0.807 2.313 0.890 1.604 0.762 0.819 0.829 1.738 0.850
Bus 29 0.782 2.335 0.811 2.767 0.662 0.773 2.069 0.867 0.873 0.775 2.538 0.702 0.774 0.927 1.967
Bus 30 0.782 2.406 0.841 2.945 0.670 0.708 2.311 0.829 0.820 0.777 2.553 0.732 0.749 0.796 2.739
M. FARHOODNEA ET AL

Figure 9 Actual and reconstructed harmonic currents at bus 15, 25, and 33 for the 5th harmonic order for the 34 bus test system.

Figure 10 Actual and reconstructed harmonic currents at bus 15, 25, and 33 for the 7th harmonic order for the 34 bus test system.

bus 7. From the results shown in Table 3, it is evident that Similar to the 30 bus test system, harmonic power flows are
harmonic source locations are located at bus 7, 16, and 30 simulated for different loading conditions by means of multi-
and this proves the accuracy of the proposed harmonic source plying the Laplace distributed random variables at all loads to
location method. create different operating conditions for each load and at each
sampling time. Again, harmonic bus voltages are computed
by solving the linear system for each harmonic order. In order
5.2 Results of the Radial 34 bus Test System to reconstruct injected harmonic currents caused by the har-
To prove the accuracy and ability of proposed method for monic sources, 500 harmonic voltage samples have been gen-
locating multiple harmonic sources in radial distribution net- erated at each bus so as to represent harmonic measurement
works, the radial IEEE 34 bus test system shown in Figure vector, X. by using preprocessing procedures, the dimension
8 is used [27]. In the test system, three harmonic sources of the harmonic measurement vectors X can be decreased so
containing the 5th , 7th , 11th , 13th , and 17th harmonic orders as to omit the redundant information of the data. Figures 9
are located at bus 15, 25, and 33. In order to use the ICA to 13 show the extracted current traces by using the Fast ICA
algorithm for reconstructing the injected harmonic currents, algorithm in comparison with the actual signals at bus 15, 25,
harmonic voltage measurements are generated using the har- and 33 for the 5th , 7th , 11th , 13th and 17th harmonic orders,
monic power flow algorithm for radial systems [28-29] which respectively, in which the red lines represent the actual har-
has been implemented using the Matlab codes. monic currents and the blue lines represent the reconstructed
 IDENTIFICATION OF MULTIPLE HARMONIC SOURCES IN POWER SYSTEMS

Figure 11 Actual and reconstructed harmonic currents at bus 15, 25, and 33 for the 11th harmonic order for the 34 bus test system.

Figure 12 Actual and reconstructed harmonic currents at bus 15, 25, and 33 for the 13th harmonic order for the 34 bus test system.

harmonic currents over 500 samples and all signals are nor-
Table 4 Correlation coefficients between the actual and the reconstructed
malized with their largest values in per-unit (p.u). signals for the 34 bus test system.
To prove high accuracy between the estimated and the actual
harmonic currents, the obtained correlation coefficients and Harmonic order Bus 15 Bus 25 Bus 33
mean square error of the signals are shown in Tables 4 and 5. 5 0.9892 0.9906 0.9945
The results of Tables 4 and 5 imply that compatibility degree 7 0.9931 0.9909 0.9978
between the reconstructed and original signals is greater than 11 0.9892 0.9907 0.9944
98% and mean square error of estimation is less than 4.33e- 13 0.9936 0.9972 0.9951
5. It is shown that the 34 bus test system has correlation 17 0.9930 0.9909 0.9944
coefficient values closer to 1, and lesser mean square error Table 5 Mean square error between the actual and reconstructed signals for
the 34 bus test system.
value in comparison with the results of the 30 bus test system.
Hence, the accuracy of the results is due to the difference in Harmonic order Bus 15 Bus 25 Bus 33
the sample size used in the two systems. 5 4.33e-5 3.74e-5 2.21e-5
To identify the location of harmonic sources, the pair-wise 7 2.24e-5 3.81e-5 8.90e-6
mutual information between the extracted current trends and 11 4.32e-5 3.73e-5 2.22e-5
bus voltage have been obtained for each harmonic frequency 13 2.57e-5 1.11e-5 1.94e-5
as shown in Table 6. In the Table, the underlined-bold mutual 17 2.79e-5 3.62e-5 2.23e-5
information values imply that the related bus has the greatest
mutual information between the estimated harmonic currents
M. FARHOODNEA ET AL

Table 6 Mutual information between estimated currents and bus voltages for the 34 bus test system.

Harmonic order
5 7 11 13 17
1 2 3 1 2 3 1 2 3 1 2 3 1 2 3
Bus 1 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
Bus 2 0.814 0.911 0.758 0.747 1.046 0.703 0.712 1.173 0.711 0.670 1.407 0.699 0.608 0.707 1.321
Bus 3 0.814 0.911 0.758 0.747 1.046 0.703 0.712 1.173 0.711 0.670 1.407 0.699 0.608 0.707 1.321
Bus 4 0.954 0.770 0.866 0.904 0.900 0.718 0.822 0.916 0.722 0.700 1.037 0.781 0.719 0.779 0.990
Bus 5 1.066 0.752 0.905 1.013 0.814 0.820 0.931 0.842 0.824 0.680 0.881 0.876 0.732 0.814 0.902
Bus 6 1.108 0.753 0.931 1.057 0.747 0.880 0.967 0.789 0.887 0.782 0.860 0.973 0.788 0.895 0.838
Bus 7 0.958 0.759 1.004 1.160 0.740 0.807 1.068 0.752 0.843 0.770 0.819 1.109 0.904 0.860 0.792
Bus 8 0.956 0.679 1.042 1.249 0.717 0.744 1.130 0.775 0.835 0.820 0.831 1.160 0.977 0.869 0.804
Bus 9 0.866 0.687 1.135 1.377 0.711 0.780 1.257 0.737 0.745 0.769 0.753 1.273 1.055 0.828 0.774
Bus 10 0.829 0.649 1.184 1.460 0.674 0.696 1.306 0.703 0.740 0.707 0.773 1.333 1.084 0.782 0.754
Bus 11 0.829 0.649 1.184 1.460 0.674 0.696 1.306 0.703 0.740 0.707 0.773 1.333 1.084 0.782 0.754
Bus 12 0.829 0.649 1.184 1.460 0.674 0.696 1.306 0.703 0.740 0.707 0.773 1.333 1.084 0.782 0.754
Bus 13 0.815 1.075 0.723 0.685 1.321 0.653 0.623 1.462 0.671 0.632 1.681 0.657 0.603 0.672 1.585
Bus 14 0.709 1.328 0.667 0.667 1.695 0.647 0.630 1.808 0.664 0.631 1.998 0.629 0.585 0.678 1.938
Bus 15 0.701 1.483 0.625 0.652 1.846 0.659 0.635 1.988 0.671 0.629 2.144 0.627 0.585 0.664 2.094
Bus 16 0.701 1.443 0.625 0.652 1.546 0.659 0.635 1.888 0.671 0.629 2.005 0.627 0.585 0.664 1.973
Bus 17 1.201 0.739 0.875 0.977 0.744 0.928 0.919 0.834 0.974 0.840 0.813 0.952 0.831 0.972 0.794
Bus 18 1.284 0.744 0.855 0.991 0.788 0.998 0.942 0.803 0.991 0.878 0.828 0.916 0.796 1.037 0.764
Bus 19 1.338 0.718 0.815 0.943 0.767 1.043 0.896 0.767 1.084 0.917 0.817 0.890 0.815 1.124 0.802
Bus 20 1.414 0.733 0.765 0.917 0.745 1.094 0.867 0.750 1.095 0.964 0.776 0.901 0.802 1.194 0.792
Bus 21 1.476 0.714 0.743 0.870 0.734 1.152 0.856 0.760 1.129 1.025 0.784 0.884 0.785 1.272 0.757
Bus 22 1.564 0.714 0.721 0.852 0.712 1.191 0.810 0.741 1.195 1.086 0.786 0.868 0.781 1.311 0.727
Bus 23 1.607 0.669 0.730 0.850 0.753 1.233 0.805 0.702 1.232 1.091 0.753 0.845 0.750 1.344 0.717
Bus 24 1.705 0.664 0.717 0.829 0.742 1.297 0.805 0.727 1.313 1.125 0.734 0.874 0.697 1.432 0.725
Bus 25 1.741 0.680 0.729 0.793 0.706 1.322 0.780 0.721 1.324 1.196 0.767 0.839 0.710 1.480 0.731
Bus 26 1.641 0.680 0.729 0.793 0.706 1.122 0.780 0.721 1.242 1.155 0.767 0.839 0.710 1.365 0.731
Bus 27 1.541 0.680 0.729 0.793 0.706 1.022 0.780 0.721 1.127 1.124 0.767 0.839 0.710 1.312 0.731
Bus 28 0.958 0.759 1.004 1.160 0.740 0.807 1.068 0.752 0.843 0.770 0.819 1.109 0.904 0.860 0.792
Bus 29 0.958 0.759 1.004 1.160 0.740 0.807 1.068 0.752 0.843 0.770 0.819 1.109 0.904 0.860 0.792
Bus 30 0.958 0.759 1.004 1.160 0.740 0.807 1.068 0.752 0.843 0.770 0.819 1.109 0.904 0.860 0.792
Bus 31 0.818 0.617 1.257 1.524 0.689 0.711 1.357 0.733 0.755 0.671 0.739 1.364 1.128 0.782 0.768
Bus 32 0.827 0.663 1.321 1.596 0.672 0.696 1.385 0.722 0.727 0.698 0.718 1.421 1.161 0.779 0.792
Bus 33 0.795 0.695 1.335 1.632 0.665 0.677 1.443 0.723 0.704 0.724 0.723 1.500 1.190 0.747 0.741
Bus 34 0.795 0.695 1.235 1.432 0.665 0.677 1.343 0.723 0.704 0.724 0.723 1.325 0.968 0.747 0.741

Table 7 Mean square errors of the mutual information values for IEEE 30 and 34 bus test system.

Harmonic orders 5th 7th 11th 13th 17th

IEEE 30 bus test system 1.29e-4 2.24e-4 2.69e-4 3.45e-4 1.26e-4
IEEE 34 bus test system 6.96e-5 8.58e-5 9.94e-5 1.21e-4 9.81e-5
 IDENTIFICATION OF MULTIPLE HARMONIC SOURCES IN POWER SYSTEMS

Figure 13 Actual and reconstructed harmonic currents at bus 15, 25, and 33 for the 17th harmonic order for the 34 bus test system.

and voltages and this bus should identify as one of the pos- REFERENCES
sible location of harmonic source. For example, the mutual
information of 1.741 between the first estimated current of the 1. K. Islam and A. Samra, "Identification of harmonic
5th harmonic order (first column) and the harmonic voltage at sources in power distribution systems," IEEE Southeast-
bus 25 (25th row) has the greatest value in the related column, [Link] new New Century., Proceedings, pp. 301
which implies that bus 25 is the location of one of the possible 303, 1997.
2. W. Xu, X. Liu, and Y. Liu, "An investigation on the validity
harmonic sources in this test system.
of power direction method for harmonic source determination,"
From the results shown in Table 6, it is evident that harmonic IEEE Power Engineering Review, vol. 22, pp. 6262, 2002.
source locations are located at bus 15, 25, and 33 for the 34 3. X. Wilsun and L. Yilu, "A method to determine customer har-
bus test system, because these buses have the highest mutual monic contributions for incentive-based harmonic control appli-
information values. Thus, the harmonic source locations have cations," in Power Engineering Society Summer Meeting, 1999.
been accurately estimated by the proposed method using ICA IEEE, 1999, pp. 361366 vol.1.
and mutual information theory. 4. X. Wilsun and [Link], "A method for determining customer and
Table 7 shows the mean square errors of the mutual infor- utility harmonic contributions at the point of common coupling,"
mation values for both 30 and 34 bus test systems. The results IEEE Transactions on Power Delivery, vol. 15, pp. 804811,
shows that the mutual information algorithm can estimate ex- 2000.
act locations of the harmonic sources with the mean square 5. C. Chaoying, L. Xiuling, D. Koval, X. Wilsun, and T. Tay-
jasanant, "Critical impedance method - a new detecting har-
errors less than 0.00035 and 0.00012 for the 30 and 34 bus
monic sources method in distribution systems," IEEE Transac-
test systems, respectively. tions on Power Delivery, vol. 19, pp. 288297, 2004.
6. N. Hamzah, A. Mohamed, and A. Hussain, "Harmonic source
location at the point of common coupling based on voltage mag-
6. CONCLUSION nitude," in TENCON 2004. 2004 IEEE Region 10 Conference,
2004, pp. 220223 Vol. 3.
In this paper, a novel technique has been presented for locating 7. S. F. de Paula Silva and J. C. de Oliveira, "The sharing of re-
multiple harmonic sources in radial and non radial distribution sponsibility between the supplier and the consumer for harmonic
systems, without prior information about the network compo- voltage distortion: A case study," Electric Power Systems Re-
nents and branch parameters. In the proposed method, ICA search, vol. 78, pp. 19591964, 2008.
which is one of the blind source separation techniques is ap- 8. T. Pfajfar, B. Blazic, and I. Papic, "Harmonic Contributions
plied for estimating the injected harmonic currents produced Evaluation With the Harmonic Current Vector Method," IEEE
by harmonic sources and the mutual information theory is ap- Transactions on Power Delivery, vol. 23, pp. 425433, 2008.
plied to identify the location of possible harmonic sources. 9. A. Cataliotti and V. Cosentino, "A New Measurement Method
ICA is also used to estimate the load profiles of the harmonic for the Detection of Harmonic Sources in Power Systems Based
on the Approach of the IEEE Std. 1459–2000," IEEE
sources with different sample sizes. The results from the case
Transactions on Power Delivery, vol. 25, pp. 332340, 2010.
studies confirm that the proposed method can accurately iden- 10. Z. P. Du, J. Arrillaga, N. R. Watson, and S. Chen, "Identification
tify the location of multiple harmonic sources in both radial of harmonic sources of power systems using state estimation,"
and non radial distribution systems. The proposed method is Generation, Transmission and Distribution, IEE Proceedings,
promising because it uses only bus voltage measurements and vol. 146, pp. 712, 1999.
there is no need for real and reactive power measurements and 11. G. DAntona, C. Muscas, and S. Sulis, "State Estimation for
other restricted power system information. the Localization of Harmonic Sources in Electric Distribution
M. FARHOODNEA ET AL

Systems," IEEE Transactions on Instrumentation and Measure- 21. E. Duggan and R. E. Morrison, "Probabilistic methods applied
ment, vol. 58, pp. 14621470, 2009. to harmonic distortion," in Probabilistic Methods Applied to
12. C. Madtharad, S. Premrudeepreechacharn, N. Watson, and R. Electric Power Systems, 1991., Third International Conference
Saeng-Udom, "An optimal measurement placement method for on, 1991, pp. 4346.
power system harmonic state estimation," IEEE Transactions 22. L. Huaiwei and D. Niebur, "Exploring independent component
on Power Delivery, vol. 20, pp. 15141521, 2005. analysis for electric load profiling," in Neural Networks, 2002.
13. C. Rakpenthai, S. Premrudeepreechacharn, S. Uatrongjit, and IJCNN 02. Proceedings of the 2002 International Joint Con-
N. Watson, "Measurement placement for power system state es- ference on, 2002, pp. 21442149.
timation using decomposition technique," Electric Power Sys- 23. L. Huaiwei and D. Niebur, "Load profile estimation in electric
tems Research, vol. 75, pp. 4149, 2005. transmission networks using independent component analysis,"
14. A. Kumar, B. Das, and J. Sharma, "Determination of location IEEE Transactions on Power Systems, vol. 18, pp. 707715,
of multiple harmonic sources in a power system," International 2003.
Journal of Electrical Power & Energy Systems, vol. 26, pp. 73 24. R. Abu-Hashim, et al., "Test systems for harmonics modeling
78, 2004. and simulation," IEEE Transactions on Power Delivery, vol. 14,
15. E. Gursoy, "Independent component analysis for harmonic pp. 579587, 1999.
source identification in electric power systems," Drexel Uni- 25. H. Saadat, "Power system analysis, WCB," ed: McGraw-Hill,
versity, 2007. Boston, 1999.
16. E. Oja and A. Hyvarinen, "Independent component analysis: al- 26. M. Grady, "PCFLO User Manual," ECE Dept, The University
gorithms and applications," Neural Networks, vol. 13, pp. 411 of Texas at Austin, . [Link]
430, 2000. 27. S. F. Mekhamer, S. A. Soliman, M. A. Moustafa, and M. E. El-
17. A. Hyvar nen, J. Karhunen, and E. Oja, "Independent Compo- Hawary, "Application of fuzzy logic for reactive-power com-
nent Analysis," 4th ed: John Wiley, 2001. pensation of radial distribution feeders," IEEE Transactions on
18. A. Hyvarinen, "Fast and robust fixed-point algorithms for in- Power Systems, vol. 18, pp. 206213, 2003.
dependent component analysis," IEEE Transactions on Neural 28. T. Jen-Hao and C. Chuo-Yean, "Backward/Forward Sweep-
Networks, vol. 10, pp. 626634, 1999. Based Harmonic Analysis Method for Distribution Systems,"
19. A. Papoulis and S. Pillai, Probability, Random Variables, and IEEE Transactions on Power Delivery, vol. 22, pp. 16651672,
Stochastic Processes: McGraw-Hill, 2002. 2007.
20. W. Li, "Mutual information functions versus correlation func- 29. J. Teng, "A network-topology-based three-phase load flow for
tions," Journal of Statistical Physics, vol. 60, pp. 823837, distribution systems," Proc Natl. Sci. Counc. ROC (A), vol. 24,
1990. pp. 259264, 2000.
View publication stats