An Improved Multioperator-Based Constrained Differential Evolution for Optimal Power Allocation in WSNs

Optimal power allocation (OPA), which can be transformed into an optimization problem with constraints, plays a key role in wireless sensor networks (WSNs). In this paper, inspired by ant colony optimization, an improved multioperator-based constrained adaptive differential evolution (namely, IMO-CADE) is proposed for the OPA. The proposed IMO-CADE can be featured as follows: (i) to adaptively select the proper operator among different operators, the feedback of operators and the status of individuals are considered simultaneously to assign the selection probability; (ii) the constrained reward assignment is used to measure the feedback of operators; (iii) the parameter adaptation is used for the parameters of differential evolution. To extensively evaluate the performance of IMO-CADE, it is used to solve the OPA for both the independent and correlated observations with different numbers of sensor nodes. Compared with other advanced methods, simulation results clearly indicate that IMO-CADE yields the best performance on the whole. Therefore, IMO-CADE can be an efficient alternative for the OPA of WSNs, especially for WSNs with a large number of sensor nodes.


Introduction
Due to their cost-effectiveness, easy deployment, and intelligence, wireless sensor networks (WSNs) have obtained considerable attention in the last few decades [1]. Nowadays, WSN is one of the promising technologies for real-world applications [2,3]. To improve their performance and lifetime, various computational intelligence techniques have been used recently for the design of WSNs [4], such as particle swarm optimization (PSO) [5], genetic algorithms (GAs), ant colony optimization (ACO) [6], machine learning algorithms [7], lion optimization [8], negatively correlated search [9], and Fibonacci tree optimization algorithm [10].
In WSNs, there are usually many sensor nodes, which have a very limited power supply. Therefore, among various design issues, optimal power allocation (OPA) plays a key role in WSNs [11][12][13][14][15]. The OPA can be transformed into an optimization problem with constraints. Due to its nonlinear and computationally expensive (for the correlated observations) nature, the optimization of the transformed problem is difficult for traditional/analytical methods; especially, there are a large number of sensor nodes with correlated observations due to the expensive computation cost [16]. In the literature, the use of nature-inspired optimization methods for solving the OPA has received more consideration [6]. For example, in [17], PSO was used for the OPA. Boussaïd et al. [18] presented a constrained bio-geography-based optimization (BBO) and differential evolution (DE)-namely, CBBO-DE-for the OPA. In [19], combined with the penalty technique, three nature-inspired methods-cat swarm optimization (CSO), cuckoo search (CS), and PSOwere compared for solving the OPA. In [20], Tsiflikiotis et al. developed a hybrid method with teaching-learning-based optimization and Jaya algorithms. Lee compared three

Background
In this section, first, the transformed optimization problem of the OPA is described. Then, the original DE is briefly introduced.

Problem Formulation
In this work, a WSN with a fusion center and K sensor nodes is considered. The optimization problem of the OPA can be formulated as follows (Since this paper mainly focuses on the development of efficient optimizer for the OPA, the data fusion problem formulation of the OPA in WSNs is not mentioned herein; the interested reader can refer to the work in [17] for the details.) [17,18]: subject to where ε is the required fusion error probability threshold, m indicates the deterministic signal, G k is the amplifier gain at node k, and e is the L-length vector with all ones. The covariance matrix is ∑ K = A T ∑ v A + ∑ w , where A = diag(H 1 G 1 , . . . , H K G K ), ∑ v , and ∑ w are the observation and receiver noise covariances, respectively. H k (k = 1, · · · , K) is the channel fading coefficient. In this work, two situations (i.i.d. and correlated observations) are considered.

Independent Observations
If the local observations and the receiver noise are both i.i.d., the probability of fusion error can be simplified to The inequality in (2) can be expressed as follows: where β = 2 m Q −1 (ε) and Q(·) is the complementary Gaussian cumulative distribution function. δ v means the variances of the observation noise and δ w represents the receiver noise. In this work, the channel fading coefficient H k follows an exponential distribution (i.e., Rayleigh fading) with unit mean [18].

Correlated Observations
If the sensor observations are spatially correlated, the observation noise covariance matrix ∑ v with the correlation degree ρ can be formulated as follows: The inequality in (2) is where d j = d(j − 1), j = 1, . . . K, which means the sensor nodes are equally spaced along a straight line. Since ∑ v is not diagonal, it is difficult to evaluate ∑ −1 K in closed form. Note that, under the correlated situation, the computational complexity for the calculation of Equation (5) is O(K 2 ), which means that the OPA with correlated observations is computationally expensive when K is large.

Differential Evolution
DE is a simple and efficient global optimizer for numerical optimization problems [23]. The flowchart of the original DE is shown in Figure 1, where there are four main proceduresi.e., population initialization, differential mutation, crossover, and selection.

Initialization
Suppose the population contains N p solutions. For the OPA, each solution x i = (x i 1 , · · · , x i K ), i = 1 · · · , N p, is a K-dimensional real-valued vector. Initially, each variable x i,k is randomly generated as where i = 1 · · · , N p and k = 1, · · · , K. The lower and upper bounds of x i k are 0 and 10, respectively. rndreal(0, 10) is a real-valued random number generator from (0, 10). Note that, for the OPA, x i k represents the amplifier gain at node k, i.e., G k .

Mutation
In DE, the core operator is the differential mutation. The commonly used "DE/rand/1" can be formulated as where v i is generated mutant vector; r1, r2, r3 ∈ {1, N p} satisfy r1 = r2 = r3 = i, and F ∈ (0, 1) is a scaling factor.

Crossover
After mutation, the crossover is performed between the target vector x i and the mutant vector v i . The binomial crossover is where u i is the trial vector, rand j ∈ {i, K} is a randomly generated integer, and Cr ∈ [0, 1] is the crossover rate.

Selection
After generating the trial vector, DE performs the one-to-one tournament selection to select the solution to survive into the next generation: The mutation, crossover, and selection will repeat until the termination criterion is met.

Related Work
As stated in [17], if the observations are i.i.d, the OPA can be derived analytically in closed form. However, if the observations are spatially correlated, the evaluation of ∑ −1 K in (6) cannot be obtained analytically in closed form. Therefore, recently, various nature-inspired methods are used for the OPA.
In [17,24], PSO was employed to obtain the OPA with the constraint of a required fusion error probability threshold. To deal with the constraint, the exterior penalty function with two additional parameters (i.e., positive penalty parameter r k and non-negative constant q) is used. In [18], Boussaïd et al. hybridized BBO with DE for the OPA, where the penalty function presented in [25] was used to handle the constraint. In [19], with the dynamically modified penalty functions, CSO, CS, and PSO were empirically compared for the OPA. The results indicated that (i) when the number of the sensor nodes is small, PSO gets better results; however, (ii) CSO outperforms CS and PSO when the network consists of a large number of sensors. Tsiflikiotis et al. [20] extended their previous work [19] and presented a hybrid TLBO and Jaya (namely, TLBO-Jaya) to solve the OPA, where the proposed TLBO-Jaya was compared with TLBO, Jaya, heat transfer search algorithm, PSO, and BBO-DE. In [21], Lee used PSO, ABC, and ACO R for optimizing the OPA, where Deb's three feasibility rules [26] were used to replace the penalty function for the constraint handling. Li et al. [22] presented a constrained DE with multiple mutation strategies, where probability matching and constrained credit assignment techniques were used.
The abovementioned methods obtained promising results for solving the OPA. However, most of them are used to solve the OPA with a small number of sensor nodes, e.g., K ≤ 50. In addition, the penalty functions usually introduce new parameters, which may influence the performance significantly. Therefore, there is still much room to develop advanced, effective, and efficient numerical methods for the OPA, especially for the OPA with a large number of sensor nodes.

Our Approach: IMO-CADE
With the main objective of finding an alternative solution for minimizing the total power consumed by the WSN, in this section, the proposed IMO-CADE is introduced in detail, including the motivations, operator pool, constrained reward assignment, improved multioperator selection technique, and parameter adaptation.

Motivations
As mentioned in Sections 1 and 3, different nature-inspired methods are presented to solve the OPA. To deal with the constraints, most methods use penalty functions. However, it is difficult to set the penalty coefficients for different problems. In [21], it adopts the feasibility rules [26]. However, the feasibility rules favor the feasible solutions, which may result in trapping into the local feasible regions.
In [27], the empirical study stated that the performance of DE may significantly be improved by the mutation operators for the OPA. It is difficult to choose the most suitable mutation operators for different cases in the OPA. In our previous work [22], probability matching was used for operator selection. However, to assign the probabilities of different operators, only operator feedback is considered. Indeed, individual information should also be considered for the probability assignment.
Based on the above considerations, in this work, we propose an improved multioperator selection technique by considering both the operator feedback and individual information. In addition, to further improve the performance, parameter adaptation of DE is also used to release the setting of parameters by the users. The main techniques used in IMO-CADE are elaborated in the following subsections.

Operator Pool
In the DE literature, various mutation operators have been proposed [28]. In this work, to efficiently solve the OPA, two mutation operators presented in [29,30] were chosen as the operator pool-"DE/current-to-pbest/1" with archive and "DE/rand-to-pbest/1" with archive: • "DE/current-to-pbest/1" with archive: • "DE/rand-to-pbest/1" with archive: where r0, r1, r2 ∈ {1, N p} and r0 = r1 = r2 = i. F i is the scaling factor of the i-th solution; x pbest is the "pbest" solution randomly chosen from the top p% solutions of the current population;x r2 is a randomly selected solution from the union of the current population P and the archive A.

Boundary Constraint Handling
After generating the mutant vector with the mutation operator, some variables may be out of their boundary constraints, i.e., v i k / ∈ [0, 10]. In this situation, the following boundary constraint handling technique (BCHT) is used:

Constrained Reward Assignment
For solving the constrained optimization problems, there are three situations of the combined parent population P and child population C.
(i) Infeasible situation: All solutions in P ∪ C are infeasible. Under this situation, the fitness F (x) of each solution x is its overall constraint violation (CV), (ii) Semifeasible situation: P ∪ C contains both the infeasible and feasible solutions. In this situation, the solutions in the parent and child populations are combined. Then, for each solution, its objective function and CV are normalized as suggested in [31]. Afterwards, the fitness is as follows: where f nor (x) and CV nor (x) are the normalized objective function and CV, respectively. The details can be found in [31]. (iii) Feasible situation: all solutions in P ∪ C are feasible. In this situation, the fitness is the objective function: According to the above fitness of each solution under situations, the relative fitness improvement is calculated as [22] where i = 1, · · · , N p, F best is the fitness of the best-so-far solution, x i is the target solution, and u i is its corresponding trial solution. Then, based on the relative fitness improvement, the reward R o of the o-th operator is calculated by

Improved Multioperator Selection
In order to adaptively select a suitable operator from the pool, the selection probability of each operator is assigned. In this work, both the operator feedback and individual information are considered simultaneously.

Probability Based on Operator Feedback
The operator feedback is measured by its constrained reward, as shown in Equation (18). Based on the reward, probability matching [32] is used to assign the probability of each operator. First, for each operator o, its quality is updated as where α ∈ (0, 1] is the adaptation rate. Then, the probability is calculated as where p min ∈ (0, 1) is the minimal probability of each operator.

Probability Based on Individual Information
Generally, for different solutions in the population, different mutation operators may be more suitable. In this work, two mutation operators, i.e., "DE/current-to-pbest/1" with archive (o = 1) and "DE/rand-to-pbest/1" with archive (o = 2), are used in the pool. "DE/current-to-pbest/1" with archive performs the local search around the target solution, and hence, is more suitable to the better solutions. "DE/rand-to-pbest/1" with archive is able to provide more diversity, and hence, is more suitable to the worse solutions. Based on this consideration, first, the solutions in P are sorted based on the feasibility rules [26] from the best to the worst. Then, according to the individual information, the probability is assigned as where p o,i is the probability for the o-th operator of the i-th solution, o = 1, 2, and i = 1, · · · , N p.

Final Probability Calculation
Finally, inspired by ACO [33], by combining operator feedback and individual information, the probability is calculated by

Parameter Adaptation
To further improve the performance of IMO-CADE for the OPA (Algorithm 1), the parameter adaptation technique proposed in [29] is also used. Initially, for each target solution x i , its associated parameters Cr i and F i are generated as and where "Gaussian(µ Cr , 0.1)" is a Gaussian random number generator with mean µ Cr and standard deviation 0.1, and "Cauchy(µ F , 0.1)" is a Cauchy random number generator with location µ F and scale 0.1. At the first generation, µ Cr = 0.5 and µ F = 0.5.

Algorithm 1 Pseudo-code of IMO-CADE
Input: Algorithm parameters: N p, µ Cr , µ F , c, p min , α, NFEs max ; WSN parameters: K, , ρ, γ 0 Output: The best solution x best = {G best 1 , · · · , G best K } 1: Initialize the population P with N p solutions; 2: Evaluate the population P and set NFEs = N p; 3: For each operator, set q o = 0 and p o = 1/2; 4: Set A = Φ; 5: while NFEs < NFEs max do 6: Sort the solutions in P based on the feasibility rules from the best to the worst; 8: For each solution in P, set p o,i according to Equation (21); 9: For each solution, calculate p o,i based on Equation (22); 10: for i = 1 to N p do 11: Select the operator o i based on the selection probability; 12: Calculate Cr i and F i based on Equations (23) and (24), respectively; 13: Generate u i and deal with the violated variables based on Equation (13); 14: Evaluate u i , NFEs = NFEs + 1, and C ← u i ; 15: end for 16: M ← P ∪ C; 17: For each solution x in M, calculate F (x) under the current situation; 18: for i = 1 to N p do 19: if F (u i ) ≤ F (x i ) then 20: A ← x i ; 21: CalculateF i based on Equation (17)  After that, if the trial solution u i is better than its target solution x i , the parameters Cr i and F i are saved in S Cr and S F , respectively. Subsequently, µ Cr and µ F are updated as follows: and where c = 0.1 is a parameter, mean A (S Cr ) calculates the arithmetic mean of S Cr , and mean L (S F ) calculates the Lehmer mean of S F [29].

Framework of IMO-CADE
By integrating the abovementioned techniques, the pseudo-code of IMO-CADE is shown in Algorithm 1, where NFEs is the current number of fitness evaluations and NFEs max is the maximal NFEs. First, the population P is initialized randomly; then, each solution is evaluated. The initial quality of each operator is set to q o = 0 and p o = 1/2. The archive A is set to be empty. In the main loop, it works as follows: • At each generation, S Cr and S F are set to be empty. • In lines 8-9, the selection probability p o,i is calculated. • In line 11, for each solution, one operator o i is selected based on the selection probability and roulette wheel selection.
• In lines 12-14, the trial solution is generated according to the selected operator and the generated parameters. The violated variables are handled based on the BCHT. • In lines 16-17, the parent population P and the child population C are combined. Then, the feasibility situation is checked, followed by the transformed fitness calculation by Equations (14)

Remarks
It is worth noting that IMO-CADE is an improved version of our previous work (PM-MDE) in [22]. However, there are several important differences between IMO-CADE and PM-MDE: (1) The core difference is that, in IMO-CADE, both the operator feedback and individual information are considered together to update the operator selection probability, whereas in PM-MDE, only the operator feedback is used.

Results and Analysis
In this section, the proposed IMO-CADE is used to solve the OPA with different sensor nodes under independent or correlated observations. Observation signal-to-noise ratio (SNR), γ 0 = 10 dB. In this work, the WSN consists of a fusion center and many spatially separated sensors, which perform amplify-and-forward local processing of their observations independently. The OPA problem was considered for the decision fusion of a deterministic signal in an inhomogeneous WSN. The channel fading coefficients follow an exponential distribution (i.e., Rayleigh fading) with unit mean. All algorithms were implemented in C++. The simulations were executed on a desktop PC with an Intel Xeon E5-2620 processor @ 2.40 GHz, 32 GB RAM, using the Windows 10 64-bit OS.

Comparison with Other Advanced DEs
To evaluate the performance of IMO-CADE, it is compared with five other DE variants for the OPA with different sensor nodes under i.i.d. observations. The five DE methods are as follows: (i) SaDE [34], (ii) JADE [29], (iii) OrSHADE [35], (iv) COLSHADE [36], and (v) CADE. The first four DE methods obtained promising results in the literature (To make a fair comparison, the parameters of the four DE variants are set as the same as those used in the original literature. All methods use NFEs max = 3000), for example, COLSHADE is one of the winners in the CEC-2020 competition for the real-world constrained single-objective optimization. CADE is a variant of IMO-CADE. The only difference between IMO-CADE and CADE is that only the operator feedback is used for operator selection in CADE.
The detailed results are reported in Table 1 for K = 10, 20, 50, 100, 150, 200, and = 0.1, 0.01, 0.001, where the "mean ± standard deviation" value is provided for each method on each case. All results are averaged over 30 runs. The best mean values are highlighted in boldface. Based on the Wilcoxon test at α = 0.05 by the KEEL software [37], in the last row of Table 1, "w/t/l" indicates that IMO-CADE performs significantly better than, similar to, or significantly worse than its competitors in w, t, or l cases, respectively. Moreover, according to the multiple-problem analysis by the Friedman test [37], the averaging rankings of the six DE variants are plotted in Figure 2. The multiple-problem analysis by the Wilcoxon test for IMO-CADE vs. other DEs is given in Table 2.   Therefore, based on the above analysis, it is clear that IMO-CADE yields the best results than other compared DE variants on the whole, especially for the OPA with a large number of sensor nodes.

Comparison with Reported Results
In this section, the results of IMO-CADE are compared with the reported results of other methods, i.e., CBBO [18], CDE [18], CBBO-DE [18], and PM-MDE [22], under both i.i.d. and correlated observations. Note that the five compared methods consume the same NFEs max = 3000 for fair comparison.

Under i.i.d. Observations
The results of the five compared methods for K = 10, 20, 50, and = 0.1, 0.01, 0.001 are shown in Table 3. In addition, the results of the multiple-problem analysis by the Wilcoxon and Friedman tests (To calculate the statistical results, "NF" is set to be 1000.) are provided in Table 4 and Figure 3, respectively.
From Table 3, we can see that IMO-CADE can consistently obtain the best results compared with the other four methods. Additionally, the results in Table 4 show that our approach can perform significantly better than CBBO, CDE, and PM-MDE. Figure 3 confirms that IMO-CADE obtains the best ranking, followed by CBBO-DE and PM-MDE.

Under Correlated Observations
As mentioned in Section 2.1.2, when the observations are correlated, i.e., ρ = 0, the optimization of OPA is more complicated and computationally expensive. In this section, IMO-CADE is compared with CBBO, CDE, CBBO-DE, and PM-MDE with K = 10, 20, 50, ρ = 0.01, 0.1, 0.5, and = 0.1, 0.01, 0.001. The detailed results are reported in Table 5. Moreover, the multiple-problem statistical results based on the Wilcoxon and Friedman tests are provided in Table 6 and Figure 4, respectively.
The results in Tables 4 and 5, and Figure 4 again confirms that IMO-CADE can get the best overall results in all cases. It is able to significantly outperform CBBO, CDE, CBBO-DE, and PM-MDE based on the Wilcoxon test. Further, it obtains the best ranking based on the Friedman test.  Based on the comparison with the reported results of other methods under both i.i.d. and correlated observations with different situations (i.e., number of sensor nodes, fusion errors, and correlation degrees), the results clearly show that the proposed IMO-CADE can provide significantly better results than the compared methods. Therefore, IMO-CADE can be an effective and efficient alternative to solve the OPA of WSNs.

Discussions
In previous subsections, the superiority of IMO-CADE is verified on both the i.i.d. and correlated observations with different numbers of sensors. Herein, the gain allocated to each sensor of different methods is discussed. Table 7 reports the results comparison between analytical and different numerical methods with K = 10, γ 0 = 10 dB, and ρ = 0. Table 8 reports the results of different numerical methods with K = 10, γ 0 = 10 dB, and ρ = 0.1.
The best results of f (x) are highlighted in boldface. In the two tables, IMO-CADE is compared with the reported results of CBBO, CDE, and CBBO-DE in [18]. Table 7. Comparison between analytical and different numerical methods to gain allocation with K = 10, γ 0 = 10 dB, and ρ = 0.  Tables 7 and 8, the following can be seen:

Sensor
• Based on the gain allocated to each sensor, the sensors with poor channels can be turned off to save system power. Based on the gain allocated to each node, we can decide that the sensors with good channel fading coefficients are assigned to more power; on the other hand, sensors with poor channels are allocated less power. • For both i.i.d. and correlated observations, IMO-CADE provides the best results of f (x) compared with other methods.

Conclusions
The OPA is one of key issues in developing WSNs. Hence, the design of effective optimization techniques for the OPA is important and obtains more attention. Based on this consideration, in this paper, we propose the IMO-CADE for the OPA. To effectively solve the OPA, in IMO-CADE, both operator feedback and individual information are considered simultaneously for multioperator selection. In addition, the parameter adaptation, efficient BCHT, and constrained reward assignment techniques are developed. IMO-CADE is extensively compared with other methods under various situations. The results clearly indicate the following: • The proposed modifications in IMO-CADE can improve its performance for the OPA under different situations. • With respect to the performance of the overall system power f (x), IMO-CADE is superior to other methods in all cases, especially for the WSN with a large number of sensor nodes. • Considering gain allocation, the numerical results of IMO-CADE agree well with the analytical results. • IMO-CADE can be an effective alternative for the OPA and other complex optimization problems of WSNs.
When the observations of sensor nodes are correlated, the OPA is computationally expensive-especially if there are a large number of sensor nodes. In the near future, combining IMO-CADE with surrogate models [38] is a promising way to efficiently solve the OPA under correlated observations. In addition, due to the similarity between the OPA and the lighting in long tunnels [39], IMO-CADE will be used for the lighting problem in long tunnels in future work. In this work, the simulations are performed to evaluate the performance of IMO-CADE. In future work, another interesting direction is using IMO-CADE in real networks.
The source code of IMO-CADE can be obtained from Dr. Gong upon request.

Acknowledgments:
The authors would like to thank the reviewers and editor for their constructive comments.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: