A Nonnegative Matrix Factorization Clustering Method for Images Based on Neurodynamics

Abstract

An image non - negative matrix factorization clustering method based on neural dynamics, comprising: 1) constructing an original data matrix; 2) constructing a sparse graph regularized non - negative low - rank matrix factorization model; 3) performing sparse graph regularized non - negative low - rank matrix factorization on the original data matrix; 4) performing clustering based on the K - means algorithm. Advantages of the present invention: 1. Sparse graph regularized non - negative low - rank matrix factorization maintains the geometric configuration of high - dimensional data represented in a low - dimensional space by adding a graph regularization term. By imposing low - rank constraints and sparse constraints, the sparse property of the obtained basis matrix is ensured; more effective information is extracted, information redundancy is reduced, and the computational difficulty is lowered; 2. The neural dynamics optimization method VF - CNO is introduced into the problem of solving the sparse graph regularized non - negative matrix low - rank matrix factorization model, improving the solving accuracy of the model. The present invention has a better dimensionality reduction effect and shows better clustering performance in image clustering tasks.

Description

Technical Field

[0001] The present invention belongs to the technical field of information processing, and specifically relates to an image non-negative matrix factorization clustering method based on neurodynamics. Background Art

[0002] With the progress of science and technology, modern medical technology has made rapid development. A variety of digital imaging technologies emerge in an endless stream, such as ultrasonic images, computed tomography, optical coherence tomography, and nuclear magnetic resonance, etc., which play a key role in pathological analysis and medical diagnosis, and thus are increasingly well-known to people. In this process, the medical image information that people can collect becomes more and more abundant. At the same time, this also causes the phenomenon of complex and diverse data types and explosive growth of data volume. These data often have the characteristics of large amount of information, high dimension, fast growth, and unstructured, etc., resulting in extremely difficult processing of medical images for us, and even the situation of "curse of dimensionality" appears. Therefore, researching how to extract effective information from intricate medical images has important practical significance and application value for realizing efficient information processing.

[0003] Non-negative Matrix Factorization (NMF) is a commonly used dimensionality reduction and feature extraction algorithm in the fields of data mining and clustering. It has been widely applied in medical image clustering, image retrieval, face recognition and other fields, and has achieved considerable success. The reason for its success in the field of images is that digital images are generally expressed in the form of matrices or tensors. In particular, digital images are stored inside a computer based on a pixel dot matrix (two-dimensional matrix). On the other hand, matrix factorization techniques can maintain the potential laws and connections existing within the data, so there is a natural advantage in dimensionality reduction. Different from traditional matrix factorization techniques, such as Principal Component Analysis (PCA), Singular Value Decomposition (SVD), Linear Discriminant Analysis (LDA) and Vector Quantization (VQ), etc., NMF guarantees the non-negativity of the resulting decomposition matrix. Generally speaking, for a given non-negative original matrix (without non-negative elements), there will definitely be no non-negative elements in its decomposition matrix, and this property exactly coincides with the real physical laws. For example, the brightness of pixel points in medical images can only be positive; there will be no negative gene sequences in DNA array analysis; there are no negative count statistics in text clustering. NMF is also a linear representation algorithm based on partial data, which coincides with the view in psychology and physiology that "the whole of an object is composed of parts" and has strong interpretability. Due to the characteristics of non-negative constraints and localized representation, NMF not only makes the data easier to be sparsely represented, but also has better robustness in its sparse representation. Inspired by the concept and principle of manifold learning algorithms, Cai D, He X, Han J et al. proposed Graph Regularized Non-negative Matrix Factorization (GNMF) in the paper "Graph regularized nonnegative matrix factorization for data representation". This algorithm uses an adjacency graph to record the potential geometric relationships of the data and embeds the manifold constraint into the matrix factorization process, maintaining the local invariance of the data in the original space and achieving good performance in image clustering.In their paper "Graph regularized non-negative low-rank matrix factorization for image clustering," Li X, Cui G, Dong Y, et al. incorporated low-rank constraints into the objective function and proposed non-negative low-rank matrix factorization (NLMF), an algorithm that can extract low-dimensional spatial information inherent in the data. While these methods improve on non-negative matrix factorization from different perspectives, none of them consider the sparsity of the resulting matrix.

[0004] In addition, the problem of solving NMF is an NP-hard (non-deterministic polynomial-time hard) problem and cannot be solved directly through mathematical deduction. Existing traditional solution methods are all based on numerical iteration for calculation, such as the Multiplicative Update Rules (MU), Alternating Non-negative Least Squares (ANLS) and Projected Gradient (PG). Their main idea is to regard NMF as an approximate decomposition process, establish the objective function, and then design appropriate iterative rules to make the objective function continuously approach the optimal position. Although these numerical methods are simple and easy to implement, they generally have problems such as low operating efficiency and unstable performance. There are two main reasons for this: (1) The objective function of the NMF model is a non-convex function. During the iteration process, the algorithm is very likely to fall into the local optimum, resulting in the optimization result not being uniquely determined; (2) Traditional iterative methods are all serial calculation methods, which have the disadvantages of long iteration time and low computational efficiency in the optimization solution process. Therefore, designing a stable, accurate, efficient and parallel NMF solution algorithm is one of the challenges that need to be urgently addressed in this field. Summary of the Invention

[0005] Based on the above technical background, this paper provides a neural dynamics-based image non-negative matrix factorization clustering method. This paper proposes a new sparse graph regularized non-negative low-rank matrix factorization algorithm model, and then optimizes and solves the proposed algorithm model using the bilateral random projection algorithm and neural dynamics methods.

[0006] The present invention is achieved through the following technical solutions.

[0007] The present invention provides a neural dynamics-based image non-negative matrix factorization clustering method, which comprises the following steps:

[0008] Step 1: Construct the original data matrix. Vectorize the images in the image dataset and construct a non - negative original data matrix V = [v1, v2,..., v n ∈ R m×n , where m is the dimension of the image samples and n is the number of image samples;

[0009] Step 2: Construct the Sparse Graph Regularized Nonnegative Low - rank Matrix Factorization (SGNLMF) model as follows:

[0010]

[0011] s.t. L = WH

[0012] L, W, H ≥ 0

[0013] ||w j ||1 = 1, j = 1,..., r

[0014] card(S) ≤ s, rank(L) ≤ l

[0015] where W ∈ R m×r is the basis matrix, H ∈ R r×n is the coefficient matrix, L ∈ R m×n is the low - rank matrix, S ∈ R m×n is the sparse matrix, l represents the range set of the rank of matrix L, s represents the sparse range set of matrix S, rank(L) represents the rank of matrix L, card(S) represents the number of non - zero elements of matrix S, α is the graph regularization parameter, L ap is the graph Laplacian matrix, β is the Tikhonov regularization parameter; Sparse Graph Regularized Nonnegative Low - rank Matrix Factorization imposes additional low - rank constraints, Tikhonov regularization constraints, sparse constraints, and graph regularization constraints during the matrix factorization process to mine sparse and effective feature information behind the data;

[0016] Step 3: Perform sparse graph regularized nonnegative low - rank matrix factorization on the original data matrix;

[0017] Step 3.1: Low - rank matrix recovery. Use the bilateral random projection algorithm to extract the potential low - rank matrix L ∈ R m×n from the original data matrix V. Mathematically, the problem can be described as:

[0018]

[0019] s.t. L = WH

[0020] L, W, H ≥ 0

[0021] card(S)≤s,rank(L)≤l

[0022] Step 3.2: Perform non-negative matrix factorization on the low-rank matrix L obtained in step 3.1 using the collaborative neural dynamics optimization method VF-CNO with variable scaling factors to obtain the basis matrix W∈R m×r and coefficient matrix H∈R r×n The basis matrix W contains the "localized" features in the original data matrix, and the coefficient matrix H represents the low-dimensional representation of the original data. At the same time, in order to enhance the decomposition performance of the model, additional Tikhonov regularization constraints, sparsity constraints, and graph regularization constraints are imposed during the matrix decomposition process. The sparsity constraints are also incorporated into the objective function. The newly formed objective function is as follows:

[0023]

[0024]

[0025] Among them, Γ = (γ1, γ2,..., γ n ) is the Lagrange multiplier of the equality constraint. ω, h are vectorized operations on the W, H matrices as follows:

[0026] ω=vec(W)=(w 11 w 21 …w m1 ;w 12 w 22 …w m2 ;…;w 1r w 2r ...w mr ) T

[0027] h=vec(H)=(h 11 h 21 …h r1 ;h 12 h 22 …h r2 ;…;h 1n h 2n …h rn ) T

[0028] Step 4: Clustering based on the K-means algorithm. Input the coefficient matrix H obtained in step 3 into the K-means classifier to achieve the clustering task of unlabeled image data.

[0029] Furthermore, step 3.1 of the present invention includes the following steps:

[0030] 3.1.1: Setting the relative reconstruction error The threshold is μ, the rank of the low-rank matrix L is set to l, and the sparse range set of the sparse matrix S is set to s;

[0031] 3.1.2: Initialize the correlation matrix and parameters, L0 = V, S0 = 0, k = 0;

[0032] 3.1.3: Perform iterative update: k = k + 1, A2 = U2,

[0033] 3.1.4: If Jump to step 3.1.2;

[0034] 3.1.5: Update matrices L and S: S k = Ω k (V - L k );

[0035] 3.1.6: If Then execute step 3.1.7, if Then execute step 3.1.3;

[0036] 3.1.7: Terminate the iteration to obtain the low-rank matrix L = L k , the sparse matrix S = S k .

[0037] Furthermore, step 3.2 of the present invention includes the following steps:

[0038] 3.2.1: Input the low-rank matrix L obtained in step 3.1 and normalize it so that its value range is [0, 1];

[0039] 3.2.2: Set the error tolerance ε, the threshold parameter ΘT, the number of particles N in the population, the maximum number of iterations k max , the regularization parameters α, β, the number of nearest neighbors ρ = 5, and the value of the decomposition dimension r;

[0040] 3.2.3: Use the K-nearest neighbor algorithm to construct a nearest neighbor graph and calculate the graph Laplacian matrix L ap = D ap - Y ap ; where Y ap is a symmetric weight matrix, and D ap is a diagonal matrix, and the diagonal elements are the column sums of Y ap ;

[0041] 3.2.4: Randomly initialize ω i (0), h i(0), Γ(0), i = 1, ..., N;

[0042] 3.2.5: Set the personal optimal solution as:

[0043] 3.2.6: Set the global optimal solution as: p * (0) = argmin[ω(0), h(0), Γ(0)], i = 1, ..., N;

[0044] 3.2.7: Conduct iterative update:

[0045] 3.2.7.1: Inspired by the design idea of FF-RNN and considering the potential time-varying characteristics of the recurrent neural network, a variant-scaling factors recurrent neural network model (VF-RNN) is proposed for the KKT point search of the constrained global optimization problem. The proposed VF-RNN model can be expressed as the following differential equation:

[0046]

[0047] where P Ω (·) is the projection operator. According to the above vector equation, its scalar form can be described as:

[0048]

[0049] where the convergence factor ∈ = (∈1, ∈2, ..., ∈ n ) T ; Ψ(∈, t) is a vector composed of the variant-scaling factor functions ψ(∈ j , t) (j ∈ {1, 2, ..., n}); for the said variant-scaling factor function ψ(∈, t), it is a function related to the convergence factor ∈ and time t and can take one of the following two forms:

[0050] (1) Exponential form:

[0051] ψ(∈, t) = ∈e t

[0052] or (2) Power-exponential form:

[0053] ψ(∈, t) = t ∈ + ∈

[0054] Compared with the traditional FF-RNN, the VF-RNN proposed in the present invention has two main advantages: (1) The VF-RNN uses a time-varying factor function ψ(∈, t) to enhance the convergence ability of the recurrent neural network, greatly reducing the convergence time required by the network; (2) By adopting a vector-type convergence parameter ∈ and setting different values of ∈ for each neural network, the exploration ability of the neural dynamic system is effectively improved, making it easier to search for the global optimal solution in the feasible domain space.

[0055] N particles perform local search according to the VF-RNN dynamic equation to obtain the personal optimal solution

[0056] 3.2.7.2: If Then p i (k + 1) = p i (k); If Then

[0057] 3.2.8: Update the global optimal position of the population: 1 ≤ i ≤ N;

[0058] 3.2.9: Calculate the initial state of the network according to the particle swarm update rule:

[0059] 3.2.9.1: Introduce the Particle Swarm Optimization (PSO). PSO performs a search operation through a group of particles and randomly creates the initial positions of the particles. Each particle represents a potential solution, and its position information is represented by the position vector x. This algorithm assumes that the particle swarm moves to the global optimum within a preset space range, and the movement speed of each particle is described by a velocity vector v. For the convenience of subsequent discussion, define x i = (x i1 ,..., x in ) T as the current position of the particle; v i = (v i1 ,..., v in ) T represents the movement speed of the particle; p i = (p i1 , …, p in ) T represents the best position of the individual during the search process; g = (g1,..., g n ) T represents the global best position among all particles during the search process. The description of the standard PSO is as follows:

[0060] v i (k + 1) = ωv i(k) + c1F1(p i -x i (k)) + c2r2(g - x i (k))

[0061] x i (k + 1) = x i (k) + v i (k + 1)

[0062] where ω ∈ [0, 1] is a fixed parameter; c1 and c2 are acceleration factors that determine the influence degree of the particle velocity by p i and g; r1 and r2 are random numbers in the range of [0, 1]; k represents the number of iterations.

[0063] First, update the initial state of the particle based on the above standard PSO update rule for updating.

[0064] 3.2.9.2: To further improve the search efficiency, wavelet mutation is used to improve the premature convergence problem of the population, thereby enhancing the population diversity. The specific rules are as follows:

[0065] First, define an evaluation function for particle diversity:

[0066]

[0067] where N represents the number of the population. Substitute into the above evaluation function;

[0068] 3.2.9.3: Compare the obtained Θ with the set threshold Θ T . If Θ ≥ Θ T , then directly execute the next step; if Θ < Θ T , then execute:

[0069]

[0070] where H i and L i represent the upper and lower bounds of x i respectively, and the parameter η can be obtained according to the following wavelet function:

[0071]

[0072] where ψ ∈ [-2.5a, 2.5a], a = exp(10(k / k max )) and k max represents the maximum number of iterations.

[0073] 3.2.10: If k ≤ k maxAnd ||p * (k + 1)-p * (k)|| ≥ ε, then jump to step 3.2.7. If k > k max And ||p * (k + 1)-p * (k)|| < ε, then the condition is satisfied, the loop terminates, and the global optimal solutions w and h are output;

[0074] 3.2.11: Re-convert the w and h obtained in 3.2.9 into matrix form to obtain the optimal matrices W and H.

[0075] Furthermore, step 4 of the present invention includes the following steps:

[0076] 4.1: Input the coefficient matrix H obtained in step 3 into the K-means classifier to implement the clustering task of unlabeled image data;

[0077] 4.2: Calculate the accuracy (Accuracy, AC) and the normalized mutual information NMI to measure the clustering performance:

[0078] 4.2.1: Accuracy AC

[0079] Accuracy is a commonly used clustering evaluation index, representing the ratio of the number of successfully predicted samples to the total number of samples. It is defined as:

[0080]

[0081] where s i and r i respectively represent the labels obtained by clustering and the true labels of the original data, n represents the total number of samples; δ(x, y) is a comparison function, which is 1 when x = y and 0 otherwise; map(x) is a label mapping function that projects the sample labels into labels equivalent to those in the dataset, and generally the Kuhn-Munkres algorithm is used to determine the best mapping; the higher the value of AC, the better the clustering effect of the algorithm.

[0082] 4.2.2: Normalized Mutual Information (NMI)

[0083] Mutual information is used to measure the similarity between the prediction result and the true result. It should be noted that NMI evaluates the clustering result based on the mutual information function and the entropy function. It is defined as:

[0084]

[0085] where c represents the number of clusters, n i represents the number of samples in the cluster of, Indicating the cluster The number of samples, n i，j Indicating the samples that belong to both cluster and cluster Generally speaking, the higher the value of NMI, the more prominent the clustering performance of the algorithm.

[0086] The beneficial effects of the present invention are as follows: The SGNLMF model proposed in the present invention has good performance in clustering analysis. From the perspective of the algorithm model, by adding a graph regularization term, SGNLMF maintains the geometric configuration of the high-dimensional data represented in the low-dimensional space. By imposing low-rank constraints and sparse constraints, the sparse properties of the decomposed basis matrix are guaranteed. Based on this, not only more effective information is extracted, but also information redundancy is reduced, and the computational difficulty is lowered. From the perspective of model optimization and solution, the neurodynamics optimization method VF-CNO is introduced into the solution problem of the SGNLMF model, improving the solution accuracy of the model. Considering the above factors, the proposed algorithm has a better dimensionality reduction effect and shows better clustering performance in image clustering tasks. BRIEF DESCRIPTION OF THE DRAWINGS

[0087] Figure 1 is the establishment process of the SGNLMF model.

[0088] Figure 2 is the solution framework diagram of the SGNLMF model.

[0089] Figure 3 is the framework diagram of the VF-CNO algorithm.

[0090] Figure 4 is the clustering framework diagram based on SGNLMF.

[0091] Figure 5 is the solution flow chart of SGNLMF. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0092] For the purpose, technical solutions and advantages of the invention to be more clearly understood, the following further detailed description of the present invention is provided in conjunction with the accompanying drawings and embodiments. The specific embodiments described herein are only for the technical solutions of the present invention and are not limited to the present invention.

[0093] The present invention provides a technical solution: An image non-negative matrix factorization clustering method based on neurodynamics.

[0094] Its process is as Figure 5 shown, including the following steps:

[0095] Step 1: Vectorize the images in the image dataset, and based on this, construct a non-negative raw data matrix V = [v1, v2,..., v n∈R m×n ;

[0096] Step 2: Construct the sparse graph regularized non - negative low - rank matrix factorization model as follows:

[0097]

[0098] s.t. L = WH

[0099] L, W, H≥0

[0100] ||w j ||1 = 1, j = 1,..., r

[0101] card(S)≤s, rank(L)≤l

[0102] Step 3: Perform sparse graph regularized non - negative matrix factorization on the original data matrix;

[0103] Step 3.1: Low - rank matrix recovery. Use the bilateral random projection algorithm to extract the latent low - rank structure L∈R in the original data matrix V m×n .

[0104] (3.1.1) Set the reconstruction relative error threshold to μ, set the rank of the low - rank matrix L to l, and set the sparse range set of the sparse matrix S to s;

[0105] (3.1.2) Initialize the relevant matrix and parameters. The low - rank matrix L0 = V, the sparse matrix S0 = 0, and the iteration number k = 0;

[0106] (3.1.3) To improve the computational efficiency of the algorithm, adopt a modification strategy based on the power - type, and use the matrix to calculate the bilateral random projection instead of V. On the one hand, and V have the same singular vectors. On the other hand, has a faster singular value decay rate. Based on this, new expressions for U1, U2, and can be obtained:

[0107]

[0108]

[0109] Then, based on the orthogonal triangular decomposition method (QR decomposition), represent the matrices U1 and U2 as:

[0110] U1 = Q1R1, U2 = Q2R2

[0111] Then L and S can be finally represented as:

[0112]

[0113] S k =Ω k (V - L k )

[0114] where Ω k (V - L k ) is the matrix element hard threshold operator, that is, project a matrix onto a set of elements and only keep the first k elements of V - L. Iteratively update based on the above update rules: k

[0115] k = k + 1, A2 = U2,

[0116] (3.1.4) Then update the low - rank matrix and the sparse matrix as follows:

[0117]

[0118] (3.1.5) Calculate the reconstruction error: If then execute step (3.1.6), if then execute step (3.1.3);

[0119] (3.1.6) Terminate the iteration to obtain the low - rank matrix L = L k , and the sparse matrix S = S k .

[0120] Step 3.2: Perform non - negative matrix factorization on the low - rank matrix L.

[0121] (3.2.1) Perform non - negative projection on the low - rank matrix L obtained in step 3.1, that is, assign zero to the non - negative elements in the low - rank matrix L;

[0122] (3.2.2) Input the low - rank matrix L obtained in (3.2.1) and normalize it so that its value range is [0, 1];

[0123] (3.2.3) Set the error tolerance ε, the threshold parameter Θ T , the number of particles N in the population, the maximum number of iterations k max , the regularization parameters α, β, the number of nearest neighbors ρ = 5, and the value of the decomposition dimension r;

[0124] (3.2.4) Use the K - nearest neighbor algorithm to construct a nearest - neighbor graph and calculate the graph Laplacian matrix L ap = D ap - Y ap ; where Y ap is the symmetric weight matrix, and D ap is a diagonal matrix, and the diagonal elements are the column sums of Y ap .​

[0125] (3.2.5) Randomly initialize the vectorized matrix ω i (0), h i (0), Γ(0), i = 1,..., N, the iteration number k = 0;

[0126] (3.2.6) Set the personal best solution as:

[0127] (3.2.7) Set the global best solution as: p * (0) = argmin[ω(0), h(0), Γ(0)], i = 1,..., N;

[0128] (3.2.8) Iteratively update w and h according to VF - RNN:

[0129]

[0130] where, P Ω (·) is a projection operator, and the rule is as follows:

[0131]

[0132] where, can be calculated as follows:

[0133]

[0134]

[0135] where, α is the graph regularization parameter, L ap is the graph Laplacian matrix, β is the Tikhonov regularization parameter. Γ = (γ1, γ2,..., γ n ) is the Lagrange multiplier of the equality constraint. can be calculated as:

[0136]

[0137] N particles perform local search according to the VF - RNN dynamic equation to obtain the personal best solution

[0138] ((3.2.9) If then p i (k + 1) = p i (k) If then

[0139] (3.2.10) Update the global optimal position of the population: 1 ≤ i ≤ N;

[0140] (3.2.11) The standard PSO update rules are as follows:

[0141] v i (k + 1) = ωv i (k) + c1r1(p i -x i (k)) + c2r2(g - x i (k))

[0142] x i (k + 1) = x i (k) + v i (k + 1)

[0143] Update the initial state of each particle according to the above rules;

[0144] (3.2.12) Calculate the evaluation function of particle diversity:

[0145]

[0146] (3.2.13) Compare the obtained Θ with the set threshold Θ T If Θ < Θ T , then execute:

[0147]

[0148] where H i and L i represent the upper and lower bounds of x i respectively, and the parameter η can be obtained according to the following wavelet function:

[0149]

[0150] where ψ ∈ [-2.5a, 2.5a], a = exp(10(k / k max )) and k max represents the maximum number of iterations.

[0151] (3.2.14) If k ≤ k max and ||p * (k + 1) - p * (k)|| ≥ ε,, then jump to step (3.2.9). If k > k max and ||p * (k + 1) - p * (k)|| < ε,, then the condition is satisfied, the loop terminates, and the global optimal solutions w and h are output;

[0152] (3.2.15) Convert the obtained \(w\) and \(h\) in (3.2.14) back into matrix form to obtain the optimal matrices \(W\) and \(H\).

[0153] Step 4: Perform clustering based on the K-means algorithm. Input the coefficient matrix \(H\) obtained in Step 3 into the K-means classifier to achieve the clustering task for unlabeled image data.

[0154] (4.1) Input the coefficient matrix \(H\) obtained in Step 3 into the K-means classifier to achieve the clustering task for unlabeled image data;

[0155] (4.2) Calculate the accuracy \(AC\) and the normalized mutual information \(NMI\) to measure the clustering performance:

[0156] 1. Accuracy (AC)

[0157] Accuracy is a commonly used clustering evaluation metric, representing the ratio of the number of successfully predicted samples to the total number of samples. Its mathematical definition is:

[0158]

[0159] where \(s\) i and \(r\) i represent the labels obtained by clustering and the true labels of the original data respectively, \(n\) represents the total number of samples; \(\delta(x, y)\) is a comparison function, which is 1 when \(x = y\) and 0 otherwise; \(map(x)\) is a label mapping function that projects the sample labels to the equivalent labels in the dataset. Generally, the Kuhn-Munkres algorithm is used to determine the optimal mapping. The higher the value of \(AC\), the better the clustering effect of the algorithm.

[0160] 2. Normalized Mutual Information (NMI)

[0161] Mutual information is used to measure the similarity between the predicted result and the true result. It should be noted that NMI is based on the mutual information function and the entropy function to evaluate the clustering result. Its mathematical definition is:

[0162]

[0163] where \(c\) represents the number of clusters, \(n\) i represents the number of samples in cluster , represents the number of samples in cluster , \(n\) i，j represents the number of samples that belong to both cluster and cluster . The higher the value of NMI, the more prominent the clustering performance of the algorithm.

[0164] The effects of the present invention can be further illustrated by the following experiments.

[0165] 1. Simulation conditions

[0166] All the experiments of the present invention were carried out on a PC equipped with an i7-6700 CPU with a main frequency of 3.40 GHz, and the simulation software was MATLAB 2019a.

[0167] 2. Experimental content

[0168] In order to verify the effectiveness of the method proposed by the present invention in the application of medical image clustering, a clustering experiment was carried out on the MedMNIST medical image dataset in this embodiment. In addition, in order to be general, supplementary experiments were also carried out on the USPS and ORL datasets in this embodiment to increase the data diversity.

[0169] (1) MedMNIST is a standard medical image dataset. In this experiment, 100 chest X-rays, breast ultrasounds, fundus photographs, retinal OCTs, and dermoscopy images were selected from it for testing. The pixel size of each image was uniformly modified to 28×28.

[0170] (2) USPS is a standard dataset from the US Postal Service handwritten digit image set. This dataset includes a total of 9,298 handwritten digit images from "0" to "9", all of which are 16×16 in size.

[0171] (3) ORL is a face image dataset created by Olivetti Research Laboratory. This dataset includes 40 groups of facial expression images of different people under different lighting conditions, with 10 images in each group. The original pixel size of the images was 112×92, and it was modified to 32×32 in this experiment.

[0172] It should be noted that the dataset was preprocessed in this embodiment, and Table 1 summarizes the information of the three datasets.

[0173] Table 1 Information of the three datasets

[0174]

[0175] In this embodiment, the proposed SGNLMF algorithm was applied to the three datasets for clustering experiments, and it was compared with the following methods:

[0176] 1. K-means: directly perform clustering analysis on the dataset without any other processing, and use it as a control experiment.

[0177] 2. PCA: Principal Component Analysis. This method can effectively extract the main components of data, thereby achieving the purpose of dimensionality reduction.

[0178] 3. NMF: Non-negative Matrix Factorization algorithm. It decomposes the data matrix into two non-negative low-dimensional sub-matrices, thereby obtaining a low-dimensional representation of the original data.

[0179] 4. GNMF: Graph Regularized Non-negative Matrix Factorization algorithm. It imposes a graph regularization constraint on the basis of NMF.

[0180] 5. SGNLMF: Sparse Graph Regularized Non-negative Low-rank Matrix Factorization algorithm proposed in the present invention. It imposes low-rank constraint, sparse constraint, Tikhonov regularization constraint and graph regularization constraint on the basis of NMF.

[0181] Among them, K-means is a widely used clustering algorithm. PCA was proposed by Wold S, Esbensen K, Geladi P, etc. in the article "Principal component analysis". NMF is introduced in detail in "Learning the parts of objects by non-negative matrix factorization". Cai D, He X, Han J, etc. proposed GNMF in the article "Graph regularized nonnegative matrix factorization for data representation".

[0182] In this experiment, p types of images are randomly selected from the dataset (p does not exceed the number of clusters in the dataset), and then they are mixed into the original data matrix V, and the number of clusters is fixed at p. In addition, the dimension r of the low-dimensional data obtained by all algorithms is set to p. Although sometimes, when r is not equal to p, the clustering effect may be better. But due to limited space, r = p is set in this experiment. p types of samples in the dataset are randomly selected for testing, and each experiment is repeated 10 times and the average value is taken.

[0183] Table 2 Clustering results on the MedMNIST dataset

[0184]

[0185] Table 3 Clustering results on the USPS dataset

[0186]

[0187] Table 4 Clustering results on the ORL dataset

[0188]

[0189] Table 2-4 gives the recognition rates and normalized mutual information of each algorithm on three different image datasets when different numbers of categories p are selected. Generally speaking, the clustering results (AC and NMI) of the SGNLMF algorithm proposed in this invention are better than those of other algorithms. In addition, as the number of categories p increases, the types and numbers of images added to the experiment also increase, and the recognition rate of the clustering algorithm shows a downward trend to a certain extent. This is because this experiment belongs to a clustering experiment without prior knowledge guidance, and the performance will decline when the experimental data increases, and this result also coincides with our intuitive feeling. Observing Table 2, it can be found that in the experiment on the MedMNIST dataset, the performances of several algorithms are relatively close, and the results of multiple experiments remain unchanged. The main reason is that when p is small, the amount of experimental data is small, and the dimensionality reduction effect is not significant. Observing the data in Table 4, when p becomes larger and the amount of data increases, the performance of the K-means algorithm will be significantly inferior to other algorithms, which shows that dimensionality reduction of high-dimensional data can effectively improve the clustering performance of the algorithm.

Claims

1. A non - negative matrix factorization clustering method for images based on neural dynamics, characterized in that According to the following steps: Step 1: Construct the original data matrix: Vectorize the images in the image dataset and construct a non - negative original data matrix \(V = [v_1, v_2, \ldots, v\) n \in\mathbb{R} m×n , where \(m\) is the dimension of the image samples and \(n\) is the number of image samples; Step 2: Construct a sparse graph regularized non - negative low - rank matrix factorization model as follows: where \(W\in\mathbb{R}\) m×r is the basis matrix, \(H\in\mathbb{R}\) r×n is the coefficient matrix, \(L\in\mathbb{R}\) m×n is the low-rank matrix, \(S\in\mathbb{R}\) m×n is the sparse matrix, \(l\) represents the range set of the rank of matrix \(L\), \(s\) represents the sparse range set of matrix \(S\), \(rank(L)\) represents the rank of matrix \(L\), \(card(S)\) represents the number of non-zero elements of matrix \(S\), \(\alpha\) is the graph regularization parameter, \(L\) ap is the graph Laplacian matrix, \(\beta\) is the Tikhonov regularization parameter; Sparse graph regularized non-negative low-rank matrix factorization imposes additional low-rank constraints, Tikhonov regularization constraints, sparse constraints and graph regularization constraints during the matrix factorization process to mine the sparse and effective feature information behind the data; Step 3: Perform sparse graph regularized non - negative low - rank matrix factorization on the original data matrix; Step 3.1: Low-rank matrix recovery: Use the bilateral random projection algorithm to extract the potential low-rank structure \(L\in\mathbb{R}\) in the original data matrix \(V\), described as: m×n , described as: Step 3.2: Perform non-negative matrix factorization on the low-rank matrix L obtained in Step 3.1 using the variable scale factor collaborative neural dynamics optimization method VF-CNO to obtain the basis matrix W ∈ R m×r and the coefficient matrix H ∈ R r×n ; where the basis matrix W contains the localized features in the original data matrix, the coefficient matrix H represents the low-dimensional representation of the original data. At the same time, in order to enhance the decomposition performance of the model, Tikhonov regularization constraints, sparse constraints, and graph regularization constraints are additionally imposed during the matrix factorization process, and the sparse constraint is incorporated into the objective function. The newly formed objective function is as follows: where Γ = (γ1, γ2, …, γ n ) are the Lagrange multipliers of the equality constraints; ω, h perform the vectorization operation on the W and H matrices as follows: ω = vec(W) = (w 11 w 21 …w m1 ; w 12 w 22 …w m2 ; …; w 1r w 2r …w mr ) T h = vec(H) = (h 11 h 21 …h r1 ; h 12 h 22 …h r2 ; …; h 1n h 2n …h rn ) T Introduce a variable - scale - factor - type recurrent neural network model VF - RNN to constrain the KKT point search for the global optimization problem. The VF - RNN model is expressed as the following differential equation: where P Ω (·) is a projection operator. According to the above vector equation, its scalar form can be described as: Among them, the convergence factor ∈ = (∈1, ∈2, …, ∈ n ) T ; Ψ(∈, t) is a vector composed of variable scale factor functions ψ(∈ j , t) (j ∈ {1, 2, …, n}); for the variable scale factor function ψ(∈, t), it is a function related to the convergence factor ∈ and time t, and takes one of the following two forms: Exponential form: ψ(∈,t) = ∈e t Or power - exponential form: ψ(∈,t) = t ∈ + ∈ Step 4: Perform clustering based on the K - means algorithm: Input the coefficient matrix H obtained in Step 3 into the K - means classifier to achieve the clustering task for unlabeled image data.

2. A non - negative matrix factorization clustering method for images based on neural dynamics according to claim 1, characterized in that The said Step 3.1 Includes the following steps: 3.1.1: Set the reconstruction relative error The threshold is μ, set the rank of the low-rank matrix L to l, and set the sparse range set of the sparse matrix S to s; 3.1.2: Initialize the correlation matrix and parameters, L0 = V, S0 = 0, k = 0; 3.1.3: Perform iterative update: k = k + 1, A2 = Y2, 3.1.4: If Jump to step 3.1.2; 3.1.5: Update matrices L and S: 3.1.6: If then execute step 3.1.

7. If then execute step 3.1.3; 3.1.7: Terminate the iteration to obtain the low-rank matrix L = L k , and the sparse matrix S = S k。 3. A non - negative matrix factorization clustering method for images based on neural dynamics according to claim 1, characterized in that The said Step 3.2 Includes the following steps: 3.2.1: Input the low - rank matrix L obtained in Step 3.1 and normalize it so that its value range is [0, 1]; 3.2.2: Set the error tolerance ε, the threshold parameter Θ T , the number of particles N in the population, the maximum number of iterations k max , the regularization parameters α, β, the number of nearest neighbors ρ = 5, and the value of the decomposition dimension r; 3.2.3: Construct a neighbor graph using the K-nearest neighbor algorithm and calculate the graph Laplacian matrix L ap = D ap - Y ap ; where Y ap is a symmetric weight matrix, and D ap is a diagonal matrix whose diagonal elements are the column sums of Y ap ; 3.2.4: Random Initialization: ω i (0), h i (0), Γ(0), i = 1, …, N; 3.2.5: Set the personal optimal solution as: 3.2.6: Set the global optimal solution as: p * (0) = argmin[ω(0), h(0), Γ(0)], i = 1, …, N; 3.2.7: Perform iterative update: 3.2.7.1: N particles perform local search according to the VF-RNN dynamic equation to obtain the personal optimal solution 3.2.7.2: If then p i (k + 1) = p i (k); If then 3.2.8: Update the global optimal position of the population: 3.2.9: Calculate the initial state of the network according to the particle swarm update rule: 3.2.9.1: Introduce the Particle Swarm Optimization algorithm PSO. Assume that the particle swarm moves to the global optimum within a preset spatial range, and the movement speed of each particle is described by a velocity vector v; define x i =(x i1 ,...,x in ) T as the current position of the particle; v i =(v i1 ,...,v in ) T represents the movement speed of the particle; p i =(p i1 ,...,p in ) T represents the best position of the individual during the search process; g = (g1,..., g n ) T represents the global best position among all particles during the search process; The description of the standard PSO is as follows: v i (k + 1) = ωv i (k) + c1r1(p i -x i (k)) + c2r2(g - x i (k)) x i (k + 1) = x i (k) + v i (k + 1) Among them, ω ∈ [0, 1] is a fixed parameter; c1 and c2 are acceleration factors that determine the degree to which the particle velocity is affected by p i and g; r1 and r2 are random numbers with a value range of [0, 1]; k represents the number of iterations; Update the initial state of the particle based on the above standard PSO update rule ; [[ID= ​ where N represents the number of populations; Substitute into the above evaluation function; 3.2.9.3: Compare the obtained Θ with the set threshold Θ T If Θ ≥ Θ T , then directly proceed to the next step; if Θ < Θ T , then execute: Among them, H i and L i represent the upper and lower bounds of x i respectively, and the parameter η is obtained according to the following wavelet function: where ψ ∈ [-2.5a, 2.5a], a = exp(10(k / k max )) and k max represents the maximum number of iterations; 3.2.10: If k ≤ k max and ‖p * (k + 1) - p * (k)‖ ≥ ε, then jump to step 3.2.

7. If k > k max and ‖p * (k + 1) - p * (k)‖ < ε, then the condition is satisfied, the loop terminates, and the global optimal solutions w and h are output; ​ 4. A non - negative matrix factorization clustering method based on neural dynamics according to claim 1, characterized in that ​ ​ ​ ​ ​ Among them, s i and r i respectively represent the labels obtained by clustering and the true labels of the original data, n represents the total number of samples; δ(x, y) is a comparison function, which is 1 when x = y and 0 otherwise; map(x) is a label mapping function that projects the sample labels into equivalent labels in the dataset, and the Kuhn-Munkres algorithm is used to determine the optimal mapping; the higher the value of AC, the better the clustering effect of the algorithm; ​ ​ where c represents the number of clusters, n i represents the number of samples in cluster , and represents the number of samples in cluster . n i,j represents the number of samples that belong to both cluster and cluster . The higher the value of NMI, the more prominent the clustering performance of the algorithm.

Images