Predictors of Efavirenz Plasma Exposure, Auto-Induction Profile, and Effect of Pharmacogenetic Variations among HIV-Infected Children in Ethiopia: A Prospective Cohort Study

(1) Background: Efavirenz plasma concentration displays wide between-patient variability partly due to pharmacogenetic variation and autoinduction. Pediatric data on efavirenz pharmacokinetics and the relevance of pharmacogenetic variation are scarce, particularly from sub-Saharan Africa, where >90% of HIV-infected children live and population genetic diversity is extensive. We prospectively investigated the short- and long-term effects of efavirenz auto-induction on plasma drug exposure and the influence of pharmacogenetics among HIV-infected Ethiopian children. (2) Method: Treatment-naïve HIV-infected children aged 3–16 years old (n = 111) were enrolled prospectively to initiate efavirenz-based combination antiretroviral therapy (cART). Plasma efavirenz concentrations were quantified at 4, 8, 12, 24, and 48 weeks of cART. Genotyping for CYP2B6, CYP3A5, UGT2B7, ABCB1, and SLCO1B1 common functional variant alleles was performed. (3) Results: The efavirenz plasma concentration reached a peak at two months, declined by the 3rd month, and stabilized thereafter, with no significant difference in geometric mean over time. On average, one-fourth of the children had plasma efavirenz concentrations ≥4 µg/mL. On multivariate analysis, CYP2B6*6 and ABCB1c.3435 C > T genotypes and low pre-treatment low-density lipoprotein (LDL) were significantly associated with higher plasma efavirenz concentration regardless of treatment duration. Duration of cART, sex, age, nutritional status, weight, and SLCO1B, CYP3A5, UGT2B7, and ABCB1 rs3842 genotypes were not significant predictors of efavirenz plasma exposure. (4) Conclusion: Pre-treatment LDL cholesterol and CYP2B6*6 and ABCB1c.3435 C > T genotypes predict efavirenz plasma exposure among HIV-infected children, but treatment-duration-dependent changes in plasma efavirenz exposure due to auto-induction are not statistically significant.


Introduction
HIV infections remain a global public health problem, particularly in sub-Saharan Africa (SSA). According to the 2020 UNAIDS estimates, around two-thirds of the 1.7 million seven HIV/AIDS treatment centers located in two adjacent administrative regions in Ethiopia: Oromia and Southern Nations Nationalities and Peoples Region (SNNPR), which have diverse ethnic groups. In the Oromia region, the study centers were Adama Hospital, Asela Referral and Teaching Hospital, and Shashemene General Hospital. The study sites in SNNPR region were Hawassa University Comprehensive Specialized Hospital, Adare General Hospital, Hosana Referral Hospital, and Arbaminch General Hospital. These regions were selected because of high HIV prevalence and population diversity, which is generalizable to the nation.

Antiretroviral Therapy and Laboratory Analysis
At enrollment, complete medical history, clinical examinations, and laboratory tests including complete and differential blood counts, platelet counts, CD4 count, HIV viral load, and liver enzymes were performed. Subsequently, study participants initiated a firstline pediatric cART regimen containing efavirenz (EFV) and lamivudine (3TC) with either zidovudine (AZT), abacavir (ABC), or tenofovir (TDF) following the Ethiopian national pediatric HIV treatment guideline [11]. All children ≥3 years old who weighed ≥10 kg, and adolescents ≥10 years older who weighed <35 kg received weight-band dosing of EFV + 3TC + AZT or ABC as described in the national HIV treatment guideline [11]. Adolescents ≥10 years old who weighed ≥35 kg received an EFV 600 mg + TDF + 3TC fixeddose combination once daily. Treatment adherence was assessed by parent/caregiver and child/adolescent interview at each study visit. The study participants were prospectively followed up to 48 weeks while on cART to monitor plasma efavirenz exposure. Whole blood samples at 12 to 16 h post efavirenz dose were collected on the 4th, 8th, 12th, 24th, and 48th weeks of starting ART. Plasma samples were prepared immediately, and aliquots were kept at −80 • C for analysis.

Chemicals and Reagents
Efavirenz reference standard was purchased from Merck (Darmstadt, Germany). Deuterated-efavirenz (rac Efavirenz-d4) was purchased from Toronto Research Chemicals (Toronto, ON, Canada). Acetone, acetonitrile, ammonium acetate, isopropanol, methanol, and acetic acid of mass spectrometry (MS) grade were purchased from Merck (Darmstadt, Germany). Blank human plasma was obtained from the blood bank of the Karolinska University Hospital Huddinge (Stockholm, Sweden).

LC-MS/MS Method and Validation
An LC-MS/MS system consisting of an Acquity Ultra Performance LC-system coupled to a Xevo TQ-S Micro (Waters, Milford, MA, USA) was used to quantify plasma efavirenz concentration. The MS/MS was operated in positive ionization mode. The data were processed using MassLynx 4.2 (Waters). In brief, plasma calibration samples were prepared by spiking blank plasma samples with efavirenz and were included in each analytical run. Quality control samples were also prepared by spiking plasma blanks to obtain low and high concentrations. The lower and upper limits of quantification were 15.8 and 15,800 ng/mL, respectively.
Plasma efavirenz concentration was determined by direct injection of the supernatant after protein precipitation using acetonitrile. In brief, a 50 µL aliquot of plasma samples was added into 200 µL acetonitrile solution containing internal standard (500 ng/mL of efavirenz-d4 in methanol), vortexed for 30 s followed by centrifugation for 5 min at 2100× g. An amount of 5 µL of the supernatant was injected onto the LC-MS/MS system. The chromatographic system consisted of an Acquity UPLC BEH C18 reversed-phase column 2.1 × 50 mm, 1.8 µm with 0.1% aqueous formic acid as mobile phase A and MeOH as mobile phase B. The initial composition of the mobile phase was 50% B, followed by a linear gradient to 99% in 1.5 min, with a flow rate of 0.4 mL/min. Total analysis/run time for each sample was approximately 3 min.
Efavirenz was monitored by the transition m/z 316 > 168 and the efavirenz-d4 by 320 > 172. The calibration curve was constructed by linear regression of the analyte/internal standard area ratios with an applied weighing of 1/x2. Accuracy was within ±10% throughout the quantification range. Precision was below 6 CV% except for LLOQ (below 10 CV%). The method was validated according to the European Medicines Agency Guideline on bioanalytical method validation [36].

Statistical Data Analysis
Study participants' sociodemographic and baseline clinical parameters are summarized as the median and interquartile range (IQR) or as frequency and percentages. For assessment of the nutritional status of the children, anthropometric measurements were converted to height for age Z score (HAZ) and weight for age Z score (WAZ) using WHO Anthro-Plus software version 1.0.4. Body mass index (BMI) for age percentiles was calculated using the age-and-sex-specific percentile for BMI as per the Centers for Disease Control calculator for children and teens aged 2-19 years.
Efavirenz plasma concentrations were transformed to log10 values before statistical analysis. The dataset was reshaped from a wide to a long format in which each individual subject had five records or less for weeks 4, 8, 12, 24, and 48. One-way analysis of variance (ANOVA) was used to identify genotype association with efavirenz plasma concentration measured at different weeks of ART. Change in efavirenz concentration over time was analyzed using a generalized linear mixed-effect model (LME). Univariate followed by multivariate analyses using LME models were conducted to determine predictors of efavirenz concentrations and to quantify the unexplained between-and within-subject variations in efavirenz exposure at different weeks of treatment. The LME modeling estimated fixed effect (predictor coefficients) and random effect parameters (between and within-subject variabilities in efavirenz concentrations). The LME models were therefore formulated as follows. EFV| = |α i × Week i | + β j × Predictor j | + , where i is a vector of weeks on treatment, in this case, these are 4, 8, 12, 24, and 48; represents random-effect parameters; j represents additional predictors that were separately added to the model for individual evaluation. The LME model was developed in 3 steps. In the first step, time on treatment (weeks) was used as the only predictor in the model. In the second step, univariate models were developed, which included time on treatment plus one additional predictor. In the last step, all additional predictors that were physiologically plausible, non-collinear, and with p-values ≤ 0.2 were included in the final LME models.
Descriptive statistics analysis for baseline characteristics was performed using Statistical Package for Social Sciences (SPSS) software, version 20.0 (IBM Corporation, Somers, NY, USA). The R software package (The R Foundation for Statistical Computing, Vienna, Austria) and GraphPad Prism version 7 (version 7, San Diego, CA, USA) were used for the analysis and presentation of the results. p values < 0.05 were considered statistically significant.

Ethical Considerations
The study received ethical approval from the SNNPR Regional Health Bureau Institutional Review Board (IRB), Addis Ababa University College of Health Sciences (Protocol No: 053/16/Pharma), and the National Research Ethics Review Committee, Ministry of Science and Higher Education-Ethiopia (Ref NO 3.10/166/2016). Prior to giving consent, study participants and their parents/guardians received written and oral information in their local language about the study. Written informed consent from parents or legal guardians was obtained. For participants ≤12 years, written informed consent was obtained from their parent or guardian. For participants >12 years of age, informed consent was obtained from the parent or guardian, and assent was obtained from the study participant.

Results
A total of 111 newly diagnosed and ART-naïve HIV-infected children were enrolled prospectively to participate in this study, of whom 63 (56.8%) were boys. Thirty-nine (35.1%) children were orphans (having lost one parent or both). Eleven (9.9%) children did not complete the study due to lost to follow-up or transfer out from the study area.

Baseline Characteristics of Study Participants
The sociodemographic, clinical, and biochemical characteristics of the study participants are presented in Table 1. The median age of the children was nine years. The dose of efavirenz received per kilogram body weight per day ranged between 12 and 24 mg/Kg with a median dose of 17 mg/kg. Study participants' body weight ranged between 8 and 49 Kg, with a median weight of 22 Kg. Following the WHO Z-score definition (<−2 ZS), 37.1% of the children were underweight; 33.9% were stunted, and 43.4% had a BMI for age percentile of less than the 5th percentile. More than half (60%) of the participants were in WHO clinical stages 2, 3, or 4. Most of the participants (53.2%) received a tenofovir TDF)/lamivudine/efavirenz regimen followed by abacavir/lamivudine/efavirenz (33%), and zidovudine/lamivudine/efavirenz regimen (13.6%).
The pattern and extent of change in plasma efavirenz concentration in each individual and the geometric mean at the 4th, 8th, 12th, 24th, and 48th weeks of ART are presented in Figure 1. Mixed effect model analysis indicated no significant within-subject variation in the geometric mean of plasma efavirenz concentrations over time (p = 0.89). Though not significant, the pattern of change in the mean log plasma efavirenz concentration over time displayed a moderate increase at the 8th week of ART, declined by the 12th week, and stabilized after that.

Effect of Genotype on Plasma Efavirenz Concentration at Each Study Time Point
Association of genotype with plasma efavirenz concentration measured at the 4th, 8th 12th, 24th, and 48th weeks of efavirenz-based antiretroviral treatment among HIV-infected children using one-way analysis of variance is presented in Supplementary Table S1. Children with the CYP2B6*6/*6 genotype had significantly higher mean log efavirenz plasma concentration than homozygous wild-type (CYP2B6*1/*1) or heterozygous (CYP2B6*1/*6) at week 4 (p = 0.004, 0.09), week 8 (p = 0.002, <0.001), week 12 (p = 0.006, 0.37), week 24 (p = 0.002, 0.01), and week 48 (p = 0.01,0.01). The CYP2B6*6 genotype was significantly associated with therapeutic plasma efavirenz concentration categories (p < 0.001). There was no significant impact of CYP3A5, UGT2B7, SLCO1B1, and ABCB1 c.4036A > G genotypes on efavirenz plasma concentration at all study time points, except for ABCB1 c.3435C > T. Carriers of at least one ABCB1 c.3435C > T allele had significantly higher mean log efavirenz plasma concentration than those with the ABCB1 c.3435 C/C genotype at weeks 4 and 12 with p-values of 0.004 and 0.03, respectively. Although not statistically significant, children with at least one variant allele of ABCB1 c.4036A > G had a higher mean log efavirenz plasma concentration at all treatment weeks. The pattern and extent of change in the geometric mean of plasma efavirenz concentration at the 4th, 8th, 12th, 24th, and 48th weeks of ART among the different CYP2B6 genotypes stratified by age group are presented in Figure 2.

Effect of Genotype on Plasma Efavirenz Concentration at Each Study Time Point
Association of genotype with plasma efavirenz concentration measured at the 4th, 8th 12th, 24th, and 48th weeks of efavirenz-based antiretroviral treatment among HIVinfected children using one-way analysis of variance is presented in Supplementary Table  S1. Children with the CYP2B6*6/*6 genotype had significantly higher mean log efavirenz plasma concentration than homozygous wild-type (CYP2B6*1/*1) or heterozygous (CYP2B6*1/*6) at week 4 (p = 0.004, 0.09), week 8 (p = 0.002, <0.001), week 12 (p = 0.006, 0.37), week 24 (p = 0.002, 0.01), and week 48 (p = 0.01,0.01). The CYP2B6*6 genotype was significantly associated with therapeutic plasma efavirenz concentration categories (p < 0.001). There was no significant impact of CYP3A5, UGT2B7, SLCO1B1, and ABCB1 c.4036A > G genotypes on efavirenz plasma concentration at all study time points, except for ABCB1 c.3435C > T. Carriers of at least one ABCB1 c.3435C > T allele had significantly higher mean log efavirenz plasma concentration than those with the ABCB1 c.3435 C/C genotype at weeks 4 and 12 with p-values of 0.004 and 0.03, respectively. Although not statistically significant, children with at least one variant allele of ABCB1 c.4036A > G had a higher mean log efavirenz plasma concentration at all treatment weeks. The pattern and extent of change in the geometric mean of plasma efavirenz concentration at the 4th, 8th, 12th, 24th, and 48th weeks of ART among the different CYP2B6 genotypes stratified by age group are presented in Figure 2.

Predictors of Plasma Efavirenz Exposure over Time
A univariate followed by a multivariate analysis was conducted to identify predictors of plasma exposure over time using the developed linear mixed effect (LME) models. The impact of sociodemographic, plausible clinical, and biochemical characteristics of the

Predictors of Plasma Efavirenz Exposure over Time
A univariate followed by a multivariate analysis was conducted to identify predictors of plasma exposure over time using the developed linear mixed effect (LME) models. The impact of sociodemographic, plausible clinical, and biochemical characteristics of the study participants presented in Table 1 were tested individually in a univariate analysis. In the following multivariate analysis, all predictors with p-values ≤ 0.2 from the univariate analysis were included in the model simultaneously, and the result is presented in Table 3. The unexplained between-subject and within-subject variability in efavirenz concentration was about 28% and 34%, respectively. CYP2B6*6 (c.516G > T), ABCB1 c. 3435 C > T genotype, and baseline LDL were significant predictors of efavirenz exposure regardless of weeks of treatment ( Table 3). The coefficients of estimates were positive for CYP2B6*6 and ABCB1 C3435T genotypes, which means carriers of these variant alleles had a higher plasma concentration than those with a wild-type allele. However, the coefficient was negative for pre-treatment LDL concentration, indicating an association of increased LDL level with a decrease in plasma efavirenz exposure. Duration of antiretroviral treatment (in weeks) as a categorical predictor was included in the model and found to have no significant effect on the plasma efavirenz exposure.

Predictors of Efavirenz Concentration over Time among CYP2B6 *1/*1 Genotype
A subgroup analysis was performed to identify predictors of plasma efavirenz concentration over time among study participants with the CYP2B6*1/*1 genotype. In a univariate analysis, the impact of sociodemographic, potential pre-treatment clinical, and biochemical characteristics of the study participants presented in Table 1 were tested individually. In the following multivariate analysis, all predictors with p-values ≤ 0.2 from the univariate analysis were included in the LME model simultaneously to predict plasma exposure over time, and the results are presented in Table 4. Among children with the CYP2B6*1/*1 genotype, the unexplained between-subject and within-subject variability in efavirenz concentration was about 22% and 26%, respectively. The effect of duration of therapy on plasma efavirenz concentration was noted at week 8 of ART. While pre-treatment higher liver enzymes (ALT and ALP), total cholesterol, and triglycerides were positive predictors, LDL was a negative predictor of plasma efavirenz exposure regardless of weeks on ART among children homozygous for the functional CYP2B6*1 variant allele.

Discussion
Previous studies reported the association of pharmacogenetic variations and plasma efavirenz concentration with antiretroviral-treatment-induced neurotoxicity and liver enzyme abnormalities among HIV-infected patients [16,34]. Therefore, identifying predictors of efavirenz plasma exposure is essential to identify patients at risk of developing treatmentassociated adverse events for personalized medicine. This study investigated the predictors of short-and long-term plasma efavirenz exposure among children on ART. The change in plasma efavirenz concentration over time was prospectively monitored as a surrogate marker to determine the time-course and extent of efavirenz autoinduction. Our main findings include (i) the presence of wide between-patient variability in efavirenz plasma exposure, but no significant time-dependent change over time among children. (ii) CYP2B6*6 (c.516G > T), ABCB1 3435C > T genotype, and low pre-treatment LDL concentration are significant predictors of high plasma efavirenz exposure regardless of the duration of therapy. (iii) Among the CYP2B6*1/*1 genotype, altered liver enzymes (ALT, ALP) and lipid profiles (total cholesterol, triglycerides, and LDL) influence plasma efavirenz exposure. (iv) There is no significant association of nutritional status (stunting and wasting), sex, age, SLCO1B1, CYP3A5, and UGT2B7 genotype on plasma efavirenz concentration or autoinduction profile. To our knowledge, this is the first study to explore the shortand long-term efavirenz plasma exposure, and autoinduction profile among HIV-infected children using efavirenz plasma level as an index and to examine any implications of pharmacogenetic variations.
Efavirenz induces its own metabolism by increasing the expression of hepatic CYP2B6 enzyme via the human constitutive androstane receptor (hCAR) [39]. Currently, treatment of HIV is lifelong, and concomitant medication for the treatment of opportunistic infection is common. Therefore, characterization of the time course and magnitude of enzyme induction is essential to avoid possible drug interactions. Studies in primary human hepatocytes have indicated that efavirenz autoinduction is non-linear, and timeand concentration-dependent [39,40]. Using an endogenous marker, we previously monitored the time course and extents of CYP3A and CYP2B6 induction among adult HIV patients from Tanzania [41][42][43], Uganda [16], and Ethiopia [31,33,35]. Results from adult HIV-infected patients indicate that efavirenz autoinduction is long-term and significantly reduces its plasma concentration over time, particularly among the CYP2B6*1/*1 genotype. This contrasts with our current findings among children, where minor changes in plasma efavirenz concentrations over time were observed regardless of the CYP2B6 genotype. Such inconsistent efavirenz autoinduction profiles between adults and children could be due to differences in age-related metabolic capacity or differences in administered dose. The observed non-significant effect of treatment duration on efavirenz plasma exposure over time may indicate that maximal induction reaches a "steady state" within the first four weeks of therapy at the given weight-band efavirenz dose for children. Nonetheless, the current efavirenz pediatric dosing based on the child bodyweight category results in no significant change on long-term plasma efavirenz exposure.
Notably, our study shows the existence of wide between-patient variability in efavirenz plasma concentration among children. In line with our previous reports from adult patients [18,19], the CYP2B6*6 genotype was significantly associated with plasma efavirenz concentration and therapeutic range categories among children. Although CYP2B6*6 allele frequency distribution is not significantly different between Ethiopian adults (29.7) [19] and children (30.1%), the proportion of patients in the various plasma efavirenz therapeutic range categories varied significantly between children and adult patients. The proportion of children with sub-therapeutic efavirenz plasma concentration levels (< 1 µg/mL) was much lower than adults [19]. However, the proportions of patients with plasma concentrations exceeding the safety threshold (≥4 µg/mL) were four-fold higher in children than adult patients [19]. The proportion of children with plasma efavirenz concentration ≥ 4 µg/mL at the 4th, 8th, 12th, 24th, and 48th weeks of ART were 19%, 30%, 29%, 22%, and 29%, respectively, while the proportion of adult patients with plasma efavirenz concentration ≥4 µg/mL at the 4th and 16th weeks of cART were 6% and 5%, respectively [19]. With the current recommended weight-band pediatric efavirenz dose, our findings indicate that on average one-fourth of Ethiopian HIV-infected children experience higher plasma drug exposure ≥4 µg/mL. This is a concern since such higher plasma efavirenz exposure is associated with treatment-induced neuropsychiatric manifestations and liver enzyme abnormalities [16,34], increased psychopathology [44], and poor long-term neurocognitive outcomes in HIV-infected children [45].
This study used linear mixed effect (LME) modeling to identify variables that predict plasma efavirenz concentration in children (Tables 3 and 4). In line with previous reports [15,18,19], CYP2B6*6 (c.516G > T) and ABCB1 c3435C > T genotypes were significant predictors of higher efavirenz plasma concentration. CYP2B6 is the main enzyme responsible for the metabolism of efavirenz. ABCB1 encodes P-glycoprotein, an ATP-dependent drug efflux pump, which is responsible for drug transport across extra-and intra-cellular membranes. Our findings reveal that sex, stunting, wasting, weight, and age had no significant role in predicting plasma efavirenz among children. The current pediatric HIV treatment guideline is only weight-based and does not consider enzyme maturation and genetic variations [46]. Hence, the current efavirenz dosing strategy for pediatric use might be inappropriate, especially for young children in sub-Saharan Africa, where the defective CYP2B6 defective variant allele occurs at a higher frequency. Consideration of a pharmacogenetics algorithm might be more useful for determining the optimal dose of efavirenz for pediatric patients [20,47].
Both HIV infection and efavirenz-based ART alter lipid profiles [48,49]. Recently we reported a significantly higher prevalence of lipid profile abnormalities (hypertriglyceridemia and reduced HDL level) among HIV-infected children and adolescents who were on cART as compared to treatment-naïve HIV-infected children [9]. A cross-sectional study from South Africa reported the association between higher plasma efavirenz concentrations with fasting plasma lipid lipoproteins among HIV-infected adult patients [50]. The current study is the first prospective study to report the association of pre-treatment lipid profile with efavirenz pharmacokinetics. VLDL/LDL acts as a drug carrier and regulates the transport and metabolism of drugs in the body, and associated drugs are eliminated faster from the body by lipoprotein apheresis, a blood purification therapy that selectively removes VLDL/LDL particles from the bloodstream [51]. In line with this, we found a significant association between low LDL-cholesterol level at baseline (pre-treatment) with higher efavirenz plasma exposure regardless of the duration of therapy. Lack of a significant change in LDL level by ART [9] means that pre-treatment LDL level can serve as a reliable predictor of efavirenz plasma exposure, as observed in this study.
This study has some limitations. We could not utilize changes in intrinsic hepatic clearance of efavirenz to study the impact of efavirenz induction on its plasma exposures, as was done in previous studies among healthy adult subjects. Intensive pharmacokinetics sampling, which requires multiple blood withdrawals at each study time point over a one-year period from vulnerable patient populations (HIV-infected children) just for study purposes, was not possible due to ethical concerns. Likewise, we could not perform fasting lipid profile analysis for ethical reasons as children did not fast for the study purpose. Instead, we characterized the time course of induction and extent of changes in efavirenz plasma concentration profile in the same individual over time, whereby each child was effectively functioning as their own control, lowering the level of unexplained variance.

Conclusions
We report the presence of wide between-patient variability in plasma efavirenz exposure among HIV-infected children, partly due to CYP2B6 and ABCB1 pharmacogenetic variations. In addition, pre-treatment plasma LDL is a significant predictor of plasma efavirenz exposure. Duration of therapy-dependent changes in plasma efavirenz concentration due to autoinduction (as observed in adult patients) is not significant with the current weight-band efavirenz dosing among children. About one-fourth of the children had plasma concentrations exceeding supratherapeutic plasma efavirenz concentration (4 µg/mL), indicating that the current weight-band efavirenz dosing maybe unnecessarily high. Weight and age, which are the basis for the current pediatric efavirenz dosing, had no significant effect on plasma efavirenz concentration. For estimating the optimum pediatric efavirenz dose, CYP2B6 genotype-based dose modifications may have more importance. Future pharmacokinetic and pharmacogenetic association studies with treatment outcome (efficacy and safety) measures are required to support a low dose efavirenz regimen among the pediatric population as confirmed in adult HIV patients.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/jpm11121303/s1, Table S1: Effect of genotype on plasma efavirenz concentration at the 4th, 8th 12th, 24th, and 48th weeks of antiretroviral treatment among HIV-infected children using one-way ANOVA.  For participants ≤12 years, written informed consent was obtained from their parent or guardian, while for participants >12 years of age, informed consent was obtained from the parent or guardian and assent obtained from the study participant.
Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.
Data Availability Statement: All data accessed and analyzed in this study are available in the article.