Significance Here we show that the spatial prevalence of human societies that believe in moralizing high gods can be predicted with a high level of accuracy (91%) from historical, social, and ecological data. Using high-resolution datasets, we systematically estimate the relative effects of resource abundance, ecological risk, cultural diffusion, shared ancestry, and political complexity on the global distribution of beliefs in moralizing high gods. The methods presented in this paper provide a blueprint for how to leverage the increasing wealth of ecological, linguistic, and historical data to understand the forces that have shaped the behavior of our own species.

Abstract Although ecological forces are known to shape the expression of sociality across a broad range of biological taxa, their role in shaping human behavior is currently disputed. Both comparative and experimental evidence indicate that beliefs in moralizing high gods promote cooperation among humans, a behavioral attribute known to correlate with environmental harshness in nonhuman animals. Here we combine fine-grained bioclimatic data with the latest statistical tools from ecology and the social sciences to evaluate the potential effects of environmental forces, language history, and culture on the global distribution of belief in moralizing high gods (n = 583 societies). After simultaneously accounting for potential nonindependence among societies because of shared ancestry and cultural diffusion, we find that these beliefs are more prevalent among societies that inhabit poorer environments and are more prone to ecological duress. In addition, we find that these beliefs are more likely in politically complex societies that recognize rights to movable property. Overall, our multimodel inference approach predicts the global distribution of beliefs in moralizing high gods with an accuracy of 91%, and estimates the relative importance of different potential mechanisms by which this spatial pattern may have arisen. The emerging picture is neither one of pure cultural transmission nor of simple ecological determinism, but rather a complex mixture of social, cultural, and environmental influences. Our methods and findings provide a blueprint for how the increasing wealth of ecological, linguistic, and historical data can be leveraged to understand the forces that have shaped the behavior of our own species.

Ecological uncertainty has long been implicated in the expression of cooperation and conflict in animal societies (1). For example, climatic variation predicts the incidence and distribution of cooperative breeding in birds (2, 3), and ecological uncertainty is associated with group living in mammals (4). These findings suggest that ecological duress can promote sociality when the benefits of cooperating during difficult periods outweigh the potential costs during benign ones (3). Although similar ecological pressures are likely to act upon humans (5, 6), the extent to which they have shaped the development of our own sociality remains disputed and is currently unclear (7⇓–9).

Human behavior is often mediated by socially transmitted conventions—that is, cultural norms—that govern expectations and behaviors during social interactions (10, 11). In particular, religious beliefs are thought to be a powerful mechanism for the enforcement of social rules (12⇓⇓–15). In support of these ideas, comparative studies have shown that a belief in moralizing high gods—defined as supernatural beings believed to have created or govern all reality, intervene in human affairs, and enforce or support human morality (7)—tends to be more prevalent among societies that recognize rights to movable property (9, 16), as well as in those that exhibit higher levels of political complexity (17), subsistence productivity (17), and norm compliance (18). In addition, psychological experiments have found that moralizing high god concepts can reduce levels of cheating (14) and increase willingness to be fair (19) and to cooperate (13). Such findings have provoked speculation that the prominence of moralizing gods at the level of cultures is a better predictor of cooperation in humans than individual religious beliefs (20). However, debates persist because statistical models assessing historical and regional dependencies are rare (e.g., refs. 7 and 9), and those that exist are limited in scope (21).

It has long been theorized that the global distribution of religious beliefs may be shaped by ecological factors (7⇓–9). Recent empirical findings indicate that beliefs in moralizing high gods not only intensify (22⇓⇓⇓⇓–27), but also promote cooperation (28⇓⇓⇓⇓–33) in situations of increased environmental risk. In addition, these findings indicate that ecological threats can strengthen mechanisms of norm enforcement in human groups (34). Given the strong correlation between cooperation and ecological uncertainty in nonhuman animals (2, 3), these findings are especially suggestive of a link between ecology and religion. In particular, the similarity between the global distribution of belief in moralizing high gods (Fig. 1), and that of social cooperation in birds (see figure 1A in ref. 2), raises the possibility that the distribution of this type of religious belief might also be shaped by ecological factors. Nevertheless, the evidence for an association between religion and ecology is currently mixed. For example, prior studies have indicated that resource scarcity is both positively (8, 9) and negatively (7) associated with a belief in moralizing high gods. These apparently contradictory findings may be the product of methodological issues. Specifically, some of the studies involved (7, 8) failed to account for the nonindependence of societies as a result of spatial proximity and common ancestry [Galton’s problem (35)], whereas others were based on a relatively small sample of cultures (9) that may have been unduly influenced by a few preclassical cultures (21). In addition, all of these studies used extremely coarse metrics of ecology [e.g., subjective ratings of resource abundance in (7, 8), or indirect measures of agricultural potential in (9)] that could have led to a distorted view of the role of environmental factors in shaping the spatial variation of religious beliefs.

Fig. 1. Global distribution of societies that exhibit beliefs in moralizing high gods (blue) or not (i.e., atheism or beliefs in nonmoralizing deities or spirits in red). The underlying map depicts the mean values of net primary productivity (i.e., the net balance of monthly consumption relative to production of carbon dioxide by living plants) in gray scale. Darker localities reflect places with greater potential for overall plant growth.

To rigorously test for an association between ecology and beliefs in moralizing high gods, we systematically modeled the effects of basic environmental variables while accounting for the effects of known covariates of religious beliefs and simultaneously adjusting for potential dependencies related to shared cultural ancestry—as measured by patterns of shared language origin (35)—and spatial diffusion—as measured by the effects of cultural traits in neighboring groups (Methods). We consider the following covariates in our models, all positively correlated with a belief in moralizing high gods: political complexity (17), practice of agriculture (17), and recognition of rights to movable property, as measured by the practice of animal husbandry (9, 16). Our analysis is based on a global cross-cultural sample of 583 human societies (36, 37) that occupy a range of different habitats and are exposed to a wide array of environmental conditions (Fig. 1). Through a comprehensive characterization of key political and cultural variables, together with a rigorous treatment of key environmental factors, we offer a systematic test for whether the global distribution of belief in moralizing high gods is associated with underlying variation in ecological parameters.

Results Environmental variables are often highly correlated (see ref. 38) and this multicollinearity is likely to lead to inaccurate and unstable estimates of statistical effects in ecological analyses (38, 39). To prevent these issues, we began by reducing the original set of ecological predictors to a smaller number of composite variables using principal components analysis (PCA) (Methods and Table 1). PC1, labeled resource abundance, captures a gradient of increasing exposure to more abundant rainfall, higher primary productivity, and greater biodiversity. This composite variable can therefore be broadly interpreted as a relative indicator of the abundance of natural resources that an environment has to offer [note that societies with similar PC1 scores may nevertheless differ in the type of resources they rely on, as well as on the extent to which they depend on them (5)]. PC2, labeled climate stability, captures a gradient of increasing exposure to more predictable annual cycles of precipitation and temperature, as well as to warmer, more stable temperatures throughout the year. As expected, both PC1 (Pearson’s ρ = −0.52, P < 0.001) and PC2 scores (ρ = −0.76, P < 0.001) decrease with latitude in our sample. Table 1. PCA of ecological variables, with varimax rotation with main contributors identified in boldface type Having derived independent ecological indicators for resource abundance and stability, we subsequently ran multiple statistical models that explored all of the possible combinations of predictors within our candidate set. These models were estimated from the same randomly chosen subset of two-thirds of societies in our dataset (n = 389 societies) (Methods). We excluded from this analysis all models that considered the effects of agriculture and animal husbandry simultaneously because these variables are highly correlated (Pearson’s χ2 = 114.3, df = 1, P < 0.001) and exhibit an almost one-to-one coding correspondence (i.e., in 82% of societies, agriculture and animal husbandry are either both absent or both present). Because no single model parameterization was clearly better supported than the others (Table 2), all inferences below are based on an average model (Table 3) computed using the multimodel inference procedures outlined in ref. 40 (Methods). Model averaging is a formal way to account for model uncertainty and reduce model selection bias (40). Briefly, the process involves estimating the overall effect of a given predictor by adding its estimated effects in the different candidate models weighted by their corresponding Akaike weights (40). Table 2. Support for alternative parameterizations of our model of religious beliefs (only the 10 best supported models are shown here) Table 3. Multimodel average for the candidate set of models of religious beliefs presented in Table 2 Our results confirm that societies that share a common language origin or coexist in close spatial proximity are likely to exhibit similar types of religious beliefs (Table 3). Similarly, they indicate that even after accounting for such historical and spatial dependencies, other predictors in our candidate set provide meaningful information on the probability of belief in moralizing high gods (Table 3). Overall, we find that belief in moralizing high gods is more prevalent in societies that practice animal husbandry (Fig. 2A), exhibit greater political complexity (Fig. 2B), and have less access to food and water (Fig. 2C). In addition, we find that the potential effects of agriculture on this type of religious beliefs are not only weakly supported (see low relative variable importance and predictive value in Table 3), but also likely to be explained by the correlation between agriculture and animal husbandry (Supporting Information and Tables S1 and S2). We also find that for the vast majority of societies (i.e., those above the 15th percentile of resource abundance), there is a further increase in the probability of belief in moralizing high gods when environments are more variable and unpredictable (closed circles and solid line in Fig. 2D). However, in societies that are chronically exposed to very low resource abundance levels (open circles and dashed line in Fig. 2D), we see an increase in the probability of belief with increasing environmental stability. The different effect of stability in extremely poor environments is also consistent with the idea that beliefs in moralizing high gods are more prevalent in societies that are more likely to be exposed to ecological duress: when environmental conditions are already as harsh as they can be, climatic changes imply that conditions may actually improve (i.e., environmental changes could bring about a sudden influx of rain or an unexpected period of milder temperatures). Fig. 2. Probability of belief in moralizing high gods among societies that vary in (A) the use of animal husbandry, (B) political complexity, (C) resource abundance, and (D) climate stability. Colors depict the belief (blue), or absence thereof (red), in moralizing high gods and lines indicate the fits of the model presented in Table 2. To facilitate the visual assessment of fits in C and D, we present here the frequency of moralizing high gods in eight equally spaced bins. In D, filled circles and the solid line depict societies above the 15th percentile of resource abundance (i.e., moderate- to high-abundance), whereas open circles and the dashed line depict societies in the lowest 15th percentile of abundance (i.e., resource-poor habitats). We now ask: How well does our model predict the distribution of beliefs in moralizing high gods? We computed the area under receiver operating characteristic curves (41), hereafter AUC (Methods), to evaluate the accuracy of our predictions for the beliefs of the 194 societies that were initially excluded from the analysis. Briefly, an AUC value of 1 indicates perfect discrimination whereas a value of 0.5 implies that predictive ability is no better than chance. The average model (Table 3) exhibits an AUC of 0.91, which indicates an excellent predictive accuracy (41) and suggests that our analysis has captured the main covariates of this aspect of religious beliefs.

Discussion Through the systematic analysis of a broad sample of human societies we have shown that the global distribution of beliefs in moralizing high gods is well explained by the combined effects of spatial proximity, political complexity, ecology, animal husbandry, agriculture, and language origin. The predictive accuracy of our model is equivalent or greater than that of a recent model of prosocial behavior in birds (AUC range = 0.67–0.92 in ref. 2). We note that in the (human) social sciences a predictive accuracy of AUC = 0.75 is already considered a large effect (42). Our results support earlier claims of a link between religion and social behavior by demonstrating that the belief in moralizing high gods is more prevalent in societies that recognize rights to movable property (9, 16), and have achieved higher levels of political complexity (9). Similarly, our finding that the recognition of rights to movable property is a better predictor of the global distribution of belief in moralizing high gods than agriculture (Table 2 and Supporting Information), is consistent with the idea that the prosocial effects of moralizing high gods are not specific to very large (“ultrasocial”) communities (43). We note as well that there appears to be no solid evidence of a quadratic effect of political complexity (see ref. 9) and that even when such an effect was considered, the nonlinear effect was minor and led to a similar fit to the one shown in Fig. 2B (thus indicating potential overfitting of the quadratic effect to the very few societies that have achieved the highest level of political complexity). In general, our findings are consistent with the notion that a shared belief in moralizing high gods can improve a group’s ability to deal with environmental duress and may therefore be ecologically adaptive (44). Philosophers have long hypothesized that the belief in morally concerned deities can reduce both anxiety (45) and existential suffering (46), and that environmental harshness can increase the appeal of a cosmic authority (47). Accordingly, we have uncovered robust evidence of a positive association between societal acceptance of moralizing high god concepts and both resource scarcity (9) and ecological risk. We emphasize that this inference applies only to variation in beliefs among social groups (as opposed to individual variation within them) because entire societies—not individual beliefs—are the unit of analysis in our study. Thus, overall, our findings indicate that the global prevalence of societies that believe in moralizing high gods is associated with higher prospects of ecological duress (i.e., conditions that are more likely to vary in reasonably productive environments or more likely to stay the same in environments that are already quite poor). In addition, the findings suggest that isolated reports of a soothing effect of ritualistic behavior in the face of randomness (48), or of the intensification of religious beliefs (22⇓⇓⇓⇓–27) and the propensity to cooperate (28⇓⇓⇓⇓–33) in riskier environments, are likely to be reflections of a general human response to threat. Other potentially important mechanisms that could help explain the global distribution of beliefs in moralizing high gods are cultural diffusion and shared ancestry. For example, the strong effects of spatial proximity in Tables 2 and 3 remind us that religious beliefs can sometimes spread by virtue of contact and proselytism. Similarly, the observed effects of language family suggest that religious beliefs could sometimes be carried along with other political and cultural traits as the descendants of ancestral cultural groups successfully colonize new sites and give rise to new societies (21). Nevertheless, although spatial proximity exhibits a high relative variable importance in our analysis, language family does not, indicating that the latter is less informative in this context (Table 3). It is interesting to note that relative variable importance [i.e., the sum of Akaike information criterion corrected for finite samples (AICc) weights of all of the models that include a given predictor] and predictive value (i.e., the predictive accuracy or AUC of the model that considers only that predictor) paint very different pictures of the role of cultural ancestry in shaping religious beliefs: although an effect of language family is poorly supported by the data (i.e., low relative variable importance), the predictive accuracy of this predictor is almost as good as the AUC of the multimodel average. A possible reason for this apparent discrepancy is that the spread of culture is likely to be shaped by both space and ecology. For example, cultural diasporas may be more successful in colonizing habitats to which they are already well adapted, and neighboring localities are likely to offer similar climates and resource levels (Fig. 3). Thus, language family may not only be a good indicator of shared cultural ancestry, but also a reasonably good proxy for ecology and spatial proximity. Fig. 3. Similarity by distance in resource abundance (black line) and climate stability (gray line). More positive values of Moran’s I indicate that societies exhibit greater similarity than expected purely by chance at a given distance, whereas more negative values indicate greater than expected negative auto-correlation. Overall, our results demonstrate the power of combining the tools of ecology and the social sciences to unravel the complex fabric of human history. Through a careful consideration of the historical, spatial, and ecological circumstances under which a wide variety of human societies exist, we have shown that the global prevalence of religious beliefs in moralizing high gods arises through a combination of social and ecological factors. A promising avenue for future research will therefore be to assess the extent to which cultural beliefs are initially transmitted through conquest, trade, proselytism, and borrowing (12, 49)—that is, the mechanisms most directly associated with spatial proximity and, to a lesser extent, with shared ancestry—and subsequently conserved or modified by the adaptive benefits they confer in their ecological settings. In moving forward we caution that the coarse estimates of shared history currently provided by language family classifications should eventually be replaced by more accurate cultural phylogenies, some of which are already available at smaller scales (50⇓–52). Similarly, we note that the belief in moralizing high gods is only one of the many aspects of religious beliefs and practices that are important for human social behavior and may be amenable to the analytical approach we have illustrated here.

Methods Data and Data Sources. All data collated for this study are included in Dataset S1. Our unit of analysis is a human society, defined as a geographically distinct group of people with a shared cultural identity at the time and place of sampling (36, 37). In 988 of the 1,267 societies included in the Ethnographic Atlas, cultural identity is reflected in the use of a unique language. In the remaining cases, groups that share a common language have been recognized as distinct societies because they exhibit pronounced differences in their regional dialects or cultural practices. The following variables were obtained from descriptions in the Ethnographic Atlas (36, 37): geographic location (v104 and v106), religious beliefs [v34, recoded as supportive of human morality (category 4) versus not (categories 1–3)], agriculture [v28, recoded as present (categories 2–6) or absent (category 1)], animal husbandry [v40, recoded as presence (categories 2–7) or absence of large domestic animals (category 1)], and political complexity (v33, number of jurisdictional hierarchy levels beyond local communities). In accordance with ref. 9, we used the presence of animal husbandry as a proxy for movable property and explored nonlinear effects of political complexity in our models. Biogeographic analyses have historically relied on latitude as a convenient proxy to characterize the series of well-known environmental changes that occur as we move away from the tropics. By taking advantage of the recent availability of high-resolution global datasets, we have replaced such imperfect surrogates with the actual ecological variables of interest to gain a better mechanistic understanding of the potential ecological factors underlying spatial variation in religious beliefs. Precipitation and temperature data were extracted from monthly global maps in the CRU-TS 3.1 Climate Database (53) at a resolution of 0.5 × 0.5° cells. For each society we measured the annual mean, variance, and predictability of climate variables from 1901 to 1950 in the corresponding cell containing its sampling locality as listed in the Ethnographic Atlas. Observations were restricted to this time frame to match the period during which the majority of the societies in our dataset were actually sampled. Predictability of climate patterns was measured via Colwell’s P (54), an information theory index that captures year-to-year variation in the onset, intensity, and duration of periodic phenomena. Colwell’s P ranges from 0 (completely unpredictable) to 1 (fully predictable). Other data sources include NASA’s Moderate Resolution Imaging Spectroradiometer data, MODIS, for net primary productivity (downloaded from neo.sci.gsfc.nasa.gov, accessed December 5, 2013), Jenkins et al. (55) for vertebrate richness, and Kreft and Jetz (56) for vascular plant richness. Of the 775 societies for which data on religious beliefs are available in the Ethnographic Atlas, we retained only the 583 societies with complete information on all other environmental and cultural variables described above. Statistics. Ecological variables were normalized [when needed (57)], centered, and scaled before PCA and the number of factors retained for subsequent analyses was informed by the Kaiser rule and parallel analysis (58). To explore the potential effects of cultural contact and spatial proximity, we included the average belief score (moralizing high gods = 1; no moralizing high gods = 0) among the 10 nearest neighbors for a given observation as a covariate in our model. In addition, we controlled for the nonindependence of societies that share a common cultural background by including language family as a random effect. Because language family classifications have changed substantially since the Ethnographic Atlas was first published, we used the updated version in ref. 59, as curated by linguists C. Bowern, Department of Linguistics, Yale University, New Haven, CT and M. Dunn, Max Planck Institute for Psycholinguistics, 6500 AH Nijmegen, The Netherlands. Although this is the best currently available estimate of dependencies because of shared cultural ancestry among human societies at a global scale, we note that cultural evolution need not map neatly onto the language family categories of the Ethnographic Atlas, and therefore caution against the possibility of measurement error when inferring conclusions about human history. All generalized linear models and generalized linear mixed models were fitted in R (60) using the lme4 (61) and MuMin (62) packages. We ran all possible models based on the alternative combinations of predictors in our candidate set and calculated the AICc for each model. Model weights were calculated from the AICc values and used to determine model average parameter estimates and unconditional SEs (40). To assess the performance of a given model, we estimated the coefficients from a randomly selected subset of two-thirds of societies in Dataset S1 and computed the AUC on our predictions for the missing one-third of our data (41). Predictive values in Table 3 were assessed by computing the AUC of models that included only the predictor or parameter of interest.

Acknowledgments We thank Q. Atkinson and K. Sterelny for comments, J. McCarter for help in the curation of the Ethnographic Atlas data, and two anonymous reviewers and our managing editor for helpful suggestions that greatly improved the clarity of our text. This work was supported in part by Grant/Cooperative Agreement G10AC00624 from the US Geological Survey (to C.A.B.), the New Zealand Marsden Fund (R.D.G.), the Templeton Foundation (28745), the Royal Society of New Zealand (11- UOA-239), and the National Evolutionary Synthesis Center (EF-0905606).

Footnotes Author contributions: C.A.B., J.B., and R.D.G. designed research; C.A.B. performed research; K.R.K. and M.C.G. contributed new reagents/analytic tools; C.A.B., B.G., K.R.K., J.B., M.C.G., and R.D.G. analyzed data; and C.A.B., J.B., and R.D.G. wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission. A.N. is a guest editor invited by the Editorial Board.

This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1408701111/-/DCSupplemental.