Currently, the development of discrete signal filtering methods that are used in computerized biometric identification systems is an urgent task.In the case of linear filtering, one of the most popular methods is the Kalman filter.A modified Kalman filtering method was proposed to improve the efficiency of digital filtering of a discrete signal. This method provides automation of parameter value determination and improves the speed and accuracy of Kalman filtering by using fewer parameters and identifying them based on immune metaheuristic methods. The proposed metaheuristic methods reduce the probability of convergence to a local extremum by using the Cauchy distribution and make parametric identification more accurate. Algorithms of immune metaheuristic methods for identifying Kalman filter parameters have been developed, which are designed for software implementation on the GPU using CUDA technology, which increases the accuracy of Kalman filtering. Further prospects of the study are to utilize the proposed immune metaheuristic methods for various general and special purpose intelligent systems.
Currently, it is an urgent task to create digital filtering methods used for noise suppression, which are used in computerized biometric identification systems as well as speech recognition and understanding systems and computer vision [1-2].
One of the popular methods of digital filtering is the use of recurrent neural networks. Such neural networks include NARMANN, ENN, JNN, GRU, LSTM. The disadvantage of such neural networks is the high computational complexity of identifying their parameters due to the lack of the ability to parallelize the learning algorithm, a large number of connections between neurons and the correct choice of activation functions and the number of neurons in the layers.
Another popular digital filtering methods is smoothing adaptive linear filtering [3-4]. For smoothing adaptive linear filtering, the identification of filter parameters plays an important role.
Approximate methods of determining parameter values based on global search do not guarantee convergence. Approximate methods of determining parameter values based on local search have a high probability of hitting a local extremum. Exact methods of determining parameter values have high computational complexity. Thus, the problem of insufficient quality of parametric identification methods arises.
Modern heuristics (or metaheuristics) are used to increase the speed of filter parameter identification and reduce the probability of hitting a local extremum [5-6]. Metaheuristics expand the capabilities of heuristics by combining heuristic methods based on a high-level strategy [7-8]. Metaheuristics often use the behavior of evolutionary and immune approaches [9-10]. Metaheuristics are approximate and, as a rule, stochastic methods [11-12]. The most effective metaheuristics use experience that accumulates during the search process and is stored in memory [13-14].
The aim of the work is to improve the quality of Kalman filtering of a discrete signal using immune metaheuristic methods.
To achieve the stated goal, it is necessary to solve the following tasks: 1. Create a modified Kalman filtering method. 2. Develop a continuous optimization method based on an artificial immune network. 3. Create a continuous optimization method based on a hybrid immune algorithm. 4. Conduct a numerical study of the proposed methods of continuous optimization and Kalman filtering.
In the pre-processing units of modern computer systems for speaker identification, multimodal interface, medical and technical diagnostics, digital filtering methods are implemented that reduce noise and allow analyzing the spectral features of the signal. Various digital filtering methods, from recurrent neural networks to linear digital filters, require the use of methods for determining parameter values. One of the popular approaches is metaheuristic.
Currently, metaheuristics are divided into evolutionary, biological, non-nature-inspired, immune, mathematical, physical, social, chemical. Most metaheuristics use stochastic search, which reduces the probability of hitting a local extremum. Metaheuristics allow solving continuous optimization problems (calculating the point at which the objective function reaches a maximum or minimum) and discrete optimization problems (e.g., clustering, knapsack problem, traveling salesman problem, assignment problem).
Modern metaheuristics have one or more of the following disadvantages: the convergence of the method may not be ensured [15-16]; the iteration number does not affect the solution finding process [17-18]; only binary potential solutions are used [19-20]; the metaheuristic method is associated with solving only one problem or there is only an abstract set of operators of this method available [21-22]; low accuracy of the method [23-24]; there is no automation of the procedure of identification of method parameters [25-26]; the method is not intended for solving problems of conditional optimization [27-28].
Based on this, the problem of constructing high quality metaheuristic optimization methods needs to be addressed. One of the most popular ones are immune metaheuristics.
In this paper, it is assumed that for each moment of time n the state matrix A ( n )= A, the control matrix B ( n )=0, the control vector u ( n )=0, the process covariance matrix Q ( n )=Q= I Q σ 2Q, the observation matrix H ( n ) are replaced by the vector h=1H σ 2H, the observed noise covariance matrix R ( n ) is replaced by the scalar σ 2R, the observation vector z ( n ) is replaced by the scalar z ( n ), the observation vector estimate ´z ( n ) is replaced by the scalar ´z ( n ), the error vector e ( n ) is replaced by the scalar e ( n ), the error vector covariance matrix D ( n ) is replaced by the scalar d ( n ), and the optimal Kalman gain matrix G ( n ) is replaced by the vector g ( n ). This reduces the amount of computation and the amount of specified data.
The modified Kalman filtering method for a discrete signal consists of the following stages: 1. Initialization.
The initial estimation of the state vector ´s ( 0 ) of length m (in this work ´s ( 0 )=0) is set. The covariance matrix of the state vector estimate C ( 0 ) is set with a size of mxm (in this work, C ( 0 ) is filled with uniformly distributed values). In this work, the state matrix A of size mxm, dispersion σ 2R, dispersion σ Q2 for calculating the matrix Q of size mxm, dispersion σ 2H for calculating the vector h of length m are set. 2. Forecast stage. 2.1. Calculating the state vector estimate In this work s ( n−1 )= A⋅ s ( n−2 ). 2.2. Calculating the covariance matrix of the state vector estimate In this work C ( n−1 )= A⋅ C ( n−2 )⋅ AT +Q, where Q= I Q σ 2Q. 3. Update stage. 3.1. Calculating the estimate of the observation vector
´s ( n−1 )= A ( n ) ´s ( n−2 )+ B ( n )u ( n ).
C ( n−1 )= A ( n ) C ( n−2 ) AT ( n )+Q ( n ).
2 In this work ´z ( n )=h⋅ ´s ( n−1 ), where h=1H σ H. 3.2. Calculating the error vector In this work e ( n )= z ( n )− ´z ( n ). 3.3. Calculating the error vector covariance matrix ´z ( n )= H ( n ) ´s ( n−1 ).
e ( n )= z ( n )− ´z ( n ).
D ( n )= H ( n ) C ( n−1 ) H T ( n )+ R ( n ).
In this work d ( n )=h⋅ C ( n−1 )⋅ hT + σ 2R, whereh=1H σ 2H. 3.4. Calculating the optimal Kalman gain matrix In this work g ( n )=
d ( n ) 3.5. Updating the state vector estimate
C ( n−1 ) hT
G ( n )=
C ( n−1 ) H T ( n )
D ( n )
.
2 , where h=1H σ H.
´s ( n )=´s ( n−1 )+G ( n ) e ( n ).
In this work ´s ( n )=´s ( n−1 )+ g ( n ) e ( n ). 3.6. Updating the covariance matrix of the state vector estimate
C ( n )=C ( n−1 )−G ( n ) H ( n ) C ( n−1 ).
2
In this work C ( n )=C ( n−1 )− g ( n )⋅ h⋅ C ( n−1 ), where h=1H σ H. 4. Selection of criteria for evaluating the effectiveness of the modified
In this paper, to evaluate the effectiveness of the modified Kalman filtering method, the accuracy
criterion is chosen, which means choosing such values of the parameters A, σ Q2, σ 2H, σ 2R, that provide a
minimum root mean square error
(
1 nmac n=1 →
min ,
A ,σQ2 ,σ2H ,σ2R
(
In accordance with the selected criterion, immune metaheuristic methods for identifying parameters A, σ Q2, σ 2H, σ 2R are proposed in this paper.
The artificial immune network was proposed by Timmis, Neal, Hunt and later modified by de Castro, von Zuben and is based on the hypothesis of representing the immune system as an idiotypic network. The disadvantage of the clonal selection theory is that it assumes that a set of cells remains unexcited when there is no antigen.
Scientist Erne proposed a hypothesis according to which the immune system is a regulated network of molecules and cells that recognize each other even in the absence of an antigen. Such structures are often called idiotypic networks, they serve as a basis for studying the behavior of the immune system. Erne's theory is interpreted as a system of differential equations describing the dynamics of the concentration of lymphocyte clones and the corresponding immunoglobulin molecules. The theory of idiotypic regulation is based on the assumption that different lymphocyte clones are not isolated from each other, but maintain communication through interactions of their receptors located on the surface of the lymphocyte.
In formulating the foundations of his theory, Jerne introduced the concepts of formal and functional networks. Formal networks serve to study issues of repertoire (recruitment), dualism and suppression. When considering functional networks, a quantitative picture of the theory is presented.
A probabilistic approach to studying idiotypic networks based on the work of Jerne was proposed by the scientist Perelson. This approach is extremely formalized and is mainly associated with the description of phase transitions. Perelson divided the plane of phase variables of the considered system of equations into a subcritical region, a transition region, and a postcritical region. Over the past 20 years, Jerne's proposed immune network theory has received considerable attention, which has led to a detailed study of many computational aspects of the corresponding mathematical models.
The modified artificial immune network method consists of the following steps: 1. Initialization. 1.1. Setting the search area: cell length M , minimum and maximum values of cell components x mji nmjax, j∈ 1 , M . Setting the maximum number of iterations N , population size K , number of clones LC. 1.2. Setting the cost function (target function) F ( x )→ min, where x – cell (real vector containing x A, σ Q2, σ 2H, σ 2 ).
R 1.3. Setting search parameters: mutation parameter α, compression threshold ε, where α >0 (the higher the α , the lower the mutation probability), ε >0. The paper it is proposed to use δ ( x mja x mjin, 0< δ <1 , instead of parameter α. 1.4. Creating an initial population P. 1.4.1. Cell number k =1, P=∅ .
1.4.2. Generating a random cell xk=( xk 1 , ... , xkM ), xkj= x mji nmjaxmjin, where U (
1.4.3. If xk∉ P, then P= P∪ { xk }, k =k +1. 1.4.4. If k ≤ K , then move to step 1.4.2. 1.5. Determine the best cell by the target function x¿=arg min F ( xk ), k ∈ 1 , K .
xk 2. Iteration number n=0. 3. Calculating the affinity of population cells P Φ ( xk )=1−
F ( xk )− min F ( xi )
i∈ 1, K max F ( xi )− min F ( xi ) i∈ 1, K i∈ 1, K 4. Order the population P by the target function, i.e. F ( xk )< F ( xk +1 ). 5. Determine the best cell by the target function
, Φ ( xk )∈ [ 0,1], k ∈ 1 , K . 6. Determining the global best cell. If F ( xk¿ )< F ( x¿ ), then x¿= xk¿. 7. Calculating the average cost value k ¿=arg min F ( xk ), k ∈ 1 , K .
k F¯ source= 1 ∑K F ( xk ).
K k=1 8. Creating the best mutated clones set H . 8.1. Set k =1, H =∅ .
~ 8.2. Creating clones set Pk={~xkl } for a population cell xk. 8.3. Creating mutated clones set P´k. 8.3.1. Clone number l=1, P´k =∅ . 8.3.2. Creating a cell (the paper proposes to use the Cauchy distribution)
x´klj=~xklj + δ ( x mja x mjin−Φ(~xkl), j∈ 1 , M ,
where Cauc h y (
x´klj=max ¿ x mji n´klj ¿, x´klj=min ¿ x mja x´klj ¿, j∈ 1 , M . 8.3.4. Calculating P´k = P´k∪ { x´kl }. 8.3.5. If l < LC, then l=l +1, go to step 8.3.2. 8.5. Calculating H = H ∪ {hk }. 8.6. If k < K , then k =k +1, go to step 8.2. 9. Calculating the average cost value 8.4. Determining the best element of set P´k by the target function hk =arg min F ( x´kl ). x´kl F¯ mutate= 1 ∑K F ( hk ).
K i=1 10. If F¯ mutate ≥ F¯ source, then move to step 8. 11. Compressing the set H and replacing population cells P with elements of the set H . 11.1. Set k =1, m=1. 11.2. Forming the ε - neighborhood of the mth element of the set H
U hm , ε={hl∨ ρ ( hm , hl ) ≤ ε , l∈ 1 , K },
where ρ – the distance between hm and hl (e.g. Euclidean distance).
(
The result is x¿.
The algorithm of the modified artificial immune network method, designed for implementation on GPU using CUDA technology, is shown in Figure 1.
This block diagram functions as follows.
Step 1 – Set the maximum number of iterations N , population size K , number of clones LC, parameter δ for generating a new solution, compression threshold ε, where 0< δ <1, ε >0.
Step 2 – Create an initial population P, using K ⋅ M threads that are grouped into K onedimensional blocks. Each thread calculates xkj= x mji nmjaxmjin. ¿ Step 3 – Determine the best cell by the target function x =arg min F ( xk ), k ∈ 1 , K . xk Step 4 – Set the iteration number n=0.
Step 5 – Calculate the minimum value of the target function a in the current population P, using K threads that are grouped into one one-dimensional block. In this block, the minimum of K elements of the form F ( xk ) is calculated based on the reduction.
Step 6 – Calculate the maximum value of the target function b in the current population P, using K threads that are grouped into one one-dimensional block. In this block, the maximum of K elements of the form F ( xk ) is calculated based on the reduction.
Step 7 – Calculate the affinity Φ ( xk ) for each cell xk of population P, using K threads that are grouped into one one-dimensional block Φ ( xk )= b− F ( xk ) .
b−a Step 8 – Order population P by the target function, i.e. F ( xk )< F ( xk+1 ). ¿ Step 9 – Determine the best cell by the target function k =arg min F ( xk ), k ∈ 1 , K . k Step 10 – Determine the global best cell. If F ( xk¿ )< F ( x¿ ), then x¿= xk¿.
Step 11 – Calculate the average cost value F¯ source, using K threads that are grouped into one onedimensional block. In this block, the sum of K elements of type F ( xk ) is calculated based on the K reduction.
Step 12 – Set the cell number k =1.
~
Step 13 – Create clones set Pk={~xkl } for population cell xk. are
Step 14 – Generate a mutant clones set P´k for a population cell xk, using LC⋅ M threads that
grouped into LC one-dimensional blocks. Each thread calculates
x´klj=~xklj+ α1 e−Φ(~xkl) Cauc h y (
Step 15 – Perform the correction of the mutant clones set P´k for the population cell xk, using LC⋅ M threads, that are grouped into LC one-dimensional blocks. Each thread calculates x´klj=max ¿ x mji n´klj ¿, x´klj=min ¿ x mja x´klj ¿.
Step 16 – Determine the best mutated clone by the target function for each population cell xk Step 17 – If k < K , then k =k +1, go to step 13.
Step 18 – Calculate the average cost value F¯mutate, using K threads that are grouped into one onedimensional block. In this block, the sum of K elements of the form F ( hk ) is calculated based on the K reduction.
Step 19 – If F¯mutate ≥ F¯ source, then go to step 12.
Step 20 – Compute a distances set { ρ ( hm , hl ) }, using K ⋅ K threads that are grouped into K onedimensional blocks. Each thread computes ρ ( hm , hl ).
Step 21 – Compress set H and replace population cells P with elements of set H . Step 22 – If k = K , then go to step 24.
Step 23 – Initialize the last K −k population cells P, using ( K −k )⋅ M threads that are grouped intoK −k one-dimensional blocks. Each thread calculates xkj= x mji nmjaxmjin.
Step 24 – If n< N −1, then n=n +1, go to step 5.
The hybrid immune algorithm was proposed by scientists Lucinska and Wierzchon and is a modification of the artificial immune network. Its distinctive feature is the use of two types of mutations.
The modified hybrid immune algorithm method consists of the following steps:
1. Initialization.
1.1. Setting the search area: cell length M , minimum and maximum values of cell components
x mji n mjax, j∈ 1 , M . Setting the maximum number of iterations N , population size K , number of
clones LC, memory size LM .
1.2. Setting the cost function (target function) F ( x ) → min, where x – cell (real vector).
x
11.3. Setting the search parameters: compression threshold ε, maximum cell age amax, memory size
excess coefficient α , where ε >0, amax – is a natural number, α >1.
1.4. Creating an initial population P.
1.4.1. Cell number k =1, P=∅ .
1.4.2. Generating a random cell xk =( xk 1 , ... , xkM ), xkj= x mji nmjaxmjin, where U (
~
3.2. Creating a set of clones Pk={(~xkl , ~akl ) } for a population cell ( xk , ak ).
3.3. Modification of a set of mutated clones P´k.
3.3.1. Clone number l=1.
3.3.2. Set λ=U (
σ´ klj={2(~xlσk´j −klj ,x´klj ) , FF (( x´x´kklljj ))<≥ FF ((~~xxkklljj )), j∈ 1 , M ,
creation of a cell (the work suggests using the Cauchy distribution)
x´klj= σ´ klj Cauchy (
The algorithm of the modified hybrid immune algorithm method, designed for implementation on GPU using CUDA technology, is shown in Figure 2.
U qi , ε={ql∨ ρ ( qi , ql ) ≤ ε , l∈ 1 ,∨Q∨¿ }¿, where ρ – the distance between qi and q j (e.g. Euclidean distance). 7.3. If ¿ U qi , ε∨¿∅ or max U qi , ε=qi, then q˘m=qi, m=m+1, Q˘ =Q˘ ∪ {q˘m }. 7.4. If i < LM , then i=i +1, go to step 7.2. 7.5. Set Q=Q˘ . 7.6. Order the set Q by the target function, i.e. F ( qi )< F ( qi+1 ). 7.7. If ¿ Q∨¿ LM, then remove from the ordered set Q the last worst ¿ Q∨− LM cells by the target function. 8. Determining the best memory cell of the set Q by the target function x¿=arg min F ( qi ).
qi
(
Step 1 – Set the maximum number of iterations N , population size K , number of clones LC, memory size LM, compression threshold ε, maximum cell age amax, memory oversize factor α , where ε >0, amax is a natural number, α >1.
Step 2 – Create an initial population P={( xk , ak ) }, using K ⋅ M threads that are grouped into K one-dimensional blocks. Each thread of each kth block calculates xkj= x mji nmjaxmjin, ak=1.
Step 3 – Set the cell number k =1.
Step 4 – Generate a set of mutated clones P´k={( x´kl , σ´ kl ) }, using L ⋅ M threads that are grouped
into LCone-dimensional blocks. Each thread calculates x´klj= x mji nmjaxmjin C
, σ´ klj=U (
Step 5 – If k < K , then k =k +1, go to step 4.
Step 6 – Initialize a set of memory cells Q=∅ .
Step 7 – Set the iteration number n=0.
Step 8 – Set the cell number k =1.
Step 9 – Create clone set ~Pk={~xkl } for population cell ( xk , ak ).
Step 10 – Calculate λl=U (
Step 11 – Modify the set of mutated clones P´k for population cell ( xk , ak ), using LC⋅ M strands that are grouped into LCone-dimensional blocks. Each strand of each lth block calculates: 2(~xlkj− x´klj ) ,
σ´ klj , If λl>0.2, then σ´ klj={
F ( x´klj )< F (~xklj ), x´klj= σ´ klj Cauc h y (
F ( x´klj )≥ F (~xklj )
If λl ≤ 0.2, then x´kl=~xkl, j=round ( 1+( M −1 )U (
Step 12 – Perform the correction of the modified set of mutated clones P´k for the population cell ( xk , ak ) using LC⋅ M strands, which are grouped intoLC one-dimensional blocks. Each thread of each lth block calculates:
If λl>0.2, then x´klj=max ¿ x mji n´klj ¿, x´klj=min ¿ x mja x´klj ¿.
Step 13 – Determine the best mutated clone by the target function for each population cell ( xk , ak ) Step 17 – If ¿ Q∨¿ α ⋅ LM, then move to step 20.
Step 18 – Compute distance set { ρ ( qi , ql ) } using Q⋅ Q threads that are grouped into Q one¿ Step 20 – Determine the best memory cell of the set Q by the target function x =arg min F ( qi ), qi dimensional blocks. Each thread calculates ρ ( qi , ql ).
Step 19 – Compress set Q. k ∈ 1 , K .
Step 21 – If n< N −1, then n=n+1, go to step 8.
The numerical study of the proposed metaheuristic methods was carried out using the Python package in the Google Colab environment. Numerical experiments were carried out using the CUDA parallel information processing technology on a GeForce 920M video card with 1025 threads in a onedimensional block.
In this work, the following parameters were used for the modified artificial immune network method: the maximum number of iterations N =100, population size K =20, number of clones LC=10, parameter for generating a new solution δ =0.1, and compression threshold ε =0.1.
In this work, the following parameters were used for the modified hybrid immune algorithm method: maximum number of iterations N =100, population size K =20, number of clones LC=10, memory size LM =10, compression threshold ε =0.1, maximum cell age amax, and memory oversize factor α =1.5.
In the work, a one-dimensional signal was generated, to which additive Gaussian noise with zero mathematical expectation and dispersion of 135 was added.
In the work, the following root mean square errors were calculated based on the formulas: RM S yz=√ RM Sx ´z=√
1 ∑200 ( y ( n )− z ( n ))2, 200 n=1
1 ∑200 ( x ( n )−´z ( n ))2, 200 n=1 (22) (23) where y ( n ) – noise-free countdown, z ( n ) – observation (noise-free countdown), ´z ( n ) – observation evaluation.
Table 1 presents the root mean square errors for immune metaheuristic methods. For all three methods RM S yz=9.25.
For example, for the modified hybrid immune algorithm method, the following parameter values were obtained: The advantage of using the proposed immune metaheuristic methods: 1. Automation of determination of Kalman filtering parameters A, σ Q2, σ 2H, σ 2 .
R 2. Immune metaheuristic methods, due to their stochastic nature, reduce the probability of convergence to a local extremum. 3. For immune metaheuristic methods, it is proposed to replace the Gaussian distribution with the Cauchy distribution, which is long-tailed, i.e. to reduce the probability of convergence to a local extremum. 4. According to Table 1, the modified hybrid immune algorithm method gives the best results in terms of the root mean square error. 5. According to Figure 3 and Figure 4, the original signal without noise and the filtered signal differ insignificantly, while there is a significant difference between the original signal without noise and the signal with noise.
Conclusions 1. A modified Kalman filtering method was developed that provides automation of parameter value determination and increases the speed and accuracy of Kalman filtering by using fewer parameters and identifying them based on immune metaheuristic methods. 2. A modified artificial immune network method was created that reduces the probability of convergence to a local extremum by using the Cauchy distribution and makes the proposed method more accurate than the existing one. 3. A modified hybrid immune algorithm method was developed, which, by using the Cauchy distribution, reduces the probability of convergence to a local extremum and makes the proposed method more accurate than the existing one. 4. Algorithms of immune metaheuristic methods for identifying Kalman filter parameters have been developed, which are intended for software implementation on GPU using CUDA technology, which increases the accuracy of Kalman filtering. The numerical studies conducted have confirmed the operability of the developed software and allow us to recommend it for practical use. 5. Further research prospects include the use of the proposed immune metaheuristic methods for various general-purpose and special-purpose intelligent systems, for example, for training neural networks.
During the preparation of this work, the authors used Grammarly in order to: Grammar and spelling check. After using this tool, the authors reviewed and edited the content as needed and take full responsibility for the publication’s content. [21] A. Slowik, Swarm Intelligence Algorithms, A Tutorial, 1st ed., Boca Raton, FL: Chapman and
Hall/CRC, New York, 2021. doi: 10.1201/9780429422614. [22] O. Bozorg Haddad, M. Solgi, H. Loaiciga, Meta-heuristic and Evolutionary Algorithms for Engineering Optimization, Hoboken, New Jersey: Wiley & Sons, 2017. doi: 10.1002/9781119387053. [23] K.-L. Du, M. N. S. Swamy, Search and Optimization by Metaheuristics. Techniques and
Algorithms Inspired by Nature, Charm: Springer, 2016. doi: 10.1007/978-3-319-41192-7. [24] A. Kaveh, T. Bakhshpoori, Metaheuristics Outlines, MATLAB Codes and Examples, Cham:
Springer, 2019. doi: 10.1007/978-3-030-04067-3. [25] A. Nayyar, D.-N. Le, N. G. Nguyen, Advances in Swarm Intelligence for Optimizing Problems in Computer Science, 1st ed., Boca Raton, FL: Chapman and Hall/CRC, New York, 2018. doi: 10.1201/9780429445927. [26] H. Emami, Anti-coronavirus optimization algorithm, Soft Computing, vol. 26 (2022) 4991-5023.
doi:10.1007/s00500-022-06903-5. [27] A. M. Khalid, K. M. Hosny, S. Mirjalili, COVIDOA: a novel evolutionary optimization algorithm based on coronavirus disease replication lifecycle, Neural Computing and Applications, vol. 34 (2022) 22465-22492. doi:10.1007/s00521-022-07639-x. [28] A. Al-Betar, Z. A. A. Alyasseri, M. A. Awadallah, I. A. Doush, Coronavirus herd immunity optimizer (CHIO), Neural Computing and Applications, vol.33 (2021) 5011-5042. doi:10.1007/s00521-020-05296-6.