Global variation in diurnal cortisol rhythms: evidence from Garisakang forager-horticulturalists of lowland Papua New Guinea

Abstract

The hypothalamic–pituitary–adrenal (HPA) axis represents an important and evolutionarily ancient biological pathway linking physical and psychological stressors with human health. Despite considerable research exploring the physiological stress response among developed populations, few studies have examined HPA activity in non-industrialized contexts, restricting understanding of variation in human stress reactivity across global socio-ecological diversity. The present study addresses this shortcoming by investigating diurnal cortisol rhythms among Garisakang forager-horticulturalists of remote, lowland Papua New Guinea. Using a large sample of repeated salivary cortisol measurements from 169 participants (age 4–70 years), multilevel growth curve models were constructed to assess Garisakang waking cortisol concentrations and diurnal cortisol slopes. As predicted, results demonstrate identifiable but substantially diminished diurnal cortisol rhythms relative to those of industrialized populations. Sample-wide, Garisakang cortisol concentrations are highest upon waking (mean = 4.86 nmol/L) and decrease throughout the day at a mean rate of only −0.18 nmol/L/h or −6.20%/h. Age and sex significantly predict evaluated cortisol parameters in ways not consistently reported among industrialized populations, suggesting that Garisakang diurnal cortisol rhythms are defined by distinct ontogenetic trajectories across the lifespan. These findings highlight cross-cultural diversity in HPA activity and have important implications for understanding basic mechanisms of the physiological stress response in contexts of chronic physical stressors such as limited nutrition, heavy burden of infectious disease, and high levels of physical activity.

Keywords:

• Stress
• HPA axis
• waking cortisol
• diurnal cortisol slope
• lifespan
• non-industrialized society

Introduction

The hypothalamic–pituitary–adrenal (HPA) axis and its primary end product, the steroid hormone cortisol, are critical components of the human physiological stress response (Dickerson & Kemeny, 2004; Sapolsky, Romero, & Munck, 2000). Triggered for release by a range of physical and psychosocial stimuli, circulating cortisol acts to maintain internal homeostasis through regulatory effects on a number of biological systems, including metabolism, reproduction, immunity, and growth (Gunnar & Quevedo, 2007; Sapolsky et al., 2000). As such, cortisol plays an important and evolutionarily ancient physiological role linking experienced stressors with human phenotypic variation and health (Dickerson & Kemeny, 2004; Weiner, 1992).

A growing body of research among industrialized populations suggests that typical, healthy human cortisol secretion follows a strong diurnal pattern (Adam & Kumari, 2009). Under basal conditions, circulating cortisol concentration is relatively high at waking, rises sharply to a peak ≍ 30–45 min post-waking, and subsequently declines over the course of the day to a nadir just prior to sleep (Gröschl, Rauh, & Dörr, 2003). Deviations from this pattern – due to, for example abnormal cortisol concentration at waking or reduced rate of diurnal cortisol decline (i.e. flattening of the diurnal cortisol slope) – are often associated with negative health outcomes, including greater incidence of depression, cardiovascular disease, and childhood growth faltering (Adam, 2006; Heim, Ehlert, & Hellhammer, 2000; Miller, Chen, & Zhou, 2007; Rosmond, 2005). On account of these relationships, the documentation of diurnal cortisol rhythms, including their ontogenetic trajectories, has become an important component of human stress research (Adam & Kumari, 2009; Heim et al., 2000; Miller et al., 2007).

Problematically, the overwhelming majority of studies investigating diurnal cortisol rhythms have been performed among industrialized populations. This shortcoming restricts understanding of variation in HPA activity across the full range of human environmental, cultural, and genetic diversity, particularly within non-industrialized, physically stressful contexts such as those experienced by most living people today and likely all humans in our recent evolutionary past.

The present study addresses this limitation by investigating diurnal cortisol rhythms across the lifespan among Garisakang forager-horticulturalists of lowland Papua New Guinea (PNG). The Garisakang practice a subsistence-based lifestyle that is characterized by chronic stressors such as limited nutrition, heavy burden of infectious disease, and high levels of physical activity. As such, research among this population provides a relatively unique opportunity to explore the impact of a non-industrialized lifestyle on HPA function and ontogeny. Using a study design involving repeated measures of waking and evening salivary cortisol, we test two hypotheses: First, that similar to many groups experiencing chronic psychosocial stress (Heim et al., 2000; Miller et al., 2007), the Garisakang, facing chronic physical stress, will have diminished diurnal cortisol rhythms (i.e. lower waking cortisol concentrations and flattened diurnal cortisol slopes) relative to healthy Western references; Second, that Garisakang diurnal cortisol rhythms – as preliminarily observed among other forager-horticulturalists (Nyberg, 2012) – will exhibit distinct relationships with age and sex, conforming to some but not all ontogenetic patterns of HPA activity found among industrialized populations with dissimilar conditions of lifetime stress (Jessop & Turner-Cobb, 2008; Saxbe, 2008).

Methods

Study population

Detailed background information for the Garisakang is provided elsewhere (Konečná & Urlacher, 2017). Briefly, the Garisakang are a subsistence-based indigenous group of ≍ 500 individuals inhabiting the remote Middle Ramu River Valley of Madang Province, PNG. Similar to many other Oceanic populations (Sahlins, 1963), Garisakang social structure is patrilineal with mixed-polygyny and established big men exercising dominant political influence. Since 1995, a large portion of the Garisakang has resided in a single community – Wanang – with a current population of ≍ 200 individuals. The Wanang region is characterized by dense evergreen tropical forest, mean daytime temperature of 26 °C, and heavy rainfall totaling ≍ 3500 mm annually (Anderson‐Teixeira et al., 2015). Government infrastructure and access to markets or modern healthcare is limited; Wanang possesses no electricity, running water, clinic, or reliable road access. To obtain basic services, community members must walk or attempt to hitchhike on a seasonal dirt logging road to Madang town located ≍ 80 km away. This combination of isolation and tropical conditions results in high rates of infectious and parasitic disease, including holoendemic malaria, as is typical of lowland PNG (Müller, Bockarie, Alpers, & Smith, 2003; Shaw, 1984). The Garisakang remain reliant on slash-and-burn horticulture, foraging, bow hunting, and fishing to supply their diet. Consumption of market foods to supplement overall nutritional constraint is uncommon, with residents in 2015 reporting eating, on average, fewer than two market items (e.g. rice, noodles, biscuits) per month. Engagement in a mixed forager-horticulturalist lifestyle fosters high levels of habitual physical activity, particularly among women who are responsible for the bulk of household domestic and gardening activities.

Data collection

All study instructions were delivered to participants in Tok Pisin, with translation into local Maghu by a field assistant when necessary. Participants provided informed consent or, in the case of young children, assent with parental consent. Research approval was obtained from village and household leaders, the Ethics Committee of the University of South Bohemia, and the Institutional Review Board of Harvard University.

Data were collected from a total of 169 Garisakang participants (aged 4.3–70.1 years) living in Wanang between 2013 and 2015. All data were collected during the annual climatic dry season from late May to early October. Women who were visibly pregnant or reported being pregnant at the time of data collection were excluded from participation. All participants reported having no known medical conditions, no current use of steroid-based medication, and no regular consumption of alcohol or caffeine. Exact estimates of age were available for 57 individuals from government-issued identification cards. Age was also estimated and cross-checked using extensive overlapping genealogies constructed from information provided by parents and other community members. In general, age estimates are considered accurate to the month for subadults and to the year for adults.

Anthropometric data facilitating description of nutritional status were collected using conventional methods (Lohman, Roche, & Martorell, 1988). Height was measured to the nearest 1.0 mm using a portable stadiometer (Seca Corporation 214, Hanover, MD). Weight was measured to the nearest 0.1 kg using an electronic scale (Tanita Corporation BF-689, Tokyo, Japan). Body mass index (BMI, kg/m²) was calculated from height and weight data.

Saliva collection

Diurnal saliva samples (n = 1654) were collected using salivette collection devices (Sarstedt, Nümbrecht, Germany), as has previously been described (Konečná & Urlacher, 2015, 2017). Sampling was devised to limit participant burden and maximize compliance. In 2013, heads of household were given salivettes the evening prior to sample collection and, following training, were directed to facilitate their use with each participating family member immediately upon waking the following morning. Heads of household were also asked to log the time of sample collection for each participant using a provided digital wristwatch. The authors retrieved the used salivettes and wristwatches from participating households early in the morning on each day of sample collection. This protocol was repeated for up to eight consecutive days per household (mean = 4.0 d, SD = 1.3 d), providing repeated waking saliva samples from each participant. In 2014 and 2015, saliva samples were collected and times were recorded directly by the lead author, with the author arriving outside households in the early morning (prior to participant waking) and in the evening to ensure compliance and accurate recording of time. Participants during these years provided saliva samples immediately upon waking and again the same evening for up to three consecutive days (2014 mean = 2.6 d, SD = 0.6 d; 2015 mean = 2.5 d, SD = 0.7 d), affording a maximum of three waking and evening saliva samples per individual. Time of collection did not significantly differ by year for waking (sample mean = 6:40; F[2, 1072] = 2.99, p = .051) or evening (sample mean = 17:56; F[1, 577] = 3.24, p = .072) samples.

Following standard protocol, waking saliva samples were not collected from any individual directly observed or reporting to have been awake for more than 5 min prior to sample collection (Stalder et al., 2016). Participants were reminded not to smoke tobacco or chew betel nut (Areca catechu) for at least 2 h prior to saliva sample collection. Samples collected from individuals reporting use of these stimulants in the preceding 2 h were discarded. Upon collection, salivettes were immediately centrifuged for 2 min at 1300 rpm using a hand-powered centrifuge (Hettich C1011, Salford, England). Extracted saliva was transferred to polystyrene vials, treated with 5% sodium azide solution preservative (Lipson & Ellison, 1989) and stored at ambient temperature until shipped to Prague for analysis and long-term storage at −20 °C.

Cortisol analysis

Saliva samples were analyzed for free cortisol concentration in the Clinical Laboratory of the Czech Institute of Endocrinology. All samples were measured in duplicate using widely available Immunotech radioimmunoassay kits (Beckman Coulter, Prague, Czech Republic) and an industrial 12-channel gamma counter (Berthold, Bad Wildbad, Germany). Manufacturer and in-house controls were included in each assay to verify expected concentration ranges. Analytical criteria agreed well with those reported by the manufacturer, with overall intra-assay CV = 8.9%, inter-assay CV = 9.2%, and analyte recovery of 82.5–107.0%.

Statistical analysis

Prior to statistical analysis, data for all continuous variables were assessed for outliers (i.e. values corresponding to z-scores > 3 or < −3) and distributional assumptions. Cortisol concentrations were positively skewed and were, therefore, log₁₀-transformed for analysis. Descriptive statistics were performed to examine sample characteristics. In some descriptive analyses, the sample was stratified by age categories for subadults (age 4.3–15.9 years, N = 83) and adults (age 16.0–70.1 years, N =  86). For comparative purposes, subadult age- and sex-specific z-scores for height-for-age (HAZ), weight-for-age (WAZ), and BMI-for-age (BAZ) were calculated using World Health Organization growth references (De Onis, 2007).

Following the analytical strategy of others examining diurnal cortisol rhythms (e.g. Adam, 2006; Hruschka, Kohrt, & Worthman, 2005), waking cortisol samples were assigned a value of +0 h for the variable Time post-waking. Evening cortisol samples were assigned a value of Time post-waking corresponding to sample-specific decimal hours since waking at time of collection (mean = 11.2 h, SD = 1.2 h). For basic descriptive analyses and population comparisons, cortisol diurnal slope was calculated as the difference in cortisol concentration between matched waking and evening saliva samples divided by Time post-waking of the evening sample.

A series of multilevel growth curve models with log₁₀-transformed cortisol concentration as the outcome variable and Time post-waking as a random factor were constructed to investigate Garisakang waking cortisol concentration (i.e. the model Intercept estimate) and cortisol diurnal slope (i.e. the model Time post-waking estimate) while accounting for variation in sample time of collection. Multilevel modeling is considered the gold standard for modeling cortisol diurnal rhythms (Adam & Kumari, 2009). Among other attributes, multilevel modeling accommodates unbalanced and incomplete datasets, can account for repeated sampling structures, and permits the simultaneous modeling of multiple diurnal cortisol parameters (Hruschka et al., 2005; Shirtcliff et al., 2012). Two-level and three-level models apportioning variance within individuals (Level 1), between individuals (Level 2) and between households or, alternatively, between years of data collection (Level 3) were assessed. To determine sample-wide mean diurnal cortisol rhythms and the partitioning of total variance across levels, a baseline, partially unconditional model was first constructed with Time post-waking as a single random predictor (at Level 1). Level-2 predictors (i.e. Age, Sex, and their interactions as well as the control variable BMI) were subsequently added to this baseline model in a stepwise fashion to test hypotheses relating to the impact of person-level variables on cortisol parameters. No predictors were included at Level 3 in any model. Akaike information criterion (AIC) and likelihood ratio tests with Satterthwaite approximations of degrees of freedom were used to assess model fit. All analyses were performed in R version 3.3.2 (http://cran.us.r-project.org/) using the lme4 package (Bates, Maechler, Bolker, & Walker, 2015) with the lmeTest extension. Results were considered statistically significant at p < .05.

Results

Descriptive analyses

Descriptive sample characteristics for subadults and adults are presented in Table 1. Relative to WHO multinational references, participants were, on average, short in stature, light in weight, and exhibited low levels of BMI. In total, 54% of subadults were classified as stunted (i.e. HAZ < 2), 23% as underweight (i.e. WAZ < 2), and 8% as wasted (i.e. BAZ < 2). For the entire sample, mean waking cortisol concentration was 4.86 nmol/L (SD = 1.42 nmol/L) with an unadjusted mean diurnal cortisol slope of −0.18 nmol/L/h (SD = 0.14 nmol/L/h). Waking cortisol concentration did not significantly differ by year of data collection (β = −0.01, SE = 0.01, p = .147). In total, 91.1% of participants had mean diurnal cortisol slopes of less than zero, reflecting an overall pattern of decreasing cortisol levels over the course of the day. Garisakang diurnal cortisol rhythms across the life course are illustrated by sex in Figure 1. Mean evening cortisol concentrations are not presented due to sample-specific variation in Time post-waking that renders unadjusted evening cortisol values largely uninformative.

Figure 1. Garisakang unadjusted mean (shaded interquartile range) waking cortisol concentration (top) and diurnal cortisol slope (bottom) across the lifespan, stratified by sex.

Table 1. Sample characteristics (mean, SD) for subadults and adults by sex.

	Subadults (age 4.3–15.9 years)		Adults (age 16.0–70.1 years)
	Females (N = 32)	Males (N = 51)	Females (N = 39)	Males (N = 47)
Age (years)	10.9 (3.6)	10.3 (3.3)	34.6 (14.6)	29.9 (14.4)
Height (cm)	131.8 (16.8)	127.4 (17.3)	152.7 (3.9)	161.6 (5.6)
Height-for-age^a	−2.0 (1.0)	−2.1 (0.8)	–	–
Weight (kg)	28.0 (10.3)	26.0 (9.1)	45.2 (6.2)	54.5 (6.6)
Weight-for-age^a	–0.6 (1.5)	−0.8 (1.4)	–	–
BMI (kg/m^2)	16.2 (1.9)	15.7 (1.5)	19.4 (2.7)	20.9 (1.9)
BMI-for-age^a	−0.9 (0.7)	−0.9 (0.8)	–	–
Saliva samples (#/participant)	8.9 (4.1)	9.9 (4.0)	10.6 (3.9)	9.6 (4.6)
Waking saliva collection time (hh:mm)	6:37 (0:30)	6:42 (0:30)	6:35 (0:29)	6:42 (0:37)
Evening saliva collection time (hh:mm)	17:48 (1:50)	17:59 (1:17)	17:55 (1:26)	17:58 (1:35)
Waking cortisol (nmol/L)	5.16 (1.42)	4.80 (1.39)	4.97 (1.47)	4.72 (1.42)
Cortisol diurnal slope (nmol/L/h)^b	−0.19 (0.14)	−0.20 (0.12)	−0.15 (0.17)	−0.16 (0.11)

^aSubadult age and sex-specific z-scores calculated using World Health Organization multinational growth references (De Onis, 2007). ^bCortisol diurnal slope calculated as slope = (evening cortisol concentration − waking cortisol concentration)/hours between waking and evening sample collection.

Multilevel modeling of diurnal cortisol rhythms

Initial three-level multilevel models revealed that 79.6% of total variation in cortisol in the sample was accounted for by within-individual factors (Level 1), 16.0% of variance by between-individual factors (Level 2) and 4.4% of variance by between household factors (Level 3). It was accordingly determined that a baseline three-level model including household did not significantly improve the fit of a simpler two-level model (χ² = 2.77, p = .096). Similar results were found for a three-level model including year of data collection at Level 3 (χ² = 2.10, p = .237). Therefore, the baseline two-level model was selected for all subsequent analyses. The inclusion of Time post-waking as a random effect into this baseline model demonstrated that 40.1% of total sample variation in cortisol concentration was attributable to differences in elapsed time since waking. Model fit was not improved by including an additional term for quadratic Time post-waking (χ² = 2.48, p = .115).

Parameter estimates for final multilevel models accounting for within- and between-individual variation in Garisakang log₁₀-transformed cortisol concentration are provided in Table 2. The baseline model (Model 1) included significant Intercept (β = 0.687, SE = 0.007, p < .001) and Time post-waking (β = −0.028, SE = 0.001, p < .001) estimates. Transformed to standard cortisol values (see Adam, 2006), these estimates indicate a sample-wide mean waking cortisol concentration of 4.86 nmol/L and a mean diurnal cortisol decline of 6.20% /h awake (Figure 2). The inclusion of person-level predictors in Model 2 identified significant effects of participant Age and Sex on diurnal cortisol rhythms. Sex was significantly related to waking log cortisol concentration (γ = −0.030, SE = 0.014, p = .039), such that being male was associated with a reduction of 0.32 nmol/L in waking cortisol. Similarly, Age was a significant predictor of cortisol diurnal slope (γ = −0.001, SE = 0.000, p = .042), indicating an increase (i.e. flattening) of 0.03% in cortisol rate of decline for every single year increase in age across the sample. Age*Sex interaction terms and the control variable BMI were not significant in any preliminary model (all p > .1) and were, therefore, removed from Model 2.

Figure 2. Garisakang all-age estimated mean cortisol concentration (shaded 95% CI) by time post-waking, stratified by sex.

Table 2. Multilevel models of Garisakang log₁₀-transformed diurnal cortisol rhythms (nmol/L).

	Model 1		Model 2
Fixed effect	Estimate (SE)	p Value	Estimate (SE)	p Value	Interpretation^a
Intercept	0.687 (0.007)	<.001	0.705 (0.015)	<.001	Mean waking cortisol of 4.86 nmol/L
Age (years)	–	–	−0.000 (0.000)	.970	N.S.
Sex (male = 1)	–	–	−0.030 (0.014)	.039	−0.32 nmol/L if male
Time post-waking	−0.028 (0.001)	<.001	−0.031 (0.002)	<.001	Mean cortisol decline of 6.20%/h
Age (years)	–	–	0.001 (0.000)	.042	+0.03%/h per year older
Sex (male = 1)	–	–	0.001 (0.002)	.645	N.S.
Random effect	Variance (SD)	p Value	Variance (SD)	p Value
Level 1 intercept	0.023 (0.151)	<.001	0.023 (0.151)	<.001	–
Level 2 intercept	0.004 (0.064)	<.001	0.004 (0.063)	<.001	–
Time post-waking	0.001 (0.010)	–	0.001 (0.010)	–	–

^a Estimates for log cortisol concentrations were transformed to standard values to enhance interpretation.

Discussion

This study is among the first to investigate diurnal cortisol rhythms across the lifespan in a non-industrialized population. Results provide insight into global diversity in the human physiological stress response.

Population-level comparisons of Garisakang diurnal cortisol rhythms

In agreement with research among industrialized populations (Adam & Kumari, 2009; Gröschl et al., 2003), Garisakang forager-horticulturalists exhibit clearly identifiable diurnal cortisol rhythms. This pattern is evident for 91% of participants, a proportion similar to that observed among other populations at various stages of economic development (Nyberg, 2012; Stone et al., 2001).

As predicted in a context of chronic socio-ecological stress (Hypothesis 1), Garisakang mean waking cortisol concentration is considerably lower and mean diurnal cortisol slope considerably flatter than reported for industrialized references (Table 3; Cohen et al., 2006; DeCaro & Worthman, 2008; Stawski, Cichy, Piazza, & Almeida, 2013). Garisakang diurnal cortisol rhythms further appear diminished relative to those of other developing populations for whom data are available (Desantis, Kuzawa, & Adam, 2015; Flinn & England, 1997; Nyberg, 2012) and, in fact, represent the lowest diurnal cortisol values on record for any population-based sample (Table 3). Indeed, the only population with reported diurnal cortisol values approaching those of the Garisakang – with waking cortisol concentrations that are 17% higher and diurnal cortisol slopes that are 21% steeper – are the Tsimane’, an Amazonian forager-horticulturalist group at a roughly similar stage of early economic development (Nyberg, 2012). It is important to note that direct comparison of diurnal cortisol values across these studies must be interpreted with caution due to slight methodological and laboratory differences. However, the consistency of suppressed HPA activity among developing populations evident in Table 3 strengthens the conclusion that Garisakang diurnal cortisol rhythms reflect the very low end of known human variation and are substantially diminished relative to those of industrialized references.

Table 3. Comparative studies reporting salivary diurnal cortisol rhythms in diverse populations. The findings of the present study are provided in bold.

Population	Age (years)	N	Samples/participant	Waking cortisol(nmol/L)	Cortisol slope(nmol/L/h)	References
Garisakang	4–15	83	9.5	4.9	−0.19	Urlacher et al. (this study)
Garisakang	16–70	86	10.0	4.8	−0.16	Urlacher et al. (this study)
Tsimane’	2–15	175	5.1	5.3	−0.21	Nyberg (2012)
Tsimane’	16–82	128	5.3	6.4	−0.26	Nyberg (2012)
Filipino	20–22	1490	2.0	7.5	−0.33	Desantis et al. (2015)
Caribbean	0–18	264	87.6	15.9	–	Flinn and England (1997)
USA (Atlanta)	5–6	28	10–14	10.4	−0.51	Decaro and Worthman (2008)
US CARDIA	18–30	781	1.9	20.3	−0.78	Cohen et al. (2006)
US MIDUS	33–84	1694	8.0	15.2	−0.73	Stawski et al. (2013)

As predicted by Hypothesis 2, Garisakang diurnal cortisol rhythms demonstrate distinct relationships with age and sex, conforming to some but not all ontogenetic patterns of HPA activity found among industrialized populations. Consistent with research in the United States and Europe (Darnall & Suarez, 2009), females exhibit significantly higher waking cortisol concentrations than males. This pattern is driven primarily by elevated female waking cortisol between puberty and mid-adulthood, presumably the result of adolescent changes in gonadotropin and androgen profiles that increase female HPA sensitivity (Romeo, 2010). However, despite elevation in waking cortisol concentration – and in opposition to findings among industrialized populations (Shirtcliff et al., 2012) – Garisakang pubertal females have substantially flatter, not steeper, diurnal cortisol slopes than males (γ = 0.071, SE = 0.010, p < .001). This pattern is the result of significantly elevated evening cortisol concentrations among females relative to males at puberty (γ = 0.985, SE = 0.049, p < .001), a phenomenon not reported among industrialized populations.

Distinct patterns of Garisakang HPA activity are also found late in life. Although, sample-wide, increasing age is related to modest flattening of the diurnal cortisol slope, the oldest Garisakang (those age 55 years or older) demonstrate steepening of the diurnal cortisol slope in association with elevated waking cortisol levels. This pattern contrasts the majority of findings among industrialized populations (Sapolsky et al., 2000), as well as among the Tsimane’ (Nyberg, 2012) which suggest that old age is typically associated with flattening diurnal cortisol slopes due to normal processes of biological senescence.

Socio-ecological stressors and determinants of Garisakang HPA activity

Socio-ecological factors – chiefly, chronic physical stressors – may largely explain diminished Garisakang diurnal cortisol rhythms. Like many non-industrialized populations, the Garisakang experience persistent physical stressors, including limited nutrition, recurrent immune activation, and physically active lifestyles. Although acute forms of these stressors are known to elicit short-term increases in HPA activity (Labsy et al., 2013; Misra et al., 2004; Rock et al., 1992; Tomiyama et al., 2010), it has been suggested that chronic exposure to these same factors may, in many cases, dampen the physiological stress response (Nyberg, 2012; Sapolsky et al., 2000). This hypothesis is supported by the observation that, among industrialized populations, chronic psychological stress is often associated with hypocortisolism and flattening of the diurnal cortisol slope (Heim et al., 2000; Miller et al., 2007; Saxbe, 2008). The biological mechanisms thought to underlie this relationship – including interference at various levels of the HPA axis (reviewed in Heim et al., 2000) – are believed to function similarly to suppress HPA activity in response to chronic physical stress (Sapolsky et al., 2000). This pathway, although requiring further testing, provides an attractive framework for understanding the diminished diurnal cortisol rhythms of the Garisakang and other non-industrialized populations globally.

As noted by others (Miller et al., 2007; Nyberg, 2012; Sapolsky et al., 2000), the suppression of HPA activity in response to chronic stress may reflect an adaptive mechanism operating to maintain long-term physiological homeostasis and avoid the deleterious effects of persistent hypercortisolism (e.g. hypertension, immunosuppression). One of the primary actions of cortisol is the initiation of catabolism to increase blood glucose (Brillon, Zheng, Campbell, & Matthews, 1995) and plasticity in diurnal cortisol rhythms should, thus, be particularly sensitive to nutritional constraint. As illustrated by widespread stunting, underweight, and low BMI, the Garisakang face chronic energy limitation. Diminished diurnal cortisol activity in this context may serve as an important mechanism preventing somatic depletion and starvation (Nyberg, 2012). Similar to findings among the Tsimane’ (Nyberg, 2012), however, we detect no evidence for a relationship between Garisakang BMI (i.e. a proxy for energy availability) and diurnal cortisol rhythms. This finding may reflect a true lack of association or, alternatively, may result from universally low BMI in the present sample. Future research in the developing world should carefully investigate the relationship between chronic nutritional stress and diurnal cortisol rhythms, particularly during energetically challenging life stages such as childhood.

Chronic immunological stress may also play an important role underlying Garisakang diminished HPA activity. Lowland PNG is characterized by high rates of infectious disease (Müller et al., 2003; Shaw, 1984), leading to persistently elevated levels of immune activity (Shinoda et al., 2012). These immune responses are energetically costly (Ulijaszek, 2000) and may suppress Garisakang diurnal cortisol rhythms indirectly via the energetic pathway described above. Importantly, persistent immune activation may also directly modify Garisakang HPA activity. Chronic inflammatory immune activity, for example, is known to dramatically down-regulate cortisol production and bioavailability over time (Silverman & Sternberg, 2012). For the Garisakang and other populations experiencing heavy burden of infectious disease, such an exposure-based suppressive effect may diminish diurnal cortisol rhythms and, ultimately, permit robust immune responses – critical for survival – while avoiding physiologically dangerous hypercortisolism. Direct testing of this hypothesis will require sensitive measures of immune function across the lifespan.

A third physical stressor that may suppress Garisakang diurnal cortisol rhythms is intensive habitual physical activity. Very little is known about the impact of habitual physical activity levels on HPA function (Rimmele et al., 2009). However, among non-industrialized populations, physical activity is costly (Dufour & Piperata, 2008) and may be predicted to suppress diurnal cortisol rhythms via the aforementioned energetic pathway. Specific patterns of division of labor may also contribute to the unique ontogeny of Garisakang HPA activity. For example, flat diurnal cortisol slopes and elevated evening cortisol concentrations among pubertal females may result from engaging in demanding chores (e.g. gardening work) that are not required of pubertal males and lead to sustained elevation of cortisol levels late in the day. This hypothesis is preliminarily supported by research demonstrating that daytime bouts of physical activity are associated with increased evening (Hansen, Blangsted, Hansen, S⊘gaard, & Sj⊘gaard, 2010; Kertes & Gunnar, 2004) but not waking (Martikainen et al., 2013) cortisol concentrations.

Although chronic physical stressors are likely the dominant factors underlying distinct patterns of Garisakang HPA activity, psychosocial stressors may also play a role. Economic development and acculturation are typically accompanied by psychological stress (Dressler, 1999; McDade, 2000). Therefore, it is possible that low levels of HPA activity among the Garisakang result from a relatively stress-free “traditional” lifestyle. This model drastically oversimplifies Garisakang way of life, however, and seems unlikely to account for the large differences in overall HPA activity from industrialized populations. One more plausible scenario is that psychosocial factors contribute to the relatively steep diurnal cortisol slopes found among elderly Garisakang; possibly due to low levels of psychological stress given a social system in which status and dominance are positively related to age (Konečná & Urlacher, 2017). Additional social factors that dynamically affect the Garisakang (e.g. alloparental care, pubertal rites, expectations for resource sharing accompanying clan-based social structure) may also contribute to a distinct pattern of diurnal cortisol rhythms. Analyses examining relationships between a range of socio-cultural variables, self-perceived psychosocial stress levels, and diurnal cortisol rhythms are currently underway to investigate such pathways.

Study limitations

This study possesses several limitations. First, given logistic difficulties of fieldwork among the Garisakang, it is limited to the collection of two saliva samples per individual per day. Although diurnal sampling to any degree is exceedingly rare among forager-horticulturalists and total sample size is relatively large, the collection of additional samples throughout the day would permit more robust modeling of diurnal cortisol rhythms, including assessment of the cortisol awakening response (Hoyt, Ehrlich, Cham, & Adam, 2016). Second, analyses did not control for incidence of periodontal disease, a condition known to increase salivary cortisol levels (Rosania, Low, McCormick, & Rosania, 2009) that may have contributed to within-population variation. Third, as previously discussed, caution must be used when making cross-population comparisons of Garisakang diurnal cortisol rhythms given some methodological differences between studies. Finally, while presenting several novel hypotheses, this study did not directly test relationships between specific stressors and diurnal cortisol rhythms, limiting the interpretation of results. Such research is needed to understand the potentially adaptive or pathological nature of global variation in the human physiological stress response.

Conclusions

This study capitalized on a rare dataset and robust multilevel modeling to investigate HPA activity among a forager-horticulturalist population. Results indicate that, relative to industrialized populations, the Garisakang have diminished diurnal cortisol rhythms and a distinct ontogeny of HPA activity. These findings highlight global variation in the human physiological stress response. Future research among the Garisakang and other non-industrialized populations should directly examine chronic stressors – including limited nutrition, heavy burden of infectious disease, and high levels of physical activity – that may sensitively calibrate HPA activity across the lifespan.

Acknowledgments

The authors thank the Garisakang for their continued research participation. We also express our gratitude to the Czech Science Foundation, the Binatang Research Center, and Nigel Baro.

Disclosure statement

The authors report no conflicts of interest.

Additional information

Funding

This work was supported by the Czech Science Foundation under grant #13-25602P.