Human social complexity was significantly lower during climate cooling events of the past 10 millennia

Human civilizations depend on the climate. Changes in climate affect the production of food and other resources that support populations and their economies. We asked whether the millennium-scale climate cooling events identified by Gerard Bond predicted social complexity in the Seshat cross-cultural database. The results show that social complexity was significantly lower during the coldest two centuries of Bond cooling events. Reductions in complexity are evident in regions north of the tropics adjacent to the Atlantic or Arctic, particularly in North Africa, Europe, and Central Eurasia.


Introduction
Like all living systems, human societies depend fundamentally on metabolism, the uptake and transformation of energy that supports a population and its activities. The food available to a population relies on the productivity of the plant and animal food species that it harvests. The productivity of these species depends on human effort, ecology, and climate, particularly in terms of temperature and precipitation (Tainter 1988;Brown et al. 2011;Burger et al. 2012;Hooper et al. 2014;Putnam et al. 2016).
Climate cooling events that occurred during the Holocene on a roughly millennial timescale were identified by Bond et al. (1997Bond et al. ( , 2001 on the basis of drift ice indices in the North Atlantic. We hypothesized that these climate events would have had a discernable impact on the extent of social complexity in human societies over time. Studies of specific regions and time periods indicate links between climate and social organization; the collapse of the Akkadian Empire and Egyptian Old Kingdom during the 4.2-kiloyear event are notable examples (Cookson et al. 2019; see Table 1). However, it has remained difficult to test relationships between climate and social variables in a global longitudinal sample. The Seshat dataset published by Turchin et al. (2015Turchin et al. ( , 2018)-a sample of 30 locations in 10 world regions spanning from 9600 BCE to 1900 CE-now allows such a test.
We asked whether the cooling events reflected in North Atlantic drift ice indices (Bond et al. 1997(Bond et al. , 2001 predict the extent of social complexity as quantified by Turchin et al. (2018). Table 1 lists the seven Bond events that overlap in time with the Seshat data. We analyzed the principal component of social complexity (PC1) calculated by Turchin et al. (2018). PC1 reflects a constellation of traits associated with social complexity, including population size, government, religious hierarchy, writing, currency, and urbanization.

Methods
We predicted that social complexity would be lowest in the last centuries of a cooling event, the period of maximum cold before productivity begins to improve again due to warming. To test this prediction, we employed a multilevel regression model estimating Turchin et al.'s (2018)  The Turchin et al. (2018) dataset includes a total of 414 observations in 30 locations (a.k.a. natural geographic areas) across 10 world regions (North America, South America, Africa, Europe, Central Eurasia, East Asia, Southwest Asia, South Asia, Southeast Asia, and Australia-Oceania). The regression model included a random effect for location to represent stable differences in complexity across locations. A 0-1 dummy variable was used to demarcate locations in the Old World (Africa, Europe, and Asia) to represent the distinctive sociopolitical trajectory of this macro-region .
To control for long-term trends in social complexity, the regression included terms for year of observation and location "age" (years from earliest observation).
Year and age were both divided by 1,000, and can therefore be interpreted in units of millennia (e.g., 1000 BCE = -1.0). Linear and quadratic terms for year and age, and the interaction between year and age were considered for inclusion in best-fit models. AIC criteria favored including a linear term for year of observation, linear and quadratic terms for age of location, and a year ⨉ age interaction term ( Table  2).
The effect of Bond cooling events was represented by a 0-1 dummy variable set to 1 for observations occurring during the last two centuries of a cooling event (i.e., the 200 years prior to and including the dates in Table 1). Understanding that the interpretation of Bond event 0 around 1500 CE may be controversial, models were fit with both the full sample of observations (N = 414) and the subsample of observations that excludes the last millennium (≤1000 CE; N = 284).
We first tested the hypothesis that cooling events are associated with lower social complexity in the global sample (Table 2), then examined estimates for specific regions. Because Bond events were described in the context of the North Atlantic, we tested whether locations north of the tropics and in one of the six regions adjacent to the Atlantic or Arctic-North America, Europe, North Africa, Central Eurasia, East Asia, and Southwest Asia-exhibit stronger signals of cooling events compared to other regions. These models included fixed effects for region and cooling ⨉ region interaction terms. We estimated the effects of cooling for northern regions as a group (Table 3), then for each region individually (Table 4). The rows corresponding to the hypothesized effects of cooling in Tables 2-4 are highlighted in gray.
Models were estimated using the lmer function in R (Bates et al. 2015; R Core Team 2021). Best-fit models were evaluated based on AIC minimization. Statistical significance was judged based on 95% confidence intervals (CI). The code and dataset used in the analysis can be downloaded from github.com/systemsscience/climate. The original data are available from the Seshat Databank (seshatdatabank.info).

Results
The analyses show that social complexity is detectably lower during the last two centuries of Bond cooling events. The model in Table 2A indicates that the social complexity measure PC1 is 0.28 (0.03-0.53 CI) units lower during cooling periods relative to other centuries. This effect is evident independent of the Little Ice Age: the model in Table 2B, which only includes observations prior to last millennium, also estimates a significant reduction in complexity (0.38 units; 0.06-0.71 CI) during cooling periods. The time controls in Table 2 suggest that-all else equal-there is a roughly linear increase in complexity over historical time, and that complexity increases at a decreasing rate with the age of the location. There is also a positive interaction between time and location age. Around the mean observation year (120 BCE) and location age (3,545 yrs), the time terms in Table 2A suggests that PC1 increases at a rate of 0.001409 units per year, or 1.409 units per millennium. This number provides a sense of scale for the other parameter estimates: the estimated effect of cooling for the full sample (0.28) equates to a loss of 197 "equivalent years" of accumulated complexity; the effect for the subsample ≤ 1000 CE (0.38) equates to a loss of 273 years. The geographic controls indicate that there is significant heterogeneity in complexity across locations, and that mean complexity is higher in Old World locations.
The models in Table 3 show that reductions in complexity occur specifically in regions north of the tropics adjacent to the Atlantic and Arctic Oceans. Grouped together, these northern regions are estimated to experience an average loss of 0.43 (0.14-0.72 CI) units of complexity during periods of cooling (Table 3A); this equates to a loss of 306 years of accumulated complexity. Prior to 1000 CE, the estimated loss associated with cooling in northern regions is 0.47 (0.11-0.83 CI), or 331 equivalent years (Table 3B). The estimate for the effect of cooling in southern regions is indistinguishable from zero in both models. Table 4 breaks down losses in complexity north of the tropics by region. In the full sample (Table 4A), the estimated effects of cooling are greatest in North Africa (0.91 ~ 642 equivalent yrs) and North America (0.89 ~ 631 yrs), followed by Europe (0.54 ~ 386 yrs), Central Eurasia (0.54 ~ 380 yrs), and East Asia (0.43 ~ 308 yrs). The estimate for the effect of cooling in North Africa is statistically significant. Unexpectedly, Southwest Asia does not follow the pattern of the other northern regions, showing a slightly positive estimate indistinguishable from zero.
In the ≤1000 CE subsample (Table 4B), the estimated losses during cooling events are greatest in Central Eurasia (1.33 ~ 941 equivalent yrs) and Europe (0.86 ~ 612 yrs), followed by North America (0.62 ~ 440 yrs), North Africa (0.50 ~ 356 yrs), and East Asia (0.45 ~ 317 yrs). The estimates for the effects of cooling in Central Eurasia and Europe are statistically significant in this earlier subsample.

Discussion
These results suggest a fundamental relationship between climate and the organization of human societies. Societies north of the tropics near the North Atlantic and Arctic exhibit lower social complexity during the most punishing years of cold events that have occurred every millennium or so throughout the Holocene. This linkage between climate and human civilization will remain salient as the climate continues to change (Brown et al. 2011;Xu et al. 2020).
These global results complement finer-scale analyses of regional trajectories and time periods, such as the population trajectories of the Sahara (Manning and Timpson 2014) and China (Wang et al. 2014) across the Holocene; the collapse of the Akkadian Empire (Cookson et al. 2019), Old Egyptian Kingdom (Krom et al. 2002) and Neolithic civilizations in China (Liu and Feng 2012;Li et al. 2018) around 2200 BCE; and the rise and fall of the Classic Maya (Haug et al. 2003;Kennett et al. 2012) and Puebloan polities (Schwindt et al. 2016;Crabtree et al. 2017) before 1500 CE.
It will be important to understand the proximate pathways by which sunlight, temperature, precipitation, seasonality, and sea levels affect energetic and cultural productivity, and consequentially social organization. Nomadic pastoralists, for example, appear to fare relatively better than agriculturalists during cooling periods, an asymmetry that may have precipitated the nomadic expansions that occurred around 200-500 CE and 1200-1500 CE (Pei and Zhang 2014;Putnam et al. 2016). Past work has also emphasized the destabilizing effects of aridity that commonly accompanies cold periods (deMenocal 2001;Li et al. 2007). In some environments, aridity may also favor more social hierarchy by increasing environmental circumscription (Kennett and Kennett 2006;Hooper et al. 2018).
This analysis represents an initial foray into understanding relationships between climate change and societal variables in the Seshat Databank. Future work with this and other datasets should leverage a variety of climate proxies to evaluate the effects of climate on multiple dimensions of social organization (Guedes et al. 2016;Peregrine 2020). While the Turchin et al. (2018) dataset is sizeable by crosscultural standards, the power of the current analysis remains limited by sample size. More comprehensive tests will be possible when high-resolution longitudinal studies have been integrated into larger comparative datasets.