Associations between Toxoplasma gondii infection and steroid hormone levels in spotted hyenas.

Toxoplasma gondii is a common parasite that infects warm-blooded animals and influences host physiology. T. gondii is known to target the host's central nervous system, affecting circulating levels of steroid hormones, fear-related behaviors, and health, although these effects appear to vary among host taxa. Here, we investigated the relationship between T. gondii infection and levels of plasma testosterone and cortisol within a wild population of spotted hyenas (Crocuta crocuta, n = 109). In our analyses, we accounted for age and sex via stratified regression analyses. We detected a negative association between circulating plasma testosterone and T. gondii infection among female cubs and subadults as well as adult male hyenas. We found no associations between T. gondii infection and cortisol in any age class or sex group of hyenas. Our work adds to a growing body of literature by characterizing the relationship between T. gondii infection and physiology in a novel host in its natural habitat. In a broader context, our findings indicate that responses to infection vary with characteristics of the host and point to a clear need for additional studies and priorities for future work that include diverse taxa and ecological settings.

Introduction

Toxoplasma gondii is a globally distributed parasite that readily infects diverse warm-blooded vertebrates, although felids are its only known definitive hosts. Experimental studies show that these infections can alter the host's physiology, health, and behavior, at least in lab animals and humans (Abdoli, 2014; Dubey, 2016). T. gondii infections originate via multiple pathways involving different parasite life stages, including ingestion of oocysts shed from felids, consumption of bradyzoites encysted within intermediate host tissues, and congenital infections originating from tachyzoites with gliding motilities (Dubey, 1998; Elmore et al., 2010). After colonizing a new host, T. gondii rapidly migrates throughout the body and preferentially encysts within the central nervous system (Di Cristina et al., 2008; Dubey, 1998). In rodent models, brain cysts are unevenly distributed and higher densities have been reported in amygdalar structures that play a key role in the regulation of fear-related behavior (Vyas et al., 2007). Still evidence of tropism is limited since T. gondii cysts have been identified throughout most regions of the rodent brain in more recent studies (Berenreiterová et al., 2011; Di Cristina et al., 2008; McConkey et al., 2013).

T. gondii infections are also associated with modulation of circulating hormones, including testosterone and glucocorticoids, albeit both positive and negative associations have been reported in the literature (Table 1). These hormones are known to target the amygdala and other interconnected regions of the brain, including the hippocampus and hypothalamus, which in turn regulate circulating steroid hormone levels that influence behavior (Heany et al., 2015; McEwen et al., 2016) and the immune response (Bereshchenko et al., 2018). Two proximate explanations for the relationship between T. gondii infection and host hormone levels are 1) a direct mechanism in which localized cysts act directly on specific regions of the brain affecting physiology and behavior, and 2) an indirect mechanism in which immune response to infection and elevated cytokine levels influence brain functions and downstream physiology and behavior (McConkey et al., 2013; Webster and McConkey, 2010). Additional hypotheses, including alterations involving dopaminergic signaling, may also explain variation in host physiology related to T. gondii infection and have been reviewed in the literature (McConkey et al., 2013; Webster and McConkey, 2010).

Table 1

A summary of studies on the relationship between T. gondii infection (positive vs negative) and steroid hormone levels (Abdulai-Saiku and Vyas, 2017, Al-Masoudi et al., 2018, Babu et al., 2007, Borráz-León et al., 2021, Colosi et al., 2015, Hamdeh et al., 2015, Hassanain et al., 2017).

In laboratory settings and in humans, the relationship between T. gondii infection and steroid hormones (i.e., testosterone, and to a lesser extent glucocorticoids), is relatively well studied. However, findings from these studies are mixed. On the one hand, T. gondii-infected humans, especially adult males, often exhibit elevated circulating levels of testosterone (Flegr et al., 2008; Shirbazou et al., 2011; Zghair et al., 2015; Zouei et al., 2018). Testes of male lab rats infected with T. gondii also produce elevated levels of testosterone, potentially due to the upregulation of genes involved in testosterone production in infected hosts (Lim et al., 2013). Similarly, positive associations between T. gondii infection and glucocorticoids have been observed in humas (Abdelazeem et al., 2015; Shirbazou et al., 2011). However, other studies of humans and lab rodents have associated T. gondii infection with both reduced plasma testosterone concentrations (Abdoli et al., 2012; Eslamirad et al., 2013; Kaňková et al., 2011) and reduced glucocorticoid concentrations (Mitra et al., 2013). Localized lesions and the inflammation that T. gondii infections can induce within the host brain may partly explain this pattern. These lesions produce meningitis in the cortex as well as perivascular cuffing in the cortex, caudoputamen thalamus, and hypothalamus (Ihara et al., 2016), which can disrupt the signal cascade that stimulates steroid hormone production. This is also observable in cases of hypogonadotropic hypogonadism, a condition that is associated with T. gondii infection (Oktenli et al., 2004; Stahl et al., 1994). Taken together, these mixed findings indicate that T. gondii infection may covary with hormone levels in wild hosts, but the nature of these relationship(s) within wild hosts cannot be predicted without empirical data.

Here, we investigate the relationship between T. gondii infections and levels of plasma testosterone and cortisol within a wild population of spotted hyenas (Crocuta crocuta). We first describe variation in infection status and hormone levels in this novel study system. We then assess associations between T. gondii infection and circulating levels of plasma testosterone and cortisol in age class and sex stratified regression analyses. We conclude with suggestions for future work that could help elucidate the mechanisms through which T. gondii influences host physiology.

Methods

Literature review

We searched Google Scholar to identify studies that measured host hormone levels with respect to T. gondii infection. First, we conducted separate searches for the following five key terms: toxoplasmosis and testosterone, toxoplasmosis and cortisol, toxoplasmosis and sex hormones, toxoplasmosis and sex hormones in animals, toxoplasmosis testosterone and animals. Here, we examined the first 100 studies of each search and retained only studies that included empirical data relating the parasite to host hormone levels. We noted an overlap in the articles found by each search phrase, and that new articles were not identified beyond the first ∼40 articles. We then supplemented our literature with some ad hoc searches to find additional studies not identified by our search terms. We summarize these results in Table 1.

Study population

We used blood samples collected between 1990 and 1999 from 109 wild spotted hyenas darted in the Masai Mara National Reserve, Kenya. All samples were collected by the Mara Hyena Project, an ongoing field study of wild spotted hyenas, from individually identified members of a single social group, the Talek clan. All data collection protocols used here have been approved by Michigan State University (MSU) Institutional Animal Care and Use Committee (IACUC) and the Kenyan Wildlife Service (KWS). We selected samples for which T. gondii infection status was known and for which steroid hormone data were available. Sex was determined based on the morphology of the erect phallus (Frank et al., 1990), and our final data set included 51 female and 58 male hyenas. Age was estimated to ±7 days based on morphological and behavioral characteristics observed shortly after the animal is born (Holekamp and Smale, 1998a). Among females, there were 28 cubs, 8 subadults and 15 adults, and among males there were 20 cubs, 9 subadults, and 29 adults. These age categories correspond to age ranges in which cubs are ≤12 months, subadults are >12 and ≤ 24 months, and adults are >24 months on the date that they were darted.

Blood sample collection, hormone assays, and Toxoplasma gondii assays

Hyenas were anesthetized using pressurized dart containing 6.5 mg/kg tiletamine-zolazepam (Telazol ®) fired from a CO2-powered rifle (Telinject Inc.). Once safely immobilized, blood was collected from the hyena's jugular vein into vacuum tubes. Following the darting procedure, plasma from sodium heparin-coated vacuum tubes was separated from red and white blood cells by centrifugation. Samples were immediately stored in liquid nitrogen, transported to Michigan State University on dry ice, and permanently stored in our −80 °C biorepository.

Steroid hormone assays

We used standardized radioimmunoassays to quantify plasma hormone concentrations. Plasma samples were assayed in duplicate using a cortisol radioimmunoassay kit (MP Biomedicals CortiCote RIA kit, #06B256440) and a total testosterone radioimmunoassay kit (ImmuChemTM Double Antibody Testosterone 125 RIA kit, #07–189,105); further details of these plasma assays are published in Holekamp and Smale, 1998b and Dloniak et al. (2004).

Toxoplasma gondii infection status assays

We used ELISA based, multi-species ID Screen® Toxoplasmosis Indirect kit (IDVET, Montpellier) to diagnose T. gondii infection in plasma samples. This assay tests for IgG reactivity to T. gondii's P-30 antigen and has been previously validated in our study population (Gering et al., 2021).

Additional covariates

In addition to sex and age class, we used long-term demographic and ecological data to measure additional variables potentially associated with circulating steroid hormone concentrations. First, we created a two-level categorical variable indicating whether or not a blood sample was collected during the annual wildebeest and zebra migration, during which hyena diet and condition undergo marked changes (Cooper et al., 1999; Holekamp et al., 1997). For adult females, we categorized reproductive state on their darting date as nulliparous, pregnant, lactating, or other (Holekamp et al., 1996), and for adult males, we categorized dispersal status as either immigrant or resident based on whether or not that male was born into the Talek clan (Holekamp and Sisk, 2003). Each time a hyena was darted we noted the time of day, the duration of time between when an animal was darted and blood was drawn, and the individual's stress level as indicated by their behavior between darting and full anesthesia. Time of day was coded as a two-level variable indicating A.M. vs P.M. because hormones, and cortisol in particular, are known to fluctuate over the course of a day (Van Meter et al., 2009). Individual behavioral stress state was scored on an ordinal scale between 1 and 5, with 1 being the lowest behavioral stress level and 5 the highest.

Statistical methods

Prior to formal analyses, we filtered and formatted the hormone data. First, we identified all samples in which the hormone concentrations (i.e., testosterone or cortisol) were below the detection limit of our assays. For these samples, we replaced zero values with a hormone concentration equivalent to one-half the value of the lowest reported concentration. We evaluated the distributions of continuous hormone levels and used a natural log transformation of both testosterone and cortisol concentrations to conform to assumptions of normality for linear regression models. Next, we conducted bivariate analyses to assess associations between each potential precision covariate with testosterone and cortisol concentrations. Here, we define precision covariates as variables potentially associated with either the explanatory or the outcome variable of interest, but not both. Unlike adjustment for confounder variables, adjustment for precision covariates does not aim to make causal inference. Such associational studies minimize model assumptions (Laubach et al., 2021) and are well suited to data sets such as ours, in which serial measurements of hormone levels and infection status were not available from focal individuals.

The main goal of this study was to explore associations between T. gondii infection status and two steroid hormones that are both involved in fear-related behavior (Agis-Balboa et al., 2009; Heany et al., 2015; Korte, 2001; Soravia et al., 2006) and immunosuppression (Wyle and Kent, 1977) and known to change following T. gondii infection within human and lab animal hosts. This approach is useful as a first step to understand the crude relationship between two key variables of interest, and to pave the way for subsequent analyses of causal inference (Laubach et al., 2021). We modeled simple associations between T. gondii infection status and steroid hormone levels within males and females of our study population. In separate models for testosterone and cortisol, T. gondii infection status was included as a two-level categorical explanatory variable and the natural log of each steroid hormone concentration was the continuous outcome. Estimates from these models were based on percentile bootstrapping using the R package ‘car’ (Fox and Weisberg, 2019). When testosterone was the outcome, we ran 4 models within sex and age class strata, including cub and subadult females (model 1), adult females (model 2), cub and subadult males (model 3), and adult males (model 4). We collapsed the cub and subadult age groups in these models due to the similarity in the hormone profiles of cubs and subadults prior to sexual maturation (Supplemental Fig. 1) and based on a priori knowledge that sex steroid hormones change at reproductive maturity (Glickman et al., 1992). When cortisol was the outcome, we conducted separate regression models for males and females, but these models were not stratified by age based on similar hormone profiles across age classes (Supplemental Fig. 2). Because our focus in this study was on baseline cortisol levels (rather than stress response dynamics), we only analyzed cortisol levels in samples that were collected within 13 min of the time an animal was darted, after which circulating cortisol rises due to the stress of the darting process (Holekamp, unpublished data). We also only included cortisol measures collected from hyenas that showed low behavioral stress between darting and full anesthesia (levels 1 or 2). In our models, statistical significance was determined from two-tailed tests with an ⍺ = 0.05.

Finally, we conducted sensitivity analyses to determine if our results were robust given that some samples were below the detection limit of the testosterone or cortisol assay. First, we classified each continuous hormone level as a two-level, binary outcome corresponding to hormone concentrations below and above the median value within sex and age class strata. Here, we used logistic regression models in which infection status was the explanatory variable of interest and the low vs. high hormone level variable was the outcome. Second, we imputed 40 data sets to generate estimated hormone values when a testosterone or cortisol level fell below the assay detection limit (Succop et al., 2004). In this step, we used the R package ‘MICE’ running the predictive mean matching method with 50 iterations to generate 40 imputed data sets (van Buuren and Groothuis-Oudshoorn, 2011). We chose to impute 40 data sets based on the percentage of the data that were below the detection limit or missing (Graham et al., 2007), and imputations were calculated using available hormone, morphometry and sociodemographic data. We then averaged the imputed testosterone and cortisol levels across all 40 data sets and replaced the samples below the assay detection limits with these average imputed values. Like previous regression models, we assessed the association between infection status and continuous hormone levels in sex and age stratified models using the imputed data. Finally, all analyses were done in program R, version 4.0.2.

Results

The measured concentrations of cortisol and testosterone for all individuals in this study (N = 109 hyenas) are given in Table 2; 63 of these animals were positively diagnosed with T. gondii infections.

Table 2

Associations of hyena sex and age with plasma steroid hormone levels among spotted hyenas.

	Testosterone (ug/dL)			Cortisol (ug/dL)
	Ntotal	Nabove detection	Mean ± SDa	Ntotal	Nabove detection	Mean ± SDa
Female age classb
Cub (<12 mos)	28	8	0.06 ± 0.04	19	15	2.72 ± 3.07
Subadults (12–24 mos)	8	1	0.05 ± NA	6	3	6.77 ± 7.39
Adult (>24 mos)	15	9	0.31 ± 0.29	10	6	3.45 ± 4.88
Male age classb
Cub (<12 mos)	20	4	0.16 ± 0.14	8	8	2.35 ± 3.23
Subadults (12–24 mos)	9	4	0.11 ± 0.09	6	5	0.49 ± 0.19
Adult (>24 mos)	29	24	1.34 ± 1.57	19	13	3.39 ± 3.81

Estimated means and standard deviation among those samples where the steroid hormone level was >0, i.e. restricted to measures not below assay detection limit.

All hyena age classes are based on the date of the hormone concentration assessment (i.e., the darting date).

Next, we conducted bivariate analyses to determine if any precision covariates (i.e., the age classes of sampled individuals, their reproductive status [females]/dispersal status [males], sample collection time of day, or whether a sample was collected during the annual wildebeest and zebra migration) were associated with hormone levels. Age class was associated with circulating testosterone but not cortisol levels. Specifically, adult hyenas of both sexes had higher plasma testosterone levels than cubs and subadults (F2, 48 = 6.17, P = 0.004, females; F2, 55 = 17.36, P < 0.001, males), while cortisol did not covary with age class in males or females (Supplemental Tables 1 and 2). Reproductive state was associated with testosterone levels among adult females (F3, 11 = 8.40, P = 0.003, Supplemental Table 1). Testosterone concentrations varied with time of day among male hyenas (F1, 56 = 6.13, P = 0.016), though it should be noted that samples collected in the PM were very rare leading to an imbalanced design. None of the other precision covariates were associated with measured plasma testosterone or cortisol levels.

In our main models, we assessed simple associations of T. gondii infection status with testosterone and cortisol concentrations within strata of sex and age class. Here, we did not control for any precision covariates given the null associations from bivariate analyses and an imbalanced design for the time of day variable. Our primary finding was that positive infection status was associated with lower testosterone among female but not male cubs and subadults (−0.74 [95%CI: -1.41, −0.10] units ln-testosterone), and similarly, T. gondii infection was associated with lower testosterone in adult males (−2.16 [95%CI: -3.52, −0.57] units ln-testosterone, Fig. 1). We observed no association between infection status and cortisol concentrations in mixed age class groups of either female or male hyenas (Fig. 2).

Fig. 1

Geometric mean ratio and 95% CI estimates from separate sex and age stratified models of the relationship between T. gondii infection and plasma testosterone. The red dashed line represents the null, and estimates are based on percentile bootstrapping (2000 simulations).

Fig. 2

Geometric mean ratio and 95% CI estimates from separate sex stratified models of the relationship between T. gondii infection and plasma cortisol. The red dashed line represents the null, and estimates are based on percentile bootstrapping (2000 simulations).

In subsequent sensitivity analyses we modeled the effect of T. gondii on hormone levels using logistic regression, in which hormone values were binned into a binary outcome, and linear regression, in which hormone values that were below assay detection limits were estimated using multiple imputation. First, logistic regression models confirmed our results from the linear regression models revealing that a positive T. gondii infection was generally associated with lower odds of having high testosterone (Supplemental Table 3). The odds of having low vs high cortisol did not differ by infection status (Supplemental Table 4), a null result that reflects the lack of association from the linear regression models in which cortisol concentration was the outcome. Second, we conducted sensitivity analyses by running models where we replaced hormone levels below the assay detection limits with imputed values. The estimates of association between infection status and testosterone levels among female cubs and subadults, as well as male adult hyenas, were in same direction as estimates from previous linear regression models, although estimates from the imputed data were attenuated towards the null (Supplemental Table 5). The estimate for the imputed testosterone data among male cubs and subadults was not in the same direction as previous models, but in all cases the confidence intervals included the null. Similarly, there were no associations between infection status and cortisol in models using imputed data (Supplemental Table 6), which agreed with previous null findings.

Discussion

In this study of wild spotted hyenas, we detected multiple associations between T. gondii infection and steroid hormone concentrations. Age class and sex stratified analyses revealed that both female cubs and subadults as well as adult male hyenas infected with T. gondii had lower plasma testosterone concentrations than uninfected animals. We found no associations between T. gondii infection and plasma cortisol concentrations. Despite relatively small sample sizes and cross-sectional data, we report novel associations between T. gondii infection and steroid hormone levels in a wild animal population that adds to an existing literature in laboratory animals and human studies (c.f., Table 1).

Our observed negative relationship between T. gondii infection and plasma testosterone has previously been found in experimentally-infected adult male lab mice (Kaňková et al., 2011), in rats (Abdoli et al., 2012), and in human hosts who naturally acquired the parasite (Eslamirad et al., 2013). Kaňková et al. (2011) speculated that reducing testosterone levels may be adaptive for these infected hosts by ameliorating testosterone's immuno-suppressive effects and thus facilitating a more robust immunological response. However, experimentally-infected rats also exhibited temporary reductions in sperm quality and quantity in parallel with reduced circulating testosterone (Abdoli et al., 2012), indicating that a heightened immune response to T. gondii infection could incur reproductive trade-offs via reduced testosterone. At a more mechanistic level, T. gondii infection corresponds with elevated pro-inflammatory cytokines, including interleukin (IL)-12 in mice (Kaňková et al., 2010) and IL-1β in humans, the latter of which negatively correlates with testosterone levels (Oktenli et al., 2004). Recent work in mice has revealed localized brain lesions associated with T. gondii infection in the cortex, caudoputamen thalamus, and hypothalamus (Ihara et al., 2016), suggesting that T. gondii-induced inflammation within the brain, particularly the hypothalamus, might be a mechanism disrupting testosterone production.

Other studies have found patterns that depart from those observed within the present study. For example, several human studies have reported positive correlations between T. gondii infection and circulating levels of testosterone (Flegr et al., 2008; Shirbazou et al., 2011; Zghair et al., 2015; Zouei et al., 2018). These positive correlations are in line with the host-manipulation hypothesis because elevated testosterone can facilitate increased boldness and thus expose infected hosts to feline predators and the parasite's definitive host. Indeed, administration of testosterone targeting the medial amygdala of rats was found to reduce fearful behavior (Singh et al., 2020; Tong et al., 2019). These results, coupled with research showing higher densities of T. gondii cysts in parts of the brain involved in fear behaviors (Vyas et al., 2007), lend support to the idea of parasite-directed, host behavioral manipulation, although additional work is needed to distill the details of a plausible mechanism (McConkey et al., 2013; Vyas and Sapolsky, 2010). At this time, however, our finding of lower plasma testosterone in infected female cubs and subadults and in adult males appears more consistent with host (vs. parasite) control over circulating hormone levels. It may therefore be the case that other modifications to the hyena nervous system (e.g., dopaminergic modulation) underlie our previous finding of increased boldness in T. gondii infected individuals (Gering et al., 2021).

Glucocorticoids also act within regions of the brain associated with fear behaviors (McEwen et al., 2016), and positive correlations between T. gondii infection and circulating cortisol levels have been reported in humans (Abdelazeem et al., 2015; Shirbazou et al., 2011). However, we did not observe associations between T. gondii infection and plasma cortisol in either male or female hyenas. This may reflect a limitation of our sample size, the influence of unmeasured variables, or true biological differences between hyenas and other T. gondii host species.

Given that positive, negative, and null associations between T. gondii infection and host steroid hormone levels have been reported in controlled lab settings, it appears likely that additional factors may influence the nature of these relationships in the animal's natural habitat. These factors might include differences in the virulence of T. gondii strains, differences in the mode and/or dosage of infection, and underlying host characteristics (e.g., species and/or individual-level differences in genotypes or condition) (Abdoli, 2014). It is also important to note that the impacts of circulating hormones on an animal's physiology and behavior are contingent on many other variables, including receptor densities and other features (e.g., methylation levels) of hormone-responsive cells and circuits; these are already known to be affected by T. gondii infection (Vyas, 2015). Lastly, because this was an associational study, we cannot exclude the possibility that T. gondii infections and plasma hormone levels are jointly influenced by unmeasured variables within hyenas – though laboratory studies of other mammals clearly support the possibility that infections and hormone levels are linked.

Strengths and limitations of the present study

A strength of our study is that we have assessed associations between T. gondii infection and circulating steroid hormone concentrations in a free-living host, the spotted hyena. The associations we observed in a natural system reflect the physiology of infected animals living in an ecological context involving trade-offs (e.g., costs and benefits of hormone production) and selection pressures that do not apply in laboratory settings. To our knowledge, this is the first study of its kind in a wild carnivore. Similar studies of other wild hosts will prove illuminating, given that T. gondii's diverse wild hosts also exhibit diverse immunological, endocrine, and nervous systems. Initially, this future work might benefit from approaches like ours - in which we have measured and report associations from simple models. Investigating simple associations requires fewer assumptions than more complex models, increasing power to detect and characterize basic relationships between measured variables (Laubach et al., 2021). Nevertheless, it is important to account for basic demographic and life-history precision covariates in any given wild host. Otherwise, these variables could mask the relationships of interest, particularly in species like hyenas that exhibit behavioral and physiological sexual dimorphism and transition through life-history stages with distinctive hormone profiles.

Our study is not without limitations. For example, our results cannot establish the direction of causality between infections and hormone levels. Confirming causation and determining the joint roles of the host and parasite in hormonal modulation falls well beyond the scope of this work. This would be an important step toward a clearer determination of whether (and how) hormonal changes can impact the fitness of the host and/or parasite in natural settings. Another limitation of our data set is the relatively small samples size for some of the groups. In particular, there was only one uninfected adult female hyena for which we also had a testosterone measurement. Finally, while a proportion of testosterone and cortisol measures for hyenas fell below our assays’ detection limits, we conducted multiple sensitivity analyses to assess the robustness of our results. Results from the sensitivity analyses generally aligned with results from our main analysis. Future studies of spotted hyenas involving larger sample sizes, alternative assays, and/or measurements of additional variables may help untangle the biological mechanisms linking T. gondii infection with host physiology, and ultimately fitness. It is also worth noting that T. gondii is not the only infectious agent to which hyenas are exposed, and the effect of co-infection on host physiology warrants additional investigation.

Conclusion

Our study corroborates numerous lab studies supporting links between T. gondii infection and steroid hormones, suggesting patterns seen in the laboratory likely also occur in nature. More work is needed to investigate the relationships between T. gondii infection and hormones spanning a wider variety of wild species: if found to be common across wild species, relationships between infection and hormone levels may have far-reaching implications given the parasite's ubiquity in terrestrial ecosystems, and the strong connections between wild animals' hormone levels, their behavior, health, and fitness (Adkins-Regan, 2005).

Funding

This material is based in part upon work supported by the under Cooperative Agreement No. DBI-0939454. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation. Sample collection and behavioral data collection were supported by grants IOS1755089 and OISE1853934. ZML was funded by the grant D19ZO-411.

Author contributions

Conceptualization, ZML, EG; data generation and curation, ZML, EG, EY, TMM, KEH; analysis, ZML; funding and resources, ZML, EG, TG, KEH; writing original draft, ZML; review & editing the manuscript, ZML, EG, EY, TMM, TG, KEH.

Data and materials availability

All data and analysis code used in this paper are available at https://zenodo.org/record/5728323#.YaARob3MKDU in a Zenodo repository DOI: 10.5281/zenodo.5728323.

Conflicts of interest

The authors have declared that no conflicts of interests exist.

Declaration of competing interest

The authors have declared that no competing interests exist.