Robust Optimal Placement and Sensitivity Analysis of Multiple Transmitting Antennas Under Uncertainty: A New Surrogate-Based Evolutionary Algorithm

: Transmitting antenna positioning or transmitter placement is a well-known NP-hard optimization problem pertinent to communication systems. Furthermore, it is of practical importance to yield an optimal location of transmitters to ensure low sensitivity with respect to potential uncertainties. Notwithstanding, incorporating uncertainties in the optimization problem can highly increase the computational expenses. This paper aims at the development of a new reduced-cost algorithm for a multi-objective robust transmitter placement under uncertainties. Toward this end, a new hybrid surrogate-metaheuristic approach is developed using Grey Wolf Optimizer (GWO) and the Kriging surrogate in a mathematical framework of a robust dual-surface model. The proposed algorithm is also able to analyze the sensitivity of the obtained optimal results. The latter is achieved by obtaining the bootstrapped confidence regions without extra simulation experiments. The paper investigates the performance of the proposed algorithm for robust optimal placing of two, three, and four transmitters, under uncertainties concerning the transmitting antenna gain. The results demonstrate the utility and the efficiency of the proposed method in rendering the robust optimal design and analyzing the sensitivity of the transmitter placement problem under practically acceptable computational efforts.


Introduction 1.1 Background of study and motivations
The antenna placement problem or cell planning problem involves locating and configuring infrastructure for cellular wireless networks. From candidate site locations, a set needs to be selected against the objectives related to issues such as financial cost and service provision [1]. Antenna placement in a multi-antenna platform involves a manual process that is challenging and time-consuming and may result in a suboptimal placement, leading to an inferior performance of communication systems. The search space becomes exponentially large with respect to the number of antennas to be placed (|search space| = , where is the number of allowable locations for each antenna, and is the number of antennas) [2]. The antenna placement problem is well-known to be NP-hard. Its solution has been attempted using heuristic approaches such as evolutionary algorithms (e.g., Genetic 1. Department of Electrical Engineering, Faculty of Engineering, Chulalongkorn University, Bangkok 10330, Thailand. 2 . Instituto de Telecomunicações Campus Universitário de Santiago P-3810-193 AVEIRO -PORTUGAL. 3. School of Engineering, University of Glasgow, Glasgow, G12 8QQ, UK. * Corresponding authors: Lunchakorn Wuttisittikulkij (Lunchakorn.W@chula.ac.th), Amir Parnianifard (amir.p@chula.ac.th). algorithm (GA) [3]), or swarm intelligence algorithms (e.g., Ant Colony Optimization (ACO) [4], Particle Swarm Optimization (PSO) [5]). However, the literature lacks studies involving new swarm intelligence optimizers such as Grey Wolf Optimizer (GWO) [6] in solving the optimal transmitter placement problem for wireless coverage systems. However, considering uncertainty as a source of variability in the computational model can significantly increase the cost of estimating the statistical figures of merit such as the mean and the standard deviation (Std) of the system response. In such cases, the application of fast surrogate models such as polynomial regression, Kriging, or radial basis function, integrated into robust design optimization flow, allows for maintaining low computational cost [7]- [11]. This paper proposes a new hybrid approach that employs Kriging interpolation surrogates and GWO as a new swarm intelligence methodology. The presented technique is combined with a dual-surface design to obtain (robust) optimal positions of the base stations in the transmitter placement problem under uncertainty. The aforementioned combination is developed to achieve the design robustness along with the accuracy in optimal allocation of the base station transmitting antennas at reduced computational cost. Also, we are interested in analyzing the sensitivity of the optimization results, in particular, computing the confidence regions caused by the randomness without extra function evaluations, which is ensured by employing the samples that are already acquired for the initial training of the model.

Related works
A transmitter placement scheme relies on a certain propagation model to assess the quality of a given transmitter allocation. For all cellular network systems, a major design step is to select the locations for the base station transmitters and to set up the optimal configurations such that the coverage of the desired area with sufficiently strong radio signals is high, and the deployment costs are low [12]. Interactive or automatic approaches have been used to produce an optimized transmitter placement based on certain propagation models [13]- [16]. When choosing the propagation model, one often needs to achieve a balance between the computational cost and the prediction accuracy [17]. Practical solution methods mainly differ in the optimization objectives and the search algorithms. Some heuristic methods have been developed for solving the transmitting placement problem. In [18], the authors adopt the ACO approach to optimize the transmitter locations to maximize the average received power. However, the benchmarking demonstrates that almost the same antenna locations are obtained with PSO and GA, whereas ACO requires more than ten times objective function evaluations (slow convergence), compared to PSO and GA. As most evolutionary algorithms, GA has proven very capable in yielding high-performance antenna designs [2], [19], [20]. Determination of the base station position for optimum signal coverage using particle swarm optimizer (PSO) has been investigated in [21]. In [22], the PSO integrated with Radial Basis Function (RBF) has been applied to obtain the optimal placement of multiple transmitters by maximizing the overall signal coverage in an objective function, and control the intersection of transmitters in a constraint. However, estimating the statistical measures due to uncertainty (the main source of antenna parameter variability) can considerably increase the computational cost. Consequently, the development of a less expensive approach, which-at the same time-offers a sufficient reliability in searching for robust optimum is a practical necessity from the perspective of real-world communication system design. To handle this, Kriging surrogate (also known as the Gaussian process) has been applied for different types of engineering design problems [23]- [25]. In the electromagnetic community, it has been applied, among others, in [26]- [28]. In [29], the Kriging surrogate is used for transmitter location optimization in both single-transmitter and three-transmitter cases. In [30], a comparative study between Kriging and GA for optimal transmitter location in an indoor environment has been performed from the fields scattered in the environment.

Main contributions
The major contributions of this study can be summarized as follows: • The development of a new hybrid Kriging/GWO approach combined with robust design optimization to estimate the Pareto front with respect to the two objectives: design robustness and its (nominal) accuracy. We also provide a sensitivity analysis of transmitting antennas placement under uncertainty.
• Conducting numerical studies concerning multiobjective constrained robust optimization through minimization of the power consumption required by signal transmission, and maximization of the signal coverage as the objective functions. The control over the signal overlap is treated as the design constraint.
• In cases with a large number of receiver test points, the computational cost of the optimization process is very high because of evaluating the intensity of signal received in all receiver points. The method developed in this work does not require exhaustive receiver evaluations in each optimization iteration.
• Carrying out extensive numerical experiments demonstrating the efficiency and reliability of the proposed algorithm in yielding robust optimal placement of multiple transmitters.
The rest of this paper is organized as follows. Section 2 provides the preliminaries required by the proposed algorithm including a definition of the free space propagation model and a formulation of the relevant multi-objective robust design optimization concepts applied in this study. The algorithmic framework of the proposed reduced-cost approach to solving and analyzing the considered optimization task is presented in Section 3. The verification examples for the robust optimal placement of two, three, and four transmitting antennas under uncertainty using the proposed approach are presented in Section 4. Section 5 concludes the paper.

Problem statement 2.1 Free space propagation model
This paper adopts the free space propagation model, which is widely used in the studies of placement problems [19], [31].
In this study, we consider a transmitting antenna and a receiving antenna separated by a large distance , . For propagation distances , much larger than the antenna size, the far field of the electromagnetic wave dominates all other components. Let and be the respective gain functions of the transmitting and receiving antenna for the direction of transmission. If the transmitted power is , the power radiated in the direction of the receiver, per unit solid angle, would be ( . If is the carrier wavelength, the receiving antenna would present a receiving cross-section ( where , 2 = {( − ) 2 + ( − ) 2 } stands for the Euclidean distance (2-norm) between locations of the receiver and the transmitter in a two-dimensional design space. The maximum operating range is determined by the signal-tonoise ratio of the detector system. If it is possible to ignore the effect of the earth on the propagation of the wave, and if is constant, it would be possible to operate the receiving system satisfactorily everywhere within the surface with radius as below: where is the minimum detectable signal for the receiver, and the surface with radius is called the "free-space coverage pattern for one-way transmission", see Figure 1.

Transmitter placement planning model
The planning model describes the environment and the mathematical model of the transmitter placement problem. The map for transmitter placement has two regions in a twodimensional area ( , ) including covered regions ( ) and placement regions ( ). The former shows the regions that need to be covered by a signal transmitted by antennas, and the latter represents the regions where the transmitters can be located. In this study, we consider the same free-space twodimensional map (ℤ 2 ) region for both covered and placement regions with grid resolution where ⊆ ℤ 2 , and ⊆ ℤ 2 . In the placement planning model, a receiver gathers signals from the transmitters. The connectivity is assessed by a signal threshold to maintain the quality of service. This paper uses a large set of receivers as test points for the coverage: a receiver ∈ has a position ( , ) ∈ with threshold . Let is the set of energy power that can be transmitted by each transmitter, then the transmitter placement problem is to create a set of transmitters = { = ( , , )| , ∈ , ∈ } and place its elements (i.e., the transmitters). In the current study, the transmitter placement problem is performed based on two objectives including maximization of percent coverage and minimization of Total Power Transmitted (TPT), when keeping the amount of coverage overlap (signal interference) under the predefined threshold. The latter is implemented as a design constraint. Percent coverage: A receiver ∈ is said to be covered by a transmitter ∈ when the signal strength is greater than the threshold, i.e., where is computed using Eq.(1). The value 1 indicates that the receiver is covered by at least one transmitter. Accordingly, the percent signal coverage of a set of transmitters can be calculated by: Note that the first objective is to maximize the percent coverage. Total power transmitted: The proper planning model needs to be designed to minimize the total power consumed for transmitting signals by all designed transmitters while achieving as high percent coverage as possible. Thus, in this paper, apart from maximizing the percent coverage (see Eq.(4)), we consider minimizing the total power transmitted, i.e., where is the power ( ) that employed in transmitter for transmitting the signal. Overlap: The coverage overlap between transmitters raises the issue of interference [12], [19]. To reduce the interference, we add a constraint to the model to keep the overlap under the predefined threshold. The relevant mathematical formulation of the overlap and the associated constraint will be explained in the next section.

Proposed Algorithm 3.1 Basic components of the optimization framework
In this study, we replace the statistical approach of Taguchi viewpoint with hybrid Kriging and GWO approach. Furthermore, we intersect two experimental designs (data sample sets). The first design is pertinent to decision variables (inner array), whereas the second one is for uncertain variables (outer array). Given the vector = (1,2, … , ) of sample points over the decision variables, and the vector = (1,2, … , ) of uncertainty scenarios, × input combinations are designed, and the true model is evaluated × times to collect relevant model responses, see Figure 2. Let be the × matrix of model responses. The mean and the standard deviation (Std) for each row of can be computed as The following Signal-to-Noise Ratio (SNR) as robustness criterion can be formulated: Since the goal is to minimize the model response, the formulation of the SNR in (8) has the opposite sign by Taguchi formulation [33].

Mathematical optimization model
Numerical optimization uses a compact mathematical model for describing the problem of concern. Here, we define the problem of multiple transmitting antenna placement under uncertainty in the framework of a robust multi-objective optimization. The goal is to obtain the set of transmitters' optimal positions on a 2D map, and the relevant optimal power for each positioned transmitter, is the set of power that can be transmitted by each transmitter. The objectives and constraints are defined as below: Subject to: In (9), the terms % and % denote the mean and the standard deviation of the coverage, estimated based on the variability of uncertain parameters in the model, using (6) and (7). This objective is formulated to enable maximization of the percent coverage as well as the robustness. The second objective function, Eq.(10), considers the minimization of the total power (expressed in mW) consumed by all transmitters. In other words, the optimization process aims at achieving the best possible placement of the transmitters, minimization of the overall power consumption required to obtain the signal coverage, as well as reducing the results sensitivity with respect to the parameter uncertainties. The terms and in Eq.(10) represent the number of transmitters in a model, and the maximum allowed power that can be allocated to each transmitter, respectively. The overlap and the intersection between the transmitters are considered in the constraint as shown in Eq. (11). The term Φ in Eq. (11) refers to the cumulative distribution function (CDF) of a standard normal distribution. The term 0 ≤ ≤ 1 can be defined by the designer and represents the allowed probability for average overlap of the transmitters in the model. The expressions . and . represent the mean and the standard deviation of the radius (see Eq. (2)). As both objective functions are expressed in the same logarithmic scale as 10 log( ), where ∈ [0,1], we can aggregate them as where ∈ [−2, +2] is the weighting factor. The optimum selection of depends on the designer's preferences and the characteristics of the optimization model. Sweeping allows for capturing the Pareto front, and for identifying trade-off designs between the considered objective functions [34]. Figure 3 provides the flow diagram of the proposed optimization approach. Its main steps are described below:

Algorithmic framework
Step 1. Initialize the model parameters.
The model parameters and constants that need to be adjusted at the beginning of the algorithm are shown in Table 1. Note that in the current study, we assume the transmitting antenna gain (̃) in Eq.(1) as an uncertain parameter that is uniformly varied in predefined range with known lower and upper bounds.

Start
Initialize the model parameters.
Construct crossed array sampling design for coverage.
Obtain resposnes for each combination of design and uncertain parameters.
Construct crossed array sampling design for power radius.
Fit one Kriging for mean and one for Std of responses.
Fit on e Kriging for mean and one for Std of responses.
Both surrogates applied in the overall objective function, Eq.  Step 2. Design of experiments using crossed array design. Design sets of experiments using the structure of crossed array design (see Figure 2) separately for two different functions, one for power absorbed in Eq.(1) and the second for the power radius in Eq. (2). For each function, also two sets of sample points are designed separately. The first set (over the design variables) in an inner array, and the second set (over an uncertain variable) in an outer array. The design variables for the absorbed power in Eq.(1) are , and , and for power radius in Eq. (2) is . The uncertain variable in both functions is ̃. Here, we apply the grid-sampling method in the class of space filling design [35] with resolution ∆≫ . In other words, the sampling grid is of considerably lower resolution than to reduce the number of required function evaluations. The information for the regions between the grid points will be obtained through interpolation using the Kriging surrogate.
Step 3. Obtain response for each combination of design and uncertain variables. Regarding the first crossed sampling design, we execute the power absorbed function using Eq.(1) and accordingly compute the "Covered(r)" using Eq.(3) with values 1 or 0. Besides, regarding the second crossed sampling design, we run the power radius function using Eq.(2) and obtain the relevant responses for each combination of design and uncertain parameters in the crossed array designs.
Step 4. Compute the statistical measures including mean and Std of responses. In this step, the mean and the standard deviation of coverage and power radius are computed for each input combination using Eq.(6) and (7) regarding a relevant crossed array design. Note that these statistical measures result from the variability of the uncertain parameter that was already considered in the previous steps.
Step 5. Construct Kriging surrogates. Four Kriging surrogates are constructed separately using the acquired input-output data pairs. The Kriging model is identified for: (i) the mean of the coverage, (ii) the standard deviation of the coverage, (iii) the mean of the power radius, and (iv) the standard deviation of the power radius.
Step 6. Validate surrogate models The surrogate model constructed in Step 5 has to be validated to ensure that its predictive power is sufficient for design optimization purposes. Here, validation is executed using the leave-one-out cross validation ( = 1) [36], [37], which works as follows. First delete the ℎ input combination and the relevant output from the complete set of the th combination ( = 1,2, … , ), i.e., to avoid the extrapolation by Kriging, we avoid dropping the sample points in the margin. The model is constructed using − 1 remaining rows and predicts the output for the left-out point ( −1 ). This is realized for all input combinations (sample points) and leads to computing predictions (̂). Finally, evaluating the standardized residuals are computed as Most of the standardized residuals should be within the interval −3 ≤ ≤ 3, and any observation outside of this interval (outlier) is potentially unacceptable with respect to its observed simulation output [36], [37].
Step 7. Set up the optimization model. The proposed robust optimization model using Eq.(9) through Eq.(12) is arranged.
Step 8. Run the GWO optimizer. The GWO optimizer is executed using the expression (12) as a fitness function, and (11) as a constraint. Here, we control the feasibility of using the death penalty for any point out of the feasible region in Eq. (11). Note that during the optimization run, all required expressions in both Eq.(11) and Eq.(12) are estimated by the relevant Kriging surrogates constructed in Step 5, so that no further evaluations of the original computational model are required. Another point is that to improve the computational efficiency of the optimization process, we do not investigate all receiver test points ⊆ ℤ 2 in each iteration. But instead, we consider a smaller set of receiver test points that randomly selected from the domain.
Step 9. Compute the two-sided BCI for an obtained robust optimal point. In stochastic simulation, each input combination X is replicated several times m ≥ 1. In the case of expensive simulations, the number possible of replications is smaller, therefore, parametric bootstrapping is unlikely to produce acceptable results (i.e., it rarely can find the exact distribution of the I/O simulation data) [24], [38], [39]. Here, to find the bootstrapped set of data, a model is resampled times ( = Table 1 The list of parameters applied in the proposed algorithm. The weighting parameter in Eq. (12) ---Initialize 1,2, … , ) (sampling with replacement). Let us assume that + is a robust optimal solution obtained from the procedure in Step 1 through Step 6. It is assumed that + is a robust optimal solution, which is obtained from the original (nonbootstrapped) surrogate. All output values at point + are estimated using all the bootstrapped surrogates. The distribution-free Bootstrapped Confidence Intervals (BCI) can be computed as below [38], [40]: Where 2 ⁄ gives two-sided BCI, whereas the Bonferroni's inequality [41] suggests that the type-I error rate of γ [38], [42] for each interval per output is divided by the number of outputs, here, the mean, Std, and the SNR. If the values of bootstrap estimate + * are sorted from low to high, then ⌊. ⌋ and ⌈. ⌉, respectively, denotes floor and ceiling function to achieve the integer part and round upwards. Here, it is assumed that the set of sample points is fixed and only old data to fit surrogate with sufficient replication is available, whereas new simulation replicating is very expensive. This augmented bootstrapping approach does not imply extra computational cost due to resampling, and the required simulation run to find a bootstrapped set of data [43]- [45]. ( = 1,2, … , ) denotes the set of sample points and each is repeated times ( = 1,2, … , ). We assume that the original set of data obtained from the original simulation model is available (size × ) when is the number of scenarios for uncertainty and is the number of input combinations. Moreover, the augmented bootstrapping procedure is sketched in Algorithm 1.

Experimental Benchmark Problems 4.1 The algorithm setup
In this section, three cases featuring a different number of transmitting antennas (e.g., two, three, and four base stations) are considered to evaluate the performance of the proposed algorithm for reduced-cost transmitting antenna placement under uncertainty. The transmitters are assumed to be homogeneous (i.e., transmitters are of only one type and the same cost). In each case, the optimization model simultaneously considers the coverage and transmit power as objective functions, and the overlap as a constraint. The freespace propagation model is used to measure the signal strength in each case. The initial parameters used in the proposed algorithm are adjusted as follows. The twodimensional map size is 120 × 120 m 2 , and the grid resolution of the map is 1 m. All radiation patterns of the transmitters are assumed to be omnidirectional. The receiving antenna gain is 1 mW, and wavelength is 0.125 m. The threshold for all receivers is 1 × 10 −7 mW. The mathematical model for each test case (i.e., with two, three, or four transmitters) is constructed according to equations (9) through (12). In each case, the design variables are , , and for each transmitter, and the transmitting antenna gain ̃ is assumed to be an uncertain parameter. The design ranges for design variables are 0 ≤ ≤ 120, 0 ≤ ≤ 120, and 0.5 ≤ ≤ 2.5, ( = 2.5 used in Eq. (10)). The uncertain parameter ̃ is also assumed to vary uniformly in the range 0.5 to 1.5. To compare the results in all cases, the parameter in (12) is set to zero. Furthermore, the sensitivity analysis is conducted with different values of in [−2, +2], separately for each case. In the optimization procedure, the standard deviation of the coverage is weighted three times of mean by considering = 3 in (15). Additionally, we allow 30 percent average intersection between all signals transmitted by multiple transmitters in a model, thus = 0.3 considered in (11). We employ the Matlab® environment for data and function analysis. The Matlab® function gridsamp is used for sample allocation for both design parameters (inner array) and uncertain parameters (outer array). The DACE [46], Matlab® toolbox has been employed to construct the Kriging surrogate with zero-order polynomial regression and Gaussian correlation functions.

Results and discussion
For each case of = 2, 3, and 4 transmitters, we run the algorithm ten times to compute the relevant statistical measures to account for the randomness in the proposed stochastic optimization algorithm. Figure 4 shows the 3D surface plot for mean and Std of coverage over input combinations in the crossed array ( , and ). As it can be seen, there is a nonlinear relationship between the input and output set of data in the crossed array design because of the existence of uncertainty in the model. Figure 5 illustrates the mean and Std of coverage over one input parameter ( , or ) while the other is fixed at a predefined value. The power radius function of Eq.(2) is applied to collect the data regarding each crossed design and uncertain parameters. Figure 6 shows the mean and Std of power radius (m) overpower transmitted (mW). Using the input-output data pairs obtained from the crossed array, two Kriging surrogates are fitted, one for the mean and the other for Std of the coverage. Two other Kriging surrogates are fitted using the input-output data pairs in the relevant crossed array design, one for the mean of power radius and another one for Std of power radius. The prediction errors for all four Kriging surrogates are computed by leave-one-out cross validation approach. As shown in Figure 7, most observed standardized residuals are within the interval [−3,3], which ensures sufficient accuracy of the surrogates. In the next step, the robust optimization model is constructed using Eq.(12) as an objective function and Eq.(11) to implement the constraint. = 0 is considered in Eq.(12) for all cases. In both functions, we employ the Kriging surrogates that are already identified for the mean and the standard deviation of the coverage and the power radius. Accordingly, the GWO is run to obtain the robust optimal solution ( , , and ) for each case (e.g., the cases with number of two, three, or four transmitters). As mentioned earlier, optimization runs are repeated 10 times for each case. Table  2 gathers the statistical results. The best result with the highest SNR among all 10 repetitions is selected as a robust optimal point. Table 3 shows the best results (corresponding to the higher SNR value) for the placement problems with two, three, and four transmitters. Figure 8 shows the coverage using the robust optimal results depicted in Table 3 for the cases with two, three, and four transmitters separately. The plots are shown for the mean of the power radius, mean+3Std, and mean-3Std. The results indicate that the mean of the coverage is 91.38% (with Std = 10.01%) by using TPT = 4.56 (mW) is obtained in the case of three transmitters. For the problem with two transmitters, the total transmitted power is 3.18 (mW) to provide 70.30% of the mean coverage with Std = 13.70%. Finally, for the    Table 3 The optimal results with the highest SNR among 10 repetitions of proposed algorithm for placement of two, three, four transmitting antennas ( = 0).  problem with four transmitters, the mean coverage 97.83% can be obtained with Std = 4.25%, when the total power is 6.41 (mW). Note that four transmitters ensure higher SNR, but it also leads to the higher average overlap. The average intersection of transmitters is kept less than = 30% due to the predefined model's constraint in Eq. (11). In general, the three transmitters case seems to be the most advantageous from the point of view of ensuring the best trade-off between the coverage, robustness, power consumption, and the overlap.

#BS
• Sensitivity analysis Here, instead of estimating a single robust optimal point using a particular response that might be inaccurate because of high variability in uncertain parameters, we derive a series of possible responses that consider the degree of uncertainty by providing confidence regions or prediction intervals. This is realized by resampling adopted to the uncertain component.
In this study, we set the bootstrapped sample size = 100, and = 0.05. For each bootstrapped sample, we randomly select 30 uncertainty scenarios from the original sample points that were already available in the crossed array design (with the same sample size of 30). Subsequently, regarding the new sample sets, the mean and the standard deviation of the coverage are computed for the obtained robust optimal points as shown in Table 3. However, 95% two-sided approximations of BCI obtained by the distribution-free bootstrapping technique for the mean and the standard deviation of the coverage in the robust optimal points and the relevant SNR are as follows:   Figure 9 illustrates the 95% confidence regions for the mean and the standard deviation of the coverage at the original robust optimal points (as shown in Table 3).
In addition to the bootstrapped confidence region, sensitivity analysis is performed with different values of the parameter in Eq. (12) to investigate the effects of on the robust optimal results of the proposed algorithm for the multiple transmitter placement problem. Toward this end, the parameter is changed from −2 to +2 with a step size of 1. The proposed hybrid algorithm and GWO are executed ten times, separately for each value of  to study the effect of randomness as well.
The relevant statistical results are given in Table 4, Table 5, Figure 8 2D Visualization of the robust optimal coverage for the problems with placement of two, three, and four transmitters. Shown are the mean, mean+3Std, and mean-3Std.

Figure 9
The 95% two-sided confidence region obtained by freedistribution bootstrapping technique for placement problems with two, three, and four transmitters. Table 4 Statistical results over 10 repetitions of proposed Kriging-GWO (hybrid surrogate-metaheuristic) approach compare to GWO (uncombined metaheuristic) in robust optimal placement of two transmitting antennas (# = 2) for different " " value in Eq. (12).  and Table 6 for the cases with two, three, and four base stations, respectively. As it can be seen, the results obtained using the hybrid of Kriging and GWO (hybrid surrogatemetaheuristic) outperforms using GWO (uncombined metaheuristic) in both the accuracy (the average objective function value) and the robustness (smaller Std). Figure 10 illustrates the mean and the standard deviation versus the total transmitted power using all samples that were already collected during the previous optimization runs. Given these results, the decision-makers select their preferred combination of the mean and standard deviation regarding the power consumed to transmit the signals in a predefined map.

Conclusion
This paper presented a novel algorithm for solving the heterogeneous transmitter placement problem under uncertainty. Our methodology employs a hybrid Kriging-GWO approach combined with robust design optimization. The proposed algorithm enables computing the statistical measures due to the source of variability (uncertainty) and evaluate the possible positions and receiver test points to obtain robust optimal placement of multiple transmitters at low computational cost. Multiple objectives including the mean and the standard deviation of the coverage, and the total required transmitted power, as well as the reliability of the model in terms of controlling the overlap, are all considered simultaneously. Compliance with Ethical Standards i) Ethical approval: This research project does not contain any studies with human participants or animals performed by any of the authors. ii) Funding: This research project is supported by Second Century Fund (C2F), Chulalongkorn University, Bangkok. iii) Conflict of interest: The authors declare that they have no conflict of interest. iv) Informed Consent: Informed consent was obtained from all individual participants included in the study.  Table 6 Statistical results over 10 repetitions of proposed Kriging-GWO (hybrid surrogate-metaheuristic) approach compare to GWO (uncombined metaheuristic) in robust optimal placement of four transmitting antennas (# = 4) for different " " value in Eq. (12).  Figure 12 The estimation of the Pareto front for the mean and standard deviation of the coverage versus the transmitted total power. The estimation is obtained regarding different values of the weighting parameter α.