A non-intrusive load identification calculation method

Abstract

This invention discloses a non-intrusive load identification calculation method. To improve the accuracy of load identification, the method first collects data from the target user, obtaining voltage and current data of the electrical equipment. An improved two-dimensional generalized likelihood ratio test is used to perform event detection on the current data. Two false detection filtering algorithms, cluster removal and derivative removal, are used to reduce the occurrence of false detections. The abrupt change points of events are identified and stored in an event list. Reactive current is used to replace the total current in the traditional V-I trajectory. Steady-state data is extracted from the event list to generate a V-reactive-I trajectory diagram. This diagram is then fused with the active and reactive power obtained from the Fourier decomposition of the measured voltage and current to form an image feature matrix. Finally, the image feature matrix is ​​input into an improved fully convolutional network model incorporating residual blocks to output the identification result.

Description

Technical Field

[0001] This invention belongs to the field of electrical load identification technology, and relates to a non-intrusive load identification method based on V-reactive I trajectory fusion of active and reactive power. Background Technology

[0002] Non-Intrusive Load Monitoring (NILM) technology is a crucial component of smart grids and plays a significant role in the rational planning and utilization of electricity. Accurate understanding of users' electricity consumption helps the grid better allocate power resources and improves users' electricity usage habits. Due to the high accuracy of supervised learning, supervised learning-based NILM methods have become the mainstream in NILM technology. While VI trajectories, as a common load feature, offer better identification results compared to other single high-frequency features, they are plotted based on normalized voltage and current, and therefore cannot fundamentally represent power characteristics. Furthermore, they struggle to distinguish between household appliances with similar VI trajectories but significant power differences. In addition, existing technologies typically use one-dimensional generalized likelihood ratios for event detection before plotting VI trajectories, which can miss slight current fluctuations, leading to missing VI trajectories and affecting the final load identification accuracy. Summary of the Invention

[0003] To address the shortcomings of existing technologies, this invention proposes a non-intrusive load identification calculation method. Based on multi-dimensional event detection, it uses a dual-channel fully convolutional network to fuse active and reactive power as load features using V-reactive-I trajectory fusion, thus solving the problem of low recognition accuracy in traditional load identification methods.

[0004] A non-intrusive load identification calculation method, the specific steps of which are as follows:

[0005] Step 1: Data Collection

[0006] The current and voltage data from the main meter and sub-meters installed on the target object are sampled. The sampled data is then denoised using median filtering, and offset points in the original data are removed to ensure data accuracy, thus obtaining sample data.

[0007] Step 2: Event Detection

[0008] Event detection identifies the moment an event occurs, thus determining the necessary steady-state data segment for the VI trajectory. To improve accuracy, this method extends the one-dimensional generalized likelihood ratio test to multi-dimensional signal detection. The generalized likelihood ratio test algorithm based on the Neyman-Pearson criterion is used for event detection. A joint window is formed by two nearest-neighbor sliding windows (pre-event and post-event windows), sliding along the current matrix. An event is considered to have occurred when the algorithm's decision function exceeds a set threshold. The specific steps are as follows:

[0009] Step 2.1: Select the 30kHz current RMS value sequence from the sample data table and arrange it sequentially according to the acquisition time to form a one-dimensional circuit sequence K, which will serve as the data points for event detection. Set the sliding window length to k sampling points and the width to 2. Multiply each value in the current sequence K by a coefficient of 0.6 to obtain a new one-dimensional current sequence G, and align it with the beginning and end of K to form a two-dimensional sequence I with a width of 2 and a length of M. m Let m be the m-th element in the two-dimensional sequence I.

[0010] Step 2.2, Set the on and off states of electrical appliances to correspond to H1 and H0 respectively. (Assumption:)

[0011] H0 assumption: The two-dimensional current sequence does not contain any inflection points, that is:

[0012] H0:I={I m ,m=M-k+1,…,M}~N(μ0,Σ)

[0013] H1 hypothesis: The two-dimensional current sequence contains a change point, that is:

[0014]

[0015] Decision function g n for:

[0016]

[0017] in:

[0018]

[0019] n c Let μ0 be the moment when the change point occurs in the two-dimensional sequence I, and let μ0 and ∑0 be the mean and covariance matrix of the Gaussian distribution under the hypothesis test H0, respectively; 1a ,∑ 1a The time n of the change point under hypothesis test H1 is respectively. c Previous mean and covariance matrices; μ 1b ,∑ 1b The time n of the change point under hypothesis test H1 is respectively.c The mean and covariance matrices follow.

[0020] Set a decision threshold h, when the decision function value g... n When the value is greater than or equal to the decision threshold h, the H1 test is valid, meaning the sequence contains a change point and an event has occurred. The time n when the change point occurs is recorded. c When the decision function value g n When the current is less than the decision threshold h, test H0 is valid, meaning there are no changing points in the sequence and no event has occurred. The decision threshold h is dynamically adjusted continuously for different electrical appliances to achieve better event detection results. Finally, the current information at times when no event has occurred is set to 0, while the current information at times when an event has occurred remains unchanged, resulting in sequence I. event .

[0021] Step 2.3: Existing event detection algorithms rely on threshold settings and struggle to balance recall and accuracy. There are no particularly effective solutions for false positives caused by long transient events resulting from equipment switching events, i.e., continuous and irregular power changes. Furthermore, fluctuations in the steady-state of some nonlinear equipment can also easily lead to false positives. Therefore, based on the existing generalized likelihood ratio transient event detection algorithm, two false positive filtering algorithms—cluster removal and derivative removal—were used to filter out false positive events.

[0022] in I event One-dimensional unsupervised clustering is performed on the sequence without setting the number of cluster centers. All consecutive non-zero values ​​are aggregated into a cluster and identified as an event, thereby removing false detections caused by small fluctuations in current under the same operating conditions, and obtaining the sequence ΔI.

[0023] Step 2.4: Cluster analysis can remove false detections of consecutive pseudo-events caused by long switching transients. These false detections occur because when the device changes its operating state, its electrical performance changes, further causing a shift in the aggregated measurement value. If the jump variable exceeds a threshold, the basic detection algorithm will falsely detect consecutive changes as events. However, for devices with long transients and a certain time interval between two detected pseudo-events, cluster analysis cannot remove them. Therefore, derivative removal is needed to further reduce false detections with time intervals. For devices with long switching transients, the derivative during the transient duration often deviates from the zero line due to the oscillation of the current waveform, and may cross the zero line several times before deviating in another direction. Based on this characteristic, the detected events can be further filtered by judging the sign of the derivative between adjacent detected events. For the current sequence ΔI obtained in Step 2.3, any current point ΔI (not an endpoint) f The first derivative can be calculated by the following formula:

[0024]

[0025] t represents the time interval between two adjacent current data points. Set the derivative threshold d. th and periodic threshold n th The absolute value of the first derivative of the sampling points corresponding to the suspected events detected by the cluster removal algorithm is compared with the derivative threshold d. th For comparison, if the absolute value of its first derivative is greater than the derivative threshold d... th If the sample point is non-zero, it is considered to be a derivative fluctuation caused by noise; otherwise, it is considered to be a derivative fluctuation caused by noise. If the number of periods with consecutive non-zero derivatives is greater than the period threshold n... th If the events before and after are not zero, then they are considered the same event. Finally, the derivatives of the detected adjacent events are compared. If the proportion of non-zero derivatives is greater than 80%, or the number of consecutive non-zero derivatives exceeds the period threshold n, then the event is considered a single event. th If the subsequent adjacent event is determined to be caused by the continuous power change of the same long transient switching event, it needs to be removed. The derivative removal algorithm can further filter out various complex long transient false detections, and this method is universal, capable of filtering out false detections regardless of power fluctuations within the long transient.

[0026] Step 3: Feature Extraction

[0027] The selected feature is the V-reactive I trajectory; a two-dimensional image matrix is ​​constructed to reflect the VI trajectory characteristics of different types of loads; the active and reactive power of each electrical appliance at different sampling points are obtained and stored, and then the reactive current sequence is obtained; then the steady-state data is extracted to obtain the voltage and current data of one cycle after final averaging.

[0028] Step 3.1: Steady-state data acquisition

[0029] Step 3.1.1: Based on the event detection results from Step 2, take two segments of steady-state data I before and after each event. before V before I after V after Each data segment lasts for 6 cycles, with the starting position of the cycle being the first positive zero crossing of the voltage data.

[0030] Step 3.1.2: Fit the discrete steady-state data obtained in step 3.1.1 using linear interpolation to obtain the corresponding continuous waveform data.

[0031] Step 3.1.3: Fit the continuous data V before V after Perform FFT decomposition to obtain the initial phase of the voltage fundamental. and offset time Where f represents the power frequency of the voltage. Based on the fundamental phase shift of the original voltage data points, the zero point of the voltage can be accurately found to ensure the integrity of each cycle.

[0032] Step 3.1.4: Add the offset time Δt obtained after phase alignment to the time of the first zero-crossing positive voltage data point in the original voltage data, and use this as the initial time t for resampling. begin , with initial time t begin Starting from the four voltage and current waveforms fitted in step 3.1.2, resample to obtain four periods of i. before v before i after v after .

[0033] Step 3.1.5: Process the voltage and current data from the four cycles according to the sampling point order: that is (i after -i before ) and (v after +v before ) / 2, the voltage and current data of the load in the event are obtained for 4 cycles by cutting off the event. The average of the sampling points at the corresponding position in each of the 4 cycles is used to obtain the voltage and current data of the final cycle.

[0034] Step 3.2, VI trajectory generation

[0035] The voltage and current waveforms of the load within one cycle after averaging are extracted. Data normalization is then applied to linearly transform the voltage and current waveform data to the range [0, 1] to obtain the normalized voltage sequence v. * Current sequence i * ;

[0036] v(n) * and i(n) * Representing sequence v respectively * i * The voltage and current values ​​at the nth sampling point are obtained. A two-dimensional matrix R of size N×N is constructed. N×N Set all elements in the matrix to 0 to obtain the initial V-reactive-I trajectory matrix. Then, map the voltage and current waveforms to the two-dimensional matrix R using the following calculations. N×N middle:

[0037]

[0038] R N×N (a,b)=1

[0039] in The floor symbol; a n Let matrix R N×Nrow index, b n Let matrix R N×N Column index; R N×N (a,b) = 1, which sets the position of the waveform within the matrix to 1. Substituting the independent variable n from 1 to the last sample into the calculation, the V-reactive-I trajectory matrix can be plotted.

[0040] Step 3.3: Create a hybrid feature matrix

[0041] The active and reactive power characteristics are processed using binary encoding, converted into binary arrays, and added to the V-reactive-I trajectory matrix. The power characteristics and the V-reactive-I trajectory matrix are then integrated to form a new hybrid feature matrix. The specific steps are as follows:

[0042] Step 3.3.1: Normalize the common features of different samples to [0, 2] according to the following formula. N Within the range of -1]

[0043]

[0044] p is the sample number, q is the feature number, and S is the feature ID. p,q s is the p-th feature value of the q-th sample after normalization; pmax s pmin These are the maximum and minimum eigenvalues ​​of the p-th feature, respectively.

[0045] Step 3.3.2: Convert the normalized eigenvalues ​​into N-bit binary numbers;

[0046] Step 3.3.3: Construct an N×2 two-dimensional feature matrix R N×2 As a supplementary feature matrix, all elements in the initial state of the supplementary feature matrix are 0; add the (N×2) order two-dimensional feature matrix F to the right side of the V-reactive I matrix M to form an N×(N+2) order mixed feature matrix, and use the following formula to fill the binary number converted in step 3.3.2 into the supplementary feature matrix F:

[0047]

[0048] c is the binary bit index, ranging from [1, N]; The supplementary feature matrix R for the q-th sample N×2 The value of the element in the c-th row and p-th column; B() is a binary function used to convert a decimal array into a binary array; B p.q (c) is the c-th bit of the binary feature value of the p-th feature of the q-th sample.

[0049] Step 4: Load Identification

[0050] A fully convolutional network (FCN) incorporating Bi-TCN residual blocks is constructed. A binary classifier and upper / lower dual channels are sequentially placed at the output position of the Bi-TCN residual block. The N×(N+2) dimensionality feature matrix obtained in step three is reduced to an N×N dimensional image feature matrix, which is then input into the Bi-TCN residual block. The binary classifier then separates the unknown appliances into the upper channel and the known appliances into the lower channel. The features of the unknown appliances are directly output after passing through a one-dimensional convolutional layer in the upper channel. The features of the known appliances are first reshaped to an N×1×N dimension in the lower channel, then passed through two convolutional blocks for encoding and decoding, and finally flattened by a Flatten layer to output a one-dimensional recognition result.

[0051] The present invention has the following beneficial effects:

[0052] 1. To further improve the load identification accuracy of features, features are fused based on the VI trajectory, and reactive current is selected to replace the total current in the traditional VI trajectory, forming a hybrid feature matrix. Selecting the V-reactive-I trajectory as the load feature can compensate for the low distinguishability of different appliance states caused by reduced sampling accuracy. Compared with traditional VI trajectory image features, the trajectory composed of reactive current has higher load identification accuracy, and even household appliances with similar VI trajectories can be better distinguished, significantly improving the load recognition rate.

[0053] 2. The network model for identifying the load employs a combination of residual blocks and a dual-channel FCN structure. In the residual blocks, dilated convolutions and batch normalization are used instead of dilated convolutions and weight normalization, giving the neural network a larger receptive field and better feature extraction capabilities. The fully convolutional structure avoids the use of fully connected layers in traditional network architectures, resulting in more efficient trainable weights and more stable gradients, and is easier to train.

[0054] 3. The generalized likelihood ratio transient event detection algorithm based on multidimensional input is selected for event detection. A false detection filtering algorithm is also added, which effectively solves the problem of false detection caused by the high-frequency disturbance of steady state caused by nonlinear components or a large number of capacitors and inductors in the load equipment, as well as the long transient state of equipment switching. It can significantly improve the accuracy of event detection. Attached Figure Description

[0055] Figure 1 This is a flowchart of an embodiment of the present invention;

[0056] Figure 2 This is a schematic diagram of the network structure according to an embodiment of the present invention;

[0057] Figure 3This is a diagram of the bidirectional TCN (Bi-tcn) residual block structure in the network architecture. Detailed Implementation

[0058] The present invention will be further explained below with reference to the accompanying drawings;

[0059] As attached Figure 1 The diagram shows the overall process of non-intrusive load identification, the neural network model, the residual block structure, and a schematic diagram of dimensionality calculation.

[0060] Step 1: Data Collection

[0061] In practice, the main meter and sub-meter data detected by the electricity meter installed on the target object, including current and voltage, can be used as samples with a sampling frequency of 30kHz. Here, the main meter and sub-meter data in the Plaid dataset are directly selected. Then, the data in the dataset is denoised using median filtering to remove offset points in the original data while ensuring that the edges of the signal are not blurred, thereby improving the accuracy of the data.

[0062] Step 2: Event Detection

[0063] Event detection identifies the moment an event occurs, thus determining the necessary steady-state data segment for the VI trajectory. To improve accuracy, this method extends the one-dimensional generalized likelihood ratio test to multi-dimensional signal detection. The generalized likelihood ratio test algorithm based on the Neyman-Pearson criterion is used for event detection. A joint window is formed by two nearest-neighbor sliding windows (pre-event and post-event windows), sliding along the current matrix. An event is considered to have occurred when the algorithm's decision function exceeds a set threshold. The specific steps are as follows:

[0064] Step 2.1: Select the 30kHz current RMS value sequence from the sample data from the sub-table, and arrange them sequentially according to the acquisition time to form a one-dimensional circuit sequence K, which serves as the data points for event detection. Set the sliding window length to k sampling points and the width to 2, with k set to 200. Multiply each value in the current sequence K by a coefficient of 0.6 to obtain a new one-dimensional current sequence G, and align it with the beginning and end of K to form a two-dimensional sequence I with a width of 2 and a length of M. m Let m be the m-th element in the two-dimensional sequence I.

[0065] Step 2.2, Set the on and off states of electrical appliances to correspond to H1 and H0 respectively. (Assumption:)

[0066] H0 assumption: The two-dimensional current sequence does not contain any inflection points, that is:

[0067] H0:I={I m ,m=M-k+1,…,M}~N(μ0,∑)

[0068] H1 hypothesis: The two-dimensional current sequence contains a change point, that is:

[0069]

[0070] Decision function g n for:

[0071]

[0072] in:

[0073]

[0074]

[0075] n c Let μ0 be the moment when the change point occurs in the two-dimensional sequence I, and let μ0 and ∑0 be the mean and covariance matrix of the Gaussian distribution under the hypothesis test H0, respectively; 1a ,∑ 1a The time n of the change point under hypothesis test H1 is respectively. c Previous mean and covariance matrices; μ 1b ,∑ 1b The time n of the change point under hypothesis test H1 is respectively. c The mean and covariance matrices follow.

[0076] Set a decision threshold h, when the decision function value g... n When the value is greater than or equal to the decision threshold h, the H1 test is valid, meaning the sequence contains a change point and an event has occurred. The time n when the change point occurs is recorded. c When the decision function value g n When the current is less than the decision threshold h, test H0 is valid, meaning there are no changing points in the sequence and no event has occurred. The decision threshold h is dynamically adjusted continuously for different electrical appliances to achieve better event detection results. Finally, the current information at times when no event has occurred is set to 0, while the current information at times when an event has occurred remains unchanged, resulting in sequence I. event .

[0077] Step 2.3: Existing event detection algorithms rely on threshold settings and struggle to balance recall and accuracy. There are no particularly effective solutions for false positives caused by long transient events resulting from equipment switching events, i.e., continuous and irregular power changes. Furthermore, fluctuations in the steady-state of some nonlinear equipment can also easily lead to false positives. Therefore, based on the existing generalized likelihood ratio transient event detection algorithm, two false positive filtering algorithms—cluster removal and derivative removal—were used to filter out false positive events.

[0078] in I eventOne-dimensional unsupervised clustering is performed on the sequence without setting the number of cluster centers. All consecutive non-zero values ​​are aggregated into a cluster and identified as an event, thereby removing false detections caused by small fluctuations in current under the same operating conditions and obtaining the sequence ΔI.

[0079] Step 2.4: Cluster analysis can remove false detections of consecutive pseudo-events caused by long switching transients. These false detections occur because when the device changes its operating state, its electrical performance changes, further causing a shift in the aggregated measurement value. If the jump variable exceeds a threshold, the basic detection algorithm will falsely detect consecutive changes as events. However, for devices with long transients and a certain time interval between two detected pseudo-events, cluster analysis cannot remove them. Therefore, derivative removal is needed to further reduce false detections with time intervals. For devices with long switching transients, the derivative during the transient duration often deviates from the zero line due to the oscillation of the current waveform, and may cross the zero line several times before deviating in another direction. Based on this characteristic, the detected events can be further filtered by judging the sign of the derivative between adjacent detected events. For the current sequence ΔI obtained in Step 2.3, any current point ΔI at a non-endpoint... f The first derivative can be calculated by the following formula:

[0080]

[0081] t represents the time interval between two adjacent current data points. Set the derivative threshold d. th and periodic threshold n th The absolute value of the first derivative of the sampling points corresponding to the suspected events detected by the cluster removal algorithm is compared with the derivative threshold d. th For comparison, if the absolute value of its first derivative is greater than the derivative threshold d... th If the sample point is non-zero, it is considered to be a derivative fluctuation caused by noise; otherwise, it is considered to be a derivative fluctuation caused by noise. If the number of periods with consecutive non-zero derivatives is greater than the period threshold n... th If the events before and after are not zero, then they are considered the same event. Finally, the derivatives of the detected adjacent events are compared. If the proportion of non-zero derivatives is greater than 80%, or the number of consecutive non-zero derivatives exceeds the period threshold n, then the event is considered a single event. th If the subsequent adjacent event is determined to be caused by the continuous power change of the same long transient switching event, it needs to be removed. The derivative removal algorithm can further filter out various complex long transient false detections, and this method is universal, capable of filtering out false detections regardless of power fluctuations within the long transient.

[0082] Step 3: Feature Extraction

[0083] The selected feature is the V-reactive I trajectory; a two-dimensional image matrix is ​​constructed to reflect the VI trajectory characteristics of different types of loads. The calculation process is simple and the restoration effect is good; the active and reactive power of each electrical appliance at different sampling points are obtained and stored, and then the reactive current sequence is obtained; then the steady-state data is extracted to obtain the voltage and current data of one cycle after final averaging.

[0084] First, for obtaining reactive current, this method is based on Fryze power theory for calculation. Based on the measured instantaneous voltage v(n) and instantaneous current i(n), discrete Fourier transform is used to obtain the current and voltage harmonics, respectively.

[0085]

[0086] s represents the harmonic order; X(s) and Y(s) are the frequency domain components of the s-th harmonic, respectively. i(n) is the n-th current sampling point, v(n) is the n-th voltage sampling point, and L is the number of sampling points. |X(s)| and |Y(s)| represent the amplitudes of the s-th current and voltage harmonics, respectively. Since the calculation of active and reactive power requires substituting the effective values, the absolute values ​​need to be taken and multiplied by a coefficient to convert them into the effective values ​​of current and voltage.

[0087] and These represent the phases of the s-th harmonics of X(s) and Y(s), respectively. This represents the phase difference.

[0088]

[0089] The active power P and reactive power Q can be obtained from the following formulas. Based on experience, s is taken as 20.

[0090]

[0091] P and Q represent active power and reactive power, respectively, V s I s These represent the RMS values ​​of voltage and current during steady-state operation, respectively. The active and reactive power of each of the η different electrical appliances at each sampling time are labeled and stored in an Excel file. The number of appliance types η can be set according to the scenario. These labels are used to match the results.

[0092] The active current i can be obtained using the following formula. a (n), and finally subtract the active current i from the instantaneous current i(n). a (n) Obtain the reactive current i b (n):

[0093]

[0094] ib (n)=i(n)-i a (n)

[0095] After acquiring the reactive current, trajectory extraction is performed:

[0096] Step 3.1: Steady-state data acquisition

[0097] Step 3.1.1: Based on the event detection results from Step 2, take two segments of steady-state data I before and after each event. before V before I after V after Each data segment lasts for 6 cycles, with the starting position of the cycle being the first positive zero crossing of the voltage data.

[0098] Step 3.1.2: Fit the discrete steady-state data obtained in step 3.1.1 using linear interpolation to obtain the corresponding continuous waveform data.

[0099] Step 3.1.3: Fit the continuous data V before V after Perform FFT decomposition to obtain the initial phase of the voltage fundamental. and offset time Where f represents the power frequency of the voltage. Based on the fundamental phase shift of the original voltage data points, the zero point of the voltage can be accurately found to ensure the integrity of each cycle.

[0100] Step 3.1.4: Add the offset time Δt obtained after phase alignment to the time of the first zero-crossing positive voltage data point in the original voltage data, and use this as the initial time t for resampling. begin , with initial time t begin Starting from the four voltage and current waveforms fitted in step 3.1.2, resample to obtain four periods of i. before v before i after v after .

[0101] Step 3.1.5: Process the voltage and current data from the four cycles according to the sampling point order: that is (i after -i before ) and (v after +v before ) / 2, the voltage and current data of the load in the event are obtained for 4 cycles by cutting off the event. The average of the sampling points at the corresponding position in each of the 4 cycles is used to obtain the voltage and current data of the final cycle.

[0102] Step 3.2, VI trajectory generation

[0103] The voltage and current waveforms of the load within one cycle after averaging are extracted. Data normalization is then applied to linearly transform the voltage and current waveform data to the range [0, 1] to obtain the normalized voltage sequence v. * Current sequence i * ;

[0104] v(n) * and i(n) * Representing sequence v respectively * i * The voltage and current values ​​at the nth sampling point are obtained. A two-dimensional matrix R of size N×N is constructed. N×N Set all elements in the matrix to 0 to obtain the initial V-reactive-I trajectory matrix. Then, map the voltage and current waveforms to the two-dimensional matrix R using the following calculations. N×N middle:

[0105]

[0106] R N×N (a,b)=1

[0107] in The floor symbol; a n Let matrix R N×N row index, b n Let matrix R N×N Column index; R N×N (a,b) = 1, which sets the position of the waveform within the matrix to 1. Substituting the independent variable n from 1 to the last sample into the calculation, the V-reactive-I trajectory matrix can be plotted.

[0108] Step 3.3: Create a hybrid feature matrix

[0109] The active and reactive power characteristics are processed using binary encoding, converted into binary arrays, and added to the V-reactive-I trajectory matrix. The power characteristics and the V-reactive-I trajectory matrix are then integrated to form a new hybrid feature matrix. The specific steps are as follows:

[0110] Step 3.3.1: Normalize the common features of different samples to [0, 2] according to the following formula. N Within the range of -1]

[0111]

[0112] p is the sample number, q is the feature number, and S is the feature ID. p,q s is the normalized eigenvalue of the p-th feature of the q-th sample; pmax s pmin These are the maximum and minimum eigenvalues ​​for the p-th feature, respectively.

[0113] Step 3.3.2: Convert the normalized eigenvalues ​​into N-bit binary numbers;

[0114] Step 3.3.3: Construct an N×2 two-dimensional feature matrix R N×2 As a supplementary feature matrix, all elements in the initial state of the supplementary feature matrix are 0; add the (N×2) order two-dimensional feature matrix F to the right side of the V-reactive I matrix M to form an N×(N+2) order mixed feature matrix, and use the following formula to fill the binary number converted in step 3.3.2 into the supplementary feature matrix F:

[0115]

[0116] c is the binary bit index, ranging from [1, N]; The supplementary feature matrix R for the q-th sample N×2 The value of the element in the c-th row and p-th column; B() is a binary function used to convert a decimal array into a binary array; B p.q (c) is the c-th bit of the binary feature value of the p-th feature of the q-th sample.

[0117] Step 4: Network Computation

[0118] The invention proposes an improved FCN network, such as Figure 2 As shown, the improved fully convolutional network combines, for example, Figure 3 The Bi-TCN residual block shown uses a fully connected structure and adds dual channels to distinguish between known and unknown electrical appliances. The improved FCN network incorporating Bi-TCN residual blocks is less affected by input size than traditional convolutional neural networks and avoids the use of fully connected layers in traditional network architectures. This new FCN-based modeling method has more efficient trainable weights and more stable gradients, and is easier to train. For computational convenience, the fused 32×(32+2) dimensional image feature matrix is ​​reduced to a 32×34 dimensional image feature matrix before being input into the entire residual block. The Bi-TCN residual block has two residual units, each of which includes, in sequence, Non-casual dilated, Batch Normalization, ReLU, and Dropout. The dilated convolutional layers used significantly increase the receptive field and provide better prediction results compared to standard convolutional layers. The BatchNormalization layer addresses vanishing and exploding gradients while also increasing learning speed. The ReLU activation function is chosen, and the Dropout layer uses a binomial function to generate a one-dimensional vector of 0s and 1s with the number of neurons as its length. A 0 represents a neuron temporarily blocked from transmitting information to subsequent neurons, thus reducing computation and preventing overfitting. The final output is a feature vector with the same size as the input.

[0119] After the residual block, a binary classifier separates the unknown and known electrical appliances. If it is an unknown electrical appliance, it enters the upper channel, passes through a one-dimensional convolutional layer, and is directly output; if it is a known electrical appliance, it enters the lower channel, which consists of two convolutional blocks: encoding and decoding. The reshape function is first used to adjust the input dimension to 32×1×32.

[0120] The obtained image matrix is ​​divided into training, testing, and prediction sets in a 6:2:2 ratio. The training set data is input into a convolutional network using the cross-entropy loss function. After multiple backpropagations to adjust the network's internal weights until they fall below a set threshold, the network is trained again using the training and validation sets to determine the final model. Finally, the test set data is input into the trained network, and the input dimensionality is continuously reduced along the convolutional blocks. The data is then flattened using a Flatten layer, transforming the multidimensional array into a one-dimensional array. The final algorithm result is then obtained.

[0121] Step 5: Compare output results with labels and evaluate performance.

[0122] After network processing, the two channels are superimposed to obtain a 14-bit binary array string, which is then compared with the tags stored in step one. If the first bit is 1 and the second bit is 0, it indicates an appliance category not stored in the tag database, and is therefore identified as an unknown appliance. If the first bit is 0 and the second bit is 1, the appliance type can be successfully identified and its on / off status obtained by matching the number of 1s in the remaining 12 bits with the tag database. The appliances represented by bits 3-12 are shown in the figure below.

[0123]

[0124] To evaluate the performance of load identification, this paper uses the recognition accuracy (Precision) to assess system performance. Precision represents the proportion of samples correctly identified as open out of all samples classified as open.

[0125]

[0126] tp represents the number of samples with the state "on" that were correctly predicted; fp represents the number of samples with the state "off" that were incorrectly predicted. Higher accuracy means higher recognition precision.

[0127] The following is a comparison of results under controlled variables:

[0128] Algorithm recognition accuracy The algorithm proposed in this invention 93.6％ Traditional FCN network model 91.9％ Traditional VI trajectory integrates active and reactive power. 90.7％

[0129] It can be seen that the improved algorithm proposed in this application has a significant improvement in the accuracy of identifying both known and unknown loads.

Claims

1. A non-intrusive load identification calculation method, characterized in that: Specifically, the following steps are included: Step 1: Data Collection The current and voltage data from the main meter and sub-meters installed on the target object are sampled. The sampled data is then denoised using median filtering, and offset points in the original data are removed to obtain the sample data. Step 2: Event Detection The generalized likelihood ratio test algorithm based on the Neyman-Pearson criterion is selected to detect events in the sample data. A joint window is formed by setting two nearest sliding windows before and after the event. The joint window slides on the current matrix. When the algorithm decision function is greater than the set threshold, it is considered that an event has occurred and the time of the event is saved. For the event detection result I obtained through generalized likelihood ratio event Further, two false detection filtering algorithms, cluster removal and derivative removal, are used for filtering, specifically: Firstly in I event One-dimensional unsupervised clustering is performed on the sequence without setting the number of cluster centers. All continuous non-zero values ​​are aggregated into a cluster and identified as an event, thereby removing false detections caused by small fluctuations in current under the same operating conditions, and obtaining the current sequence ΔI. Then calculate the current point ΔI at any non-endpoint in the current sequence ΔI. f First derivative ΔI f ': t represents the time interval between two adjacent current data points; Set the derivative threshold d th and periodic threshold n th The absolute value of the first derivative of the sampling points corresponding to the suspected events detected by the cluster removal algorithm is compared with the derivative threshold d. th For comparison, if the absolute value of its first derivative is greater than the derivative threshold d... th If the value is zero, the sampling point is considered non-zero; otherwise, it is considered a derivative fluctuation caused by noise. If the number of periods in which the continuous derivative is not zero is greater than the period threshold n th If the events before and after are considered as the same event, then the derivatives of the detected adjacent events are compared. If the proportion of non-zero derivatives is greater than 80% or the number of consecutive non-zero derivatives exceeds the period threshold n, then the event is considered as the same event. th If the next adjacent event is determined to be caused by the continuous change in power of the same long transient switching event, it needs to be removed. Step 3: Feature Extraction Construct a two-dimensional image matrix to reflect the trajectory characteristics of different types of loads (VI); The active and reactive power of each electrical appliance at different sampling points are acquired and stored, and then the reactive current sequence is obtained. Subsequently, steady-state data is extracted to obtain the voltage and current data of one cycle after final averaging. Step 4: Load Identification An FCN network incorporating Bi-TCN residual blocks is constructed. Each Bi-TCN residual block has two residual units, and each residual unit includes Non-casual dilated, Batch Normalization, ReLU, and Dropout in sequence. A binary classifier and upper and lower dual channels are sequentially set at the output position of the Bi-TCN residual block. The image feature matrix obtained in step three is dimensionality reduced and then input into the Bi-TCN residual block. The binary classifier then splits the unknown appliances into the upper channel and the known appliances into the lower channel. The features of the unknown appliances are directly output after passing through a one-dimensional convolutional layer in the upper channel. The features of the known appliances are first adjusted to N×1×N by a reshape function in the lower channel, and then passed through two convolutional blocks to complete the encoding and decoding process. Finally, the features are flattened by a Flatten layer to output a one-dimensional recognition result.

2. The non-invasive load identification calculation method as described in claim 1, characterized in that: The specific process of step 2 is as follows: Step 2.1: Select the 30kHz current RMS value sequence from the sample data table and arrange it sequentially according to the acquisition time to form a one-dimensional circuit sequence K, which serves as the data points for event detection; set the sliding window length to k sampling points and the width to 2; multiply each value in the current sequence K by a coefficient of 0.6 to obtain a new one-dimensional current sequence G, and align it with the beginning and end of K to form a two-dimensional sequence I with a width of 2 and a length of M. m Let m be the m-th element in the two-dimensional sequence I; Step 2.2, Set the on and off states of electrical appliances to correspond to H1 and H0 respectively. (Assumption:) H0 assumption: The two-dimensional current sequence does not contain any inflection points, that is: H0:I={I m ,m=M-k+1,…,M}~N(μ0,∑) H1 hypothesis: The two-dimensional current sequence contains a change point, that is: Decision function g n for: in: n c Let μ0 be the moment when the change point occurs in the two-dimensional sequence I, and let μ0 and ∑0 be the mean and covariance matrix of the Gaussian distribution under the hypothesis test H0, respectively; 1a ,∑ 1a The time n of the change point under hypothesis test H1 is respectively. c Previous mean and covariance matrices; μ 1b ,∑ 1b The time n of the change point under hypothesis test H1 is respectively. c The mean and covariance matrices follow; Set a decision threshold h, when the decision function value g... n When the value is greater than or equal to the decision threshold h, the H1 test is valid, meaning the sequence contains a change point and an event has occurred. The time n when the change point occurs is recorded. c When the decision function value g n When the current is less than the decision threshold h, test H0 is valid, meaning there are no changing points in the sequence and no event has occurred. The decision threshold h is dynamically adjusted continuously for different electrical appliances to achieve better event detection results. Finally, the current information at times when no event has occurred is set to 0, while the current information at times when an event has occurred remains unchanged, resulting in sequence I. event .

3. The non-invasive load identification calculation method as described in claim 1, characterized in that: The specific process of step 3 is as follows: Step 3.1: Steady-state data acquisition Step 3.1.1: Based on the event detection results from Step 2, take two segments of steady-state data I before and after each event. before V before I after V after Each data segment lasts for 6 cycles, with the starting position of the cycle being the first positive zero crossing of the voltage data. Step 3.1.2: Fit the discrete steady-state data obtained in Step 3.1.1 using linear interpolation to obtain the corresponding continuous waveform data; Step 3.1.3: Fit the continuous data V before V after Perform FFT decomposition to obtain the initial phase of the voltage fundamental. and offset time Where f represents the power frequency of the voltage; based on the fundamental phase shift of the original voltage data points, the zero point of the voltage can be accurately found to ensure the integrity of each cycle; Step 3.1.4: Add the offset time Δt obtained after phase alignment to the time of the first zero-crossing positive voltage data point in the original voltage data, and use this as the initial time t for resampling. begin , with initial time t begin Starting from the four voltage and current waveforms fitted in step 3.1.2, resample to obtain four periods of i. before v before i after v after ; Step 3.1.5: Process the voltage and current data from the four cycles according to the sampling point order: that is (i after -i before ) and (v after +v before ) / 2, the voltage and current data of the load in the event are obtained for 4 cycles by cutting off the event. The average of the sampling points at the corresponding position in each of the 4 cycles is used to obtain the voltage and current data of the final cycle. Step 3.2, VI trajectory generation The voltage and current waveforms of the load within one cycle after averaging are extracted. Data normalization is then applied to linearly transform the voltage and current waveform data to the range [0, 1] to obtain the normalized voltage sequence v. * Current sequence i * ; v(n) * and i(n) * Representing sequence v respectively * i * The voltage and current values ​​at the nth sampling point are used to construct a two-dimensional matrix R of size N×N. N×N Set all elements in the matrix to 0 to obtain the initial V-reactive-I trajectory matrix. Then, map the voltage and current waveforms to the two-dimensional matrix R using the following calculations. N×N middle: R N×N (a,b)=1 in The floor symbol; a n Let matrix R N×N row index, b n Let matrix R N×N Column index; R N×N (a,b)=1 means setting the position passed by the waveform in the matrix to 1; substitute the independent variable n from 1 to the last sample into the calculation to draw the V-reactive-I trajectory matrix; Step 3.3: Create a hybrid feature matrix The active and reactive power characteristics are processed using binary encoding, converted into binary arrays, and added to the V-reactive-I trajectory matrix. The power characteristics and the V-reactive-I trajectory matrix are then integrated to form a new hybrid feature matrix. The specific steps are as follows: Step 3.3.1: Normalize the common features of different samples to [0, 2] according to the following formula. N Within the range of -1] p is the sample number, q is the feature number, and S is the feature ID. p,q s is the p-th feature value of the q-th sample after normalization; pmax s pmin These are the maximum and minimum eigenvalues ​​of the p-th feature, respectively; Step 3.3.2: Convert the normalized eigenvalues ​​into N-bit binary numbers; Step 3.3.3: Construct an N×2 two-dimensional feature matrix R N×2 As a supplementary feature matrix, all elements in the initial state of the supplementary feature matrix are 0; add the (N×2) order two-dimensional feature matrix F to the right side of the V-reactive I matrix M to form an N×(N+2) order mixed feature matrix, and use the following formula to fill the binary number converted in step 3.3.2 into the supplementary feature matrix F: c is the binary bit index, ranging from [1, N]; The supplementary feature matrix R for the q-th sample N×2 The value of the element in the c-th row and p-th column; B() is a binary function used to convert a decimal array into a binary array; B p.q (c) is the c-th bit of the binary feature value of the p-th feature of the q-th sample.

4. The non-intrusive load identification calculation method as described in claim 1, characterized in that: Reactive current is calculated based on Fryze power theory.

5. A non-invasive load identification calculation method as described in claim 1 or 4, characterized in that: Discrete Fourier transform is applied to instantaneous voltage v(n) and instantaneous current i(n) to obtain current and voltage harmonics, respectively: Where s is the harmonic order, X(s) and Y(s) are the frequency domain components of the s-th harmonic, respectively; i(n) is the n-th current sampling point, v(n) is the n-th voltage sampling point, and L is the number of sampling points; |X(s)| and |Y(s)| represent the amplitudes of the s-th current and voltage harmonics, respectively. These represent the phases of the s-th harmonics of X(s) and Y(s), respectively. Phase difference; The calculation methods for active power P and reactive power Q are as follows: V s I s These represent the effective values ​​of voltage and current during steady-state operation, respectively; reactive current i b The formula for calculating (n) is: i b (n)=i(n)-i a (n) Where i a i(n) represents the active current, and i(n) represents the instantaneous current.

Images