Grey Wolf Optimizer with Behavior Considerations and Dimensional Learning in Three-Dimensional Tooth Model Reconstruction

Three-dimensional registration with the affine transform is one of the most important steps in 3D reconstruction. In this paper, the modified grey wolf optimizer with behavior considerations and dimensional learning (BCDL-GWO) algorithm as a registration method is introduced. To refine the 3D registration result, we incorporate the iterative closet point (ICP). The BCDL-GWO with ICP method is implemented on the scanned commercial orthodontic tooth and regular tooth models. Since this is a registration from multi-views of optical images, the hierarchical structure is implemented. According to the results for both models, the proposed algorithm produces high-quality 3D visualization images with the smallest mean squared error of about 7.2186 and 7.3999 μm2, respectively. Our results are compared with the statistical randomization-based particle swarm optimization (SR-PSO). The results show that the BCDL-GWO with ICP is better than those from the SR-PSO. However, the computational complexities of both methods are similar.


Introduction
Three-dimensional reconstruction from multi-view images has been used in many applications, including in orthodontics.It is also used in diagnostic and treatment planning processes in adults and children with dental caries [1,2], especially dental caries that are a cause of chronic diseases in children [3].With recent technology, e.g., laser or CT, there are several 3D reconstructions from multimodal images research works [4][5][6][7][8][9] from those technologies.However, in rural areas, there is limited access to these sophisticated devices and also insufficient oral healthcare [10].In Thailand, the Dental Innovation Foundation, under royal patronage, has provided dental care access in rural communities for a long time.Due to limited access to sophisticated devices in those communities, taking multi-view teeth images inside children's mouths is a very difficult task.Hence, multiview teeth optical images are collected to be utilized in the 3D reconstruction system.In 3D reconstruction from multi-view images, there are several processes, including image registration, which is the transforming process of different sets of data into one coordinate system.There are existing works on 2D medical image registration [11][12][13][14][15][16][17][18][19][20][21][22][23][24][25][26].Since point cloud coordinates are used in our 3D teeth reconstruction, the 3D registration is more proper.In the literature, there are several works on 3D registration [27][28][29][30][31] utilizing several features in the registration process, including point cloud coordinates representing the 3D shapes of objects [32][33][34][35][36][37].These mentioned works used a variation swarm optimization (PSO) in the location matching between the source and target images.In the hope of improving the registration accuracy, there was study conducted on using the grey wolf optimizer (GWO) instead of PSO in 3D registration [38].However, only rotation and/or translation were used in those matching locations.Hence, in our previous work [39], the statistical randomization-based particle swarm optimization (SR-PSO) algorithm with the iterative closet point (ICP) method was used to find the optimal affine transform (translation, scaling, rotation, and shearing (shortened from a shearing mapping that displaces each point in a fixed direction by an amount proportional to its signed distance from a given line parallel to that direction)) between teeth optical images.
In particular, a system with the 3D registration using a modified grey wolf optimizer that can reconstruct a 3D image from teeth optical images is developed.However, due to a research ethical approval requirement, we will not use real images taken from children.Hence, we postulate scanned images from two commercial tooth models and then create point cloud images [39].To avoid premature convergence and to balance between exploration and exploitation, we modify the grey wolf optimization algorithm [40] with behavior considerations and dimensional learning strategies [41][42][43][44], called BCDL-GWO, to find the suitable affine transform between the source and target images.This can also enhance global and local searching and improve an ability to escape from a local optima.It has been shown in [41][42][43][44] that the BCDL-GWO is a good optimization tool when used to find suitable parameters in several applications, including engineering design problems (pressure vessel design, tension/compression spring design, and welded beam design problems), biomedical real-life problems (breast cancer and heart disease detection), and 14 real-world problems from the 2011 IEEE Congress on Evolutionary Computation.Furthermore, to refine the resulting registration, the iterative closet point (ICP) method [45,46] is used because of its ability to refine registered results [11,34,39].In the final step, we reconstruct the 3D tooth models.

Registration Method
In this section, we will briefly review the 3D registration method used in this paper.The registration between two point cloud images (target (P = [p i ] M × 4 , M is the number of target point cloud points) and source (Q = [q j ] N × 4 , N is the number of source point cloud points) point cloud images) can be found by the following transformation: where H is the geometry transform estimated by finding the nearest neighbor between a set of point pairs (p j and q j ) [46,47], and f (•) is an objective function (minimum distance error between two corresponding points).Hence, the mean squared error (MSE) can be used as f (•) to find a suitable H* as In this case, there are 15 unknown parameters, i.e., 3, 3, 3, and 6 parameters for scaling (S), translation (T), rotation (R), and shearing (SH), respectively [48].To give a simpler equation, let cox = cos(ϕ x ), coy = cos(ϕ y ), coz = cos(ϕ z ), six = sin(ϕ x ), siy = sin(ϕ y ), and siz=sin(ϕ z ); then, the 3D transformation matrix H is computed as Hence, where It is worthwhile noting that a through j are non-rigid transformations resulting from the combination of scaling, shearing, and rotation properties.
To find the optimal H, the proposed grey wolf optimization algorithm with behavior considerations and dimensional learning strategies (BCDL-GWO) algorithm described in the following section is utilized.Table 1 shows the defined search space with 15-dimensional individuals in the swarm.

Overview of Grey Wolf Optimizer Algorithm (GWO)
The gray wolf optimizer (GWO) [40] algorithm is divided into five mathematical models, i.e., (1) social hierarchy, (2) encircling prey, (3) hunting prey, (4) attacking prey (exploitation), and (5) seeking prey (exploration).The wolves are first generated as a set of candidate solutions (search agents) by randomization.At each generation, the wolves, called omega (ω), are guided by their three leaders, named alpha (α), beta (β), and delta (δ), to find more favorable regions in search spaces while searching or hunting for prey.Let X = {x l |l = 1 . . .K }; x l = x lj |j = 1 . . .d be a set of K search agents (individuals) with d-dimensional feature space.The encircling behavior of the lth grey wolf (x l ) around the pth prey (x p ) in the jth dimension at iteration (t) is where D t p is the distance between x l and x p at iteration t.The A t p and C t p are defined as where the components of a decrease linearly from 2 to 0 over the course of iterations.r 1 and r 2 are random vectors in [0, 1].Therefore, each element in A t p will be a random value in [−a, a], whereas that in C t p will be a random value in [0, 2].The position update equation of each individual will follow the 3 leaders, i.e., α, β, and δ represented by x t α , x t β , and x t β , respectively. then 2.2.The Modified GWO Algorithm with Behavior Considerations and Dimensional Learning (BCDL-GWO) The modified GWO with behavior considerations and dimensional learning is based on the idea of [41][42][43][44].Firstly, we incorporate the Sine Cosine Algorithm (SCA) [41] in the alpha grey wolf to alleviate the unbalancing between exploration and exploitation and to help with the premature convergence by Each jth element of C t α , C t β , and C t δ is modified following the method in [49] as where r 3 , r 4 , and r 5 are the uniformly distributed random numbers in [0, 1].c is decreasing linearly from 1 to 0 over the course of iterations as follows: Therefore, each element in vector C is stochastically generated in [0, 2] in the first iteration and decreases to 1 at the final iteration.This process helps to provide a better exploration capability.The control vectors A t α , A t β , and A t δ are calculated as in Equation ( 10) with a calculated by [50] where T is the total number of iterations, and 0 < µ ≤ 2.
Our hypothesis is that the alpha grey wolf is the most important leader; hence, the updated position of each omega grey wolf is modified to [51] where From Equation ( 20), w 1 ≥ w 2 ≥ w 3 .w 1 is close to 1, and w 2 and w 3 are close to 0 at the beginning.These values will finally be close to 1/3 in the last iteration.Now, we are ready to incorporate real-life behavioral considerations into the algorithm [52] by discarding the wolves during the migration (prey searching) with low fitness values and allowing mating (crossover and mutation as in the genetic algorithm) to improve the pack's diversity.However, in our GWO, the first half of the iteration is focused on the exploration behavior (when |A| > 1), whereas the remaining half is transformed into exploitation behavior (when |A| < 1).Hence, we applied the Lévy flights (LF) [42] and Random Opposition Learning (ROL) [43] to improve the pack's diversity and to enhance the capability of the global and local search.The LF is also applied to each element j of the three leaders as The LF is defined as where and Γ is a standard grammar function.In the experiment, we set β to 1.5, and α decreases over time [53] as follows: with lb j and ub j are the lower and upper bound, respectively, of the search space in the jth dimension.To add in the ROL in the exploitation behavior [43], suppose xl = xl,j |j = 1 . . .d ; xl,j ∈ lb j , ub j be a d-dimensional vector with where x b l is the individual best [54] of the lth grey wolf.Now, we introduce a new candidate solution (x BC l ) as It is worth noting that the fitness function in the BCDL-GWO (f (•)) is the aforementioned MSE.
The new update position is For the dimensional learning part [55], we suppose x DL l = x DL l,j |j = 1 . . .d .The distance between the current and the next position of x t l is The neighborhood of x t l is defined by Hence, each jth element of x DL l is calculated by where x t n and x t r are randomly selected from N t l and search agents, respectively.Then, the updated position will be Finally, the position of each grey wolf will only change if the next fitness value is better than the current one.Hence, the final update position equation will be The BCDL-GWO algorithm is summarized as shown in Algorithm 1.
The BCDL-GWO optimal solution is selected from the global best in the last population.The iterative closest point algorithm (ICP) method, as in [45,55] with the Nelder-Mead simplex method [56], is utilized to fine-tune the registration results.Table 2 shows the parameters used in our experiment.

Experimental Results
The diagram of the 3D reconstruction with the optimal transformation matrix H −1 (transform from source point cloud to target point cloud) found by the BCDL-GWO is shown in Figure 1.To fine-tune the resultant H −1 , the ICP method is used.Finally, the 3D tooth models are reconstructed based on the registered source and target point clouds.

Experimental Results
The diagram of the 3D reconstruction with the optimal transformation matrix H −1 (transform from source point cloud to target point cloud) found by the BCDL-GWO is shown in Figure 1.To fine-tune the resultant H −1 , the ICP method is used.Finally, the 3D tooth models are reconstructed based on the registered source and target point clouds.We first test our system on synthetic cylindrical and pyramid shapes as shown in Figure 2. We first test our system on synthetic cylindrical and pyramid shapes as shown in Figure 2.
The transformation matrix used to generate from source to target point clouds for both shapes is Hence, the H −1 that the BCDL-GWO needs to find is We compare the results from the BCDL-GWO with µ varied from 0.5 to 2.0 with a step size of 0.5 with those from our previous work (SR-PSO) [39].Tables 3 and 4 show the best registration MSEs results using the BCDL-GWO algorithm with and without refining the ICP method on the synthetic cylindrical and pyramid shapes, respectively.The best result on the synthetic cylindrical shape with MSE of 2.71 × 10 −27 from BCDL-GWO is at µ = 0.5.Whereas that on the synthetic pyramid shape with MSE of 7.79 × 10 −20 is at the same µ.The final best H −1 from both synthetic shapes is the same.Because the MSE is extremely small, both final best H −1 are the same as H −1 shown in Equation (34).The resulting registration of two cylindrical images and two pyramid images are shown in Figure 3a,b, respectively.It can be said that the best result from the BCDL-GWO is comparable with the best one from SR-PSO (α = 1.5).And when we look at the results from the BCDL-GWO with the other µ, those are better than that from the SR-PSO with the other α.To confirm this result, we also report the average ± standard deviation from several experiments of this algorithm on the same data set shown in Tables 5 and 6 for the synthetic cylindrical and pyramid shapes, respectively.Table 3.The mean squared error (MSE) in pixels 2 of the registration on the synthetic cylindrical shape.[39] α   We also compare the results from several experiments of both synthetic data sets shown in Tables 7 and 8 with those achieved by the butterfly optimization algorithm (BOA) [57], Harris hawks optimization (HHO) [58], slime mold algorithm (SMA) [59], and whale optimization algorithm (WOA) [60], whereas each method has the best parameter setting.It can be seen that the result from BCDL-GWO without ICP is better than that from all compared methods with ICP.Hence, we can assume that our BCDL-GWO can escape local minima.

Research work in
We also provide an indirect registration comparison of our BCDL-GWO with other methods without the utilization of swarm intelligence, i.e., Zhan et al. [33], Li et al. [61], and Du et al. [62], and as shown in Table 9.The results also confirm that our BCDL-GWO provides better results than its counterparts.Table 6.The average registration mean squared error (MSE) ± standard deviation in pixels 2 on the synthetic pyramid shape.[39] α     From the synthetic data set results, we are certain that the BCDL-GWO can be used in the tooth model 3D reconstruction.The regular tooth model and orthodontic tooth model from [39] were used in the experiment.For each model, six consecutive point cloud coordinate (x, y, z) views with an interval of 30 degrees are used in the experiment.Table 10 shows the information on the tooth point cloud data.

SR-PSO
In this experiment, the size of the original image in each view was randomly sampled to 60% with the assumption that there was an overlap between each consecutive view.The voxel hull method [63][64][65] was used to select representative points inside the overlapping area.After that, the registration process with the parameter setting shown in Table 2 was implemented.Since there were six consecutive views, the hierarchical registration with F = 6 was used to increase the registration performance shown in Figure 4.The survival at each level was the best final registration result (BCDL-GWO algorithm with the ICP method), and that result proceeded to the next level of the hierarchical registration.
Table 11 shows the registration MSE results from the BCDL-GWO without the ICP of the regular tooth model at hierarchical level 1, whereas those with ICP are shown in Table 12.Figures 5 and 6 show the best registration of each consecutive pair without and with ICP, respectively.In this experiment, the size of the original image in each view was randomly sampled to 60% with the assumption that there was an overlap between each consecutive view.The voxel hull method [63][64][65] was used to select representative points inside the overlapping area.After that, the registration process with the parameter setting shown in Table 2 was implemented.Since there were six consecutive views, the hierarchical registration with F = 6 was used to increase the registration performance shown in Figure 4.The survival at each level was the best final registration result (BCDL-GWO algorithm with the ICP method), and that result proceeded to the next level of the hierarchical registration.Table 11 shows the registration MSE results from the BCDL-GWO without the ICP of the regular tooth model at hierarchical level 1, whereas those with ICP are shown in Table 12.Figures 5 and 6 show the best registration of each consecutive pair without and with ICP, respectively.We compare our results with those from SR-PSO [39] as well.The MSEs of the regular tooth model at hierarchical level 1 are shown in Table 13.From the results, we can see that both methods provide comparable results.However, when we look at the MSEs of the final 3D reconstruction of the regular tooth model from six consecutive views shown in Table 14, we can see that the results from the best BCDL-GWO with ICP (7.2186 µm 2 when µ = 0.5) are better than SR-PSO with ICP (7.3666 µm 2 ).The final registration of the regular tooth model is shown in Figure 7.We can see that the reconstruction result provides a good visualization.

SR-PSO with ICP BCDL-GWO with ICP
Bioengineering 2024, 11, x FOR PEER REVIEW 16 of 25 We compare our results with those from SR-PSO [39] as well.The MSEs of the regular tooth model at hierarchical level 1 are shown in Table 13.From the results, we can see that both methods provide comparable results.However, when we look at the MSEs of the final 3D reconstruction of the regular tooth model from six consecutive views shown in Table 14, we can see that the results from the best BCDL-GWO with ICP (7.2186 µm 2 when µ = 0.5) are better than SR-PSO with ICP (7.3666 µm 2 ).The final registration of the regular tooth model is shown in Figure 7.We can see that the reconstruction result provides a good visualization.Finally, we implement the BCDL-GWO on the orthodontic tooth model to observe more experiments.The MSE registration results at hierarchical level 1 from the BCDL-GWO and BCDL-GWO with the ICP are shown in Tables 15 and 16, respectively.Figures 8 and 9 show the best registration results for each consecutive pair.Finally, we implement the BCDL-GWO on the orthodontic tooth model to observe more experiments.The MSE registration results at hierarchical level 1 from the BCDL-GWO and BCDL-GWO with the ICP are shown in Tables 15 and 16, respectively.Figures 8 and 9 show the best registration results for each consecutive pair.Again, we compare our MSEs on the registration results with those from the SR-PSO, as shown in Table 17.We can see that our proposed without the ICP method provides better registration results than the SR-PSO without the ICP method, except for a pair of 2 vs. 3.However, for the fine registration, our BCDL-GWO with ICP is comparable with those from the SR-PSO with ICP.But when we look at the final registration 3D orthodontic reconstruction model, as shown in Table 18, the results from the BCDL-GWO with ICP are better than the best results from the SR-PSO with ICP (7.4130 µm 2 ).While the best result from the BCDL-GWO with ICP at µ = 1.5 is 7.3999 µm 2 .The final 3D reconstruction of the orthodontic tooth model is shown in Figure 10.We can see that the reconstruction result can still provide a good visualization.

SR-PSO with ICP BCDL-GWO with ICP
* K denotes population size, and D indicates the number of dimensions in search spaces.To confirm that our method is good enough, we indirectly compare our results with those existing methods in the literature on the dental 3D registration data sets.However, those methods were performed on different data sets.The comparison results are shown in Table 19.Again, the results from our BCDL-GWO are better.One might wonder what the computational complexity of BCDL-GWO is compared with the SR-PSO, BOA, HHO, SMA, and WOA shown in Table 20.To compute the complexity of the BCDL-GWO, we start with the population initialization step.Since there are K grey wolves in the population and each grey wolf is represented by a D-dimensional vector, the computational complexity in this step is O(K × D).For the control parameter step, the first operation is the GWO-SCA, which needs O(K × D).The next operation in this step is Equation ( 15), which will need O(D).Then, both the position update and fitness evaluation will need O(K × D).Hence, in this step, it is O(K × D).However, in the fitness comparison step, it will need O(K).The next step the behavior consideration procedure.In this step, the new candidate solution calculation from Equation (26) will need O(K × D).The fitness calculation in this step will need O(K × D).The update position using Equation ( 27) will be O(K × D).Hence, in this step, the complexity will be O(K × D).Finally, in the dimensional learning procedure step, the distance calculation will need O(K × D).Equation ( 29) used in the search agents will need O(K × D 2 ).Again, the position update needs O(K × D).However, in this step, the total complexity is O(K × D 2 ).Since there are T iterations, the total complexity of the BCDL-GWO will be O(T × K × D 2 ).For other algorithms, the complexities are calculated similarly.Table 20 shows the Big O of each step in each algorithm.The complexities of all algorithms are very similar.Even though the BCDL-GWO has a slightly higher complexity than the others, the tradeoff with the performance of our BCDL-GWO is still good.

Conclusions
To help in dental diagnostic and treatment planning in rural areas with limited access to sophisticated devices, a 3D reconstruction from multi-view optical images is needed.To provide a good 3D reconstruction, a good 3D registration process is required.In this paper, we developed the grey wolf optimization algorithm with behavior considerations and dimensional learning strategies (BCDL-GWO) with iterative closet point (ICP) to find the optimal affine transform in the 3D registration process.We compare the results with those from the statistical randomization-based particle swarm optimization (SR-PSO).We found that the final best result of BCDL-GWO with the ICP yields a mean squared error (MSE) of 7.2186 µm 2 for 3D reconstruction from six consecutive views of the regular tooth model, whereas that of SR-PSO with the ICP method is 7.3666 µm 2 .The MSE of the BCDL-GWO with the ICP method is 7.3999 µm 2 for the orthodontic tooth model, while the SR-PSO with the ICP provides 7.4130 µm 2 .We can say that the 3D reconstruction of the regular and orthodontic tooth models from the BCDL-GWO with ICP is better than the SR-PSO with ICP.
We also estimate the computational complexity of both the BCDL-GWO and the SR-PSO.We could say that they are comparable.However, from the nature of the BCDL-GWO, we can also say that it can cope with a premature convergence, an unbalance between exploration and exploitation, and finally, it increases a pack's diversity.
Currently, there is only one research work involving the 3D model to assess dental caries [70].This shows that there is a need for a 3D model for dental caries assessment.

Figure 2 .
Figure 2. The original shape (target point cloud) and the 3D transformation (source point cloud) of (a) synthetic cylindrical and (b) synthetic pyramid shapes.

Figure 3 .
Figure 3.The final best registration using BCDL-GWO with µ = 0.5 of synthetic (a) cylindrical and (b) pyramid images.

Figure 4 .
Figure 4. Hierarchical structure for multiple-views registration.Table 11.MSE of BCDL-GWO on the regular tooth model at hierarchical level 1.

Figure 7 .
Figure 7.The final registration result of six consecutive views of the regular tooth model from the BCDL-GWO with the ICP.

Figure 7 .
Figure 7.The final registration result of six consecutive views of the regular tooth model from the BCDL-GWO with the ICP.

Figure 10 .
Figure 10.The final registration result of six consecutive views of the orthodontic tooth model from the BCDL-GWO with the ICP.

Figure 10 .
Figure 10.The final registration result of six consecutive views of the orthodontic tooth model from the BCDL-GWO with the ICP.

Table 20 .*
Computational time complexities of BCDL-GWO and SR-PSO.K denotes population size, and D indicates the number of dimensions in search spaces.

Table 1 .
Parameters boundaries in optimization process.

Table 4 .
The mean squared error (MSE) in pixels2of the registration on the synthetic pyramid shape.

Table 5 .
The average registration mean squared error (MSE) ± standard deviation in pixels 2 on the synthetic cylindrical shape.

Table 7 .
The best average registration mean squared error (MSE) ± standard deviation in pixels 2 on the synthetic cylindrical shapes.

Table 8 .
The best average registration mean squared error (MSE) ± standard deviation in pixels 2 on the synthetic pyramid shapes.

Table 9 .
The indirect comparison of several registration data sets with other existing methods.

Table 10 .
Tooth data set information.

Table 11 .
MSE of BCDL-GWO on the regular tooth model at hierarchical level 1.

Table 12 .
MSE of the BCDL-GWO with ICP on the regular tooth model at hierarchical level 1.

Table 12 .
MSE of the BCDL-GWO with ICP on the regular tooth model at hierarchical level 1.

Table 13 .
MSE comparison on regular tooth model at level 1.

Table 14 .
MSE of the final registration of six consecutive views (µm 2 ) for the regular tooth model (the best value is in bold).

Table 13 .
MSE comparison on regular tooth model at level 1.

Table 14 .
MSE of the final registration of six consecutive views (µm 2 ) for the regular tooth model (the best value is in bold).

Table 15 .
MSE of BCDL-GWO on the orthodontic tooth model at hierarchical level 1.

Table 16 .
MSE of BCDL-GWO with ICP on the orthodontic tooth model at hierarchical level 1.

Table 17 .
MSE comparison on orthodontic tooth model at hierarchical level 1.

Table 18 .
MSE of the final registration of six consecutive views (micrometer 2 ) for the orthodontic-tooth model (the best value is in bold).

Table 19 .
Indirect comparison results on dental 3D registration data sets.