CN112821559A - Non-invasive household appliance load depth re-identification method - Google Patents

Non-invasive household appliance load depth re-identification method Download PDF

Info

Publication number
CN112821559A
CN112821559A CN202110089840.XA CN202110089840A CN112821559A CN 112821559 A CN112821559 A CN 112821559A CN 202110089840 A CN202110089840 A CN 202110089840A CN 112821559 A CN112821559 A CN 112821559A
Authority
CN
China
Prior art keywords
sample
formula
identification
value
load
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110089840.XA
Other languages
Chinese (zh)
Other versions
CN112821559B (en
Inventor
张志禹
周咪
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shenzhen Wanzhida Technology Co ltd
Wuxing Technology Shenzhen Co ltd
Original Assignee
Xian University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xian University of Technology filed Critical Xian University of Technology
Priority to CN202110089840.XA priority Critical patent/CN112821559B/en
Publication of CN112821559A publication Critical patent/CN112821559A/en
Application granted granted Critical
Publication of CN112821559B publication Critical patent/CN112821559B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • H—ELECTRICITY
    • H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J13/00—Circuit arrangements for providing remote indication of network conditions, e.g. an instantaneous record of the open or closed condition of each circuitbreaker in the network; Circuit arrangements for providing remote control of switching means in a power distribution network, e.g. switching in and out of current consumers by using a pulse code signal carried by the network
    • H02J13/00002—Circuit arrangements for providing remote indication of network conditions, e.g. an instantaneous record of the open or closed condition of each circuitbreaker in the network; Circuit arrangements for providing remote control of switching means in a power distribution network, e.g. switching in and out of current consumers by using a pulse code signal carried by the network characterised by monitoring
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00—Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/04—Architecture, e.g. interconnection topology
    • G06N3/045—Combinations of networks
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/08—Learning methods
    • H—ELECTRICITY
    • H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00—Circuit arrangements for ac mains or ac distribution networks
    • H—ELECTRICITY
    • H02—GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02J—CIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2310/00—The network for supplying or distributing electric power characterised by its spatial reach or by the load
    • H02J2310/10—The network having a local or delimited stationary reach
    • H02J2310/12—The local stationary network supplying a household or a building
    • H02J2310/14—The load or loads being home appliances
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02B—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO BUILDINGS, e.g. HOUSING, HOUSE APPLIANCES OR RELATED END-USER APPLICATIONS
    • Y02B70/00—Technologies for an efficient end-user side electric power management and consumption
    • Y02B70/30—Systems integrating technologies related to power network operation and communication or information technologies for improving the carbon footprint of the management of residential or tertiary loads, i.e. smart grids as climate change mitigation technology in the buildings sector, including also the last stages of power distribution and the control, monitoring or operating management systems at local level
    • Y02B70/3225—Demand response systems, e.g. load shedding, peak shaving
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y04—INFORMATION OR COMMUNICATION TECHNOLOGIES HAVING AN IMPACT ON OTHER TECHNOLOGY AREAS
    • Y04S—SYSTEMS INTEGRATING TECHNOLOGIES RELATED TO POWER NETWORK OPERATION, COMMUNICATION OR INFORMATION TECHNOLOGIES FOR IMPROVING THE ELECTRICAL POWER GENERATION, TRANSMISSION, DISTRIBUTION, MANAGEMENT OR USAGE, i.e. SMART GRIDS
    • Y04S20/00—Management or operation of end-user stationary applications or the last stages of power distribution; Controlling, monitoring or operating thereof
    • Y04S20/20—End-user application control systems
    • Y04S20/222—Demand response systems, e.g. load shedding, peak shaving

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Physics & Mathematics (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Computing Systems (AREA)
  • Mathematical Physics (AREA)
  • Computational Linguistics (AREA)
  • General Health & Medical Sciences (AREA)
  • Molecular Biology (AREA)
  • Biophysics (AREA)
  • Biomedical Technology (AREA)
  • Health & Medical Sciences (AREA)
  • Software Systems (AREA)
  • Power Engineering (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Biology (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

The invention discloses a non-invasive household appliance load depth re-identification method, which is implemented according to the following steps: denoising the acquired high-frequency data of the user; event detection is carried out on the data; extracting multi-dimensional load characteristics for the detected event change points; screening relevant features aiming at the extracted multi-dimensional features; and taking the obtained features as load marks, identifying the working state of the household appliance load in the user by using the maximum membership degree through establishing a load feature library, training a parameter smoothing factor of the GRNN, and performing deep re-identification on a sample to be detected by using the trained GRNN to finish the final identification of the non-invasive household appliance load. The problem of load identification accuracy rate low caused by a single algorithm in the prior art is solved.

Description

Non-invasive household appliance load depth re-identification method
Technical Field
The invention belongs to the technical field of household appliance load identification, and relates to a non-invasive household appliance load depth re-identification method.
Background
The load monitoring technology is one of basic technologies for building ubiquitous power internet of things and transparent power grids, and through knowing the power utilization law, a user can be helped to save the electric energy expenditure by 14%, and the non-intrusive load monitoring (NILM) technology becomes an attention object of power grid service providers and users due to the characteristic of low cost.
The non-invasive load monitoring means that monitoring equipment is installed at an electric power inlet, the types and the running conditions of single loads in a load cluster are obtained by monitoring signals such as voltage, current and the like at the position and analyzing, the non-invasive household appliance load identification is a user-side-oriented non-invasive load monitoring technology, and the process is concentrated in three steps: event detection, feature extraction and load identification. In the aspect of event detection, the generalized likelihood estimation (GLR) mainly considers the mean value, so that misjudgment is easy to occur, and the F test can reflect the sample variance; in the aspect of feature extraction, feature types are generally determined subjectively according to experience, existing information is underutilized, semi-supervised Relief-F and filtering feature selection combining maximum correlation and minimum redundancy (mRmR) can fully utilize known label information and quantize features; in addition, the prior art has the problem that the load identification accuracy is low due to power overlapping and a single algorithm.
Disclosure of Invention
The invention aims to provide a non-invasive household appliance load depth re-identification method, which solves the problem of low load identification accuracy caused by a single algorithm in the prior art.
The technical scheme adopted by the invention is that the non-invasive household appliance load depth re-identification method is implemented according to the following steps:
step 1, acquiring user high-frequency data, and performing denoising processing on the acquired user high-frequency data;
step 2, carrying out event detection on the data in the step 1 through improved generalized likelihood ratio test, if an event is detected, executing the step 3, otherwise, returning to the step 1;
step 3, extracting multi-dimensional load characteristics for the detected event change points;
step 4, aiming at the multi-dimensional features extracted in the step 3, screening relevant features by using a semi-supervised algorithm combining a Relief-F and a mRMR;
step 5, taking the characteristics obtained in the step 4 as load marks, establishing a load characteristic library through a self-adaptive FCM algorithm, identifying the load working state of the household appliances in the user by using the maximum membership degree, finishing the identification if the identification results of the two FCM algorithms are consistent, and executing the step 6 if the identification results are inconsistent;
and 6, training the parameter smoothing factor of the GRNN by adopting an SA-BAS (simulated annealing-Tianniu whisker) algorithm, and performing deep re-identification on the sample to be detected by utilizing the trained GRNN to finish the final identification of the non-invasive household appliance load.
The invention is also characterized in that:
the step 1 is implemented according to the following steps:
step 1.1, acquiring high-frequency household appliance load data containing electrical parameters, specifically representing voltage, current and corresponding power;
step 1.2, denoising processing of the power signal, wherein because an isolated noise point is easily identified as an event by an event detection algorithm, a median filtering method is selected to process the original power signal, so that the noise is eliminated without changing edge information: suppose that there is a sequence x of digital signalsj(-∞<j<+ ∞) is filtered, a window with an odd number L is first defined, where L is 2N +1, N is a positive integer, and it is assumed that at a certain time i, the signal sample in the window is xi-N,…,xi,…,xi+NWherein x isiIs the sample value of the signal located in the center of the window, and after rearranging the L signals from small to large, the value is defined as the output value of the median filter.
The step 2 is implemented according to the following steps:
step 2.1, calculating the active power P of the fundamental wave according to the formula (1)1Taking the combined active power P as a two-dimensional power time sequence
Figure BDA0002912022900000031
A binary hypothesis test is proposed according to the formula (2);
Figure BDA0002912022900000032
Figure BDA0002912022900000033
in the formula, V1Is a fundamental voltage, I1Is the current of the fundamental wave,
Figure BDA0002912022900000034
is the phase difference between the two. n iscFor the time of occurrence of the change point, k is the total length of the window, n is the last sample time in the window, μ0,∑0Testing for hypotheses H0Mean of Gaussian distribution, covariance matrix, μ under the conditionsa,∑aIs H1Multi-dimensional signal mean, multi-dimensional covariance matrix, mu, before occurrence of change point under conditionb,∑bIs H1A multidimensional signal mean and a multidimensional covariance matrix after the variable point occurs under the condition;
step 2.2, defining two consecutive windows W within this time sequenceaAnd WbThe sample in the two windows is Xn＝{xmWhere m is n-k +1, …, n, and the window lengths are both k/2, the μ sum Σ in each of the two windows is calculated from equations (3) and (4), and then the decision function g is calculated from equation (5)n；
Figure BDA0002912022900000035
Figure BDA0002912022900000037
Figure BDA0002912022900000036
Step 2.3, mixing gnAnd a threshold value h1And comparing and searching suspicious points of event occurrence: when the decision function value is larger than h1When it is, refuse H0The data distribution in the two windows is not consistent, and at the point-changing time ncThere is a possibility of an event occurring; when the decision function is less than h1When it is, refuse H1And the data distribution of the two windows is consistent, and no event occurs. And when the sampling points of the detected multiple events are within 3 continuously or at intervals, the sampling points are regarded as having the maximum gnA corresponding event occurs;
and 2.4, taking the event point as a base point, and carrying out F-test screening on the false detection events. First, it is assumed that the variance values of the window data before and after the base point are equal:
Figure BDA0002912022900000041
calculating the value of statistic F according to equation (6), giving significance level α, rejecting hypothesis H if the value of F satisfies equation (7)0Judging that the two maternal variances have obvious difference, and judging that an event occurs at the point;
Figure BDA0002912022900000042
Figure BDA0002912022900000043
or
Figure BDA0002912022900000044
In the formula (I), the compound is shown in the specification,
Figure BDA0002912022900000045
respectively the sample volume and variance of the window before and after the base point.
Step 3 is specifically carried out as follows: the method for extracting the power characteristics of the variable points specifically comprises the following steps: active power, fundamental active power, reactive power, fundamental reactive power, apparent power, distortion power, power factor angle, fundamental power factor; extracting harmonic characteristics at a variable point, specifically including voltage, amplitude of each harmonic of one to nine times of voltage, content rate of each harmonic, difference of content rates of each harmonic and total harmonic distortion rate; the current waveform characteristics comprise the wave peak value, the average value and the wave crest coefficient; the method for extracting the V-I track characteristics at the variable points specifically comprises the following steps: symmetry, surrounding direction, surrounding area, number of intersection points, Y-axis intercept, Y-axis span, midline curvature, trace middle part peak value, left and right part area, middle part shape, and instantaneous admittance standard deviation.
Step 4 is specifically implemented according to the following steps:
step 4.1, labeling part of multi-dimensional characteristic samples, and recording as marked samples S1And the remainder is marked as the unmarked sample S2The label is C, and the multidimensional feature set is A ═ A1,A2,…,AN)；
Step 4.2, Slave S1In which a sample s of class C is randomly drawnq(CqE.g. C); from S2D neighboring samples are selected and recorded as
Figure BDA0002912022900000051
At S1From C to CqEach class C of the otherpRespectively solving a nearest neighbor sample x of s in the epsilon C (p is not equal to q); and at S2D neighbor samples of medium solution x, noted
Figure BDA0002912022900000052
Wherein the neighbor formula is shown as formula (8); updating all feature weights based on equation (9);
Figure BDA0002912022900000053
Figure BDA0002912022900000054
Figure BDA0002912022900000055
in the formula, M represents the iteration number, d is the number of the adjacent samples,
Figure BDA0002912022900000056
for the t-th neighbor sample in the q class to which the sample s belongs,
Figure BDA0002912022900000057
representing the t-th neighbouring sample in the P-th class, which is different from the class of the sample s, P (C)p) Representing the probability of the p-th class of objects, AkAs a result of the k-th feature,
Figure BDA0002912022900000058
representing a sample s and a sample
Figure BDA0002912022900000059
With respect to feature AkThe distance of (d);
step 4.3, circularly executing the step 4.2 for M times, and obtaining the characteristic weight omega finally output by the semi-supervised Relief-Fk'(k is 1,2, …, w'), removing features with weight coefficient less than theta to obtain a candidate feature subset, and then performing semi-supervised mRmR feature selection on the candidate feature subset;
step 4.4, marking the sample S according to the formula (11)1Calculating the correlation degree between each feature and the sample label, and obtaining the unmarked sample S according to the formula (12)2Calculating the redundancy of each characteristic;
Figure BDA00029120229000000510
Figure BDA00029120229000000511
in the formula (I), the compound is shown in the specification,
Figure BDA00029120229000000512
is shown in the labeled sample S1In (1) calculation ofkMutual information with label C, R (A)k,Sm-1) Representing a subset S of existing featuresm-1Which contains m-1 features, not Sm-1Features A in the subsetkRedundancy with selected features;
step 4.5, establishing a characteristic candidate set H, and selecting the maximum correlation degree DmaxCorresponding features as candidate set leader H1Sequentially selecting the kth feature A according to equation (13)kPutting the obtained product into H until the designated characteristic number w is selected;
Figure BDA0002912022900000061
and 4.6, calculating each characteristic weight according to the formula (14).
Figure BDA0002912022900000062
Step 5 is specifically implemented according to the following steps:
step 5.1, setting the initial clustering number b, determining an initial clustering center based on a maximum-minimum distance algorithm, and simultaneously, enabling a counter k to return to zero;
a. calculating the mean of the data set
Figure BDA0002912022900000063
The sample point farthest from the mean is V1j(j ═ 1,2, …, w); b. calculating the minimum distance D between each data point and the selected cluster center according to the formula (15)xSelection of DxThe maximum value point is used as a new clustering center; c. repeating the step b until b initial clustering centers are selected;
Dx＝min d(xi,Z'k) k'＝1,…,kselected (15)
Figure BDA0002912022900000064
step 5.2, constructing a loss function L of the FCM1 algorithm1As shown in formula (17), for uij，λjAnd vikCalculating and combining the partial derivatives, and performing u by using the derived formula (18) and formula (19)ijAnd vikThe iterative update of (2);
Figure BDA0002912022900000065
Figure BDA0002912022900000071
Figure BDA0002912022900000072
wherein n is the number of samples, b is the number of clusters, w is the number of features,
Figure BDA0002912022900000073
is a sample weight, uijThe membership degree of the jth sample belonging to the ith class, v is a clustering center matrix, TuE (0, infinity) replaces the fuzzy index m for controlling the entropy value, and the maximum entropy is introduced as regularization;
step 5.3, calculating J1A value of, if
Figure BDA0002912022900000074
Finishing clustering of the initial clustering number b, terminating iteration, otherwise, turning to the step 5.2 to continue iteration until the requirement is met;
step 5.4, calculating the clustering effectiveness index of the initial clustering number b according to the formula (20), namely the contour coefficient SbIf it satisfies Sb-2＜Sb-1And Sb-1＞SbIf yes, the self-adaptive FCM1 classification is ended, otherwise, whether the clustering number b reaches the maximum value or not is judged
Figure BDA0002912022900000075
If the maximum value is reached, all S are takenbIf b, u and v corresponding to the medium maximum value do not reach the maximum value, b is equal to b +1, and the step 5.1 is executed again;
Figure BDA0002912022900000076
wherein n is the number of samples, a (i) is the average distance between the sample point i and the remaining sample objects in the same cluster, b (i) is the minimum value of the average distance between the sample point i and the sample objects in each of the remaining clusters, SbIn [ -1,1 [)]The larger the S value is, the better the clustering effect is;
step 5.5, setting the optimal clustering number b of the FCM1 as the clustering number of the FCM2, similarly adopting a maximum-minimum distance algorithm to determine an initial clustering center, and simultaneously, enabling a counter k to return to zero;
step 5.6, constructing a loss function J of the FCM2 algorithm2As shown in formula (21), for vikCalculating a partial derivative, and performing v using the derived formula (22)ikIterative updating of (1), sampling projection gradient descent method pair uijCarrying out iteration:
Figure BDA0002912022900000081
Figure BDA0002912022900000082
in the formula, lambda is more than or equal to 0 and is an orthogonal constraint parameter, the second term of the function is to add approximate orthogonal constraint to the membership matrix to balance the influence of unbalanced data, wherein uipThe membership degree of the ith sample point belonging to the class cluster p is represented by an original membership function uip' pass through
Figure BDA0002912022900000083
Transforming to obtain;
step 5.7, calculate J2A value of, if
Figure BDA0002912022900000084
FCM2 clustering is complete and iteration terminates, otherwise k is k +1, go to step 5.6 and continue iteration until it is satisfied
Figure BDA0002912022900000085
Load(s)Determining the clustering number and the clustering center of the FCM algorithm in the feature library;
and 5.8, inputting a sample to be detected, determining an initial clustering center by the load feature library, obtaining a membership function through FCM clustering, taking and outputting membership vectors corresponding to the sample to be identified, sequencing the membership, taking a category corresponding to the maximum membership, and recording as the belonged load category of the sample to be identified. If the sample identification results of the FCM1 and the FCM2 are consistent, the algorithm identification is finished, the category corresponding to the maximum membership degree is the final sample identification result, and if the sample identification results are not consistent, the step 6 is executed.
Step 6 is implemented according to the following steps:
6.1, combining samples with consistent FCM identification results to serve as a training set to train the improved GRNN neural network;
step 6.2, setting parameters: the number of neurons in the input layer is a characteristic dimension w; the number n of samples in a sample set of which the pattern layer neurons are a training set, and the transfer function of the ith neuron is shown as a formula (23); the summation layer comprises two types of neurons, one is to perform arithmetic summation on the outputs of all the mode layer neurons, and the transfer function of the sum is shown as a formula (24), and the other is to perform weighted summation on the outputs of all the mode layer neurons, and the transfer function of the sum is shown as a formula (25); the number of neurons in the output layer is equal to the dimension of an output vector in a training sample, namely the number b of clustered electric appliance categories, and the network output of the electric appliance is the division of two neurons, as shown in a formula (26);
Figure BDA0002912022900000091
Figure BDA0002912022900000092
Figure BDA0002912022900000093
Figure BDA0002912022900000094
wherein X is an input sample, XiFor the ith sample data of the training set, σ is the smoothing factor, YijFor the ith output sample YiThe jth element in (a);
6.3, the network parameters of GRNN only have smooth factors, so that an SA-BAS algorithm is used for searching the optimal smooth factor sigma; firstly, inputting training data, and initializing parameters: longicorn initial position X 01, the temperature T is 200, the step factor alpha is 0.95, the maximum iteration number N is 200, the annealing cycle number L is 100, the distance d from the center of mass to the whisker is 1.5, and the counter h is 1;
step 6.4, setting the step length S of the longicorn to be T, and meanwhile, enabling a counter T to be 1;
6.5, updating the left and right beard positions of the longicorn according to a formula (27), and updating the position X of the longicorn according to a formula (28)t+1Calculating the probability according to Metropolis criterion, and judging whether to accept Xt+1As a new solution, the step S is updated as shown in equation (29)t＝αSt-1Updating the distance from the center of mass to the tentacle according to the formula (30), and judging whether the counter t is more than or equal to the annealing cycle number L, if so, executing the step 6.6, otherwise, if t is t +1, and executing the step 6.5 in a circulating manner;
Figure BDA0002912022900000101
Figure BDA0002912022900000102
Figure BDA0002912022900000103
Figure BDA0002912022900000104
in the formula (I), the compound is shown in the specification,
Figure BDA0002912022900000105
is a random unit vector, the direction of the right hair pointing to the left hair is taken as the orientation of the longicorn in the space, sign () is a sign function, the value is greater than zero, 1 is taken, is less than zero, 1 is taken, is equal to zero, 0 is taken, T is the current temperature, and delta T is f (X)t+1)-f(Xt)，
Figure BDA0002912022900000106
As an expression for the fitness function, ωjAs a feature weight, xijIn order to train the sample network output values,
Figure BDA0002912022900000107
expected output values for the training samples;
and 6.6, updating the adaptive factor beta according to the formula (31), wherein T is beta T, judging whether the counter h is more than or equal to the iteration number N, if so, outputting the optimal longicorn position as the optimal smooth factor, and finishing the training, otherwise, h is h +1, and executing the step 6.4.
Figure BDA0002912022900000108
In the formula (f)hFor the current fitness value, fminThe historical optimal fitness value is obtained;
and (3) building a GRNN network by combining the trained smooth factors with the conditions of the step 6.2, so that deep re-identification of the sample to be detected is performed, and final identification of the non-invasive household appliance load is completed.
The invention has the beneficial effects that: the invention discloses a non-invasive household appliance load depth re-identification method, which solves the problem of low load identification accuracy caused by a single algorithm in the prior art. The event detection algorithm combining GLR and F detection is used for locating the occurrence time of the event, the algorithm combining semi-supervised Relief-F and mRMR is used for screening out the features with high correlation with the electric label, the combination of self-adaptive FCM1 and FCM2 is used for carrying out primary judgment, the improved GRNN is used for carrying out secondary judgment, the misjudgment and feature redundancy are reduced, and the method has the advantages of high identification accuracy and high identification speed.
Drawings
FIG. 1 is a flow chart of a non-invasive appliance load depth re-identification method of the present invention;
FIG. 2 is a flow chart of event detection based on GLR and F tests for a non-invasive method for deep re-identification of appliance load according to the present invention;
FIG. 3 is a flow chart of feature selection based on semi-supervised Relief-F and mRmR of the non-invasive household appliance load depth re-identification method of the present invention;
fig. 4 is a load identification flow chart of a non-invasive household appliance load deep re-identification method based on one-time identification of an adaptive FCM1 algorithm according to the present invention;
FIG. 5 is a flowchart of load identification based on the adaptive FCM2 algorithm for one-time identification of a non-invasive household appliance load deep re-identification method according to the present invention;
fig. 6 is a load identification flow chart based on the improved secondary discrimination of the GRNN neural network of the non-invasive household electrical appliance load depth re-identification method of the present invention.
Detailed Description
The present invention will be described in detail below with reference to the accompanying drawings and specific embodiments.
The invention discloses a non-invasive household appliance load depth re-identification method, which is implemented according to the following steps as shown in figure 1:
step 1, acquiring user high-frequency data, and performing denoising processing on the acquired user high-frequency data;
the step 1 is implemented according to the following steps:
step 1.1, acquiring high-frequency household appliance load data containing electrical parameters, specifically representing voltage, current and corresponding power;
step 1.2, denoising processing of the power signal, wherein an isolated noise point is easily identified as an event by an event detection algorithm, so a median filtering method is selected to process the original power signal, and the noise is eliminated without changingEdge information: suppose that there is a sequence x of digital signalsj(-∞<j<+ ∞) is filtered, a window with an odd number L is first defined, where L is 2N +1, N is a positive integer, and it is assumed that at a certain time i, the signal sample in the window is xi-N,…,xi,…,xi+NWherein x isiIs the sample value of the signal located in the center of the window, and after rearranging the L signals from small to large, the value is defined as the output value of the median filter.
Step 2, carrying out event detection on the data in the step 1 through improved generalized likelihood ratio test, if an event is detected, executing the step 3, otherwise, returning to the step 1;
as shown in fig. 2, step 2 is specifically implemented according to the following steps:
step 2.1, calculating the active power P of the fundamental wave according to the formula (1)1Taking the combined active power P as a two-dimensional power time sequence
Figure BDA0002912022900000121
A binary hypothesis test is proposed according to the formula (2);
Figure BDA0002912022900000122
Figure BDA0002912022900000123
in the formula, V1Is a fundamental voltage, I1Is the current of the fundamental wave,
Figure BDA0002912022900000124
is the phase difference between the two. n iscFor the time of occurrence of the change point, k is the total length of the window, n is the last sample time in the window, μ0,∑0Testing for hypotheses H0Mean of Gaussian distribution, covariance matrix, μ under the conditionsa,∑aIs H1Multi-dimensional signal mean, multi-dimensional covariance matrix, mu, before occurrence of change point under conditionb,∑bIs H1A multidimensional signal mean and a multidimensional covariance matrix after the variable point occurs under the condition;
step 2.2, defining two consecutive windows W within this time sequenceaAnd WbThe sample in the two windows is Xn＝{xmWhere m is n-k +1, …, n, and the window lengths are both k/2, the μ sum Σ in each of the two windows is calculated from equations (3) and (4), and then the decision function g is calculated from equation (5)n；
Figure BDA0002912022900000131
Figure BDA0002912022900000132
Figure BDA0002912022900000133
Step 2.3, mixing gnAnd a threshold value h1And comparing and searching suspicious points of event occurrence: when the decision function value is larger than h1When it is, refuse H0The data distribution in the two windows is not consistent, and at the point-changing time ncThere is a possibility of an event occurring; when the decision function is less than h1When it is, refuse H1And the data distribution of the two windows is consistent, and no event occurs. And when the sampling points of the detected multiple events are within 3 continuously or at intervals, the sampling points are regarded as having the maximum gnA corresponding event occurs;
and 2.4, taking the event point as a base point, and carrying out F-test screening on the false detection events. First, it is assumed that the variance values of the window data before and after the base point are equal:
Figure BDA0002912022900000134
calculating the value of statistic F according to equation (6), giving significance level α, rejecting hypothesis H if the value of F satisfies equation (7)0Judging that the two maternal variances have obvious difference, and judging that an event occurs at the point;
Figure BDA0002912022900000135
Figure BDA0002912022900000136
or
Figure BDA0002912022900000137
In the formula (I), the compound is shown in the specification,
Figure BDA0002912022900000138
respectively the sample volume and variance of the window before and after the base point.
Step 3, extracting multi-dimensional load characteristics for the detected event change points;
step 3 is specifically carried out as follows: the method for extracting the power characteristics of the variable points specifically comprises the following steps: active power, fundamental active power, reactive power, fundamental reactive power, apparent power, distortion power, power factor angle, fundamental power factor; extracting harmonic characteristics at a variable point, specifically including voltage, amplitude of each harmonic of one to nine times of voltage, content rate of each harmonic, difference of content rates of each harmonic and total harmonic distortion rate; the current waveform characteristics comprise the wave peak value, the average value and the wave crest coefficient; the method for extracting the V-I track characteristics at the variable points specifically comprises the following steps: symmetry, surrounding direction, surrounding area, number of intersection points, Y-axis intercept, Y-axis span, midline curvature, trace middle part peak value, left and right part area, middle part shape, and instantaneous admittance standard deviation. Step 4, aiming at the multi-dimensional features extracted in the step 3, screening relevant features by using a semi-supervised algorithm combining a Relief-F and a mRMR;
as shown in fig. 3, step 4 is specifically implemented according to the following steps:
step 4.1, labeling part of multi-dimensional characteristic samples, and recording as marked samples S1And the remainder is marked as the unmarked sample S2The label is C, and the multidimensional feature set is A ═ A1,A2,…,AN)；
Step 4.2, Slave S1In which a sample s of class C is randomly drawnq(CqE.g. C); from S2D neighboring samples are selected and recorded as
Figure BDA0002912022900000141
At S1From C to CqEach class C of the otherpRespectively solving a nearest neighbor sample x of s in the epsilon C (p is not equal to q); and at S2D neighbor samples of medium solution x, noted
Figure BDA0002912022900000142
Wherein the neighbor formula is shown as formula (8); updating all feature weights based on equation (9);
Figure BDA0002912022900000143
Figure BDA0002912022900000144
Figure BDA0002912022900000145
in the formula, M represents the iteration number, d is the number of the adjacent samples,
Figure BDA0002912022900000146
for the t-th neighbor sample in the q class to which the sample s belongs,
Figure BDA0002912022900000147
representing the t-th neighbouring sample in the P-th class, which is different from the class of the sample s, P (C)p) Representing the probability of the p-th class of objects, AkAs a result of the k-th feature,
Figure BDA0002912022900000148
representing a sample s and a sample
Figure BDA0002912022900000149
With respect to feature AkThe distance of (d);
step 4.3, circularly executing the step 4.2 for M times, and obtaining the characteristic weight omega finally output by the semi-supervised Relief-Fk'(k is 1,2, …, w'), removing features with weight coefficient less than theta to obtain a candidate feature subset, and then performing semi-supervised mRmR feature selection on the candidate feature subset;
step 4.4, marking the sample S according to the formula (11)1Calculating the correlation degree between each feature and the sample label, and obtaining the unmarked sample S according to the formula (12)2Calculating the redundancy of each characteristic;
Figure BDA0002912022900000151
Figure BDA0002912022900000152
in the formula (I), the compound is shown in the specification,
Figure BDA0002912022900000153
is shown in the labeled sample S1In (1) calculation ofkMutual information with label C, R (A)k,Sm-1) Representing a subset S of existing featuresm-1Which contains m-1 features, not Sm-1Features A in the subsetkRedundancy with selected features;
step 4.5, establishing a characteristic candidate set H, and selecting the maximum correlation degree DmaxCorresponding features as candidate set leader H1Sequentially selecting the kth feature A according to equation (13)kPutting the obtained product into H until the designated characteristic number w is selected;
Figure BDA0002912022900000154
and 4.6, calculating each characteristic weight according to the formula (14).
Figure BDA0002912022900000155
Step 5, taking the characteristics obtained in the step 4 as load marks, establishing a load characteristic library through a self-adaptive FCM algorithm, identifying the load working state of the household appliances in the user by using the maximum membership degree, finishing the identification if the identification results of the two FCM algorithms are consistent, and executing the step 6 if the identification results are inconsistent;
as shown in fig. 4 and 5, step 5 is specifically implemented according to the following steps:
step 5.1, setting the initial clustering number b, determining an initial clustering center based on a maximum-minimum distance algorithm, and simultaneously, enabling a counter k to return to zero;
a. calculating the mean of the data set
Figure BDA0002912022900000161
The sample point farthest from the mean is V1j(j ═ 1,2, …, w); b. calculating the minimum distance D between each data point and the selected cluster center according to the formula (15)xSelection of DxThe maximum value point is used as a new clustering center; c. repeating the step b until b initial clustering centers are selected;
Dx＝min d(xi,Z'k) k'＝1,…,kselected (15)
Figure BDA0002912022900000162
step 5.2, constructing a loss function L of the FCM1 algorithm1As shown in formula (17), for uij，λjAnd vikCalculating and combining the partial derivatives, and performing u by using the derived formula (18) and formula (19)ijAnd vikThe iterative update of (2);
Figure BDA0002912022900000163
Figure BDA0002912022900000164
Figure BDA0002912022900000165
wherein n is the number of samples, b is the number of clusters, w is the number of features,
Figure BDA0002912022900000166
is a sample weight, uijThe membership degree of the jth sample belonging to the ith class, v is a clustering center matrix, TuE (0, infinity) replaces the fuzzy index m for controlling the entropy value, and the maximum entropy is introduced as regularization;
step 5.3, calculating J1A value of, if
Figure BDA0002912022900000167
Finishing clustering of the initial clustering number b, terminating iteration, otherwise, turning to the step 5.2 to continue iteration until the requirement is met;
step 5.4, calculating the clustering effectiveness index of the initial clustering number b according to the formula (20), namely the contour coefficient SbIf it satisfies Sb-2＜Sb-1And Sb-1＞SbIf yes, the self-adaptive FCM1 classification is ended, otherwise, whether the clustering number b reaches the maximum value or not is judged
Figure BDA0002912022900000168
If the maximum value is reached, all S are takenbIf b, u and v corresponding to the medium maximum value do not reach the maximum value, b is equal to b +1, and the step 5.1 is executed again;
Figure BDA0002912022900000171
wherein n is the number of samples, a (i) is the average distance between the sample point i and the remaining sample objects in the same cluster, b (i) is the minimum value of the average distance between the sample point i and the sample objects in each of the remaining clusters, SbIn [ -1,1 [)]In between, the larger the S value, the more polymerizedThe better the class effect;
step 5.5, setting the optimal clustering number b of the FCM1 as the clustering number of the FCM2, similarly adopting a maximum-minimum distance algorithm to determine an initial clustering center, and simultaneously, enabling a counter k to return to zero;
step 5.6, constructing a loss function J of the FCM2 algorithm2As shown in formula (21), for vikCalculating a partial derivative, and performing v using the derived formula (22)ikIterative updating of (1), sampling projection gradient descent method pair uijCarrying out iteration:
Figure BDA0002912022900000172
Figure BDA0002912022900000173
in the formula, lambda is more than or equal to 0 and is an orthogonal constraint parameter, the second term of the function is to add approximate orthogonal constraint to the membership matrix to balance the influence of unbalanced data, wherein uipThe membership degree of the ith sample point belonging to the class cluster p is represented by an original membership function uip' pass through
Figure BDA0002912022900000174
Transforming to obtain;
step 5.7, calculate J2A value of, if
Figure BDA0002912022900000175
FCM2 clustering is complete and iteration terminates, otherwise k is k +1, go to step 5.6 and continue iteration until it is satisfied
Figure BDA0002912022900000176
Determining the clustering number and the clustering center of the FCM algorithm in the load characteristic library;
and 5.8, inputting a sample to be detected, determining an initial clustering center by the load feature library, obtaining a membership function through FCM clustering, taking and outputting membership vectors corresponding to the sample to be identified, sequencing the membership, taking a category corresponding to the maximum membership, and recording as the belonged load category of the sample to be identified. If the sample identification results of the FCM1 and the FCM2 are consistent, the algorithm identification is finished, the category corresponding to the maximum membership degree is the final sample identification result, and if the sample identification results are not consistent, the step 6 is executed.
And 6, training the parameter smoothing factor of the GRNN by adopting an SA-BAS algorithm, and performing secondary identification on the sample to be detected by using the trained GRNN to finish secondary identification of the non-invasive household appliance load.
As shown in fig. 6, step 6 is specifically implemented according to the following steps:
6.1, combining samples with consistent FCM identification results to serve as a training set to train the improved GRNN neural network;
step 6.2, setting parameters: the number of neurons in the input layer is a characteristic dimension w; the number n of samples in a sample set of which the pattern layer neurons are a training set, and the transfer function of the ith neuron is shown as a formula (23); the summation layer comprises two types of neurons, one is to perform arithmetic summation on the outputs of all the mode layer neurons, and the transfer function of the sum is shown as a formula (24), and the other is to perform weighted summation on the outputs of all the mode layer neurons, and the transfer function of the sum is shown as a formula (25); the number of neurons in the output layer is equal to the dimension of an output vector in a training sample, namely the number b of clustered electric appliance categories, and the network output of the electric appliance is the division of two neurons, as shown in a formula (26);
Figure BDA0002912022900000181
Figure BDA0002912022900000182
Figure BDA0002912022900000183
Figure BDA0002912022900000184
wherein X is an input sample, XiFor the ith sample data of the training set, σ is the smoothing factor, YijFor the ith output sample YiThe jth element in (a);
6.3, the network parameters of GRNN only have smooth factors, so that an SA-BAS algorithm is used for searching the optimal smooth factor sigma; firstly, inputting training data, and initializing parameters: longicorn initial position X 01, the temperature T is 200, the step factor alpha is 0.95, the maximum iteration number N is 200, the annealing cycle number L is 100, the distance d from the center of mass to the whisker is 1.5, and the counter h is 1;
step 6.4, setting the step length S of the longicorn to be T, and meanwhile, enabling a counter T to be 1;
6.5, updating the left and right beard positions of the longicorn according to a formula (27), and updating the position X of the longicorn according to a formula (28)t+1Calculating the probability according to Metropolis criterion, and judging whether to accept Xt+1As a new solution, the step S is updated as shown in equation (29)t＝αSt-1Updating the distance from the center of mass to the tentacle according to the formula (30), and judging whether the counter t is more than or equal to the annealing cycle number L, if so, executing the step 6.6, otherwise, if t is t +1, and executing the step 6.5 in a circulating manner;
Figure BDA0002912022900000191
Figure BDA0002912022900000192
Figure BDA0002912022900000193
Figure BDA0002912022900000194
in the formula (I), the compound is shown in the specification,
Figure BDA0002912022900000195
is a random unit vector, the direction of the right hair pointing to the left hair is taken as the orientation of the longicorn in the space, sign () is a sign function, the value is greater than zero, 1 is taken, is less than zero, 1 is taken, is equal to zero, 0 is taken, T is the current temperature, and delta T is f (X)t+1)-f(Xt)，
Figure BDA0002912022900000196
As an expression for the fitness function, ωjAs a feature weight, xijIn order to train the sample network output values,
Figure BDA0002912022900000197
expected output values for the training samples;
and 6.6, updating the adaptive factor beta according to the formula (31), wherein T is beta T, judging whether the counter h is more than or equal to the iteration number N, if so, outputting the optimal longicorn position as the optimal smooth factor, and finishing the training, otherwise, h is h +1, and executing the step 6.4.
Figure BDA0002912022900000201
In the formula (f)hFor the current fitness value, fminThe historical optimal fitness value is obtained;
and (3) building a GRNN network by combining the trained smooth factors with the conditions of the step 6.2, so that the depth re-identification of the sample to be detected is carried out, and the non-invasive secondary identification of the load of the household appliance is completed.
The invention relates to a non-invasive household appliance load depth re-identification method, which is characterized in that the occurrence time of an event is positioned through an event detection algorithm combining GLR and F detection, a semi-supervised algorithm combining Relief-F and mRmR is used for screening out characteristics with high correlation with an electric label, primary judgment is carried out through the combination of self-adaptive FCM1 and FCM2, final judgment is carried out through improved GRNN, misjudgment and characteristic redundancy are reduced, and compared with other literature methods, the result is high in identification accuracy and high in identification speed.

Claims (7)

1. A non-invasive household appliance load depth re-identification method is characterized by comprising the following steps:
step 1, acquiring user high-frequency data, and performing denoising processing on the acquired user high-frequency data;
step 2, carrying out event detection on the data in the step 1 through improved generalized likelihood ratio test, if an event is detected, executing the step 3, otherwise, returning to the step 1;
step 3, extracting multi-dimensional load characteristics for the detected event change points;
step 4, aiming at the multi-dimensional features extracted in the step 3, screening relevant features by using a semi-supervised algorithm combining a Relief-F and a mRMR;
step 5, taking the characteristics obtained in the step 4 as load marks, establishing a load characteristic library through a self-adaptive FCM algorithm, identifying the load working state of the household appliances in the user by using the maximum membership degree, finishing the identification if the identification results of the two FCM algorithms are consistent, and executing the step 6 if the identification results are inconsistent;
and 6, training the parameter smoothing factor of the GRNN by adopting an SA-BAS algorithm, and performing secondary identification on a sample to be detected by using the trained GRNN to complete non-invasive household appliance load depth re-identification.
2. The method for non-invasive home appliance load depth re-identification according to claim 1, wherein the step 1 is implemented by the following steps:
step 1.1, acquiring high-frequency household appliance load data containing electrical parameters, specifically representing voltage, current and corresponding power;
step 1.2, denoising processing of the power signal, wherein because an isolated noise point is easily identified as an event by an event detection algorithm, a median filtering method is selected to process the original power signal, so that the noise is eliminated without changing edge information: suppose that there is a sequence x of digital signalsj(-∞<j<+ ∞) is filtered, a window with the length of an odd number L is firstly defined,l is 2N +1, N is a positive integer, and it is assumed that at a certain time i, the signal sample in the window is xi-N,…,xi,…,xi+NWherein x isiIs the sample value of the signal located in the center of the window, and after rearranging the L signals from small to large, the value is defined as the output value of the median filter.
3. The method for non-invasive depth re-identification of appliance load according to claim 1, wherein the step 2 is implemented by the following steps:
step 2.1, calculating the active power P of the fundamental wave according to the formula (1)1Taking the combined active power P as a two-dimensional power time sequence
Figure FDA0002912022890000021
A binary hypothesis test is proposed according to the formula (2);
Figure FDA0002912022890000022
Figure FDA0002912022890000023
in the formula, V1Is a fundamental voltage, I1Is the current of the fundamental wave,
Figure FDA0002912022890000024
is the phase difference between the two. n iscFor the time of occurrence of the change point, k is the total length of the window, n is the last sample time in the window, μ0,∑0Testing for hypotheses H0Mean of Gaussian distribution, covariance matrix, μ under the conditionsa,∑aIs H1Multi-dimensional signal mean, multi-dimensional covariance matrix, mu, before occurrence of change point under conditionb,∑bIs H1A multidimensional signal mean and a multidimensional covariance matrix after the variable point occurs under the condition;
step 2.2,Two successive windows W are defined in the time seriesaAnd WbThe sample in the two windows is Xn＝{xmWhere m is n-k +1, …, n, and the window lengths are both k/2, the μ sum Σ in each of the two windows is calculated from equations (3) and (4), and then the decision function g is calculated from equation (5)n；
Figure FDA0002912022890000025
Figure FDA0002912022890000031
Figure FDA0002912022890000032
Step 2.3, mixing gnAnd a threshold value h1And comparing and searching suspicious points of event occurrence: when the decision function value is larger than h1When it is, refuse H0The data distribution in the two windows is not consistent, and at the point-changing time ncThere is a possibility of an event occurring; when the decision function is less than h1When it is, refuse H1And the data distribution of the two windows is consistent, and no event occurs. And when the sampling points of the detected multiple events are within 3 continuously or at intervals, the sampling points are regarded as having the maximum gnA corresponding event occurs;
and 2.4, taking the event point as a base point, and carrying out F-test screening on the false detection events. First, it is assumed that the variance values of the window data before and after the base point are equal:
Figure FDA0002912022890000033
calculating the value of statistic F according to equation (6), giving significance level α, rejecting hypothesis H if the value of F satisfies equation (7)0Judging that the two maternal variances have obvious difference, and judging that an event occurs at the point;
Figure FDA0002912022890000034
Figure FDA0002912022890000035
in the formula, n1,
Figure FDA0002912022890000036
n2,
Figure FDA0002912022890000037
Respectively the sample volume and variance of the window before and after the base point.
4. The method for non-invasive depth re-identification of appliance load according to claim 1, wherein the step 3 is implemented as follows: the method for extracting the power characteristics of the variable points specifically comprises the following steps: active power, fundamental active power, reactive power, fundamental reactive power, apparent power, distortion power, power factor angle, fundamental power factor; extracting harmonic characteristics at a variable point, specifically including voltage, amplitude of each harmonic of one to nine times of voltage, content rate of each harmonic, difference of content rates of each harmonic and total harmonic distortion rate; the current waveform characteristics comprise the wave peak value, the average value and the wave crest coefficient; the method for extracting the V-I track characteristics at the variable points specifically comprises the following steps: symmetry, surrounding direction, surrounding area, number of intersection points, Y-axis intercept, Y-axis span, midline curvature, trace middle part peak value, left and right part area, middle part shape, and instantaneous admittance standard deviation.
5. The method for non-invasive depth re-identification of appliance load according to claim 1, wherein the step 4 is implemented by the following steps:
step 4.1, labeling part of multi-dimensional characteristic samples, and recording as marked samples S1And the remainder is marked as the unmarked sample S2The label is C, multidimensionalThe syndrome is A ═ A1,A2,…,AN)；
Step 4.2, Slave S1In which a sample s of class C is randomly drawnq(CqE.g. C); from S2D neighboring samples are selected and recorded as
Figure FDA0002912022890000041
At S1From C to CqEach class C of the otherpRespectively solving a nearest neighbor sample x of s in the epsilon C (p is not equal to q); and at S2D neighbor samples of medium solution x, noted
Figure FDA0002912022890000042
Wherein the neighbor formula is shown as formula (8); updating all feature weights based on equation (9);
Figure FDA0002912022890000043
Figure FDA0002912022890000044
Figure FDA0002912022890000045
in the formula, M represents the iteration number, d is the number of the adjacent samples,
Figure FDA0002912022890000046
for the t-th neighbor sample in the q class to which the sample s belongs,
Figure FDA0002912022890000047
representing the t-th neighbouring sample in the P-th class, which is different from the class of the sample s, P (C)p) Representing the probability of the p-th class of objects, AkAs a result of the k-th feature,
Figure FDA0002912022890000048
representing a sample s and a sample
Figure FDA0002912022890000049
With respect to feature AkThe distance of (d);
step 4.3, circularly executing the step 4.2 for M times, and obtaining the characteristic weight omega finally output by the semi-supervised Relief-Fk'(k is 1,2, …, w'), removing features with weight coefficient less than theta to obtain a candidate feature subset, and then performing semi-supervised mRmR feature selection on the candidate feature subset;
step 4.4, marking the sample S according to the formula (11)1Calculating the correlation degree between each feature and the sample label, and obtaining the unmarked sample S according to the formula (12)2Calculating the redundancy of each characteristic;
Figure FDA0002912022890000051
Figure FDA0002912022890000052
in the formula (I), the compound is shown in the specification,
Figure FDA0002912022890000053
is shown in the labeled sample S1In (1) calculation ofkMutual information with label C, R (A)k,Sm-1) Representing a subset S of existing featuresm-1Which contains m-1 features, not Sm-1Features A in the subsetkRedundancy with selected features;
step 4.5, establishing a characteristic candidate set H, and selecting the maximum correlation degree DmaxCorresponding features as candidate set leader H1Sequentially selecting the kth feature A according to equation (13)kPutting the obtained product into H until the designated characteristic number w is selected;
Figure FDA0002912022890000054
and 4.6, calculating each characteristic weight according to the formula (14).
Figure FDA0002912022890000055
6. The method for non-invasive depth re-identification of appliance load according to claim 1, wherein the step 5 is implemented by the following steps:
step 5.1, setting the initial clustering number b, determining an initial clustering center based on a maximum-minimum distance algorithm, and simultaneously, enabling a counter k to return to zero;
a. calculating the mean x of the data set, the sample point farthest from the mean being V1j(j＝1,2,…,w)；
b. Calculating the minimum distance D between each data point and the selected cluster center according to the formula (15)xSelection of DxThe maximum value point is used as a new clustering center; c. repeating the step b until b initial clustering centers are selected;
Dx＝min d(xi,Z'k)k'＝1,…,kselected (15)
Figure FDA0002912022890000061
step 5.2, constructing a loss function L of the FCM1 algorithm1As shown in formula (17), for uij，λjAnd vikCalculating and combining the partial derivatives, and performing u by using the derived formula (18) and formula (19)ijAnd vikThe iterative update of (2);
Figure FDA0002912022890000062
Figure FDA0002912022890000063
Figure FDA0002912022890000064
wherein n is the number of samples, b is the number of clusters, w is the number of features,
Figure FDA0002912022890000065
is a sample weight, uijThe membership degree of the jth sample belonging to the ith class, v is a clustering center matrix, TuE (0, infinity) replaces the fuzzy index m for controlling the entropy value, and the maximum entropy is introduced as regularization;
step 5.3, calculating J1A value of, if
Figure FDA0002912022890000066
Finishing clustering of the initial clustering number b, terminating iteration, otherwise, turning to the step 5.2 to continue iteration until the requirement is met;
step 5.4, calculating the clustering effectiveness index of the initial clustering number b according to the formula (20), namely the contour coefficient SbIf it satisfies Sb-2＜Sb-1And Sb-1＞SbIf yes, the self-adaptive FCM1 classification is ended, otherwise, whether the clustering number b reaches the maximum value or not is judged
Figure FDA0002912022890000067
If the maximum value is reached, all S are takenbIf b, u and v corresponding to the medium maximum value do not reach the maximum value, b is equal to b +1, and the step 5.1 is executed again;
Figure FDA0002912022890000068
wherein n is the number of samples, and a (i) is the sample point i and the remaining sample pairs in the same clusterAverage distance of the image, b (i) minimum of average distance of sample point i from sample object of each remaining cluster, SbIn [ -1,1 [)]The larger the S value is, the better the clustering effect is;
step 5.5, setting the optimal clustering number b of the FCM1 as the clustering number of the FCM2, similarly adopting a maximum-minimum distance algorithm to determine an initial clustering center, and simultaneously, enabling a counter k to return to zero;
step 5.6, constructing a loss function J of the FCM2 algorithm2As shown in formula (21), for vikCalculating a partial derivative, and performing v using the derived formula (22)ikIterative updating of (1), sampling projection gradient descent method pair uijCarrying out iteration:
Figure FDA0002912022890000071
Figure FDA0002912022890000072
in the formula, lambda is more than or equal to 0 and is an orthogonal constraint parameter, the second term of the function is to add approximate orthogonal constraint to the membership matrix to balance the influence of unbalanced data, wherein uipThe membership degree of the ith sample point belonging to the class cluster p is represented by an original membership function uip' pass through
Figure FDA0002912022890000073
Transforming to obtain;
step 5.7, calculate J2A value of, if
Figure FDA0002912022890000074
FCM2 clustering is complete and iteration terminates, otherwise k is k +1, go to step 5.6 and continue iteration until it is satisfied
Figure FDA0002912022890000075
Determining the clustering number and the clustering center of the FCM algorithm in the load characteristic library;
and 5.8, inputting a sample to be detected, determining an initial clustering center by the load feature library, obtaining a membership function through FCM clustering, taking and outputting membership vectors corresponding to the sample to be identified, sequencing the membership, taking a category corresponding to the maximum membership, and recording as the belonged load category of the sample to be identified. If the sample identification results of the FCM1 and the FCM2 are consistent, the algorithm identification is finished, the category corresponding to the maximum membership degree is the final sample identification result, and if the sample identification results are not consistent, the step 6 is executed.
7. The method for non-invasive depth re-identification of appliance load according to claim 1, wherein the step 6 is implemented by the following steps:
6.1, combining samples with consistent FCM identification results to serve as a training set to train the improved GRNN neural network;
step 6.2, setting parameters: the number of neurons in the input layer is a characteristic dimension w; the number n of samples in a sample set of which the pattern layer neurons are a training set, and the transfer function of the ith neuron is shown as a formula (23); the summation layer comprises two types of neurons, one is to perform arithmetic summation on the outputs of all the mode layer neurons, and the transfer function of the sum is shown as a formula (24), and the other is to perform weighted summation on the outputs of all the mode layer neurons, and the transfer function of the sum is shown as a formula (25); the number of neurons in the output layer is equal to the dimension of an output vector in a training sample, namely the number b of clustered electric appliance categories, and the network output of the electric appliance is the division of two neurons, as shown in a formula (26);
Figure FDA0002912022890000081
Figure FDA0002912022890000082
Figure FDA0002912022890000083
Figure FDA0002912022890000084
wherein X is an input sample, XiFor the ith sample data of the training set, σ is the smoothing factor, YijFor the ith output sample YiThe jth element in (a);
6.3, the network parameters of GRNN only have smooth factors, so that an SA-BAS algorithm is used for searching the optimal smooth factor sigma; firstly, inputting training data, and initializing parameters: longicorn initial position X01, the temperature T is 200, the step factor alpha is 0.95, the maximum iteration number N is 200, the annealing cycle number L is 100, the distance d from the center of mass to the whisker is 1.5, and the counter h is 1;
step 6.4, setting the step length S of the longicorn to be T, and meanwhile, enabling a counter T to be 1;
6.5, updating the left and right beard positions of the longicorn according to a formula (27), and updating the position X of the longicorn according to a formula (28)t+1Calculating the probability according to Metropolis criterion, and judging whether to accept Xt+1As a new solution, the step S is updated as shown in equation (29)t＝αSt-1Updating the distance from the center of mass to the tentacle according to the formula (30), and judging whether the counter t is more than or equal to the annealing cycle number L, if so, executing the step 6.6, otherwise, if t is t +1, and executing the step 6.5 in a circulating manner;
Figure FDA0002912022890000091
Figure FDA0002912022890000092
Figure FDA0002912022890000093
Figure FDA0002912022890000094
in the formula (I), the compound is shown in the specification,
Figure FDA0002912022890000095
is a random unit vector, the direction of the right hair pointing to the left hair is taken as the orientation of the longicorn in the space, sign () is a sign function, the value is greater than zero, 1 is taken, is less than zero, 1 is taken, is equal to zero, 0 is taken, T is the current temperature, and delta T is f (X)t+1)-f(Xt)，
Figure FDA0002912022890000096
As an expression for the fitness function, ωjAs a feature weight, xijIn order to train the sample network output values,
Figure FDA0002912022890000097
expected output values for the training samples;
and 6.6, updating the adaptive factor beta according to the formula (31), wherein T is beta T, judging whether the counter h is more than or equal to the iteration number N, if so, outputting the optimal longicorn position as the optimal smooth factor, and finishing the training, otherwise, h is h +1, and executing the step 6.4.
Figure FDA0002912022890000098
In the formula (f)hFor the current fitness value, fminThe historical optimal fitness value is obtained;
and (3) building a GRNN network by combining the trained smooth factors with the conditions of the step 6.2, so as to carry out deep re-identification on the sample to be detected and complete the final identification of the non-invasive household appliance load.
CN202110089840.XA 2021-01-22 2021-01-22 Non-invasive household appliance load depth re-identification method Active CN112821559B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110089840.XA CN112821559B (en) 2021-01-22 2021-01-22 Non-invasive household appliance load depth re-identification method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110089840.XA CN112821559B (en) 2021-01-22 2021-01-22 Non-invasive household appliance load depth re-identification method

Publications (2)

Publication Number Publication Date
CN112821559A true CN112821559A (en) 2021-05-18
CN112821559B CN112821559B (en) 2023-08-01

Family

ID=75858877

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110089840.XA Active CN112821559B (en) 2021-01-22 2021-01-22 Non-invasive household appliance load depth re-identification method

Country Status (1)

Country Link
CN (1) CN112821559B (en)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113595242A (en) * 2021-06-28 2021-11-02 湖南大学 Non-invasive load identification method based on deep CNN-HMM
CN113723479A (en) * 2021-08-18 2021-11-30 南京工程学院 Non-invasive load identification method based on GRNN and mean shift algorithm
CN114325081A (en) * 2021-12-29 2022-04-12 润建股份有限公司 Non-invasive load identification method based on multi-modal characteristics
CN114384999A (en) * 2021-11-19 2022-04-22 福州大学 User irrelevant myoelectricity gesture recognition system based on self-adaptive learning
CN116893314A (en) * 2023-09-04 2023-10-17 国网浙江省电力有限公司余姚市供电公司 Non-invasive power load monitoring method, device, equipment and storage medium

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20160141032A (en) * 2015-05-27 2016-12-08 전자부품연구원 Non-Intrusive Appliance Load Monitoring Method using a Switching Factorial Hidden Markov Model and System applying the same
CN108846530A (en) * 2018-09-28 2018-11-20 国网上海市电力公司 One kind being based on the short-term load forecasting method of " cluster-recurrence " model
CN108898154A (en) * 2018-09-29 2018-11-27 华北电力大学 A kind of electric load SOM-FCM Hierarchical clustering methods
CN110555369A (en) * 2019-07-16 2019-12-10 浙江工业大学 MLCDTL-based non-intrusive load identification method
CN110866841A (en) * 2019-11-20 2020-03-06 江苏方天电力技术有限公司 Power consumer industry dimension power consumption pattern identification analysis method and system based on double clustering method
CN110991818A (en) * 2019-11-14 2020-04-10 广西电网有限责任公司电力科学研究院 Load identification method integrating event detection and neural network
CN112101110A (en) * 2020-08-13 2020-12-18 西安理工大学 Non-invasive load identification method for user side of power system

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20160141032A (en) * 2015-05-27 2016-12-08 전자부품연구원 Non-Intrusive Appliance Load Monitoring Method using a Switching Factorial Hidden Markov Model and System applying the same
CN108846530A (en) * 2018-09-28 2018-11-20 国网上海市电力公司 One kind being based on the short-term load forecasting method of " cluster-recurrence " model
CN108898154A (en) * 2018-09-29 2018-11-27 华北电力大学 A kind of electric load SOM-FCM Hierarchical clustering methods
CN110555369A (en) * 2019-07-16 2019-12-10 浙江工业大学 MLCDTL-based non-intrusive load identification method
CN110991818A (en) * 2019-11-14 2020-04-10 广西电网有限责任公司电力科学研究院 Load identification method integrating event detection and neural network
CN110866841A (en) * 2019-11-20 2020-03-06 江苏方天电力技术有限公司 Power consumer industry dimension power consumption pattern identification analysis method and system based on double clustering method
CN112101110A (en) * 2020-08-13 2020-12-18 西安理工大学 Non-invasive load identification method for user side of power system

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
PAULO RICARDO ZAMBELLI TAVEIRA ET AL.: "Non-Intrusive Identification of Loads by Random Forest and Fireworks Optimization", 《IEEE》 *
江帆 等: "基于广义回归神经网络的非侵入式负荷识别方法", 《电测与仪表》 *
汪繁荣 等: "基于聚类特征及seq2seq深度CNN的家电负荷识别方法研究", 《电测与仪表》 *
满蔚仕 等: "利用快速S变换的电能质量扰动识别方法", 《西安交通大学学报》 *

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113595242A (en) * 2021-06-28 2021-11-02 湖南大学 Non-invasive load identification method based on deep CNN-HMM
CN113595242B (en) * 2021-06-28 2024-02-20 湖南大学 Non-invasive load identification method based on depth CNN-HMM
CN113723479A (en) * 2021-08-18 2021-11-30 南京工程学院 Non-invasive load identification method based on GRNN and mean shift algorithm
CN114384999A (en) * 2021-11-19 2022-04-22 福州大学 User irrelevant myoelectricity gesture recognition system based on self-adaptive learning
CN114384999B (en) * 2021-11-19 2023-07-21 福州大学 User-independent myoelectric gesture recognition system based on self-adaptive learning
CN114325081A (en) * 2021-12-29 2022-04-12 润建股份有限公司 Non-invasive load identification method based on multi-modal characteristics
CN116893314A (en) * 2023-09-04 2023-10-17 国网浙江省电力有限公司余姚市供电公司 Non-invasive power load monitoring method, device, equipment and storage medium

Also Published As

Publication number Publication date
CN112821559B (en) 2023-08-01

Similar Documents

Publication Publication Date Title
CN112821559B (en) Non-invasive household appliance load depth re-identification method
Ahmed et al. Feature selection–based detection of covert cyber deception assaults in smart grid communications networks using machine learning
CN107101813B (en) A kind of frame-type circuit breaker mechanical breakdown degree assessment method based on vibration signal
CN109145516B (en) Analog circuit fault identification method based on improved extreme learning machine
CN109886464B (en) Low-information-loss short-term wind speed prediction method based on optimized singular value decomposition generated feature set
CN112732748B (en) Non-invasive household appliance load identification method based on self-adaptive feature selection
CN112560596B (en) Radar interference category identification method and system
CN113780684A (en) Intelligent building user energy consumption behavior prediction method based on LSTM neural network
CN111563827A (en) Load decomposition method based on electrical appliance physical characteristics and residential electricity consumption behaviors
CN112085111A (en) Load identification method and device
CN113641906A (en) System, method, device, processor and medium for realizing similar target person identification processing based on fund transaction relation data
CN115375921A (en) Two-stage non-intrusive load identification method and terminal
CN115079052A (en) Transformer fault diagnosis method and system
Wu et al. Weak NAS Predictor Is All You Need
CN115114484A (en) Abnormal event detection method and device, computer equipment and storage medium
Abbas et al. Volterra-system identification using adaptive real-coded genetic algorithm
Herath et al. Comprehensive analysis of convolutional neural network models for non-instructive load monitoring
CN117278314A (en) DDoS attack detection method
CN112422546A (en) Network anomaly detection method based on variable neighborhood algorithm and fuzzy clustering
Al-Khamees et al. Data Stream Clustering Using Fuzzy-based Evolving Cauchy Algorithm.
CN116527346A (en) Threat node perception method based on deep learning graph neural network theory
CN111797979A (en) Vibration transmission system based on LSTM model
CN116400168A (en) Power grid fault diagnosis method and system based on depth feature clustering
CN116318925A (en) Multi-CNN fusion intrusion detection method, system, medium, equipment and terminal
CN114124437B (en) Encrypted flow identification method based on prototype convolutional network

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
TA01 Transfer of patent application right
TA01 Transfer of patent application right

Effective date of registration: 20230706

Address after: 518000 1002, Building A, Zhiyun Industrial Park, No. 13, Huaxing Road, Henglang Community, Longhua District, Shenzhen, Guangdong Province

Applicant after: Shenzhen Wanzhida Technology Co.,Ltd.

Address before: 710048 Shaanxi province Xi'an Beilin District Jinhua Road No. 5

Applicant before: XI'AN University OF TECHNOLOGY

Effective date of registration: 20230706

Address after: 628, building 118, No.72 Guowei Road, Xianhu community, Liantang street, Luohu District, Shenzhen, Guangdong 518000

Applicant after: Wuxing Technology (Shenzhen) Co.,Ltd.

Address before: 518000 1002, Building A, Zhiyun Industrial Park, No. 13, Huaxing Road, Henglang Community, Longhua District, Shenzhen, Guangdong Province

Applicant before: Shenzhen Wanzhida Technology Co.,Ltd.

GR01 Patent grant
GR01 Patent grant