Stochastic Leader Gravitational Search Algorithm for Enhanced Adaptive Beamforming Technique

In this paper, stochastic leader gravitational search algorithm (SL-GSA) based on randomized k is proposed. Standard GSA (SGSA) utilizes the best agents without any randomization, thus it is more prone to converge at suboptimal results. Initially, the new approach randomly choses k agents from the set of all agents to improve the global search ability. Gradually, the set of agents is reduced by eliminating the agents with the poorest performances to allow rapid convergence. The performance of the SL-GSA was analyzed for six well-known benchmark functions, and the results are compared with SGSA and some of its variants. Furthermore, the SL-GSA is applied to minimum variance distortionless response (MVDR) beamforming technique to ensure compatibility with real world optimization problems. The proposed algorithm demonstrates superior convergence rate and quality of solution for both real world problems and benchmark functions compared to original algorithm and other recent variants of SGSA.


Introduction
Increasing interference due to multiple users and other signal sources is one of the fundamental problems in wireless communication and has been extensively studied for many years. Smart antenna systems decrease interference by adaptive beamforming techniques like minimum variance distortionless response (MVDR). It is one of the commonly utilized adaptive array beamforming techniques [1], but it is often not able to form nulls towards any nearby interference sources satisfactorily. Consequently, MVDR may lead to significant performance degradation in the case of unexpected interfering signals [2]. It is difficult and time consuming to solve these problems through conventional empirical approach, and sometimes, in the applied cases, is impractical. Recently, the employment of meta-heuristics algorithm has been growing instead of exhaustive and exact procedures in similar applications [3][4][5][6][7].
Consequently, meta-heuristics and exploratory methods need to provide mathematically reliable solution for this complicated class of optimization problems. However, the performance of these algorithms is often unsatisfactory for cases with three or more interference sources due to issues such as premature convergence and lack of sufficient exploration. Several methods are suggested for increasing the search diversity of SGSA, such as increasing the initial number of leaders (initial kbest). However, this significantly increases computational complexity of the force equation in GSA without properly addressing the key issue of dominant agents, with large masses, causing premature convergence. Primary reason for the search pattern domination is because agents are allowed to exert a force proportional to their performance and most SGSA variants allow the best agents to consistently influence all agents. Therefore, this paper, suggests a stochastic leader gravitational search algorithm (SL-GSA) to enhance MVDR beamforming performance by preventing premature convergence and improving overall exploration.
Standard gravitational search algorithm (SGSA) [8] was proposed as a global optimization method for computationally complex real world problems. In SGSA, the particles, called agents, move based on Newton's law of universal gravitation. The search space is represented as an 'n' dimensional space and the position of each agent is represented by a coordinate vector of length n. The mass of these agents are determined based on their fitness. The performance of each agent is calculated using the fitness function and their positions are updated accordingly. All the SGSA search agents (individuals) globally move toward the agents with heavier masses due to their gravitational force. Hence, superior solutions of the problems are represented by the heavier masses. The global search ability and high performance of SGSA in solving several nonlinear functions have been confirmed previously [8].
The balance between exploration and exploitation is critical for heuristic algorithms to achieve robust and reliable performance. In SGSA, this balance is achieved using the time variant linearly decreasing kbest parameter, which determines the number of agents that are allowed to exert force on the others in a given iteration. Thus, the parameter kbest is initially large and linearly reduced to provide some protection from premature convergence. This technique still allows the optimization process to be heavily influenced by agents with superior fitness resulting in poor exploration properties. As kbest agents are chosen based on their current fitness, it allows agents with superior fitness to attract the others towards optimal solutions. Thus, the algorithm is highly dependent on the best performing agents. However, if the kbest agents stagnate at a local optimum, the other agents become practically helpless to prevent premature convergence. The SGSA agents gravitate towards 'kbest' optimum agents. This allows convergence towards superior solutions but also allows the search to stagnate at local optima.
In this paper, SL-GSA randomly selects agents from a gradually reducing set that removes agents with inferior performance based on the adaptive parameter, γ. This directly prevents the domination of the search pattern by any individual agent. Thus, SL-GSA is far less likely to stagnate in a local optimum because it randomly ignores the best particles sometimes. This allows more efficient exploration before final convergence. The proposed new parameter, γ, prevents selection of the agents with the worse fitness in the later part of the optimization. This, in conjunction with the linear decrease of the parameter k, allows SL-GSA to converge faster than SGSA. This is verified by applying the proposed algorithm to six benchmark functions and two case studies of MVDR beamforming technique. High performance of convergence and quality of final solution compared to original algorithm is achieved as discussed in simulation results. The rest of this paper is organized as follows: Section 2 introduces the brief review of SGSA. The proposed SL-GSA is presented in section 3. The basics of adaptive beamforming and the conventional MVDR technique are explained in section 4 and 5, respectively. The testing of the proposed SL-GSA via benchmark functions and the simulation results obtained via SGSA and its variants are reported in section 6. Section 7 shows the incorporation of MVDR in SL-GSA. The efficiency of SL-GSA for different interferences in two case studies is also reported in this section. Finally, Section 8 concludes this investigation.

Standard Gravitational Search Algorithm
Standard gravitational search algorithm (SGSA) was presented as one of the recent heuristic population search algorithms [8] based on the mass interactions and laws of gravity. This method uses the Newtonian laws of gravitation to govern the motion of virtual agents through a search space. The i th GSA agent is defined as . .,x i n ), for i = 1,2,. . .,N that will be evaluated through their masses. In the definition of X i , N is the total number of agents, n is the total number of dimensions in the search space. These objects move together toward the agents with heavier masses due to gravitational force. The new position and velocity of i th agent at iteration t along the d th dimension will be upgraded as follows: where x i d and v i d is position and velocity of i th agent at dimension d respectively. rand i augments a randomized characteristic to the search pattern. The acceleration of agent i in iteration number t by the law of motion, is calculated as below: where M i is the normalized inertial mass of i th agent that evaluated by fitness function. More efficient agents have heavier masses that are calculated using fitness values. The updated inertia masses are as follows: where fit i (t) is the fitness value of the i th agent at iteration t. worst (t) and best (t) are the worst and best fitness of population, respectively that are defined as below (for a minimization problem): One way to increase the performance of SGSA is to make a suitable tradeoff between exploitation and exploration. In order to avoid local optimum stagnation, the algorithm needs to emphasize exploration at beginning. However, in later iterations, exploration must fade out and exploitation needs to increase. kbest is a time function that is initialized as the total number of agents and linearly decreased to one finally. Thus, at the beginning, all agents apply the force and, at the end, only one agent applies force to the others [8]. Furthermore, this reduces the number of overall forces that need to be summed to calculate the net force. Hence, Eq (8) is modified to the following: where rand j is the random value in the interval [0,1] and F ij d is force acting from mass j on mass i at dimension t as follows: where ε is the small constant value to prevent division by zero, R i,j (t) and G(t) are the Euclidean distance between agents i and j and the gravitational constant respectively as below: GðtÞ where t is the current iteration, t max is the maximum iteration number, β is a constant value. The gravitational constant determines how much influence the agents have on each other. High gravitational constant basically allows larger acceleration. It is exponentially reduced to allow convergence at the end of the optimization process. The acceleration of the i th agent at dimension d and iteration t is calculated by the Newtonian motion laws as following: Stochastic Leader Gravitational Search Algorithm The method of selecting best agents to apply the force on others in SGSA using a linear time based function, kbest, is an elegant way to improve the performance of the algorithm by efficiently compromising between the exploration and exploitation properties. However, selecting best agents to be the leaders may cause the algorithm to get trapped in local optima soon after reaching a near optimal region due to adherence to the leaders. In this paper we propose a modification to overcome this drawback by allowing other agents to lead randomly according to their performances. This is done by modifying Eq (9) as follows: where LAI is a vector representing the indices of leaders to be selected from a vector of agents sorted in descending order of fitness, at iteration t and can be described as follows: LAI ðtÞ ¼ ½z 1 ðtÞ; z 2 ðtÞ; . . . ; z i ðtÞ; . . . ; z K ðtÞ ð15Þ where z i (t) is the i th leader index number in the vector of agents which are sorted from the best to the worst at iteration t. γ is an adaptive factor between 0 and 1 determining the range of selection and is linearly adjusted during the search process from an initial value to a final value based on Eq (17); rand i produces a uniform random number between 0 and 1. In each iteration, the agents are sorted in a vector based on their performances from heavier masses to lighter ones. Then K leaders are selected randomly from the vector according to Eq (16). The selection mechanism is designed to allow selection of weaker agents in early iterations to lead others, resulting in improved exploration. As the iterations proceed, the chance of selecting the weaker agents is reduced by Eq (17) and superior agents attract others, which results in enhanced exploitation and convergence behaviour. By this modification the agents will be more active and the chance of getting trapped in local optima is reduced during first few iterations, while not degrading the overall convergence performance.
The flowchart for SL-GSA, illustrated in Fig 1, explains the algorithm in detail. It shows that the algorithm is rather convenient and robust as it requires calibration of only one additional parameter, γ.
This modification is designed to improve the early segment of the search process, while allowing the proposed algorithm to behave similarly to SGSA by discriminating against agents with poor performance.

Adaptive Beamforming
The aim of adaptive beamforming is to steer radiation beam towards the signal of interest (SOI) and produce null towards the interference sources. Beamforming, by signal processing technique, automatically recognizes and adjusts incoming SOIs and interference from received data. Signals from individual antennas were combined linearly after being scaled with the matching weights. This procedure improves the antenna array to obtain maximum gain in the desired signal radiation and nulls in the interference sources. An adaptive algorithm, like MVDR, uses the SOI direction (direction of main beam) and interference directions to produce a linear combination of the signals from each antenna, as illustrated in Fig 2, to achieve optimal signal to interference and noise ratio (SINR) performance. Furthermore, the algorithm must ensure a distortion-less response towards the SOI.
The output y is given by a linear combination of the data from the M sensors, with signal from each element x(n) and weight w(n) where Ã symbolizes the complex conjugate as shown below: where W is length N weights vector, X is length N received signals vector, and Superscript H is the conjugate transpose (Hermitian vector).

Conventional MVDR Beamforming
The minimum output energy (MOE) uses as a beamforming technique that implements several optimality standards. MOE has a fixed array gain on desired signal while simultaneously minimizing the output energy. Therefore, when gain is constant, any decline in the energy is because of suppressing interference as below: where h 0 is the steering vector. Lagrange multipliers method can solve this minimization. By solving the gain constraint and using the constraint on the weight, the parameter of Lagrange λ can be acquired. The weights can be written as: where R is the covariance matrix, with size of (M, 1), for the received signal x. Determining arbitrary constant, c = 1, creates the minimum variance distortionless response (MVDR). One of the beamforming techniques presented in literature [1] is known as MVDR due to minimum energy (variance) required for output signal and unity signal gain (desired signal is not distorted). Moreover, the output power that is subjected to a unity gain constraint in the direction of wanted signal must be minimized. The output power can be obtained as below: MVDR beamforming keeps a distortionless main lobe response towards the desired signal and minimize the array output power simultaneously. The Lagrange multiplier method does not explore the entire search space. Instead, the search is limited to locations identified by Eq (22), where the response towards the SOI is 0 dB. Furthermore, Lagrange multiplier method is reliant on the 'smoothness' of the fitness function to achieve satisfactory performance. Thus, this method may have unsatisfactorily low nulling in multiple interference scenarios, which may have significantly adverse effect on any system performance [9]. Hence, the introduction of a meta-heuristic technique with improved global exploration capability is proposed in this paper for beamforming application.

Model Verification
In order to analyse, compare and validate the effectiveness and efficiency of the proposed SL-GSA, the proposed algorithm is implemented on six standard benchmark functions and the results are compared with SGSA [8], OBGSA [10] and MGSA [11]. These functions and parameters used in this study are summarized in Tables 1 and 2, respectively. The first four functions are unimodal and the other two are multimodal. The results shown in Tables 3, 4 and 5 are recorded from 30 simulation of each algorithm. The worst, mean, median, best and standard deviation are calculated from aforementioned 30 runs and compared.
; uðx i ; a; k; mÞ ¼ kðx i À aÞ m ::::::; x i > a 0:::::::::::::; Àa  These benchmark functions provide an efficient method for comparison between heuristic algorithms along with providing insight into specific aspects such as rate of convergence and stability of the algorithm. The parameters presented in Table 2 are chosen to conduct a fair comparison between different variants of SGSA.
Tables 3 and 4 compare the optimization performance of the proposed algorithm to the results reported in literature for 500 and 1000 iterations, respectively.
According to the results of Tables 3 and 4, proposed SL-GSA provides considerably lower mean fitness values in 500 iterations than SGSA in 1000 iterations for all functions. For the functions F 1 , F 2 , and F 3 , proposed algorithm also provides better results in 500 iterations compared with MGSA [11] at equal number of iterations. Furthermore, for the same number of iterations, results achieved by proposed SL-GSA are better than those calculated by OBGSA [10] for the functions F 4 and F 5 . The comparison for 1000 iterations of SGSA and SL-GSA on F 6 concludes the comparisons with a prevalent trend of superior SL-GSA performance. Moreover, the standard deviations for 30 runs of proposed algorithm are generally smaller than those of SGSA, OBGSA and MGSA for both 500 and 1000 iterations, which strongly indicate that SL-GSA is significantly more stable and can produce better solutions more consistently. The results indicate, for all tested functions, the mean fitness values achieved by proposed algorithm are significantly lower than SGSA, OBGSA and MGSA, in same number of iterations. Therefore, the results of Tables 3 and 4 establish the proposed SL-GSA as more robust and capable than SGSA and its variants. The convergence histories of the mean best fitness of SGSA and SL-GSA in the 30 independent runs for six benchmark functions are shown in Figs 3-8. The Figs show that the convergence rate of the proposed SL-GSA is better than SGSA in 500 iterations.  show that the SL-GSA has faster varying curves of mean best fitness than SGSA and achieves lower and better solution compared with SGSA. Furthermore, Figs 6 and 7 indicate that SGSA suffers from premature convergence during optimization of F 4 and F 5 respectively. They stopped improving approximately before 50 iterations and hence produced an inferior solution. The overall results show that SL-GSA performs significantly better than SGSA for both unimodal and multimodal optimization problems consistently. The difference in exploration performance in the early phase between SGSA and proposed SL-GSA is vividly visible in most cases as the SL-GSA curve falls away from the SGSA before the hundredth iteration. This early advantage translates into further superior performance in the later phase of exploitation as the SL-GSA graph continues to separate further from the SGSA curve in case of F 1 , F 4 , F 5 and F 6 .

Sensitivity Analysis of γ Parameters
The γ parameters provide the distinctive characteristics to the proposed algorithm. Thus, in this section, the effects of the γ parameters on the proposed algorithm's performance are analyzed in detail. The γ initial determines the performance cutoff for eligible agents at the beginning of the optimization. The value 0.6 roughly corresponds to the top 60 percentile in terms of performance. This cutoff criterion becomes stricter with progressive iterations. As γ decreases linearly towards the γ final value of 0.1, the agents must be within the approximately top 10 percentile to be eligible for selection. The algorithm randomly chooses the leaders from the eligible set of agents. This process implements equal probability of selection for each eligible agent. Furthermore, an eligible agent may be selected multiple times. The proposed algorithm is more dependent on the choice of value for γ final than γ initial . Higher values for γ initial are likely to yield more exploration in the early phase, but degrade the overall convergence rate. γ final value of zero indicates that only the best particle is eligible for selection. This strategy is similar to the common practice of reducing kbest to 1 in SGSA and similar algorithms. Implementation of the less stringent γ final value of 0.1 in this study allows the proposed algorithm to outperform the SGSA and some of its variants. However, varying γ from 1 to 0 is detrimental to the search process. This causes all agents, including agents with the worst performance, to be eligible for selection in the early phase. The algorithm can be misguided if the poorest performing agents are repeatedly selected, resulting in poor stability. Table 5 illustrates a performance comparison for different pairs of γ initial and γ final . The algorithm was applied 30 times for each configuration to the multimodal function F6 to study the overall effects of γ. The worst performance is recorded for the extreme values of γ initial = 1 and γ final = 0. It is the least stable configuration as it exhibits the highest standard deviation. γ initial = 0.9 and γ final = 0.5 shows slightly better results. However the accuracy and stability of convergence are still unsatisfactory as indicated by the considerably high mean and standard deviation. The third case of γ initial = 0.9 and γ final = 0.1, in Table 5, can be considered adjustment from the first case results in. The algorithm is extremely stable and accurate in this configuration. The proposed algorithm exhibits similar behavior for the last two cases  Table 5, which illustrates that it is not too sensitive to the γ parameters provided they are within the general guidelines mentioned above.
An overall sensitivity analysis for γ initial = 0.6 ± 0.1 and γ final = 0.1 ± 0.1 on unimodal function F 3 , presented in Table 6, further illustrates the stable performance of the proposed algorithm. The last two rows in Table 6 show the cumulative effect of parameter value alterations. 'Up' row represents the results obtained from increasing both γ parameters. Conversely, 'Down' considers the scenario with both parameters reduced. Each configuration was implemented 30 times and the deviations are recorded as percentage of the results presented in Table 6, as shown in Eq (24).
The results in Table 6 are attained using the recommended configuration of γ initial = 0.6 and Y final = 0.1. The highest average performance deviation in Table 6 is 2.352% for increasing γ initial by 0.1. This does not accumulate further as the 'Up' scenario shows deviation of 2.2403%. Thus, the performance degradation caused by γ initial is partially mitigated by increasing γ final . It is evident from the p-values of Wilcoxon Rank Sum test that the recommended configuration achieves a statically significant performance improvement over only one of the six scenarios tested in Table 6. The 'Down' configuration resulted in p-value less than 0.05, which means the null hypothesis (no statistical difference) must be rejected and alternate hypothesis (recommended configuration outperforms 'Down' configuration) must be accepted. Thus, there is no statistically The results presented in Tables 5 and 6 indicate that the proposed algorithm is applicable to both multimodal and unimodal functions. Overall, the analyses show that the proposed algorithm is sufficiently stable for application in real world optimization problems.

Model Application
In the previous section, the robustness and effectiveness of the suggested algorithm (SL-GSA) have been tested through six benchmark functions. The performance of the proposed method to enhance signal to interference and noise ratio (SINR) of MVDR beamforming technique will be investigated and compared with conventional MVDR and MVDR-SGSA in this section. Two interference scenarios will be considered in this study. The first case is with two interference sources located at 30°and 50°while the second case is with four interference sources located at 30°, 50°, 25°and 60°. Both cases have one user at 0°and number of elements is four. The proposed algorithm was implemented by using MATLAB 1 .

Signal to Interference and Noise Ratio Calculation
In this paper, the SL-GSA was employed to improve the MVDR beamforming performance by increasing the SINR value of four element antenna array. The SL-GSA and SGSA were applied on MVDR technique to maximize the SINR by varying the complex weights. In these algorithms, w mvdr (MVDR weight vector) will be replaced as the first agent of initial population. The rest of the population is initialized randomly, as shown in Eq (26). This system will create a population of N agents to represent weight vectors. The weight vectors in every agent will contain M number of weights, where M is the number of elements in the array. The weights considered in the optimization process, thus, are represented as the matrix W NM .
where dimension d is analogous to number of sensor m. Agent index, i, in SL-GSA and SGSA is set to n, as shown in Eq (26). The population weight vectors can be defined in the matrix format as follows: where W NM is the weight vectors of N agents with M sensors and W mvdr is the weight vectors from conventional MVDR beamformering. The fitness function employed in this study is the SINR as shown in Eq (27). Therefore, the optimization processes will try to maximize the fitness function by finding the corresponding optimal weight vectors.
where P u is the target user power, P i is the power of interference at i interference, I is the number of interference sources and N is the noise.

Simulation
In this section, the validation of the proposed approach for real world applications with two interference scenarios is studied. The parameters of SGSA are chosen according to the recommendations and guidelines presented in literature [8]. These configurations of SGSA have also been used extensively after the development of SGSA [10][11][12][13]. The algorithms are simulated 20 times and the best results are recorded. The maximum iteration number and population size set to 100 and 50 respectively.
7.2.1. Case 1: One user two interferences. Two interference sources at 30°, 50°and user at 0°have been assumed in the first case study.
According to the weights in Table 7, the performance of power response for MVDR, SGSA-MVDR and SL-GSA-MVDR beamforming is plotted in Fig 9 to explain the target user and interference by using different values. Fig 9 clearly illustrates that the proposed algorithm achieves significantly deeper nulls than SGSA-MVDR and conventional MVDR.
In Table 8, 67.79% and 105.84% improvements of SINR are achieved by using SGSA and SL-GSA, respectively, over conventional MVDR. The superior performance of SL-GSA in this problem, along with the results of the benchmark problems, establishes SL-GSA's as a suitable interference mitigation algorithm. 7.2.2. Case 2: One user four interferences. Four interference sources at 30°, 50°, 25°and 60°and user at 0°have been assumed as the second case study.
The power response for MVDR, SGSA-MVDR and SL-GSA-MVDR is plotted in Fig 10  according to the weights in Table 9 to explain the target user and interference by using different values. Conventional MVDR rarely manages to create sufficiently deep nulls towards more       than 3 interference sources as shown in Fig 10. The two interferences at 25°and 30°produce one shallow null for MVDR and SGSA-MVDR. Furthermore, the nulls towards the other two interference sources are similarly shallow. Two of the three nulls created by MVDR are significantly improved by the implementation of the proposed algorithm. This clearly shows that SL-GSA can outperform the beamforming techniques in complex scenarios involving more than 3 interfering sources. The deep nulls correspond to superior SINR performance as shown in Table 10. Table 10 illustrates the improvement of SINR that are 0.32% and 4.84% by SGSA and SL-GSA, respectively, as compared to conventional MVDR method. This case represents the most complex interference mitigation optimization problem in this paper. The results are consistent with the previous case, indicating SL-GSA to be generally superior for real world optimization problems compared to SGSA.
Referring to the both simulation cases, when the interference sources become closer, it is more difficult for both conventional MVDR and artificial intelligence based technique to enhance null level at the interference direction. From the simulation result, it is clearly observed that, SGSA performs inferior than SL-GSA. The superior exploratory and exploitive properties of SL-GSA over SGSA have resulted in better beam steering and interference mitigation performances in both cases.

Conclusions
In this paper, a stochastic leader gravitational search algorithm via randomizing the choice of kbest utilized in total force calculation of the original algorithm was proposed. The use of kbest in SGSA that allows the heaviest few agents to apply force on other agents may often causes premature convergence due to one or a group of dominant agents. The proposed SL-GSA alleviates this problem while simultaneously enhancing the exploration in the early phase of optimization by randomly selecting the leading agents. This stochastic process prevents the aforementioned domination and results in superior performance for both benchmark functions and real world optimization problems. The parameter γ provides an elegant method for balancing the exploration and exploitation properties of the proposed algorithm. Overall, SL-GSA exhibits higher robustness, convergence rate, accuracy, and stability without increasing computational difficulty in comparison with SGSA. The achievement of satisfactorily deep nulls in complex interference environments indicates that the proposed algorithm is ready to be implemented in real world applications.