A method for identifying cancer cooperative driving pathways based on a weighted mutation matrix

Abstract

The application discloses a method for identifying cancer collaborative driving pathways based on a weighted mutation matrix, and comprises the following steps: 1) setting a nonlinear CMDP model; 2) setting an MPGA: 2.1) chromosome coding and initial population; 2.2) setting a fitness function; 2.3) selecting an operator; 2.4) a mutation operator; and 2.5) an MPGA algorithm specific process. The method can identify more biological information, has strong expansibility and practicability, and can identify genes enriched on multiple important driving pathways.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of bioinformatics, and particularly to a method for identifying cancer cooperatively driven pathways based on a weighted mutation matrix. BACKGROUND

[0002] The pathogenesis of cancer is very complex, which seriously threatens human health. With the rapid development of high-throughput sequencing, some large-scale cancer genomics projects, such as The Cancer Genome Atlas (TCGA) and the International Cancer Genome Consortium (ICGC), have provided a large amount of data for the study of the formation and development of cancer.

[0003] Distinguishing driver genes (whose mutations have a major impact on the proliferation and spread of cancer cells) from passenger genes (whose mutations have no contribution to the development of cancer) is a major challenge in cancer research. Previous studies have shown that since driver genes are enriched in many important cellular signaling or regulatory functional pathways, any abnormality of a driver gene in a certain pathway will lead to dysfunction of the pathway, thereby forming a tumor. Currently, there are three main directions for identifying driver pathways: identifying single driver pathways, identifying cooperatively driven pathways, and identifying public driver pathways and specific driver pathways for pan-cancer. The present method mainly identifies cooperatively driven pathways.

[0004] In the past decade, some progress has been made in this research field. Leiserson et al. proposed the Multi-Dendrix method by extending the Dendrix method for identifying single driver pathways. According to the characteristics of high coverage (i.e., the genes in the driver pathway cover a large number of samples) and high exclusivity (i.e., the mutation of a gene in a single pathway can lead to cancer formation) of cancer driver pathways, a summation model of multiple pathways was designed to find cooperatively driven pathways. Zhang et al. proposed the CoMDP model based on somatic mutation data using the co-occurrence between a pair of driver pathways and the exclusivity within a single driver pathway. They also used integer linear programming to solve the model. Cai et al. proposed the ICMDP model, which uses gene expression data to adjust the co-occurrence between a pair of driver pathways, and designed a genetic algorithm to solve the model. Yang et al. proposed the CDPath method, which identifies core driver modules from mutation data and clusters them in a pathway interaction network to discover cooperatively driven pathways. Wang et al. proposed the CDPathway method, which first constructs a heterogeneous network of genes, miRNAs, and pathways, and then performs matrix decomposition on the heterogeneous network to discover cooperatively driven pathways.

[0005] In the above research methods, most of the methods only use somatic mutation data and identify the model based on the high coverage and high exclusivity of the driver pathway, or use existing pathway prior knowledge for identification, but some driver gene combinations are not highly exclusive, and these special cases cannot be well identified, and limited pathway prior knowledge may limit the discovery of new driver gene combinations. Therefore, the MPGA-CMDP method for identifying non-linear model collaborative driver pathways based on numerical matrix is proposed to solve the above problems. SUMMARY

[0006] The purpose of the present application is to overcome the shortcomings of the prior art and provide a non-linear intelligent optimization method for identifying cancer collaborative driver pathways. This method can provide more useful biological information, has strong expansibility and strong practicability, and can identify genes enriched in multiple important driver pathways.

[0007] The technical scheme of the present application is:

[0008] A method for identifying cancer collaborative driver pathways based on a weighted mutation matrix, comprising the following steps:

[0009] 1) Set a non-linear CMDP model:

[0010] Suppose there is a somatic mutation matrix A copy number variation matrix A cancer sample gene expression matrix And a normal sample gene expression matrix All rows of the first three matrices represent cancer samples, and the rows of the fourth matrix represent normal samples. Their columns represent candidate genes G S , G C and G E , for each element in the matrix s ij =1 indicates that the jth gene has a mutation in the ith sample, otherwise no mutation occurs; for each element in the matrix , if the copy number of the jth gene in the ith sample is amplified or deleted, then c ij ∈{-2,-1,1,2}, where i=1,2,...,|P|, j=1,2,...,|G C |, otherwise c ij =0; for each element in the matrices and , the size of the numerical value represents the expression level of the jth gene in the ith sample.

[0011] First, a binary mutation matrix is constructed, where G A= G S ∪ G C , if s ij and c ij are not equal to 0, then a ij = 1, i = 1, 2,..., |P|, j = 1, 2,..., |G A |, otherwise 0, then a difference fold change matrix D is constructed using the cancer sample gene expression matrix and the normal sample gene expression matrix Specifically, d ij = |log2e ij - m j |, where m j is the median of the logarithm of all expression values of the jth gene in the normal sample gene expression matrix, and finally a numerical mutation matrix A is generated using the quartile method |P|×|G| The specific formula is as follows:

[0012]

[0013] Where Q3(d -j ) is the third quartile of the jth column value of the difference fold change matrix D.

[0014] Suppose M and N are two gene non-repeated sub-matrices of A |P|×|G| , three formulas are designed based on the numerical mutation matrix A |P|×|G| to calculate the three features of the collaboration driven pathway, namely common coverage, non-common coverage and mutual exclusivity within a single pathway.

[0015] The specific calculation of common coverage is as follows, let λ M (i) = max{a ij |a -j ∈ M}, λ N (i) = max{a ij |a -j ∈ N}, then the common coverage of the ith sample of the two matrices M and N is:

[0016]

[0017] The common coverage of the whole two matrices M and N is

[0018] The specific calculation of non-common coverage is as follows: suppose Γ(g) = {a i- |a ig ≠ 0, i = 1, 2,..., |P|} represents the sample set of gene a -g mutation, and represent the sample sets covered by matrices M and N respectively, then the non-common coverage of the ith sample of the two matrices M and N is:

[0019]

[0020] The common coverage of the two matrices M and N is:

[0021]

[0022] where |U| = |Γ(M)∪Γ(N)|, |Q| = |Γ(M)∩Γ(N)|;

[0023] The specific calculation of the intra-path mutual exclusion is as follows: a measurement method v(i, M) (resp. v(i, N)) based on the coefficient of variation is designed to estimate the dispersion degree of the weight in the ith row, i = 1, 2,..., |P|, and the calculation formula is: where σ(i, M) (resp. σ(i, N)) and μ(i, M) (resp. μ(i, N)) represent the standard deviation and the average weight of the ith row of the matrix M and the matrix N, respectively, and the mutual exclusion of the two matrices M and N as a whole is:

[0024]

[0025]

[0026] Based on the above definitions of common coverage CO, non-common coverage DIF, and intra-path mutual exclusion ME, a nonlinear CMDP model is constructed: given a numerical mutation matrix A |P|×|G| and a parameter K, 0 < K < |G|, determine two sub-matrices and k1+k2=K, to maximize the weight function W(M, N):

[0027] W(M, N) = CO(M, N) - DIF(M, N) + ME(M) + ME(N);

[0028] 2) Set the MPGA algorithm:

[0029] 2.1) Chromosome coding and initial population:

[0030] A chromosome is a solution to the problem, in the MPGA algorithm, a chromosome is represented by , which is a solution made by two matrices M and N, where The specific coding method of a chromosome is as follows: 1. Calculate the wa value of each gene in the matrix A |P|×|G| and arrange the genes in descending order according to the calculated value, the wa calculation formula is 2. Based on the arranged order, construct an initial chromosome for every K genes. For an initial chromosome, the first... One gene is Part, after One gene is Part, such a chromosome is divided into For each combination, repeat this step. This is the first step to obtain the initial population;

[0031] 2.2) Define the fitness function:

[0032] The fitness function is used to evaluate the quality of a solution and plays a guiding role in the evolution of the population, given a chromosome. and Representing matrix A respectively |P|×|G| Two submatrices of size |P|×k1 and |P|×k2, The fitness function is used to evaluate the quality of a solution. The higher the fitness value, the better the solution. The fitness function is defined as follows:

[0033]

[0034] 2.3) Selection Operator:

[0035] In the MPGA algorithm, different selection operators are used at different stages of population evolution. In the initial stage of evolution, the exponential ranking selection algorithm is used to increase the diversity of the population by increasing the probability of selecting poor chromosomes. In the remaining stages of evolution, roulette wheel selection is used to ensure that offspring evolve rapidly in the optimal direction.

[0036] 2.4) Mutation Operator:

[0037] 2.4.1) Single-point exchange mutation:

[0038] Given a chromosome exist and In the process, a pair of interchange sites s1 and s2 are randomly selected, where 1 ≤ s1 ≤ k1 and 1 ≤ s2 ≤ k2. By exchanging the genes at sites s1 and s2, a new offspring X′ is obtained. If Fitness(X′) is greater than Fitness(X), then X is replaced with X′; otherwise, it remains unchanged.

[0039] 2.4.2) Two-point crossover mutation:

[0040] Given a chromosome exist and Two pairs of interchange sites s were randomly selected from each of the two groups.11 ,s 12 , and s 21 ,s 22 , 1≤s 11 ≤s 12 ≤k1, 1≤s 21 ≤s 22 ≤k2, by swapping s 11 and s 21 Genes at the locus and s 12 and s 22 For each gene at a given locus, a new offspring X″ is generated. If Fitness(X″) is greater than Fitness(X), then X is replaced with X″; otherwise, it remains unchanged.

[0041] 2.4.3) Greedy Mutation:

[0042] Given a chromosome Construct a candidate gene set, and sequentially in and Greedy mutation is performed on top, from Randomly delete a gene g from the process to generate an initial offspring individual, i.e. Then, based on a greedy strategy, from the candidate gene set... A gene is randomly selected from a subset of the size and inserted into... To obtain better offspring individuals, that is The same process will then occur. The process is performed on the above and a new offspring X″′ is generated. If Fitness(X″′) is greater than Fitness(X), then X is replaced with X″′; otherwise, it remains unchanged.

[0043] 2.5) Specific process of MPGA algorithm:

[0044] In Algorithm 1, the MPGA algorithm is summarized. Step 1 sets the parameters used in the MPGA algorithm. Step 2 generates multiple initial populations. Steps 3 to 22 allow these populations to evolve independently. Step 5 performs selection on the parent populations to generate offspring. Steps 6 to 12 attempt to improve the worst individuals in the populations. Steps 13 to 17 perform mutation operations to construct new populations. The entire evolutionary process is controlled by two parameters: the maximum number of generations (maxg) and the threshold (maxt) that keeps the optimal solution unchanged. Finally, the best individual among the multiple populations is output.

[0045]

[0046]

[0047] The technical solution studies the identification method from the perspective of the characteristics of the cooperation-driven pathway, and the main contributions are as follows:

[0048] (1) A numerical mutation matrix without parameter control is constructed by using multi-omics data;

[0049] (2) A non-parametric nonlinear identification model is designed based on three characteristics of the cooperation-driven pathway, which has higher identification accuracy and practical value;

[0050] (3) A multi-population evolutionary algorithm is proposed to solve the model, and the new identification algorithm is more suitable for solving complex combined data;

[0051] (4) The single cancer cooperation-driven pathway found by the whole technical method contains more important driving pathways and genes.

[0052] This method provides more useful information, strong scalability, strong practicality, and can identify more important genes enriched in different driving pathways. BRIEF DESCRIPTION OF DRAWINGS

[0053] Figure 1 An example diagram for constructing a mutation matrix in the embodiment is shown;

[0054] Figure 2 An example diagram for constructing a numerical matrix in the embodiment is shown;

[0055] Figure 3 An example diagram for running the cooperation-driven pathway for glioblastoma GBM with a driving pathway size of 4 is shown;

[0056] Figure 4 An example diagram for running the cooperation-driven pathway for glioblastoma GBM with a driving pathway size of 8 is shown. DETAILED DESCRIPTION

[0057] The application will be further described in detail below in conjunction with the drawings and specific embodiments, but is not limited to the application.

[0058] The example experiment was performed on a computer (Intel(R) Core(TM) i7-6700 3.40GHz CPU, 32G memory), the operating system was Windows 10, and the compilation and running tool was R4.1.2.

[0059] This example is described for the cooperation-driven pathway identification problem.

[0060] Embodiment:

[0061] A nonlinear intelligent optimization method for identifying cancer cooperation-driven pathways based on a numerical matrix, comprising the following steps

[0062] 1) Set the nonlinear CMDP model:

[0063] Two cancer datasets, namely glioblastoma (GBM) and ovarian cancer (OVCA), were used. Mutational data of the two cancers were from Zhao et al., and mutation and normal expression data of the two cancers were from TCGA (http: / / tcgadata.nci.nih.gov / tcga / ). The specific data are shown in Table 1:

[0064] Table 1: Data of two cancers (sample x gene)

[0065]

[0066] In this example, the somatic mutation matrix of glioblastoma (GBM) is used Copy number variation matrix Cancer sample gene expression matrix and normal sample gene expression matrix All rows of the first three matrices represent cancer samples, and the rows of the fourth matrix represent normal samples. Their columns represent the entire candidate gene G S , G C and G E , respectively. For each element in matrix , s ij = 1 indicates that the jth gene is mutated in the ith sample, otherwise it is not mutated. For each element in matrix , if the jth gene has significant copy number amplification or deletion in the ith sample, then c ij ∈ {-2, -1, 1, 2} (i = 1, 2,..., |P|, j = 1, 2,..., |G C |), otherwise c ij = 0. For each element in matrices and , the size of the numerical value represents the expression level of the jth gene in the ith sample.

[0067] First, a binary mutation matrix is constructed, where G A = G S ∪ G C . If s ij and c ij are not equal to 0, then a ij = 1 (i = 1, 2,..., |P|, j = 1, 2,..., |G A |), otherwise 0. After the construction, genes with a mutation rate less than 1% in the sample are filtered out, as shown in Figure 1 , where |P| = 88, |G A| = 919, then use the cancer sample gene expression matrix and the normal sample gene expression matrix to construct the difference fold matrix Specifically, d ij = |log2e ij -m j |, where m j is the median of all expression values of the jth gene in the normal sample gene expression matrix, and finally a numerical mutation matrix A is generated using the quartile method |P|×|G| , and the specific formula is as follows:

[0068]

[0069] Where Q3(d -j ) is the third quartile of the jth column value of the difference fold matrix D, and the numerical mutation matrix is as shown in Figure 2 ;

[0070] Suppose M and N are two gene non-repeating A |P|×|G| sub-matrices, and based on the numerical mutation matrix A |P|×|G| Three formulas are designed to calculate the three features of the collaborative driving path (common coverage; non-common coverage; mutual exclusivity within a single path);

[0071] The specific calculation method of common coverage is as follows, let λ M (i) = max{a ij |a -j ∈M}, λ N (i) = max{a ij |a -j ∈N}, then the common coverage of the ith sample of the two matrices M and N is:

[0072]

[0073] The common coverage of the two matrices M and N as a whole is

[0074] The specific calculation method of non-common coverage is as follows: suppose Γ(g) = {a i- |a ig ≠0, i = 1, 2,..., |P|} represents the sample set of gene a -g mutation, and represent the sample sets covered by matrices M and N respectively, then the non-common coverage of the ith sample of the two matrices M and N is:

[0075]

[0076] The common coverage of the two matrices M and N as a whole is:

[0077]

[0078] where |U| = |Γ(M)∪Γ(N)|, |Q| = |Γ(M)∩Γ(N)|;

[0079] The specific calculation method of the intra-path mutual exclusion is as follows: a measurement method v(i, M) (resp. v(i, N)) based on the coefficient of variation is designed to estimate the dispersion degree of the weight in the ith (i = 1, 2,..., |P|) row, and the calculation formula is: where σ(i, M) (resp. σ(i, N)) and μ(i, M) (resp. μ(i, N)) represent the standard deviation and the average weight of the ith row of the matrix M and the matrix N respectively, and the mutual exclusion of the two matrices M and N as a whole is:

[0080]

[0081]

[0082] Based on the above definitions of the common coverage CO, the non-common coverage DIF and the mutual exclusion ME within the single path, a nonlinear CMDP model is constructed: given a numerical mutation matrix A |P|×|G| and a parameter K (0 < K < |G|), two sub-matrices and are determined to maximize the weight function W(M, N):

[0083] W(M, N) = CO(M, N) - DIF(M, N) + ME(M) + ME(N);

[0084] 2) Set the MPGA algorithm:

[0085] 2.1) Chromosome coding and initial population:

[0086] A chromosome is a solution to the problem, in the MPGA algorithm, a chromosome is represented by , which is a solution made by two matrices, where The specific coding method of a chromosome is as follows: 1. Calculate the wa value of each gene in the matrix A |P|×|G| , and arrange the genes in descending order according to the calculated value, the wa calculation formula is 2. According to the order after the arrangement, an initial chromosome is constructed every K genes, for an initial chromosome, the first genes are part, the last genes are part, such a chromosome is divided into This step is repeated for each combination times to obtain an initial population;

[0087] 2.2) Set fitness function:

[0088] Fitness function is used to evaluate the quality of a solution, which plays a guiding role in the evolution of the population. Given a chromosome and represent two sub-matrices of size |P| x k1 and |P| x k2 of matrix A |P|×|G| respectively, which is used to evaluate the quality of a solution, the greater the fitness value, the better the corresponding solution, and the fitness function is defined as follows:

[0089]

[0090] 2.3) Selection operator:

[0091] In the MPGA algorithm, different selection operators are used in different stages of population evolution. In the initial stage of evolution, the exponential ranking selection algorithm is used to improve the diversity of the population by increasing the probability of selecting poor chromosomes. In the remaining stages of evolution, roulette wheel selection is used to ensure that the offspring quickly evolves in the optimal direction;

[0092] 2.4) Mutation operator:

[0093] 2.4.1) Single-point crossover mutation:

[0094] Given a chromosome In and , a pair of crossover sites s1 (1≤s1≤k1) and s2 (1≤s2≤k2) are randomly selected, and a new offspring X' is obtained by exchanging the genes at sites s1 and s2. If Fitness(X') is greater than Fitness(X), replace X with X', otherwise keep it unchanged;

[0095] 2.4.2) Two-point crossover mutation:

[0096] Given a chromosome In and , two pairs of crossover sites s 11 ,s 12 (1≤s 11 ≤s 12 ≤k1) and s 21 ,s 22 (1≤s21 ≤s 22 ≤k2), by exchanging s 11 and s 21 genes at s 12 and s 22 sites, a new offspring X" is obtained, if Fitness(X") is greater than Fitness(X), then X is replaced by X", otherwise it remains unchanged;

[0097] 2.4.3) Greedy mutation:

[0098] Given a chromosome A candidate gene set is constructed, and greedy mutation is performed on and respectively, a gene g is randomly deleted from to generate an initial offspring individual, i.e. Then, according to the greedy strategy, a gene is randomly selected from the candidate gene set of size to be inserted into to obtain a better offspring individual, i.e. The same process is then performed on and a new offspring X'" is generated, if Fitness(X"') is greater than Fitness(X), then X is replaced by X"', otherwise it remains unchanged;

[0099]

[0100] In Algorithm 1, the MPGA algorithm is summarized, Step 1 sets some parameters used in the MPGA algorithm, Step 2 generates the initial multiple populations, Steps 3 to 22 evolve these populations independently, Step 5 performs the selection operator on the parent population to generate the offspring population, Steps 6 to 12 try to improve the worst individual in the population, Steps 13 to 17 perform mutation operations to construct a new population, the entire evolution process is controlled by two parameters, i.e., the maximum evolution generation maxg and the threshold maxt for keeping the optimal solution unchanged, finally, the best individual in the multiple populations is output, the parameter settings of MPGA are as follows: the population size N is set to the maximum evolution generation maxg is equal to 300 and the threshold maxt for keeping the optimal solution unchanged is equal to 100;

[0101]

[0102]

[0103]

[0104] AsFigure 3 and Figure 4 The results of the identification of the driving pathways with a scale of 4 and 8 are shown in Table 2. The maximum fitness values of the running results indicate that the identified genes are driving genes, and they are enriched in two different driving pathways, respectively.

Claims

1. A method for identifying collaborative driver pathways in cancer based on a weighted mutation matrix, characterized in that, Comprising the following steps: 1) Set up the nonlinear CMDP model: Suppose there is a somatic mutation matrix , a copy number variation matrix , a cancer sample gene expression matrix and a normal sample gene expression matrix , all rows of the first three matrices represent cancer samples, rows of the fourth matrix represent normal samples, and their columns represent candidate genes respectively , and , for each element in matrix , indicates that the th gene has a mutation in the th sample, otherwise no mutation occurs; for each element in matrix , if the th gene has copy number amplification or deletion in the th sample, then , where , , otherwise for each element in matrices and , the magnitude of its value represents the expression level of the th gene in the th sample; Firstly, a binary mutation matrix is constructed where if and are not equal to 0, then , , otherwise 0, then the differential fold matrix is constructed using the cancer sample gene expression matrix and the normal sample gene expression matrix , specifically where is the median of all expression values of the th gene of the normal sample gene expression matrix after taking the logarithm, and finally the numerical mutation matrix is generated using the quartile method , the specific formula is as follows: ， wherein refers to the third quartile of the jth column values of the fold change matrix refers to the third quartile of the jth column values of the fold change matrix Assume and are two gene non-redundant sub-matrices, based on numerical mutation matrix Three formulas are designed to calculate the three features of the collaboration-driven pathway, namely, common coverage, non-common coverage, and mutual exclusivity within a single pathway; The specific calculation of the common coverage is as follows, let , , and The common coverage of the first samples of the two matrices is: ， and The common cover of both matrices as a whole is ; The non-common coverage is calculated as follows: Assume represents the set of genes mutated in the sample set and represents the set of genes covered by the matrix and represents the set of samples covered by the matrix respectively, then the non-common coverage of the th sample of the two matrices is given by: ， and The common cover of both matrices as a whole is: ， wherein , , ; The calculation of the intra-row exclusivity for the single-path case is as follows: a measure based on the coefficient of variation is designed to estimate the dispersion of the weights in the row, given by where and denote the standard deviation and the mean weight of the first row of the matrix and the matrix respectively, and the intra-row exclusivity of the two matrices as a whole is given by and ​ ， ， Based on the above definitions of co-coverage , non-co-coverage , and mutual exclusivity within a single pathway , a nonlinear CMDP model is constructed: given a numerical mutation matrix and a parameter , , determine two sub-matrices and , to maximize a weight function : ； 2) Set up the MPGA algorithm: 2.1) Chromosome coding and initial population: A chromosome represents a solution to the problem. In the MPGA algorithm, a chromosome is represented by... Let represent a solution formed by two matrices, where , , , , The specific encoding method for a chromosome is as follows:

1. Calculate the matrix Each gene in The genes are then sorted in descending order of their calculated values. The calculation formula is ; 2. According to the arranged order, each One gene constructs an initial chromosome. For an initial chromosome, the first... One gene is Part, after One gene is Part, such a chromosome is divided into For each combination, repeat this step. This is the first step to obtain the initial population; 2.2) Set up the fitness function: The fitness function is used to evaluate the quality of a solution, and plays a guiding role in the evolution of the population. Given a chromosome X=( , ), and represent two sub-matrices of size and of the matrix , is used to evaluate the quality of the solution, the greater the fitness value, the better the corresponding solution, and the fitness function is defined as follows: ； 2.3) Selection operator: In the MPGA algorithm, different selection operators are used in different stages of population evolution; in the initial stage of evolution, the exponential ranking selection algorithm is used to improve the diversity of the population by increasing the probability of selecting poor chromosomes, and in the remaining stages of evolution, roulette wheel selection is used to ensure that offspring evolve quickly in the optimal direction; 2.4) Mutation operator: 2.4.1) Single-point crossover mutation: Given a chromosome , , , , , , in and , a pair of crossover sites are randomly chosen and , , , a new offspring is obtained by swapping the genes at the and sites , if is greater than , replace with , otherwise keep it unchanged; 2.4.2) Two-point crossover mutation: Given a chromosome , , , , , , , randomly select two pairs of crossover sites and in and , , , respectively, by exchanging the genes at sites and and the genes at sites and , a new offspring is obtained, if is greater than , replace with , otherwise keep it unchanged; 2.4.3) Greedy mutation: Given a chromosome , , , , , , , a candidate gene set is constructed, and greedy mutation is performed on and in turn, a gene is randomly deleted from to generate an initial offspring, i.e. , , , then a gene is randomly selected from a subset of size from the candidate gene set , and the gene is inserted into to obtain a better offspring, that is , the same process is performed on , and a new offspring is generated, if is greater than , replace with , otherwise keep it unchanged; 2.5) MPGA algorithm specific process: In Algorithm 1, the MPGA algorithm is summarized, Step 1 sets the parameters used in the MPGA algorithm, Step 2 generates the initial multiple populations, Steps 3 to 22 evolve these populations independently, Step 5 executes the selection operator on the parent population to generate the offspring population, Steps 6 to 12 attempt to improve the worst individual in the population, Steps 13 to 17 perform mutation operations to construct new populations, the entire evolution process is controlled by two parameters, namely the maximum evolution generation maxg and the threshold value maxt that keeps the optimal solution unchanged, finally, the best individual in multiple populations is output;

Images