# Simultaneous Load Identification Method Based on Hybrid Features and Genetic Algorithm for Nonintrusive Load Monitoring

Authors: Shu-Hui Yi; Jian Wang; Jun-Jie Liu
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

Nonintrusive load monitoring (NILM) is a widely accepted technology to conduct load monitoring. Many effective methods have been established to make NILM more practical. However, the focus of current methods is mainly on the identification accuracy and efficiency of single load under the individual appliance operated independently, which have limited support for the identification problem under multiple appliances operated simultaneously. Therefore, a simultaneous identification method is proposed to efficiently identify the total load under multiple appliances operated simultaneously in this paper. The proposed identification method mainly consists of three parts: hybrid features extraction, simultaneous identification optimization model construction, and frequency-weighting-factor-based genetic algorithm (FWF-GA). Firstly, the hybrid feature model, which integrates the features of active power, reactive power, and harmonic magnitude, is constructed by hybrid features extraction. Secondly, the simultaneous identification optimization model is constructed by employing the features of active and reactive power. Thirdly, the developed FWF-GA is used to solve the simultaneous identification optimization problem. In FWF-GA, the relative errors of active power, reactive power, and the frequency-weighting factor of harmonic magnitude are used to evaluate the fitness of an individual. Finally, a NILM practice to identify household appliances is used to demonstrate the validity of the proposed method.


---

## Body

## 1. Introduction

As a high-grade energy carrier, electricity plays a critical role in economic development [1]. With the development of the economy and society, the demand for electrical power will grow continually. Since electricity production relies on other energies as well as its feature is unsuited to mass storage, energy-saving becomes an urgent need for modern power systems [2]. One best way to optimize energy savings is to conduct load monitoring, which has been widely used in smart grids. Many investigations reveal that both electricity consumers and providers can benefit from load monitoring [3–5]. Customers can smartly use their appliances with the detailed information of load monitoring, such as making a more rational strategy to reduce energy consumption and identifying the failure appliances and inefficient appliances. For electricity providers, some notable benefits derived from load monitoring include making public policy on electricity consumption, segmenting customers, predicting energy demand, and even detecting illegal electricity behaviors [6].Traditional load monitoring is conducted by analyzing the collected data that are obtained from the sensors installed on each appliance. This monitoring manner is also called intrusive load monitoring (ILM). The technique advantage of ILM is that it can directly obtain the relatively accurate and reliable monitoring data of the individual appliance, while the key disadvantage lies in the hardware cost, the complicated maintenance, as well as the large number of installation work. In contrast, a cost-effective technique named nonintrusive load monitoring (NILM) has been developed to make load monitoring more practical. NILM only needs a single minimal set of plug-panel sensors to acquire the aggregate electric consumption data, and then, the acquired data is treated as a disaggregation problem. NILM, therefore, also refers to load disaggregation [7]. Benefiting from the advantage of technology, NILM has gained tremendous attention and is becoming more and more popular.NILM aims to identify the operation state of the individual appliance from the information of aggregate load. There are typically four steps summarized in general NILM procedures, including data acquiring, event detecting, feature extracting, and load identifying. Since feature extracting and load identifying dominate more to the efficiency and accuracy, the current researches around NILM mainly focus on these two aspects. The task of feature extraction is to obtain relevant features so as to discriminate different types of appliance loads. It is the inherent requirement to extract features accurately and effectively. Researchers therefore pay especially attention to the feature representation. Many feature representations are developed based on current, voltage, and power. Among them, real power and reactive power are most widely used for their stable performance [8]. Besides, harmonic-based features were proposed to conduct feature extraction as a supplement to real and reactive power features [9–11]. To better detect the transient events, time-frequency tools were introduced to model the nonstationary features. Lin and Tsai [12] proposed a transient feature extraction, where the transient response was modeled based on S-transform. Shaw et al. [13] adopted the short-time Fourier transforms to compute the spectral envelope coefficients during the NILM implementation. Gillis and Morsi [14] applied discrete wavelet transform to represent transient features in event-efficient identification. Another new feature representation method is shape representation, which can directly obtain efficient feature information from the raw current and voltage signals [15–17]. The shape representation was also extended to transient features modeling, such as Stockwell transform-based time-frequency feature [18]. Moreover, He et al. [19] adopted the probability density function to model multiple operating states, where both the power measurements and current signals were considered. As different feature representation indicates the different characteristics of appliances, the synthetic feature is increasingly used in NILM for accurate and efficient identification [20, 21].The second research focus of the current NILM is load identification. It is a significant goal of NILM that the operating state of individual electrical appliances can be accurately discriminated from the aggregate consumption signal. To achieve this, load identification as one crucial step in NILM is used to distinguish the different types of load. Large quantities of methods were developed in the load identification area, which can be roughly grouped into two categories. The first category is called machine learning-based methods. This type of methods can be implemented without event detection and usually has high efficiency. During its implementation, the discrimination of different type loads is generally treated as a multilabel classification problem, and therefore, the machine learning methods, including unsupervised, semisupervised, and supervised methods, all can be employed to deal with the identification problems. Tabatabaei et al. [22] developed two supervised multilabel classifiers based on a k-nearest neighbor. Lai et al. [23] proposed a hybrid model in which the support vector machine is integrated with the Gaussian mixture model. The hidden Markov model and its improved models were used as the unsupervised classification to implement the load identification, providing a well performance of data clustering [24–26]. Yang et al. [27] proposed a semisupervised multilabel classification method based on deep learning to further improve the efficiency of the deep learning-based identification methods. The second category is optimization-based methods. In the optimization-based method, the load identification is treated as a mathematical optimization problem, where the total consumption signal is considered as the combination of the individual appliance to be identified. For a long time, the optimization-based method was adopted as another identification method, which is always used to identify one load at a time. Therefore, efficiency and accuracy were especially concerned. Some well-known methods, such as mixed-integer linear programming [28] and quadratic programming (QP) [29, 30], were already successfully applied in optimization-based load identification.Recently, the attention to load identification is also drawn to simultaneous identification problems for multiple electrical loads exist. Hua et al. [31] proposed an event-based NILM method, in which the mixed linear integer programming model was developed to achieve the simultaneous disaggregation of multiple appliances. Liu et al. [32] proposed a deep dictionary learning model and an adaptive-window-based detection approach to conducting the NILM with simultaneous switching operations. Sun et al. [33] proposed the dynamic adaptive particle swarm optimization algorithm (DAPSO) to realize nonintrusive household load identification. Fang et al. [34] proposed a NILM method that leverages advances in statistical learning. Yang et al. [35] proposed an intelligent event-driven nonintrusive load monitoring method based on a convolutional neural network, which extends the NILM application range to the handling process in terms of various and combined characteristics. Furthermore, Zhou et al. [36] proposed a NILM method based on the CNN-LSTM hybrid model to improve the performance of the whole network system. Since different users have different habits, various cases of power consumption are possible. The simultaneous identification therefore is a direction worthy of development and research [37].The above feature representation and load identification methods can efficiently improve the quality of NILM. There is a clear trend that more and more feature representations use the synthetic feature deriving from the different aspects of appliance load. This is because the synthetic feature can better distinguish different appliances compared with a single feature. Furthermore, the good feature representation is always combined with an excellent load identification method to strengthen its practicability. Although many load identification methods have already been developed and well used in practical NILM, their focus is mainly on the single load identification at a time. The support for the identification problem of the total load under multiple appliances operated simultaneously is still very limited, and both the identification model and the efficient solution algorithm are needed to be further developed. Therefore, to achieve the identification of the total load under multiple appliances operated simultaneously, a simultaneous load identification method based on hybrid features and a genetic algorithm is proposed in this paper, in which the frequency weighting factor-based genetic algorithm (FWF-GA) is developed to efficiently solve the simultaneous identification problem.The main contributions are highlighted as follows:(1)
This paper proposes for the first time a fitness function based on the frequency weighting factor that is used to evaluate the fitness of individuals powerfully in FWF-GA. This proposed fitness function not only takes the influence of the harmonic feature into account but also considers the relative error of active power and reactive power, and more importantly, during the evaluation of the fitness function, the influence of negligence due to the large gap between active power and reactive power will be avoided.(2)
The test results show that the proposed method has better accuracy and efficiency compared with the method of QP, DAPSO, and CNN. Meanwhile, the proposed method does not need a large amount of data for training like the method based on CNN.The rest of the paper is organized as follows. Section2 introduces related preliminaries about NILM. Section 3 proposes the simultaneous identification method, including hybrid features, a simultaneous identification optimization model, and the procedure of FWF-GA. Section 4 mainly uses the load data of five household appliances to test the effectiveness of the proposed method. Finally, Section 5 concludes this work and suggests future work.


## 2. Preliminaries

### 2.1. Discrete Fourier Transform

The signal of the time domain can be transformed into the signal of the frequency domain by using discrete Fourier transform (DFT) [38]. The DFT consists in decomposing a signal into a sum of elementary signals, which have the property of being easy to implement and observe [39]. On the basis of DFT, the time signal A=ak=a0,a1,…,aK−1 is transformed into the frequency sequence B=bm=b0,b1,…,bS−1. The DFT is defined by(1)B=bm=DFTak,bm=∑k=0K−1ake−j2π/Kmk,0≤m≤K−1,where ak is the kth value of the signal A, bm is the mth DFT value of the signal A, K is the size of A, and S is the size of B.For each value ofm in equation (2), the number of calculations to obtain each value is K, and therefore, a total calculative number K2 is required to get the frequency sequence by DFT. When the number K becomes large, the computational cost of common DFT will be extremely large. Consequently, the fast Fourier transform (FFT) is proposed to improve computational efficiency.According to the symmetry ofe−j2π/Kmk, the time signal A can be decomposed into the sum of two sequences, which is defined by(2)A=A1+A2,where A1 and A2 are the time sequences of size K/2 and A1 is the even sequence, while A2 is the odd sequence. Let WK=e−j2π/K, then equation (2) can be transformed into the following form:(3)B=bm=B1+WKmB2,where B1=bm1 and B2=bm2 are, respectively, the frequency sequence after the DFT of the time signal A1 and A2. bm can be obtained by(4)bm=∑k=0K/2−1a1kWK2km+∑k=0K/2−1a2kWK2k+1m=∑k=0K/2−1a1kWK2km+WKm∑k=0K/2−1a2kWK2km=bm1+WKmbm2,0≤m1,m2≤K2−1,0≤m≤K−1.


### 2.2. Load Feature

The function and internal structure of electrical equipment are different, and the load features are also different. Therefore, different load features can describe the operation of different electrical equipment. In the current study on NILM, the common load features are the current waveform features, current harmonic features, power features, and so on.


#### 2.2.1. Current Waveform Features

From the current waveform, three features can be extracted, which are given by(5)IS=1K∑k=0Kak2,IP=maxak,IF=IPIS,where IS, IP, and IF are the root mean square value, peak magnitude, and crest factor of the current waveform, respectively.


#### 2.2.2. Harmonic Features

The harmonic features can reflect the structure of electrical equipment. The harmonic features can be calculated by DFT, which can be given by the following equations:(6)M=bm=Rebm2+Imbm2,P=argbm=arctanImbmRebm,kTHD=∑i=2LGiG12,where M, P, and kTHD are the magnitude, phase, and total harmonic distortion (THD) of the frequency signal, respectively; Im and Re represent the imaginary and the real part of the DFT, respectively; Gi is the harmonic component of the ith order; G1 is the fundamental component; and L is the maximum of effective harmonic order.


#### 2.2.3. Power Features

The power features are the most common indicators to describe the electrical behavior of electrical equipment. The power features contain active power and reactive power, which can be calculated as follows:(7)P=∑k=0∞Pk=∑k=0∞UkIkcosφk,Q=∑k=0∞Qk=∑k=0∞UkIksinφk,where P is the active power, Q is the reactive power, U is the voltage, I is the current, k is the harmonic order, and φ is the phase difference between voltage and current.


## 2.1. Discrete Fourier Transform

The signal of the time domain can be transformed into the signal of the frequency domain by using discrete Fourier transform (DFT) [38]. The DFT consists in decomposing a signal into a sum of elementary signals, which have the property of being easy to implement and observe [39]. On the basis of DFT, the time signal A=ak=a0,a1,…,aK−1 is transformed into the frequency sequence B=bm=b0,b1,…,bS−1. The DFT is defined by(1)B=bm=DFTak,bm=∑k=0K−1ake−j2π/Kmk,0≤m≤K−1,where ak is the kth value of the signal A, bm is the mth DFT value of the signal A, K is the size of A, and S is the size of B.For each value ofm in equation (2), the number of calculations to obtain each value is K, and therefore, a total calculative number K2 is required to get the frequency sequence by DFT. When the number K becomes large, the computational cost of common DFT will be extremely large. Consequently, the fast Fourier transform (FFT) is proposed to improve computational efficiency.According to the symmetry ofe−j2π/Kmk, the time signal A can be decomposed into the sum of two sequences, which is defined by(2)A=A1+A2,where A1 and A2 are the time sequences of size K/2 and A1 is the even sequence, while A2 is the odd sequence. Let WK=e−j2π/K, then equation (2) can be transformed into the following form:(3)B=bm=B1+WKmB2,where B1=bm1 and B2=bm2 are, respectively, the frequency sequence after the DFT of the time signal A1 and A2. bm can be obtained by(4)bm=∑k=0K/2−1a1kWK2km+∑k=0K/2−1a2kWK2k+1m=∑k=0K/2−1a1kWK2km+WKm∑k=0K/2−1a2kWK2km=bm1+WKmbm2,0≤m1,m2≤K2−1,0≤m≤K−1.


## 2.2. Load Feature

The function and internal structure of electrical equipment are different, and the load features are also different. Therefore, different load features can describe the operation of different electrical equipment. In the current study on NILM, the common load features are the current waveform features, current harmonic features, power features, and so on.


### 2.2.1. Current Waveform Features

From the current waveform, three features can be extracted, which are given by(5)IS=1K∑k=0Kak2,IP=maxak,IF=IPIS,where IS, IP, and IF are the root mean square value, peak magnitude, and crest factor of the current waveform, respectively.


### 2.2.2. Harmonic Features

The harmonic features can reflect the structure of electrical equipment. The harmonic features can be calculated by DFT, which can be given by the following equations:(6)M=bm=Rebm2+Imbm2,P=argbm=arctanImbmRebm,kTHD=∑i=2LGiG12,where M, P, and kTHD are the magnitude, phase, and total harmonic distortion (THD) of the frequency signal, respectively; Im and Re represent the imaginary and the real part of the DFT, respectively; Gi is the harmonic component of the ith order; G1 is the fundamental component; and L is the maximum of effective harmonic order.


### 2.2.3. Power Features

The power features are the most common indicators to describe the electrical behavior of electrical equipment. The power features contain active power and reactive power, which can be calculated as follows:(7)P=∑k=0∞Pk=∑k=0∞UkIkcosφk,Q=∑k=0∞Qk=∑k=0∞UkIksinφk,where P is the active power, Q is the reactive power, U is the voltage, I is the current, k is the harmonic order, and φ is the phase difference between voltage and current.


## 2.2.1. Current Waveform Features

From the current waveform, three features can be extracted, which are given by(5)IS=1K∑k=0Kak2,IP=maxak,IF=IPIS,where IS, IP, and IF are the root mean square value, peak magnitude, and crest factor of the current waveform, respectively.


## 2.2.2. Harmonic Features

The harmonic features can reflect the structure of electrical equipment. The harmonic features can be calculated by DFT, which can be given by the following equations:(6)M=bm=Rebm2+Imbm2,P=argbm=arctanImbmRebm,kTHD=∑i=2LGiG12,where M, P, and kTHD are the magnitude, phase, and total harmonic distortion (THD) of the frequency signal, respectively; Im and Re represent the imaginary and the real part of the DFT, respectively; Gi is the harmonic component of the ith order; G1 is the fundamental component; and L is the maximum of effective harmonic order.


## 2.2.3. Power Features

The power features are the most common indicators to describe the electrical behavior of electrical equipment. The power features contain active power and reactive power, which can be calculated as follows:(7)P=∑k=0∞Pk=∑k=0∞UkIkcosφk,Q=∑k=0∞Qk=∑k=0∞UkIksinφk,where P is the active power, Q is the reactive power, U is the voltage, I is the current, k is the harmonic order, and φ is the phase difference between voltage and current.


## 3. Simultaneous Identification Method

The proposed simultaneous identification method is mainly used to deal with the problem involving simultaneously identifying multiple loads in NILM. The proposed method is developed on the basis of hybrid features and genetic algorithm, which mainly consists of three parts: hybrid features extraction, simultaneous identification optimization model construction, and optimization model solution by using frequency weighting factor-based genetic algorithm (FWF-GA). The frame of the proposed simultaneous identification method is shown in Figure1.Figure 1
Frame of the proposed simultaneous identification method.


### 3.1. Hybrid Features

In this paper, a hybrid feature model is proposed for NILM prediction. The proposed model is a set of features that is composed of active power, reactive power, and harmonic magnitude. The hybrid features of the appliances are expressed in the following matrix form:(8)H=P1Q1M11…M1j…M1LP2Q2M21…M2j…M2L⋮⋮⋮⋮⋮⋮⋮PiQiMi1⋯Mij…MiL⋮⋮⋮⋮⋮⋮⋮PNQNMN1…MNj…MNL,where Pi and Qi are the active and reactive power of ith appliance, respectively. Mij is the jth harmonic magnitude of ith appliance, and N is the total number of appliances. For clarity of presentation, the matrix H also can be regarded as the combination of matrix H=HP,HQ,HM. HP, HQ, and HM are given as follows:(9)HP=P1P2⋯Pi⋯PNT,HQ=Q1Q2⋯Qi⋯QNT,HM=M11…M1j…M1LM21…M2j…M2L⋮⋮⋮⋮⋮Mi1…Mij…MiL⋮⋮⋮⋮⋮MN1…MNj…MNL,where the superscript T denotes the transposition of vectors or matrices. It is noted that when performing the simultaneous identification, both the hybrid feature model of the individual appliance operated independently and multiple appliances operated simultaneously need to be given out.The single load under the individual appliance operated independently and the total load under multiple appliances operated simultaneously are obtained by real-time monitoring. Through the feature extracting and DFT, the hybrid feature models of single load and total load can be given by(10)HS=HS−P,HS−Q,HS−MHT=HT−P,HT−Q,HT−M,where the subscript S denotes the single load that is a load of the individual electrical appliance operated independently and the subscript T denotes the total load that is a load of multiple appliances operated simultaneously. It is noted that HT−P and HT−Q are both constant values, which are defined as HT−P=P and HT−Q=Q, respectively, while HT−M is one-dimensional vector, that is, HT−M=MT1,MT2,…,MTi,…,MTL, where MTi denotes the ith harmonic magnitude of the total load.


### 3.2. Simultaneous Identification Optimization Model

The hybrid feature modelH is divided into three components due to different uses. The vector HP and vector HQ are used to construct the optimization model of simultaneous identification, while the matrix HM is introduced to the algorithm developed for solving the optimization model. The object of load identification is to find the best combination of the appliance load to approach the total power consumption, and thus, the most possible appliances switched on can be inferred. The objective function therefore is given as follows:(11)minfx=HT−P−xHS−P+HT−Q−xHS−Q,where x is the vector of design variables corresponding to each appliance, that is, x=x1,x2,…,xi,…,xN, and xi is the design variable corresponding to the ith appliance.Considering the limit case of simultaneous appliances switched on, the simultaneous identification problem can be established by an optimization problem, which is given by(12)minfx=HT−P−xHS−P+HT−Q−xHS−Qs.t.∑i=1Nxi≤Nxi=0or1.The value ofxi is 0 or 1, and 0 denotes that the ith appliance is not running, while 1 denotes that the ith appliance is running. Therefore, the optimization problem is a 0-1 integer optimization problem. The genetic algorithm [40, 41] is useful for the integer optimization problem. Therefore, in order to improve the identification accuracy and efficiency, a genetic algorithm based on frequency weighting factors (FWF-GA) is developed to solve the simultaneous identification problem of equation (12).


### 3.3. FWF-GA

The procedure of simultaneous identification with FWF-GA is shown in Figure2.Figure 2
Procedures of simultaneous identification with FWF-GA.The procedure can be summarized in the following steps:Step 1: Calculating frequency weighting factorsSince the frequency domain can directly reflect the inherent characteristics of appliance load from a frequency perspective, the frequency weighting factor is defined on the basis of frequency domain information. The frequency weighting factor is defined as follows:(13)ξi=1LMi1MT1+Mi2MT2⋅⋅⋅+MijMTj+⋅⋅⋅+MiLMTL,whereξi is the frequency weighting factor of ith appliance.As shown by equation (13), the frequency weighting factor is actually the average value of L terms weight sum. To further clearly illustrate the defined frequency weighting factor, the preceding six terms weight are taken as an example, which is shown in Figure 3. It is worth noting that the value of L is not absolutely fixed, that is to say, the value of L still depends on an individual decision considering which ones play an important role.Step 2: Constructing fitness functionThe fitness evaluation is one crucial step in the GA process, whose fitness function directly determines the accuracy and efficiency of the algorithm. To further enhance the practicability of simultaneous identification problem solving, a new fitness function is developed based on the frequency-weighting factor, which is used to evaluate the fitness of an individual more efficiently.The new fitness function is defined as follows:(14)Ft=1−1P∑i=1nxiPi+1−1Q∑i=1nxiQi+1−∑i=1nxiξi∑i=1LMTi.Step 3: Generating populationThe binary coding is adopted in the proposed genetic algorithm. Each gene stands for the state of one appliance, and therefore, the number of genes that compose an individual equals the number of total appliances. Besides, the orders of gene on an individual correspond to different appliances, for example, four appliances (that are, respectively, marked as 1, 2, 3, and 4) with the state ON, OFF, ON, and ON, are encoded as 1011, in which the value 0 in the second place stands for the state of appliance marked as 2 is switched off. For the first step, the initial population can be obtained randomly.Step 4: Calculating fitness valueAccording to equation (14), the fitness value under the current population is calculated. The smaller the value of individual fitness, the better the individual, that is, the combination of the appliance corresponding to the individual is closer to the realistic state.Step 5: Judging convergenceThe convergence criterion can be given by(15)FtCycle−FtCycle−1FtCycle≤εorCycle≤Cyclemax,where Cycle denotes iteration time;FtCycle and FtCycle−1 are the fitness value corresponding to the(Cycle)th and (Cycle − 1)th best individual, respectively; ε is the convergence accuracy; and Cyclemax is the maximum number of iterations set. If the convergence is achieved, stop the solution and go to step 9. Otherwise, go to step 6.Step 6: Selecting and copyingIn this step, the superior individuals are selected to generate the next population. The roulette wheel selection approach [42] is adopted to select the superior individuals. The selective probability of an individual is given as follows:(16)Prs1=Fts1∑j=1nFtsj,where Pr stands for the probability,s1 stands for the individual, and n stands for the total number of individuals.Step 7: Performing crossover and mutationThe crossover operations are first performed for the copied individuals, which are based on the crossover rate and the gene bit to be exchanged. Then, the mutation operations can be performed according to the given mutation rate.Step 8: Updating populationThrough step 7, the new population can be obtained. Then go to step 3 and continue with the next iteration.Step 9: Obtaining a best individualWhen the solution is stopped, give out the best individual, that is, give out the simultaneous identification results.Figure 3
Schematic diagram of the frequency-weighting factor.


## 3.1. Hybrid Features

In this paper, a hybrid feature model is proposed for NILM prediction. The proposed model is a set of features that is composed of active power, reactive power, and harmonic magnitude. The hybrid features of the appliances are expressed in the following matrix form:(8)H=P1Q1M11…M1j…M1LP2Q2M21…M2j…M2L⋮⋮⋮⋮⋮⋮⋮PiQiMi1⋯Mij…MiL⋮⋮⋮⋮⋮⋮⋮PNQNMN1…MNj…MNL,where Pi and Qi are the active and reactive power of ith appliance, respectively. Mij is the jth harmonic magnitude of ith appliance, and N is the total number of appliances. For clarity of presentation, the matrix H also can be regarded as the combination of matrix H=HP,HQ,HM. HP, HQ, and HM are given as follows:(9)HP=P1P2⋯Pi⋯PNT,HQ=Q1Q2⋯Qi⋯QNT,HM=M11…M1j…M1LM21…M2j…M2L⋮⋮⋮⋮⋮Mi1…Mij…MiL⋮⋮⋮⋮⋮MN1…MNj…MNL,where the superscript T denotes the transposition of vectors or matrices. It is noted that when performing the simultaneous identification, both the hybrid feature model of the individual appliance operated independently and multiple appliances operated simultaneously need to be given out.The single load under the individual appliance operated independently and the total load under multiple appliances operated simultaneously are obtained by real-time monitoring. Through the feature extracting and DFT, the hybrid feature models of single load and total load can be given by(10)HS=HS−P,HS−Q,HS−MHT=HT−P,HT−Q,HT−M,where the subscript S denotes the single load that is a load of the individual electrical appliance operated independently and the subscript T denotes the total load that is a load of multiple appliances operated simultaneously. It is noted that HT−P and HT−Q are both constant values, which are defined as HT−P=P and HT−Q=Q, respectively, while HT−M is one-dimensional vector, that is, HT−M=MT1,MT2,…,MTi,…,MTL, where MTi denotes the ith harmonic magnitude of the total load.


## 3.2. Simultaneous Identification Optimization Model

The hybrid feature modelH is divided into three components due to different uses. The vector HP and vector HQ are used to construct the optimization model of simultaneous identification, while the matrix HM is introduced to the algorithm developed for solving the optimization model. The object of load identification is to find the best combination of the appliance load to approach the total power consumption, and thus, the most possible appliances switched on can be inferred. The objective function therefore is given as follows:(11)minfx=HT−P−xHS−P+HT−Q−xHS−Q,where x is the vector of design variables corresponding to each appliance, that is, x=x1,x2,…,xi,…,xN, and xi is the design variable corresponding to the ith appliance.Considering the limit case of simultaneous appliances switched on, the simultaneous identification problem can be established by an optimization problem, which is given by(12)minfx=HT−P−xHS−P+HT−Q−xHS−Qs.t.∑i=1Nxi≤Nxi=0or1.The value ofxi is 0 or 1, and 0 denotes that the ith appliance is not running, while 1 denotes that the ith appliance is running. Therefore, the optimization problem is a 0-1 integer optimization problem. The genetic algorithm [40, 41] is useful for the integer optimization problem. Therefore, in order to improve the identification accuracy and efficiency, a genetic algorithm based on frequency weighting factors (FWF-GA) is developed to solve the simultaneous identification problem of equation (12).


## 3.3. FWF-GA

The procedure of simultaneous identification with FWF-GA is shown in Figure2.Figure 2
Procedures of simultaneous identification with FWF-GA.The procedure can be summarized in the following steps:Step 1: Calculating frequency weighting factorsSince the frequency domain can directly reflect the inherent characteristics of appliance load from a frequency perspective, the frequency weighting factor is defined on the basis of frequency domain information. The frequency weighting factor is defined as follows:(13)ξi=1LMi1MT1+Mi2MT2⋅⋅⋅+MijMTj+⋅⋅⋅+MiLMTL,whereξi is the frequency weighting factor of ith appliance.As shown by equation (13), the frequency weighting factor is actually the average value of L terms weight sum. To further clearly illustrate the defined frequency weighting factor, the preceding six terms weight are taken as an example, which is shown in Figure 3. It is worth noting that the value of L is not absolutely fixed, that is to say, the value of L still depends on an individual decision considering which ones play an important role.Step 2: Constructing fitness functionThe fitness evaluation is one crucial step in the GA process, whose fitness function directly determines the accuracy and efficiency of the algorithm. To further enhance the practicability of simultaneous identification problem solving, a new fitness function is developed based on the frequency-weighting factor, which is used to evaluate the fitness of an individual more efficiently.The new fitness function is defined as follows:(14)Ft=1−1P∑i=1nxiPi+1−1Q∑i=1nxiQi+1−∑i=1nxiξi∑i=1LMTi.Step 3: Generating populationThe binary coding is adopted in the proposed genetic algorithm. Each gene stands for the state of one appliance, and therefore, the number of genes that compose an individual equals the number of total appliances. Besides, the orders of gene on an individual correspond to different appliances, for example, four appliances (that are, respectively, marked as 1, 2, 3, and 4) with the state ON, OFF, ON, and ON, are encoded as 1011, in which the value 0 in the second place stands for the state of appliance marked as 2 is switched off. For the first step, the initial population can be obtained randomly.Step 4: Calculating fitness valueAccording to equation (14), the fitness value under the current population is calculated. The smaller the value of individual fitness, the better the individual, that is, the combination of the appliance corresponding to the individual is closer to the realistic state.Step 5: Judging convergenceThe convergence criterion can be given by(15)FtCycle−FtCycle−1FtCycle≤εorCycle≤Cyclemax,where Cycle denotes iteration time;FtCycle and FtCycle−1 are the fitness value corresponding to the(Cycle)th and (Cycle − 1)th best individual, respectively; ε is the convergence accuracy; and Cyclemax is the maximum number of iterations set. If the convergence is achieved, stop the solution and go to step 9. Otherwise, go to step 6.Step 6: Selecting and copyingIn this step, the superior individuals are selected to generate the next population. The roulette wheel selection approach [42] is adopted to select the superior individuals. The selective probability of an individual is given as follows:(16)Prs1=Fts1∑j=1nFtsj,where Pr stands for the probability,s1 stands for the individual, and n stands for the total number of individuals.Step 7: Performing crossover and mutationThe crossover operations are first performed for the copied individuals, which are based on the crossover rate and the gene bit to be exchanged. Then, the mutation operations can be performed according to the given mutation rate.Step 8: Updating populationThrough step 7, the new population can be obtained. Then go to step 3 and continue with the next iteration.Step 9: Obtaining a best individualWhen the solution is stopped, give out the best individual, that is, give out the simultaneous identification results.Figure 3
Schematic diagram of the frequency-weighting factor.


## 4. Example Verification

In this section, a load identification problem that contains five household appliances are used to demonstrate the proposed method. To verify the effectiveness of the proposed method, two aspects, efficiency and accuracy are mainly considered.


### 4.1. Example Description

Five different appliances including a notebook computer (NC), a fluorescent light (FL), a water fountain (WF), a hairdryer (HD), and a television (TV) are used to verify the proposed method. The single load and total load are obtained by real-time monitoring. The single load is the load when five appliances operate independently, but the total load is the mixed load when five appliances operate under different states.The single load of five appliances is shown in Figure4. As shown in Figure 4, there is only one operating state for the appliance of NC, FL, and TV, while there are many operating states for the appliance of WF and HD. For the WF, there are three states: heating, cooling, and keeping warm. For the HD, there are four states: 1st hot air, 1st cold air, 2nd hot air, and 2nd cold air. The active power and reactive power of each appliance are listed in Table 1. By using the FFT, the current harmonic characteristics of each appliance are shown in Figure 5.Figure 4
Single load of five appliances.Table 1
Features of each appliance.
ApplianceActive power (W)Reactive power (Var)NC24.206.79FL40.740.62WFHeating: 398.22Heating: 0.56Cooling: 64.96Cooling: 12.82Keep warm: 11.62Keep warm: 0.03HD1st hot air: 253.441st hot air: 241.931st cold air: 155.871st cold air: 25.202nd hot air: 1298.122nd hot air: 8.812nd cold air: 313.232nd cold air: 7.54TV108.9535.36Figure 5
Current harmonic characteristics of each appliance.Figure6 shows the total load. It is noted that Figure 6 only gives the partial waveforms of mixed operation states. Eight group data marked in Figure 6 are used to illustrate the effectiveness of the proposed method. The mixed operating states of five appliances corresponding to eight group data are listed in Table 2.Figure 6
Current and voltage waveforms of multiple appliances operated simultaneously.Table 2
Mixed operating states of five appliances corresponding to eight group data.
No.Total active power (W)Total reactive power (Var)Mixed operating states168.61313.983[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]2484.69314.753[1, 1, 0, 0, 1, 0, 0, 0, 1, 0]3517.540.329[1, 0, 1, 0, 0, 0, 0, 0, 0, 1]4729.702102.524[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]5783.11412.995[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]6881.37142.043[1, 1, 1, 0, 0, 0, 0, 0, 1, 1]71,485.59642.08[1, 1, 0, 0, 1, 0, 0, 1, 0, 1]81,883.86239.746[0, 1, 1, 0, 0, 0, 0, 1, 0, 1]In Table2, the active power and reactive power are obtained by feature extraction from the total load. The mixed operating states are the vector of on or off state corresponding to the appliance of NC, FL, heating of WF, cooling of WF and keep warm of WF, 1st hot air of HD, 1st cold air of HD, 2nd hot air of HD, 2nd cold air of HD, and TV. For example, the mixed operating states [1, 1, 1, 0, 0, 0, 1, 0, 0, 1] of fourth group data denote that the current operating state is NC, FL, heating of WF, 1st cold air of HD, and TV.


### 4.2. Solution Analysis of FWF-GA

The FWF-GA is used to perform the load identification. Considering the different operating states of various electrical appliances, it needs to set 10 design variables, which isx=x1,x2,x3,x4,x5,x6,x7,x8,x9,x10, where x1, x2, …, x10, respectively, stand for the on or off state corresponding to the appliance of NC, FL, heating of WF, cooling of WF, keep warm of WF, 1st hot air of HD, 1st cold air of HD, 2nd hot air of HD, 2nd cold air of HD, and TV. In addition, not only the sum of all the design variables cannot be greater than the total number of appliances, but also the sum of the design variables corresponding to WF and HD cannot be greater than the maximum number of running states. Consequently, the simultaneous identification optimization model of the example can be established by(17)minfx=P−xHP+Q−xHQs.t.∑i=110xi≤N;∑i=35xi≤1;∑i=69xi≤1,xi=0or1.In the proposed FWF-GA, some system parameters significantly influence the evaluation performance, such as the crossover rate and population size. The sensitivity of these parameters is tested by the data set of mixed load exacted from Figure6. The length of the data set is 5,000. Table 3 shows the sensitivity of the crossover rate and population size.Table 3
Sensitivity of crossover rate and population size.
ParameterParameter valueSIP (%)AIP (%)NCFLWFHDTVCrossover rate0.692.494.790.890.592.592.60.794.695.192.993.294.793.90.896.497.494.394.795.896.40.995.796.893.794.195.295.6Population size1090.191.589.789.390.290.82092.893.690.491.692.892.73094.595.292.193.794.694.44096.597.894.695.396.396.85096.798.195.095.796.496.9In Table3, SIP is the single identification probability, that is, the probability of identifying each appliance in the total load identification, while AIP is the average identification probability, that is, the probability of identifying all the appliances in the total load identification. As can be seen from Table 3, when the crossover rate is 0.8, the SIP and AIP reach the highest. For the sensitivity of population size, when the population size is 50, the SIP and AIP reach the highest, but the values are very close to those when the population size is 40. Consequently, 0.8 and 50 are chosen as the crossover rate and population size, respectively. Then, the mixed load corresponding to the eight group data can be identified by FWF-GA, and the identification results are shown in Table 4.Table 4
Results by FWF-GA.
No.Design variablesxAccuracy1[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]True2[1, 1, 0, 0, 1, 0, 0, 0, 1, 0]True3[1, 0, 1, 0, 0, 0, 0, 0, 0, 1]True4[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]True5[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]True6[1, 1, 1, 0, 0, 0, 0, 0, 0, 1]False7[1, 1, 0, 0, 1, 0, 0, 1, 0, 1]True8[0, 1, 1, 0, 0, 0, 0, 1, 0, 1]TrueFrom Tables2 and 4, the design variable values of the other seven groups’ data are the same as the values of the actual mixed operating states, except for the 6th group data. For the 6th group data, the actual mixed operating states is NC, FL, heating of WF, 2nd cold air of HD, and TV, but the identification result is NC, FL, heating of WF, and TV. That is, the running of 2nd cold air of HD is not identified for the 6th group data. This demonstrates that there is one false identification result of eight group data for FWF-GA.


### 4.3. Different Fitness Function Evaluation

In order to better illustrate the effectiveness of FWF-GA, the algorithms by different fitness functions are used to perform the identification for the eight group data. The algorithms are called Re-GA and Co-GA, and the fitness function is defined by(18)RE−GA;Ft=1−1P∑i=1nxiPi+1−1Q∑i=1nxiQiCo−GA;Ft=P−xHP+Q−xHQ.The results of different fitness functions are listed in Tables5 and 6.Table 5
Results by different fitness function.
No.Re-GACo-GADesign variables,xAccuracyDesign variables,xAccuracy1[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]True[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]True2[1, 1, 0, 0, 1, 0, 0, 0, 1, 0]True[1, 1, 1, 0, 0, 0, 0, 0, 0, 0]False3[1, 0, 1, 0, 0, 0, 0, 0, 0, 1]True[0, 0, 1, 0, 0, 0, 0, 0, 0, 1]False4[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]True[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]True5[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]True[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]True6[1, 1, 1, 0, 0, 0, 0, 0, 0, 1]False[1, 1, 1, 0, 0, 0, 0, 0, 1, 1]True7[0, 1, 0, 0, 1, 0, 0, 1, 0, 1]False[1, 1, 0, 0, 1, 0, 0, 1, 0, 1]True8[0, 1, 1, 0, 0, 0, 0, 1, 0, 1]True[1, 1, 1, 0, 0, 0, 0, 1, 0, 1]FalseTable 6
Power and relative error by different fitness function.
No.FWF-GARe-GACo-GAPower (W, Var)Relative errorPower (W, Var)Relative errorPower (W, Var)Relative error1[64.96, 12.82][0.0533, 0.0832][64.96, 12.82][0.0533, 0.0832][64.96, 12.82][0.0533, 0.0832]2[389.79, 14.98][0.1958, 0.0154][389.79, 14.98][0.1958, 0.0154][463.16, 4.97][0.0444, 0.4598]3[531.37, 42.71][0.0268, 0.059][531.37, 42.71][0.0268, 0.059][507.17, 35.92][0.0200, 0.1093]4[727.98, 68.53][0.0024, 0.3316][727.98, 68.53][0.0024, 0.3316][727.98, 68.53][0.0024, 0.3316]5[784.81, 284.64][0.0022, 0.3108][784.81, 284.64][0.0022, 0.3108][784.81, 284.64][0.0022, 0.3108]6[572.11, 43.33][0.3509, 0.0306][572.11, 43.33][0.3509, 0.0306][885.34, 50.87][0.0045, 0.2100]7[1483.63, 51.61][0.0013, 0.2265][1459.43, 44.82][0.0176, 0.0651][1483.63, 51.61][0.0013, 0.2265]8[1846.03, 45.35][0.0201, 0.1410][1846.03, 45.35][0.0201, 0.1410][1870.23, 52.14][0.0072, 0.3118]In Table6, the power is the vector of active power and reactive power, and the relative error is the vector of the relative error of active power and reactive power, which can be calculated by(19)ε1=PI−PCPC,ε2=QI−QCQC,where PI and QI are the active power and reactive power by identification, respectively, and PC and QC are the active power and reactive power by actual measurement.As shown in Table5, the identification of two group data (the 6th and 7th groups) are false for Re-GA, and the identification of three group data (the 2nd, 3rd, and 8th groups) are false for Co-GA, while only the 6th group is false for FWF-GA which is shown in Table 4. These data indicate that the proposed FWF-GA has better accuracy than Re-GA and Co-GA.From Table6, comparing the results by Co-GA with that by FWF-GA, the relative error of active power by Co-GA is relatively smaller than that of reactive power in eight group data, which is caused by the numerical gap between active power and reactive power. This is because the fitness function of Co-GA is the absolute error of active and reactive power. The numerical value of active power is relatively larger than reactive power, and the absolute error of active power has a great impact on the optimization results. Consequently, the influence of active power is mainly considered and the influence of reactive power is ignored in order to ensure the minimum target value in the optimization process, which maybe make the identification results incorrect. Additionally, comparing the results by Re-GA with that by FWF-GA, the results of the 6th group are both false by FWF-GA and RE-GA, but the results of the 7th group are also false by RE-GA. FWF-GA and RE-GA both use the relative error of active and reactive power, while FWF-GA also uses the frequency weighting factor of harmonic feature as the adjustment. For further description of the influence of the harmonic feature, the frequency-weighting factor corresponding to the 7th group data is given in Table 7.Table 7
Frequency weighting factor corresponding to the 7th group.
iFrequency-weighting factor10.987620.072430.897140.951350.574260.653370.175380.923590.6791100.1048From Tables5 and 2, the mistake of the results of the 7th group data is that the design variable x1 corresponding to NC is 0, which means that NC is wrongly judged as not running. This is because RE-GA adopts the relative error of active and reactive power as the optimization objective. For the 7th group data, the relative error of reactive power has a greater impact on the objective value than that of active power. From Table 1, the active power of NC is small, but the reactive power is relatively large. The optimization objective under the running state of NC will be enlarger than that under the closed state. Consequently, RE-GA gives false results. However, for FWF-GA, not only the relative errors of active and reactive power are adopted as the optimization objective but also the magnitude of harmonic feature are used to adjust the influence of active and reactive power. From Table 7, the frequency weighting factor ξ1 corresponding to NC is larger than other factors. Through the adjustment of the harmonic feature, the influence on the optimization result of the open/closed state of NC is changed. Consequently, FWF-GA gives true results.In order to clearly describe the identification efficiency of different fitness functions, taking the results of the 4th group to draw the iterative history as shown in Figure 7, and the convergence times of 8 group data are shown in Table 8.Figure 7
Iterative history of the 4th group data.Table 8
Convergence times of eight mixed loads.
No.Convergence timesFWF-GARE-GACo-GA11071421210153101013415102051211166161518715917810712From Figure7, the three algorithms of FWF-GA, RE-GA, and Co-GA tend to be stable after 15 iterations, 10 iterations, and 20 iterations, respectively. As shown in Table 8, the number of convergence times of FWF-GA (which is range from 10 to 16) is less than that of Co-GA (which is range from 12 to 20) while larger than that of RE-GA (which is range from 7 to 15). These show that the FWF-GA has better computational efficiency than Co-GA. However, due to increasing the adjustment of harmonic features for FWF-GA, the computational efficiency of FWF-GA is not as good as that of RE-GA.


### 4.4. Comparison of Different Methods

In order to verify the effectiveness of the proposed method, the existing method such as DAPSO, CNN, and QP are also used to perform the identification of the data set of mixed load. Table9 gives the identification accuracy and test time.Table 9
Comparison of each method.
MethodSIP (%)AIP (%)Costtime (s)NCFLWFHDTVDAPSO94.395.692.191.395.393.7211.82CNN95.898.293.595.794.895.4314.25QP96.498.594.795.295.596.6675.67Proposedmethod96.798.195.095.796.496.9260.73As shown in Table9, the minimum of SIP for five appliances (NC is 94.3%, FL is 95.6%, WF is 92.1%, HD is 91.3%, and TV is 94.8%) is obtained by DAPSO and CNN, respectively, while the maximum (NC is 96.7%, FL is 98.1%, WF is 95.0%, HD is 95.7%, and TV is 96.4%) is obtained by the proposed method. Furthermore, the minimum of AIP is 93.7% by DAPSO, and the maximum is 96.9 by the proposed method. It can be also derived from Table 9 that DAPSO takes the least time (211.82 s in total), and QP takes the longest time (675.67 s in total). The proposed method takes 260.73 s, which is larger than DAPSO and lower than CNN and QP. These demonstrate that the proposed method has a better accuracy and calculation efficiency.


### 4.5. Influence Analysis of Noise

In load identification, noise inevitably exists. For the purpose of test, the influence of the noise on the identification result of FWF-GA, Gaussian white noise with various signal-to-noise ratios (SNR) is added to the test data. The SNR is designed to be 50 dB, 40 dB, 30 dB, and 20 dB respectively. The identification accuracy with various SNRs is shown in Table10.Table 10
Identification accuracy with various SNR by FWF-GA.
SNR (dB)SIP (%)AIP (%)NCFLWFHDTV2094.794.692.893.095.294.53095.295.393.293.795.495.44095.696.193.994.695.796.15096.397.594.895.396.196.7As shown in Table10, the SIP and AIP with different SNR are reduced compared with those without noise, but all the values of identification accuracy are still larger than 90%, which demonstrates that the proposed identification method still maintains a high accuracy under various SNR.


## 4.1. Example Description

Five different appliances including a notebook computer (NC), a fluorescent light (FL), a water fountain (WF), a hairdryer (HD), and a television (TV) are used to verify the proposed method. The single load and total load are obtained by real-time monitoring. The single load is the load when five appliances operate independently, but the total load is the mixed load when five appliances operate under different states.The single load of five appliances is shown in Figure4. As shown in Figure 4, there is only one operating state for the appliance of NC, FL, and TV, while there are many operating states for the appliance of WF and HD. For the WF, there are three states: heating, cooling, and keeping warm. For the HD, there are four states: 1st hot air, 1st cold air, 2nd hot air, and 2nd cold air. The active power and reactive power of each appliance are listed in Table 1. By using the FFT, the current harmonic characteristics of each appliance are shown in Figure 5.Figure 4
Single load of five appliances.Table 1
Features of each appliance.
ApplianceActive power (W)Reactive power (Var)NC24.206.79FL40.740.62WFHeating: 398.22Heating: 0.56Cooling: 64.96Cooling: 12.82Keep warm: 11.62Keep warm: 0.03HD1st hot air: 253.441st hot air: 241.931st cold air: 155.871st cold air: 25.202nd hot air: 1298.122nd hot air: 8.812nd cold air: 313.232nd cold air: 7.54TV108.9535.36Figure 5
Current harmonic characteristics of each appliance.Figure6 shows the total load. It is noted that Figure 6 only gives the partial waveforms of mixed operation states. Eight group data marked in Figure 6 are used to illustrate the effectiveness of the proposed method. The mixed operating states of five appliances corresponding to eight group data are listed in Table 2.Figure 6
Current and voltage waveforms of multiple appliances operated simultaneously.Table 2
Mixed operating states of five appliances corresponding to eight group data.
No.Total active power (W)Total reactive power (Var)Mixed operating states168.61313.983[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]2484.69314.753[1, 1, 0, 0, 1, 0, 0, 0, 1, 0]3517.540.329[1, 0, 1, 0, 0, 0, 0, 0, 0, 1]4729.702102.524[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]5783.11412.995[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]6881.37142.043[1, 1, 1, 0, 0, 0, 0, 0, 1, 1]71,485.59642.08[1, 1, 0, 0, 1, 0, 0, 1, 0, 1]81,883.86239.746[0, 1, 1, 0, 0, 0, 0, 1, 0, 1]In Table2, the active power and reactive power are obtained by feature extraction from the total load. The mixed operating states are the vector of on or off state corresponding to the appliance of NC, FL, heating of WF, cooling of WF and keep warm of WF, 1st hot air of HD, 1st cold air of HD, 2nd hot air of HD, 2nd cold air of HD, and TV. For example, the mixed operating states [1, 1, 1, 0, 0, 0, 1, 0, 0, 1] of fourth group data denote that the current operating state is NC, FL, heating of WF, 1st cold air of HD, and TV.


## 4.2. Solution Analysis of FWF-GA

The FWF-GA is used to perform the load identification. Considering the different operating states of various electrical appliances, it needs to set 10 design variables, which isx=x1,x2,x3,x4,x5,x6,x7,x8,x9,x10, where x1, x2, …, x10, respectively, stand for the on or off state corresponding to the appliance of NC, FL, heating of WF, cooling of WF, keep warm of WF, 1st hot air of HD, 1st cold air of HD, 2nd hot air of HD, 2nd cold air of HD, and TV. In addition, not only the sum of all the design variables cannot be greater than the total number of appliances, but also the sum of the design variables corresponding to WF and HD cannot be greater than the maximum number of running states. Consequently, the simultaneous identification optimization model of the example can be established by(17)minfx=P−xHP+Q−xHQs.t.∑i=110xi≤N;∑i=35xi≤1;∑i=69xi≤1,xi=0or1.In the proposed FWF-GA, some system parameters significantly influence the evaluation performance, such as the crossover rate and population size. The sensitivity of these parameters is tested by the data set of mixed load exacted from Figure6. The length of the data set is 5,000. Table 3 shows the sensitivity of the crossover rate and population size.Table 3
Sensitivity of crossover rate and population size.
ParameterParameter valueSIP (%)AIP (%)NCFLWFHDTVCrossover rate0.692.494.790.890.592.592.60.794.695.192.993.294.793.90.896.497.494.394.795.896.40.995.796.893.794.195.295.6Population size1090.191.589.789.390.290.82092.893.690.491.692.892.73094.595.292.193.794.694.44096.597.894.695.396.396.85096.798.195.095.796.496.9In Table3, SIP is the single identification probability, that is, the probability of identifying each appliance in the total load identification, while AIP is the average identification probability, that is, the probability of identifying all the appliances in the total load identification. As can be seen from Table 3, when the crossover rate is 0.8, the SIP and AIP reach the highest. For the sensitivity of population size, when the population size is 50, the SIP and AIP reach the highest, but the values are very close to those when the population size is 40. Consequently, 0.8 and 50 are chosen as the crossover rate and population size, respectively. Then, the mixed load corresponding to the eight group data can be identified by FWF-GA, and the identification results are shown in Table 4.Table 4
Results by FWF-GA.
No.Design variablesxAccuracy1[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]True2[1, 1, 0, 0, 1, 0, 0, 0, 1, 0]True3[1, 0, 1, 0, 0, 0, 0, 0, 0, 1]True4[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]True5[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]True6[1, 1, 1, 0, 0, 0, 0, 0, 0, 1]False7[1, 1, 0, 0, 1, 0, 0, 1, 0, 1]True8[0, 1, 1, 0, 0, 0, 0, 1, 0, 1]TrueFrom Tables2 and 4, the design variable values of the other seven groups’ data are the same as the values of the actual mixed operating states, except for the 6th group data. For the 6th group data, the actual mixed operating states is NC, FL, heating of WF, 2nd cold air of HD, and TV, but the identification result is NC, FL, heating of WF, and TV. That is, the running of 2nd cold air of HD is not identified for the 6th group data. This demonstrates that there is one false identification result of eight group data for FWF-GA.


## 4.3. Different Fitness Function Evaluation

In order to better illustrate the effectiveness of FWF-GA, the algorithms by different fitness functions are used to perform the identification for the eight group data. The algorithms are called Re-GA and Co-GA, and the fitness function is defined by(18)RE−GA;Ft=1−1P∑i=1nxiPi+1−1Q∑i=1nxiQiCo−GA;Ft=P−xHP+Q−xHQ.The results of different fitness functions are listed in Tables5 and 6.Table 5
Results by different fitness function.
No.Re-GACo-GADesign variables,xAccuracyDesign variables,xAccuracy1[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]True[0, 0, 0, 1, 0, 0, 0, 0, 0, 0]True2[1, 1, 0, 0, 1, 0, 0, 0, 1, 0]True[1, 1, 1, 0, 0, 0, 0, 0, 0, 0]False3[1, 0, 1, 0, 0, 0, 0, 0, 0, 1]True[0, 0, 1, 0, 0, 0, 0, 0, 0, 1]False4[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]True[1, 1, 1, 0, 0, 0, 1, 0, 0, 1]True5[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]True[1, 0, 1, 0, 0, 1, 0, 0, 0, 1]True6[1, 1, 1, 0, 0, 0, 0, 0, 0, 1]False[1, 1, 1, 0, 0, 0, 0, 0, 1, 1]True7[0, 1, 0, 0, 1, 0, 0, 1, 0, 1]False[1, 1, 0, 0, 1, 0, 0, 1, 0, 1]True8[0, 1, 1, 0, 0, 0, 0, 1, 0, 1]True[1, 1, 1, 0, 0, 0, 0, 1, 0, 1]FalseTable 6
Power and relative error by different fitness function.
No.FWF-GARe-GACo-GAPower (W, Var)Relative errorPower (W, Var)Relative errorPower (W, Var)Relative error1[64.96, 12.82][0.0533, 0.0832][64.96, 12.82][0.0533, 0.0832][64.96, 12.82][0.0533, 0.0832]2[389.79, 14.98][0.1958, 0.0154][389.79, 14.98][0.1958, 0.0154][463.16, 4.97][0.0444, 0.4598]3[531.37, 42.71][0.0268, 0.059][531.37, 42.71][0.0268, 0.059][507.17, 35.92][0.0200, 0.1093]4[727.98, 68.53][0.0024, 0.3316][727.98, 68.53][0.0024, 0.3316][727.98, 68.53][0.0024, 0.3316]5[784.81, 284.64][0.0022, 0.3108][784.81, 284.64][0.0022, 0.3108][784.81, 284.64][0.0022, 0.3108]6[572.11, 43.33][0.3509, 0.0306][572.11, 43.33][0.3509, 0.0306][885.34, 50.87][0.0045, 0.2100]7[1483.63, 51.61][0.0013, 0.2265][1459.43, 44.82][0.0176, 0.0651][1483.63, 51.61][0.0013, 0.2265]8[1846.03, 45.35][0.0201, 0.1410][1846.03, 45.35][0.0201, 0.1410][1870.23, 52.14][0.0072, 0.3118]In Table6, the power is the vector of active power and reactive power, and the relative error is the vector of the relative error of active power and reactive power, which can be calculated by(19)ε1=PI−PCPC,ε2=QI−QCQC,where PI and QI are the active power and reactive power by identification, respectively, and PC and QC are the active power and reactive power by actual measurement.As shown in Table5, the identification of two group data (the 6th and 7th groups) are false for Re-GA, and the identification of three group data (the 2nd, 3rd, and 8th groups) are false for Co-GA, while only the 6th group is false for FWF-GA which is shown in Table 4. These data indicate that the proposed FWF-GA has better accuracy than Re-GA and Co-GA.From Table6, comparing the results by Co-GA with that by FWF-GA, the relative error of active power by Co-GA is relatively smaller than that of reactive power in eight group data, which is caused by the numerical gap between active power and reactive power. This is because the fitness function of Co-GA is the absolute error of active and reactive power. The numerical value of active power is relatively larger than reactive power, and the absolute error of active power has a great impact on the optimization results. Consequently, the influence of active power is mainly considered and the influence of reactive power is ignored in order to ensure the minimum target value in the optimization process, which maybe make the identification results incorrect. Additionally, comparing the results by Re-GA with that by FWF-GA, the results of the 6th group are both false by FWF-GA and RE-GA, but the results of the 7th group are also false by RE-GA. FWF-GA and RE-GA both use the relative error of active and reactive power, while FWF-GA also uses the frequency weighting factor of harmonic feature as the adjustment. For further description of the influence of the harmonic feature, the frequency-weighting factor corresponding to the 7th group data is given in Table 7.Table 7
Frequency weighting factor corresponding to the 7th group.
iFrequency-weighting factor10.987620.072430.897140.951350.574260.653370.175380.923590.6791100.1048From Tables5 and 2, the mistake of the results of the 7th group data is that the design variable x1 corresponding to NC is 0, which means that NC is wrongly judged as not running. This is because RE-GA adopts the relative error of active and reactive power as the optimization objective. For the 7th group data, the relative error of reactive power has a greater impact on the objective value than that of active power. From Table 1, the active power of NC is small, but the reactive power is relatively large. The optimization objective under the running state of NC will be enlarger than that under the closed state. Consequently, RE-GA gives false results. However, for FWF-GA, not only the relative errors of active and reactive power are adopted as the optimization objective but also the magnitude of harmonic feature are used to adjust the influence of active and reactive power. From Table 7, the frequency weighting factor ξ1 corresponding to NC is larger than other factors. Through the adjustment of the harmonic feature, the influence on the optimization result of the open/closed state of NC is changed. Consequently, FWF-GA gives true results.In order to clearly describe the identification efficiency of different fitness functions, taking the results of the 4th group to draw the iterative history as shown in Figure 7, and the convergence times of 8 group data are shown in Table 8.Figure 7
Iterative history of the 4th group data.Table 8
Convergence times of eight mixed loads.
No.Convergence timesFWF-GARE-GACo-GA11071421210153101013415102051211166161518715917810712From Figure7, the three algorithms of FWF-GA, RE-GA, and Co-GA tend to be stable after 15 iterations, 10 iterations, and 20 iterations, respectively. As shown in Table 8, the number of convergence times of FWF-GA (which is range from 10 to 16) is less than that of Co-GA (which is range from 12 to 20) while larger than that of RE-GA (which is range from 7 to 15). These show that the FWF-GA has better computational efficiency than Co-GA. However, due to increasing the adjustment of harmonic features for FWF-GA, the computational efficiency of FWF-GA is not as good as that of RE-GA.


## 4.4. Comparison of Different Methods

In order to verify the effectiveness of the proposed method, the existing method such as DAPSO, CNN, and QP are also used to perform the identification of the data set of mixed load. Table9 gives the identification accuracy and test time.Table 9
Comparison of each method.
MethodSIP (%)AIP (%)Costtime (s)NCFLWFHDTVDAPSO94.395.692.191.395.393.7211.82CNN95.898.293.595.794.895.4314.25QP96.498.594.795.295.596.6675.67Proposedmethod96.798.195.095.796.496.9260.73As shown in Table9, the minimum of SIP for five appliances (NC is 94.3%, FL is 95.6%, WF is 92.1%, HD is 91.3%, and TV is 94.8%) is obtained by DAPSO and CNN, respectively, while the maximum (NC is 96.7%, FL is 98.1%, WF is 95.0%, HD is 95.7%, and TV is 96.4%) is obtained by the proposed method. Furthermore, the minimum of AIP is 93.7% by DAPSO, and the maximum is 96.9 by the proposed method. It can be also derived from Table 9 that DAPSO takes the least time (211.82 s in total), and QP takes the longest time (675.67 s in total). The proposed method takes 260.73 s, which is larger than DAPSO and lower than CNN and QP. These demonstrate that the proposed method has a better accuracy and calculation efficiency.


## 4.5. Influence Analysis of Noise

In load identification, noise inevitably exists. For the purpose of test, the influence of the noise on the identification result of FWF-GA, Gaussian white noise with various signal-to-noise ratios (SNR) is added to the test data. The SNR is designed to be 50 dB, 40 dB, 30 dB, and 20 dB respectively. The identification accuracy with various SNRs is shown in Table10.Table 10
Identification accuracy with various SNR by FWF-GA.
SNR (dB)SIP (%)AIP (%)NCFLWFHDTV2094.794.692.893.095.294.53095.295.393.293.795.495.44095.696.193.994.695.796.15096.397.594.895.396.196.7As shown in Table10, the SIP and AIP with different SNR are reduced compared with those without noise, but all the values of identification accuracy are still larger than 90%, which demonstrates that the proposed identification method still maintains a high accuracy under various SNR.


## 5. Conclusions

In this paper, a simultaneous load identification method is proposed to achieve the identification of the total load under multiple appliances operated simultaneously. Three parts including hybrid features extraction, simultaneous identification optimization model, and FWF-GA are mainly involved in the proposed method. The simultaneous identification optimization model is constructed by using the hybrid feature model. FWF-GA is developed based on the GA framework to efficiently solve the identification optimization model. In FWF-GA, the fitness function is constructed on the basis of hybrid features, and the relative error of active power, reactive power, and the frequency weighting factor are used to evaluate the fitness of individuals. The practicability of the developed fitness function is examined by comparing with Re-GA and Co-GA, and the accuracy and efficiency of the proposed FWF-GA are compared with the existing DAPSO, CNN, and QP. The test results all demonstrated that FWF-GA has better accuracy and efficiency than the other methods. In addition, the impact of noise on the identification results is analyzed by adding Gaussian white noise with various SNRs to the test data. The result also demonstrates that FWF-GA still has a high accuracy under various SNRs.For the shortcomings of the method proposed in this paper, FWF-GA is proposed to identify the stable operation state of all electrical appliances, while instantaneous state on or off with on or off is not considered. In addition, the uncertainty of the appliance operation has an important impact on the identification result. Another important problem is that, in order to get better accuracy, the adjustment of system parameters needs a large calculation cost. Consequently, it would be beneficial to develop an adaptive identification method considering the instantaneous state and the operating uncertainty in future work.



---

*Source: 7876380-2022-05-19.xml*