Study sample
This study was approved by the Research Ethics Board of the Stomatological Hospital of Chongqing Medical University (No. 2020-013). All patients gave informed consent prior to participation. Patients for this study were recruited from consecutive adult patients visiting the Department of Orthodontics, Stomatological Hospital of Chongqing Medical University, Chongqing, China. A total of 18 females aged 18–26 years were recruited for this study. The mean age of the patients were 22.6 years. Of the 18 patients, 15 received extraction of four first or second premolars prior to the start of orthodontic treatment while 3 patients received orthodontic treatment without tooth extraction.
Patients eligible for this study should satisfy all of the following criteria: adult females of Chinese ethnicity, ANB angle between 0 to 4 degrees, and Body Mass Index within the range of 18.5 to 25 \(kg/{m}^{2}\), which represented individuals of normal weight (18). Patients with obvious facial asymmetry, craniofacial anomalies, previous history of orthodontic treatment, and defective dentitions were excluded.
All patients were treated with the same fixed appliances (0.022 × 0.028-inch bracket slot). Mandibular bracket was bonded one month after treatment began. The order in which nickel-titanium archwire sequence used were 0.012 inch, 0.014 inch, 0.016 inch, 0.016\(\times\)0.022 inch, and 0.018\(\times\)0.025 inch.
Facial surface imaging
Digital facial stereophotogrammetry (Morpheus 3D, Korea) was used to capture 3D facial surfaces for each individual. Patients were imaged following standard facial image acquisition protocol (19). Patients were asked to gently close mouth, maintain neutral facial expression, and assume natural head position during imaging. Images were taken at three treatment stages: baseline (T0), three months after treatment initiation (T1), and six months after treatment initiation (T2).
Spatially dense facial quasi-landmarking
3D facial images obtained from the Morpheus 3D systems was stored in .M3D format, which was converted to the OBJ format by the company. Facial images in OBJ format was further converted to ASCII PLY format. Each facial image in PLY format was imported into the IDAV Landmark Editor v.3.0.0.6 to digitize five anchoring points (right exocanthus, left exocanthus, pronasale, right cheilion, left cheilion) in a fixed order (20).
Facial images in OBJ format, together with coordinates of the five anchoring points, were imported into the MeshMonk toolbox of MATLAB (R2018b) for spatially dense facial quasi-landmarking (21). Based on the five anchoring points, an anthropometric mask was mapped to each facial image through rigid and non-rigid registration algorithms. This resulted in 7160 3D facial quasi-landmarks that capture facial region of interest while removing irrelevant structures such as hair, ears, and any dissociated polygons (22).
Generalized Procrustes Analysis
Human face is internally symmetric around the midsagittal plane (23). Each quasi-landmark on the right side had a homologous quasi-landmark on the left side. We reflected and relabeled the quasi-landmark configuration of each patient’s face. All configurations and their relabeled reflections were superimposed through Generalized Procrustes Analysis (GPA). This removed among-configuration variation in size, location, and orientation (24) and resulted in Procrustes shape coordinates that characterized facial shape (25).
Although orthodontic treatment may impact left and right side of the face differentially, this is not the focus of the present study. Therefore, all GM facial analysis in this study was performed based on the symmetric component of facial shape.
Changes in mean facial shape across treatment stages
Statistical significance of changes of facial shape from T0 through T1 to T2 was evaluated using permutational multivariate analysis of variance (MANOVA) of distance matrices comprising pairwise Procrustes distance (PD) between configurations (26). Procrustes shape coordinates were used as response variable and treatment stage, treatment modality (tooth extraction vs non-extraction), and their interactions were used as explanatory variables. Unlike traditional MANOVA in which groups under comparison were independent, facial shape at T0, T1, and T2 were correlated because they represented repeated measurements on the same group of patients. Therefore, we constrained permutation to be within each participant, as recommended by Anderson and Braak (2003), to preclude the confounding effect of among-patient variation on evaluations of facial shape changes during treatment (13). Empirical p values were calculated as the probability that the permuted pseudo F-statistic was larger or equal to the observed F-statistic (26), based on 10000 permutations. In addition, coefficient of partial determination (\(partial {R}^{2}\)) was calculated for each explanatory variable following the formula by Kleinbaum et al. (27).
To gain an understanding of facial regions where shape changed significantly with treatment, we determined statistical significance and relative magnitude of positional changes of each facial quasi-landmark. The analyses were performed based on the same model used above while replacing Procrustes shape coordinates for the entire face with 3D coordinates for each facial quasi-landmark. Similar approach has been adopted by Zaidi et al. (28) and Claes et al. (29).
To ascertain the specific treatment stages between which facial shape changed significantly, post-hoc pairwise comparisons were performed. The magnitude of shape change was quantified by PD between mean facial shape of the two treatment stages under comparison. The permutational MANOVA was performed using the Adonis function of the vegan package version 2.5-6 in R version 4.0.0 (30, 31). 3D Facial heatmaps were used for visualization of study findings.
Changes in variance of facial shape across treatment stages
Variance of the facial shape was estimated as Procrustes variance (32). The morphol.disparity function in geomorph package version 3.2.1 in R was used to quantify pairwise differences in variance of facial shape among all treatment stage-by-modality groups (29,33).
In additional to the above analyses of facial shape, facial form was obtained by multiplying Procrustes shape coordinates of each facial image with its corresponding centroid size. Statistical significance of changes in mean and variance of facial form was performed following the same procedures described above. The level of statistical significance was set at 0.05 for all analyses. Quasi-landmark repeatability error of MeshMonk has been reported to be as low as 0.002 in unit of Procrustes distance, which is robust to manual landmark digitization errors associated with digitization of the five anchoring points (34).