SCIENCE ADVANCES | RESEARCH ARTICLE ANTHROPOLOGY Copyright © 2020 The Authors, some The life history of human foraging: Cross-cultural rights reserved; exclusive licensee and individual variation American Association for the Advancement of Science. No claim to Jeremy Koster1,2*, Richard McElreath2,3, Kim Hill4, Douglas Yu5,6,7, Glenn Shepard Jr.8, original U.S. Government Nathalie van Vliet9, Michael Gurven10, Benjamin Trumble4,11, Rebecca Bliege Bird12, Douglas Bird12, Works. Distributed Brian Codding13, Lauren Coad9,14, Luis Pacheco-Cobos15, Bruce Winterhalder3, Karen Lupo16, under a Creative Dave Schmitt16, Paul Sillitoe17, Margaret Franzen18, Michael Alvard19, Vivek Venkataraman20, Commons Attribution Thomas Kraft10, Kirk Endicott21, Stephen Beckerman13, Stuart A. Marks22,23, Thomas Headland24, NonCommercial License 4.0 (CC BY-NC). Margaretha Pangau-Adam25,26, Anders Siren27, Karen Kramer13, Russell Greaves13, Victoria Reyes-García28,29, Maximilien Guèze29, Romain Duda29,30, Álvaro Fernández-Llamazares31, Sandrine Gallois32, Lucentezza Napitupulu29, Roy Ellen33, John Ziker34, Martin R. Nielsen35, Elspeth Ready2,36, Christopher Healey37, Cody Ross2 Human adaptation depends on the integration of slow life history, complex production skills, and extensive sociality. Refining and testing models of the evolution of human life history and cultural learning benefit from increasingly accurate measurement of knowledge, skills, and rates of production with age. We pursue this goal by inferring hunters’ increases and declines of skill from approximately 23,000 hunting records generated by more than 1800 indi- viduals at 40 locations. The data reveal an average age of peak productivity between 30 and 35 years of age, although high skill is maintained throughout much of adulthood. In addition, there is substantial variation both among individuals and sites. Within study sites, variation among individuals depends more on heterogeneity in rates of decline than in rates of increase. This analysis sharpens questions about the coevolution of human life history and cultural adaptation. INTRODUCTION Progress is made in these debates via models that integrate growth, Among hominoids, humans are distinguished by a suite of life his- reproduction, cognitive development, skill development, sociality, tory traits that includes a prolonged juvenile and adolescent period, and cultural evolution. The most advanced attempt that we know is short interbirth intervals, and an extended postreproductive life span the optimal control model of González-Forero and Gardner (6). (1). Multiple conceptual models have been advanced to explain the Drawing on observed rates of brain and somatic growth, their model evolution of these traits, focusing on distinctive human behaviors estimates the relative importance of different ecological and social such as pair bonding and alloparental care from grandparents and challenges to the evolution of intelligence. In comparison to alter- others (2–4). Any satisfactory model of human life history evolution natives, the ecological challenge of acquiring food emerges as the must simultaneously account for the large brains that characterize strongest predictor of the observed pattern of human growth. In this our species. Foraging complexity and competitive social challenges model, brains develop first followed by the body because this se- have been alternately championed as the evolutionary prime mover quence allows a longer period of learning and ultimately higher adult of encephalization, while others combine perspectives by citing the productivity. This finding complements recent comparative work advantages of flexible cultural learning among juveniles as a solution on the prominence of foraging complexity as a predictor of primate to social and ecological challenges (5). brain sizes (7). These findings and predictions direct our attention 1 Department of Anthropology, University of Cincinnati, Cincinnati, OH, USA. 2Department of Human Behavior, Evolution, and Culture, Max Planck Institute for Evolutionary Anthropology, Leipzig, Germany. 3Department of Anthropology and Graduate Group in Ecology, University of California, Davis, Davis, CA, USA. 4School of Human Evolution and Social Change, Arizona State University, Tempe, AZ, USA. 5State Key Laboratory of Genetic Resources and Evolution, Kunming Institute of Zoology, Kunming, China. 6 School of Biological Sciences, University of East Anglia, Norwich, UK. 7Center for Excellence in Animal Evolution and Genetics, Chinese Academy of Sciences, Kunming, China. 8Division of Human Sciences, Museu Paraense Emílio Goeldi, Belém, Brazil. 9Center for International Forestry Research, Bogor, Indonesia. 10Department of Anthropology, University of California, Santa Barbara, Santa Barbara, CA, USA. 11Center for Evolution and Medicine, Arizona State University, Tempe, AZ, USA. 12Department of Anthro- pology, Pennsylvania State University, University Park, PA, USA. 13Department of Anthropology, University of Utah, Salt Lake City, UT, USA. 14Interdisciplinary Centre for Conservation Science, University of Oxford, Oxford, UK. 15Facultad de Biologia, Universidad Veracruzana, Xalapa, Veracruz, Mexico. 16Department of Anthropology, Southern Methodist University, Dallas, TX, USA. 17Anthropology Department, Durham University, Durham, UK. 18Unaffiliated. 19Department of Anthropology, Texas A&M University, College Station, TX, USA. 20Institute for Advanced Study in Toulouse, Toulouse, France. 21Department of Anthropology, Dartmouth College, Hanover, NH, USA. 22Depart- ment of Biology, University of Richmond, Richmond, VA, USA. 23Department of Anthropology, University of the Free State, 205 Nelson Mandela Drive, Bloemfontein 9300, South Africa. 24Department of Anthropology, SIL International, Dallas, TX, USA. 25Department of Conservation Biology, University of Göttingen, Göttingen, Germany. 26 Biology Department–FMIPA, Cenderawasih University, Papua, Indonesia. 27Department of Geography and Geology, University of Turku, Turku, Finland. 28Institució Catalana de Recerca i Estudis Avançats (ICREA), Barcelona, Spain. 29Institut de Ciència i Tecnologia Ambientals, Universitat Autònoma de Barcelona, Barcelona, Spain. 30 Unit Medical Anthropology and Ecology of Disease Emergence, Institut Pasteur, Paris, France. 31Helsinki Institute of Sustainability Science (HELSUS), Faculty of Biological and Environmental Sciences, University of Helsinki, Helsinki, Finland. 32Faculty of Archaeology, Leiden University, Leiden, Netherlands. 33School of Anthropology and Conservation, University of Kent, Canterbury, UK. 34Department of Anthropology, Boise State University, Boise, ID, USA. 35Department of Food and Resource Economics, University of Copenhagen, Frederiksberg, Denmark. 36Department of Anthropology, University of Florida, Gainesville, FL, USA. 37School of Culture, History and Language, Australian National University, Canberra, Australian Capital Territory, Australia. *Corresponding author. Email: jeremy.koster@uc.edu Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 1 of 7 SCIENCE ADVANCES | RESEARCH ARTICLE to age-related variation in foraging skill in human societies. To the 18, this fictional average hunter has 89% of maximum skill. In addi- extent that foraging complexity underlies the evolution of human tion, skill declines slowly, such that skill falls below 89% of maximum life history traits, we anticipate protracted mastery of foraging tasks only after age 56. The blue shading around the posterior mean in across the life span. this plot shows the entire posterior distribution, fading out to trans- Here, we focus on age as a predictor of harvests in the largest yet parent as probability declines. There is correspondingly a lot of in- assembled database of hunting returns (Fig. 1). We focus on sub- formation in this sample about the global mean. sistence hunting for multiple reasons. First, because hunting returns While the overall pattern is clear, not every site nor individual are evident at the conclusion of an excursion; the modeling of hunt- forager exhibits the same pattern. The site-level plots in Fig. 2 illus- ing inputs and outputs generally requires fewer assumptions than trate this variation. Each site displays the mean skill function for each the analogous analysis of agricultural or pastoralist returns. Second, hunter in the sample. While there is substantial uncertainty about increased hunting is often cited as an adaptive shift that coincides with individual skill curves, there is good evidence of individual-level other distinguishing features of the genus Homo (8). However, we variation at some sites, such as the Matsigenka (9 MTS), the Colombian note that none of the societies in our sample rely exclusively on hunt- site (11 CLB), the Aché (16 ACH), and the Martu (35 MRT). Differ- ing (or foraging more generally). Contemporary hunting is not a prim- ences among individuals can be quite large. Some individuals have itive economy, but rather a recurring one that integrates with other half the adult skill of others in the same community. means of production. The data should be evaluated accordingly. For each site, the figure also shows the age of peak skill for a sta- tistically average hunter at that site, as indicated by the vertical dashed lines. While these peak ages cluster around 30 years of age, there is RESULTS noteworthy variation. On the low end, the Matsigenka (9 MTS) and There are many ways to summarize the model inferences. We focus Wola (40 WOL) peak early, near 24 years of age. Note that the best on three foundational issues that motivated the analysis. hunters at these sites tend to peak even earlier, a trend that is also 1) What is the overall pattern of skill development? evident among the Barí (7 BAR). On the higher end, the Aché (16 ACH) 2) How variable is this pattern within and between societies? and Valley Bisa (28 BIS) peak at 37 and 45, respectively, but with 3) Which components of the model—increases early in life or relatively slow declines. declines later in life—describe variation? These skill functions are inputs into site-specific production func- At the highest level of pooling, the model provides a statistical tions that also include labor and technology inputs that vary in im- answer to the question, “What is a typical human life history of for- portance across sites. This means that the relationship between aging skill?” This is very much an abstraction, one that attempts to production and age can be different from the relationship between factor away all the variation in production functions and associated skill and age. In the Supplementary Materials, we therefore present elasticities to reveal an underlying, dimensionless life history. It alternative versions of the plot that incorporate the other components cannot say much about absolute levels of production, either within of production (figs. S5 to S7). One feature of the production compo- or between societies. However, it can inform comparisons of rela- nents is that variation can arise from different sources, which, in turn, tive skill at different stages of life. have different implications for age-related variation in harvests at The statistically average hunter in this sample peaks at 33 years the individual level. Furthermore, at some sites (e.g., the Dolgan, of age (top left plot, Fig. 2). However, this peak is not sharp. At age site 30), the skill functions for individual hunters cluster around a Fig. 1. Distribution of study sites. For the key, see Table 1. Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 2 of 7 SCIENCE ADVANCES | RESEARCH ARTICLE central mean. However, this does not necessarily support inferences Table 1. Study sites and their numerical and text codes. See the help that hunters at these sites have equal skill because there may be in- file of the cchunts package for related citations. sufficient evidence to distinguish them. For an alternative perspective Dataset in cchunts on the anticipated variation among hunters within sites, we simu- Number Code Country Group package late variation from the posterior samples of the model (fig. S8). 1 CRE Canada Cree Winterhalder Last, skill functions vary both within sites and between sites. 2 MYA Belize Maya Pacheco Which components of skill contribute to this variation? To address 3 MYN Nicaragua Mayangna Koster this question, we examine the model parameters that measure vari- ation in the components k (rate of increase) and m (rate of decline) 4 QUI Ecuador Quichua Siren of the skill function. Since this is a nonlinear model, we cannot exact- 5 ECH Colombia Embera Chami Ross ly partition total variance. The impact of variation in a component 6 WAO Ecuador Waorani Franzen of skill depends on the values of all the other components. We can, 7 BAR Venezuela Bari Beckerman however, consider relative sizes of components of variation on the 8 INU Canada Inuit Ready latent scale. 9 MTS Peru Matsigenka Yu_et_al First, we find moderately greater variation in m than k within 10 PIR Peru Piro Alvard sites (note the cyan curve in the top left plot of Fig. 3). By contrast, between-site heterogeneity in age-related skill is divided roughly Van_Vliet_et_al_ 11 CLB Colombia South_America_sites evenly between variation in m and k (see the orange curve in the top left plot in Fig. 3). Some caution is necessary in these comparisons 12 PME Venezuela Pume Kramer_Greaves because the relationship between m and k is not additive. However, Fernandez_ 13 TS1 Bolivia Tsimane the implication is that among individuals within a given ecology, Llamazares skill varies more later in life than earlier in life. 14 TS2 Bolivia Tsimane Reyes-Garcia We also find a modest positive correlation between k and m (top 15 TS3 Bolivia Tsimane Trumble_Gurven right plot), suggesting that hunters who develop skill relatively 16 ACH Paraguay Aché Hill_Kintigh quickly also show advanced skill later in life. Each cyan density in 17 GB1 Gabon Coad the top right plot is the posterior correlation between hunters’ k and 18 GB2 Gabon Van_Vliet_et_al_Gabon m parameters at a given site. This correlation is particularly pro- 19 GB3 Gabon Van_Vliet_et_al_Ovan nounced for the Aché and modest otherwise. This may reflect the lack of longitudinal data on individual hunters at most study sites, Van_Vliet_et_al_ 20 CN1 DR Congo limiting what can be learned about this correlation. In contrast, the Phalanga Van_Vliet_et_al_ Aché site contains enough longitudinal data on individuals to make 21 GB4 Gabon stronger inferences about the correlation. Djoutou 22 BK1 Cameroon Baka Gallois Last, relative variation in m and k can be decomposed within and between sites. We show the posterior distributions of the SDs in 23 BK2 Cameroon Baka Duda both k (bottom left) and m (bottom right) in Fig. 3. The cyan densi- 24 CN2 Congo Van_Vliet_et_al_Ingolo ties are the SDs within sites. This corresponds to the plausible values Van_Vliet_et_al_ for variation among individuals. The orange densities are the SDs 25 CN3 Congo Ngombe among sites, corresponding to the plausible values for variation in Central African site-level averages. The dashed curves in both plots show the prior 26 BFA Bofi and Aka Lupo_Schmitt Republic distributions, which were the same for both within and between 27 CN4 DR Congo Van_Vliet_et_al_Baego components. For both k and m, there is relatively less information 28 BIS Zambia Valley Bisa Marks about variation among sites. As a result, the orange curves remain 29 HEH Tanzania Nielsen flatter than the cyan curves. There is substantially more informa- 30 DLG Russia Dolgan Ziker tion about variation within sites, and so the cyan curves are rather 31 BTK Malaysia Batek Venkataraman_et_al peaked in both cases. While there is a hint that variation between 32 PN1 Indonesia Punan Gueze sites contributes more to variation in k, while variation within sites contributes more to variation in m, strong inferences cannot be drawn 33 PN2 Indonesia Punan Napitupulu until more information is available for inferring the between-site 34 AGT Philippines Agta Headland variance. 35 MRT Australia Martu Bird_Bird_Codding 36 NUA Indonesia Nuaulu Ellen 37 NIM Indonesia Nimboran Pangau_Adam DISCUSSION Papua New Overall, these results provide an empirical counterpoint to compu- 38 NEN Nen Guinea Healey_Nen_PNG tational models of life history evolution [e.g., (6, 9)]. On the one hand, Papua New there is agreement among models about the central tendencies for 39 MAR Maring Healey Guinea the ontogeny of skill, which accelerates most rapidly during child- Papua New hood and adolescence before reaching a plateau during adulthood. 40 WOL Wola Sillitoe Guinea In all study sites, skill peaks after physical and reproductive maturity. This result is largely consistent with predictions of embodied capital Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 3 of 7 SCIENCE ADVANCES | RESEARCH ARTICLE Global mean 0 --KEY-- Ids (obs) 1 CRE 16 (127) 2 MYA 59 (464) 3 MYN 52 (359) 4 QUI 32 (189) Obs range Peak 26 31 35 33 30 0 18 33 56 80 0 40 80 0 40 80 40 80 40 80 0 40 80 5 ECH 2 (6) 6 WAO 48 (373) 7 BAR 18 (233) 8 INU 15 (29) 9 MTS 69 (1441) 10 PIR 42 (274) 31 33 30 28 24 30 0 40 80 40 40 80 11 CLB 12 PME 23 (172) 13 TS1 14 TS2 37 (139) 15 TS3 16 ACH 147 (14364) 26 28 32 36 32 37 0 40 80 0 40 80 40 80 40 80 40 80 17 GB1 18 GB2 19 (114) 19 GB3 15 (69) 20 CN1 14 (116) 21 GB4 12 (44) 22 BK1 80 (249) 35 33 34 29 31 29 40 40 80 0 40 80 40 80 0 40 80 40 80 23 BK2 24 CN2 16 (70) 25 CN3 6 (76) 26 BFA 59 (433) 27 CN4 3 (33) 28 BIS 24 (231) 26 30 33 28 29 45 40 80 0 40 80 0 40 80 40 80 0 40 80 40 80 29 HEH 45 (45) 30 DLG 26 (76) 31 BTK 27 (268) 32 PN1 35 (119) 33 PN2 23 (125) 34 AGT 44 (211) 31 33 34 28 32 29 40 40 80 40 80 0 40 80 0 40 80 40 80 35 MRT 77 (758) 36 NUA 36 (140) 37 NIM 26 (180) 38 NEN 7 (7) 39 MAR 6 (28) 40 WOL 27 (410) 30 32 39 32 33 24 40 40 80 40 80 0 40 80 0 40 80 0 40 80 Fig. 2. Skill functions. The figure depicts the global average of skill (top left plot) and skill at the respective study sites. Within study sites, each curve is the posterior mean skill for an individual hunter, standardized to the maximum within each site. In the header of each plot, the site number and three-letter code are shown along with the number of individual hunters in each sample, followed by the number of observed harvests in parentheses. The orange span of ages corresponds to ages observed within each site, while the gray ranges were unobserved and are instead implied by the underlying model. The vertical dashed lines show the average ages at peak within sites. theory (4), although, on average, the empirical data do not reveal sponse to variables such as experience, motivation, opportunities for the distinct peaks in skill during middle age that are evident in previous social learning, and the physical and cognitive demands of hunting studies of the Aché (10, 11). Instead, the empirical model suggests in different socioecological environments. As opposed to a canalized that skill starts to plateau by early adulthood and that only moderate human life history strategy, this study suggests potential develop- increases are evident subsequently. Declines are typically slow, such mental plasticity in traits associated with hunting skill, which manifest that an 80-year-old may retain two-thirds of maximum skill. not only in contemporary settings but also potentially in ancestral Another noteworthy result is the extent of variability in skill, both settings. These results further imply that singular study sites can rarely among and within sites. Cross-cultural variation is evident in the be viewed as straightforward analogs for evolutionarily relevant en- rate at which hunters develop peak skill. Within sites, the rate at vironments (13). which hunters develop skill is relatively homogeneous compared to Within the respective study sites, the model brings new attention the variation that distinguishes young hunters in different study sites. to the variation in skill among individual hunters. What explains To explain cross-cultural variation in the development of foraging this variation? Data and theory suggest that physical strength and proficiency, it is common and reasonable for anthropologists to stamina, accumulated knowledge, and motivation all plausibly con- emphasize ecological predictors, such as extrinsic mortality risks tribute to age-related variation (14, 15). In most empirical datasets, [e.g., (12)]. However, varying rates of skill development may stem however, only data on the hunters’ ages are available as predictors, as well from mediating social factors that relate only indirectly to not other attributes of the individuals. The peak in the average skill ecological differences. Additional theorizing is needed to generate function at approximately 33 years old is also near the age when hypotheses about the cross-cultural ontogeny of foraging skill in re- physical strength and ecological knowledge plateau (15–18). From a Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 4 of 7 SCIENCE ADVANCES | RESEARCH ARTICLE Relative variation k/m Correlation k,m among individuals 2.0 Within 2.0 Global average [Aché] 1.5 1.5 Between Density Density 1.0 1.0 0.5 0.5 Sites 0.0 0.0 0.0 0.5 1.0 1.5 2.0 −1.0 −0.5 0.0 0.5 1.0 (SD k)/(SD m) Correlation Variation in k (increase rate) Variation in m (decline rate) 6 6 5 5 Between 4 4 Within Within Density Density 3 3 Between 2 2 1 1 Prior 0 0 0.0 0.2 0.4 0.6 0.8 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 SD SD Fig. 3. Variation in components of skill. (Top left) Relative variation in k and m. The horizontal axis is the ratio of the SD of k to the SD of m. The vertical dashed line at 1 indicates equality of variances. The orange density is between-site variation. The cyan density is within-site variation. There is more variation in declines (m) than increases (k) in skill within sites, whereas the ratio is roughly equivalent across sites. (Top right) Correlation between k and m among individuals within sites. The orange density is the global average. Each cyan density represents a single site. The Aché stand out and are labeled separately. (Bottom left and bottom right) Variation in k (left) and m (right) comparing variation within and between sites. theoretical perspective, an optimal life history should develop these has received relatively little consideration (22). Our analysis suggests components together, with the important caveat that brain growth that this variability typifies communities of human foragers while may need to precede body growth, to enable learning (9). To distin- concomitantly altering the effectiveness of food sharing for buffer- guish the relative importance of phenotypic traits to variation in ing risk. That is, when hunters vary substantially in their skill and hunting skill across the life span, data are needed that link measure- productivity, there are asymmetric benefits to participation in risk-­ ments of these attributes to individual-level hunting returns, ideally pooling distribution systems. To the extent that prosociality and other longitudinally. In the interim, the currently available data suggest that traits in the human lineage stem from the cooperative challenges individual hunters develop physical and cognitive abilities in con- posed by this asymmetry, the high variation in hunting skill across cert, resulting in high hunting success by their late 20s and early 30s. the life span merits further attention. Estimates from empirical studies provide inputs to the parame- terization of computational life history models. Our analysis of hunting by age provides refined estimates of the average skill function, a target MATERIALS AND METHODS of inference for recent theoretical work (6, 9). However, the added The total sample contains 23,747 observations of 1821 individual clarity about average skill is belied by the substantial heterogeneity foragers across 40 study sites (Fig. 1). There is substantial imbalance that is evident among individual hunters. For future theoretical de- in sample size across units. One site contributes only six trips from velopments about the unique human life history pattern, this varia- two individuals. Another contributes more than 14,000 trips from tion in skill merits careful attention. Prevailing hypotheses about the 147 individuals. Some individuals contribute only a single outcome, adaptive shift to hunting by human ancestors assert that reciprocal while others contribute dozens. The majority of the sample comprises food sharing in small bands was necessary to smooth variance in male hunters, with too little data on female hunters to infer general- consumption, given that daily harvests by hunters are unpredictable izable sex differences. Most sites contribute primarily cross-sectional (19–21). In this literature, variability in the skill of individual hunters data, while a few others exhibit impressive time series. Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 5 of 7 SCIENCE ADVANCES | RESEARCH ARTICLE Because skill cannot be directly observed, what is required is a model ​​S(x ) = M(x ) K ​(x)​​  b​​ ​ (4) with latent age-varying skill, which informs a production function for observable foraging returns. The model is described in detail in where b is represents the elasticity for knowledge in the skill function. the Supplementary Materials. Building on earlier research (10), our As detailed in the Supplementary Materials, the von Bertalanffy modeling framework was developed in a grant proposal and reviewed model permits diverse functional forms of skill across the life span, before seeing the assembled sample. Using a Cobb-Douglas produc- ranging from approximately sigmoidal to roughly quadratic shapes. tion function common to economic research, we model hunting returns Although we assume that the growth and senescence components (harvest) on excursions as a standard log-linear function of skill, labor of skill relate to proximate mechanisms, such as age-related varia- inputs, and auxiliary inputs from technology and cooperation tion in ecological knowledge and physical abilities, the available data do not allow us to examine those proximate mechanisms [cf. (15, 18)]. ​harvest = ​S​​  η​ ​L​​  β​  α​ (1) As a result, age-related variation in skill must be inferred from its effects on the observed productivity of hunters of heterogeneous ages. where S represents the hunter’s skill with its elasticity, , the labor The statistical model allows the von Bertalanffy parameters to vary input and its elasticity are represented by L, and  is a linear model across individuals, reflecting different rates of increasing skill or se- for covariates such as group size, the number of assistants, the use of nescence among hunters. The parameters also vary across study dogs, and the use of firearms. These latter variables have been shown sites, allowing for different rates of increase and decline in heteroge- to influence hunting returns (23–26). Note that production requires neous environments. As noted previously, we model hunting returns both skill and labor; if either is zero, then there is no harvest. The using a zero-augmented gamma model, and the respective skill func- respective elasticities reflect the proportional effects that skill and tions of hunters and societies are estimated jointly from the Bernoulli labor have on harvests. That is, increasing either skill or labor results and gamma functions. It would have been possible to estimate sep- in increased harvests, but the scale of the increase is reflected in the arate, correlated individual-level and site-level skill functions for the elasticities. Our parameterization of the function does not impose Bernoulli and gamma functions, respectively. Our modeling approach constant returns to scale [cf. (27)]. instead assumes that increases in skill have comparable effects in Data on hunting returns pose the particular challenge of including terms of reducing the probability of unsuccessful excursions and many zeroes (for hunting trips in which nothing was acquired), and increases in the amount of meat acquired on successful outings. This the harvests on successful hunting trips exhibit positive skew. We there- assumption receives support from an earlier analysis of the Aché fore adapt the Cobb-Douglas function in Eq. 1 to a zero-augmented dataset, which suggests that there is a positive correlation (and no model in which the zeroes and nonzero harvests are modeled via evident tradeoff) between hunters’ probabilities of acquiring some- separate functions, as detailed in the Supplementary Materials. As thing and the amounts that are harvested on successful trips (10). in previous work (10), we use a Bernoulli distribution to model the In addition to the varying effects on the skill parameters, the model probability of success versus failure, and the distribution of nonzero also allows the parameters for labor inputs and covariates in the harvests is assumed to follow a gamma distribution. linear model to vary across study sites. Missing data, particularly for To model latent skill across the life span, we adapt the von Bertalanffy hunt duration and technological variables, are common and addressed growth model (28). The benefit of the model is that individuals’ skill using Bayesian imputation and averaging methods. is assumed to be lowest at birth with eventual declines due to senes- At some study sites, hunters work cooperatively to harvest prey, cence. Within these constraints, the functional form of the model and the data on these excursions assign the hunting returns to the potentially exhibits considerable diversity depending on the empirical group, not individual hunters. In those cases, we replace individual data. Age-related variation in skill is determined by a rate of growth, hunter skill in the production equation with the weighted average of k, and a rate of decline, m. The growth component, k, potentially the skill of the group members. The statistical model follows the includes ecological knowledge, strength, cognitive function, and other principles of a multiple membership model (29, 30). When hunters traits that underlie foraging success but that exhibit reduced accel- are observed in different combinations of groups, it is possible to eration with age. For simplicity, the composite of these attributes can distinguish differences in skill between them. be dubbed knowledge. For a hunter of age x, the growth function, K, We cannot rule out selection biases that complicate inferences. is represented by For instance, if there were a study site where highly skilled hunters are active regardless of environmental conditions and the relatively ​​K(x ) = 1 − exp(− kx)​​ (2) unskilled hunters are active only when returns are expected to be particularly favorable, then the estimated variation in hunters’ skill where k is a parameter greater than zero that reflects the rate of would likely be lower than a site where hunting activity occurs inde- increase. pendently of skill. The declining component, m, reflects the senescence of traits re- We fit our model using a Hamiltonian Monte Carlo sampling lated to hunting skill. For a forager of age x, the decline in produc- algorithm in the Rstan package, version 2.16.1 (31). We implemented tive ability, M, decreases at a constant rate given by the model both as a forward simulation and as a statistical model. The forward simulation validates that the statistical model can re- ​​M(x ) = exp(− mx ) #​​ (3) cover parameters from data with known values. The statistical model, estimated from 10 chains of 500 iterations, exhibits efficient mixing where m is a parameter greater than zero that reflects the rate of and adequate diagnostics. The data and convenience functions are decline. included as part of the cchunts R package, which is available along- Age-specific skill is then represented by the weighted product of side the model posterior, coding scripts, and accompanying infor- the two preceding functions mation about the Open Science Framework (https://osf.io/2kzb6/). Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 6 of 7 SCIENCE ADVANCES | RESEARCH ARTICLE SUPPLEMENTARY MATERIALS 21. H. S. Kaplan, E. Schniter, V. L. Smith, B. J. Wilson, Risk and the evolution of human Supplementary material for this article is available at http://advances.sciencemag.org/cgi/ exchange. Proc. Biol. Sci. 279, 2930–2935 (2012). content/full/6/26/eaax9070/DC1 22. R. Boyd, The evolution of reciprocity when conditions vary, in Coalitions and Alliances in Humans and Other Animals, A. Harcourt, F. DeWaal, Eds. (Oxford Univ. Press, 1992), pp. 473–489. REFERENCES AND NOTES 23. R. B. Hames, A comparison of the efficiencies of the shotgun and the bow in neotropical 1. J. H. Jones, Primates and the evolution of long, slow life histories. Curr. Biol. 21, forest hunting. Hum. Ecol. 7, 219–252 (1979). R708–R717 (2011). 24. E. A. Smith, Inuit foraging groups: Some simple models incorporating conflicts of interest, 2. C. O. Lovejoy, The origin of man. Science 211, 341–350 (1981). relatedness, and central-place sharing. Ethol. Sociobiol. 6, 27–47 (1985). 3. K. Hawkes, J. F. O’Connell, N. B. Jones, H. Alvarez, E. L. Charnov, Grandmothering, 25. M. Alvard, Shotguns and sustainable hunting in the Neotropics. Oryx 29, 58–66 menopause, and the evolution of human life histories. Proc. Natl. Acad. Sci. U.S.A. 95, (1995). 1336–1339 (1998). 26. J. M. Koster, Hunting with dogs in Nicaragua: An optimal foraging approach. Curr. 4. H. Kaplan, K. Hill, J. Lancaster, A. M. Hurtado, A theory of human life history evolution: Anthropol. 49, 935–944 (2008). Diet, intelligence, and longevity. Evol. Anthropol. 9, 156–185 (2000). 27. K. J. Arrow, H. B. Chenery, B. S. Minhas, R. M. Solow, Capital-labor substitution 5. E. Herrmann, J. Call, M. V. Hernández-Lloreda, B. Hare, M. Tomasello, Humans have and economic efficiency. Rev. Econ. Stat. 43, 225–250 (1961). evolved specialized skills of social cognition: The cultural intelligence hypothesis. Science 28. L. von Bertalanffy, Untersuchungen Über Die Gesetzlichkeit Des Wachstums : I. Teil: 317, 1360–1366 (2007). Allgemeine Grundlagen Der Theorie; Mathematische Und Physiologische 6. M. González-Forero, A. Gardner, Inference of ecological and social drivers of human Gesetzlichkeiten Des Wachstums Bei Wassertieren. Wilhelm Roux Arch. Entwickl. Mech. brain-size evolution. Nature 557, 554–557 (2018). Org. 131, 613–652 (1934). 7. A. R. DeCasien, S. A. Williams, J. P. Higham, Primate brain size is predicted by diet but not 29. W. J. Browne, H. Goldstein, J. Rasbash, Multiple membership multiple classification sociality. Nat. Ecol. Evol. 1, 0112 (2017). (MMMC) models. Stat. Model. 1, 103–124 (2001). 8. R. A. Foley, The evolutionary consequences of increased carnivory in hominids, in The 30. B. Aven, H. Hillmann, Structural role complementarity in entrepreneurial teams. Early Human Diet: The Role of Meat; Meat-Eating and Human Evolution, C. B. Stanford, Manage. Sci. 64, 5688–5704 (2018). H. T. Bunn, Eds. (Oxford Univ. Press, 2001), pp. 305–331. 9. M. González-Forero, T. Faulwasser, L. Lehmann, A model for brain life history evolution. 31. Stan Development Team, RStan: The R interface to Stan (2016). R package version PLOS Comput. Biol. 13, e1005380 (2017). 2.14.1. 10. R. McElreath, J. Koster, Using multilevel models to estimate variation in foraging returns. Effects of failure rate, harvest size, age, and individual heterogeneity. Hum. Nat. 25, Acknowledgments: Audiences in London, Lausanne, Nijmegen, Aarhus, UCLA, and the 100–120 (2014). University of Utah contributed useful feedback on draft analyses and interpretations, as did 11. R. Walker, K. Hill, H. Kaplan, G. McMillan, Age-dependency in hunting ability among the anonymous reviewers. The map was created was the multimedia department at the Max the Ache of Eastern Paraguay. J. Hum. Evol. 42, 639–657 (2002). Planck Institute for Evolutionary Anthropology. Funding: This study was funded by the NSF 12. N. Blurton Jones, K. Hawkes, P. Draper, Foraging returns of !Kung adults and children: (#1534548). Author contributions: J.K. and R.M. conceived the study, analyzed the data, and Why didn’t !Kung children forage? J. Anthropol. Res. 50, 217–248 (1994). wrote the paper. All other authors contributed hard-earned field data. Competing interests: 13. W. Irons, Adaptively relevant environments versus the environment of evolutionary The authors declare that they have no competing interests. Data and materials availability: adaptedness. Evol. Anthropol. 6, 194–204 (1998). All data needed to evaluate the conclusions in the paper are present in the paper and/or the 14. N. Blurton Jones, F. W. Marlowe, Selection for delayed maturity : Does it take 20 years Supplementary Materials. Supplementary data and code are available at https://osf.io/2kzb6/. to learn to hunt and gather? Hum. Nat. 13, 199–238 (2002). Additional data related to this paper may be requested from the authors. 15. M. Gurven, H. Kaplan, M. Gutierrez, How long does it take to become a proficient hunter? Implications for the evolution of extended development and long life span. J. Hum. Evol. Submitted 3 May 2019 51, 454–470 (2006). Accepted 17 March 2020 16. S. Zent, E. López-Zent, Ethnobotanical convergence, divergence, and change among Published 24 June 2020 the Hoti of the Venezuelan Guayana. Adv. Econ. Bot. 15, 37–78 (2004). 10.1126/sciadv.aax9070 17. K. Demps, F. Zorondo-Rodŕıguez, C. Garćıa, V. Reyes-Garćıa, Social learning across the life cycle: Cultural knowledge acquisition for honey collection among the Jenu Kuruba, India. Citation: J. Koster, R. McElreath, K. Hill, D. Yu, G. Shepard Jr., N. van Vliet, M. Gurven, B. Trumble, Evol. Hum. Behav. 33, 460–470 (2012). R. B. Bird, D. Bird, B. Codding, L. Coad, L. Pacheco-Cobos, B. Winterhalder, K. Lupo, D. Schmitt, 18. J. Koster, O. Bruno, J. L. Burns, Wisdom of the elders? Ethnobiological knowledge across P. Sillitoe, M. Franzen, M. Alvard, V. Venkataraman, T. Kraft, K. Endicott, S. Beckerman, S. A. Marks, the lifespan. Curr. Anthropol. 57, 113–121 (2016). T. Headland, M. Pangau-Adam, A. Siren, K. Kramer, R. Greaves, V. Reyes-García, M. Guèze, R. Duda, 19. G. Isaac, The food-sharing behavior of protohuman hominids. Sci. Am. 238, 90–108 (1978). Á. Fernández-Llamazares, S. Gallois, L. Napitupulu, R. Ellen, J. Ziker, M. R. Nielsen, E. Ready, C. Healey, 20. B. Winterhalder, Diet choice, risk, and food sharing in a stochastic environment. C. Ross, The life history of human foraging: Cross-cultural and individual variation. Sci. Adv. 6, J. Anthropol. Archaeol. 5, 369–392 (1986). eaax9070 (2020). Koster et al., Sci. Adv. 2020; 6 : eaax9070 24 June 2020 7 of 7