Different environmental variables predict body and brain size evolution in Homo

Body and brain size database

The fossil dataset consists of the hitherto largest collection of body (n = 204) and brain size estimates (n = 166) from Homo in the past ~1.0 Ma (Fig. 1). The data on hominin body size estimates are derived from our own previous study⁶ plus additional estimates⁶⁰ and updated chronometric ages from more recent literature. Individual body size estimates are provided by specimen in Supplementary Data 1 with data sources. The bulk of data on hominin brain sizes (endocranial volume, in cm³) is derived from recent meta-analyses^{7,12,13,61,62} and updated chronometric information. Specific sources of these data are indicated in Supplementary Data 2, with some assessments bearing larger errors due to the incomplete state of the crania on which they are based (e.g., Arago 21, Vértesszőlős, Zuttiyeh). Each body and brain size estimate is associated with information on estimated chronometric age (dating method and data source), geographical location (longitude and latitude), and taxonomic attribution. For the exact locations per specimen, see interactive map in Supplementary Note 1. We divided the dataset into three taxonomic units: Pleistocene H. sapiens, Neanderthals, and Mid-Pleistocene Homo. Whereas hypodigms of H. sapiens and Neanderthal remains are generally agreed upon, we use “Mid-Pleistocene Homo” as a strictly analytical unit to denote African and European Middle Pleistocene hominins that predate Neanderthals and are not assigned to Homo naledi, between ~800 and 130 ka. We refrained from further division of this group due to the often fragmentary nature of fossils, unclear alpha taxonomy, and small sample size. Analyses performed within these taxonomic units minimize phylogenetic effects of, e.g., significantly different brain sizes (e.g., ref. ²). Specimens from H. naledi and Homo floresiensis had to be excluded from this analysis as for each taxon they derive from a single location and age bracket, precluding assessment of paleoclimatic variation. Limitations to the fossil datasets (see e.g., refs. ^2,3,4,6) include imprecision of brain and body size estimates due to methodical and taphonomic problems, uncertainties of absolute ages that translate into uncertainties of the associated climate, and unequal sampling of hominin fossils across time and space. These limitations were incorporated into the construction of the synthetic dataset to assess the extent of their effects on the overall results for the actual fossil dataset. For all further analyses, brain and body size values were log-transformed as they increase multiplicatively.

Climate reconstructions

Each body and brain size estimate required corresponding estimates of relevant climatic variables. Our climate records are numerical model estimates based on global climate reconstructions for the past 1 Ma using the global climate model emulator GCMET²⁷. The main idea behind GCMET is that global climate model (GCM) simulations of the past 120,000 years contain sufficient information about long-term climatic changes on time scales of ≥1000 years. Given that we know the external boundary conditions, we can reconstruct previous glacial–interglacial climatic changes. The Quaternary climate is largely determined by dynamics of the Northern Hemisphere ice sheets, which, in turn, are affected by orbital variations of the Earth around the Sun and variations of atmospheric CO₂. Using these factors as external boundary conditions, GCMET can emulate the climate of the Quaternary in a similar way as a state-of-the-art GCM²⁷.

The atmospheric CO₂ record of the past 1 Ma that we use in this study is a composite of the EPICA CO₂ record from an Antarctic ice core⁶³ and of output from a carbon cycle model (CYCLOPS)⁶⁴. The EPICA record covers the past ~800 ka, whereas we use the CYCLOPS model output to cover the time up until 1.0 Ma. Orbital variations are based on calculations by Berger and Loutre⁶⁵. Ice-sheet extents for the past 800 ka are based on numerical ice-sheet model output⁶⁶. For the period before 800 ka, we assumed present-day ice-sheet configurations. This is an appropriate assumption given that all but one specimen of the fossil record before 800 ka in our datasets are within Africa or southeast Asia and thus far away from ice-sheet margins, with the local GCMET climate reconstructions not affected by this simplification.

For each fossil site location from the body and brain size database, we extracted a time series of the relevant climate variables, see Table 1 (also Supplementary Figs. 10 and 11). The time series were used to attach the value of each climate variable to the fossil record, both for the actual fossil data as well as for the synthetic fossil datasets.

Linear models

The null and two alternative linear models used throughout this manuscript are defined as follows. The null model simply estimates the mean for each taxonomic group, and we refer to this model as LM-T (linear model with taxa):

Here Y corresponds to either body or brain size (or the log-transformed thereof), whereas β₀ is the intercept, which is equivalent to the mean size of the reference taxon, and β₁ is a factor that reflects the deviation from this mean size for a taxonomic group (thus giving independent intercepts for Mid-Pleistocene Homo, Neanderthals, or Pleistocene H. sapiens).

The first alternative model contains the effect of the climate variable X (across all taxa):

Here β₀ and β₁ are the intercept terms, giving taxon-specific values, and β₂ is the slope, which is the same across all taxa. We refer to this model as LM-TC (linear model with taxonomic differences plus a climate effect).

The second alternative model takes taxonomic differences for the slope of the climate effect into account. This is done via an interaction term, β₃, which acts as a modifier for the slopes (i.e., different intercepts, given by β₀ and β₁, and slopes, given by β₂ and β₃, for each taxonomic group):

We refer to this model as LM-T*C (linear model with taxonomic differences plus a taxon-specific climate effect). The slopes, β₂ and β₃ in Eqs. (2) and (3), respectively, are presented in the main text in Tables 2 and 3 (for the log-transformed sizes) and in Supplementary Tables 1 and 2 (for the natural units of the sizes).

Synthetic datasets and power analysis

Apart from determining the smallest sample size suitable to detect the effect of a given test at the desired level of significance, power analysis can also be used as a formal way to test whether a relationship between dependent and independent variables can be detected with the available data and proposed methods (i.e., linear models in our case) assuming that such a relationship exists. Before testing for any true association between local climate and the fossil record, we use such a power analysis to assess our power to detect relationships of different effect sizes given the uncertainties, for example, in body/brain sizes, dating, and climate reconstructions. We generated 1000 synthetic datasets for each of the ten climate variable associations (MAP, MAT, NPP, mean temperature of coldest quarter, mean precipitation of driest quarter, and the logarithm of their running standard deviation over a 10,000-year window) with body and brain size. For each association, we assumed a strong, a medium, and a weak relationship between size and climate.

By strong, we refer to 1/4 of the maximum possible slope given by the range of the climate and size. Subsequently, medium is half the slope of the strong relationship, (1/8 maximum possible slope), and weak is half the slope of the medium relationship (1/16 maximum possible slope). For example, the strong association between MAT and body size (Bergmann’s rule) is −0.34 kg/°C, based on the above defined rule. This is close to the estimated association between temperature and body size of about −0.4 kg/°C found for modern humans in a recent study (ref. ³¹, their Fig. 5A). Unfortunately, there are no empirical data about other climatic relationships and body (or brain) size. For simplicity, we therefore applied the same rule of strong, medium, and weak associations for all other climate variables and for brain size.

Before generating a synthetic dataset, we estimated the intercepts β₀ and β₁ and the slope β₂ for the LM-TC model, Eq. (2). However, for the real fossil analysis we used the model with the interaction term, LM-T*C, Eq. (3). First, we looked up the climate record for each location and time from the climate time series and attached it to the respective empirical fossil records. We calculate the maximum slope from the X and Y ranges as β₁ = range(Y)/range(X). Assigning an actual relationship factor, e.g., strong (=1/4), the intercept β₀ can be calculated using the X– and Y-midpoints, β₀ = Y_midpoint − 1/4β₁X_midpoint.

For the synthetic fossil datasets, we assume an age uncertainty range of 10% (±5%) for radiocarbon-dated fossils, i.e., younger than 50 ka cal BP (e.g. ref. ⁶⁷), and 20% (±10%) for fossils older than 50 ka coming from other dating methods with higher uncertainty such as luminescence, U-series, or ESR (e.g., ref. ⁶⁸). Furthermore, we assume a standard error of 2 K for mean annual and mean temperature of coldest quarter. For all other climate variables, we assume a 20% error range (±10%). The 2 K and the 20% are in line with climate model biases as estimated in a recent study⁶⁹. Within a taxonomic unit of the genus Homo, we assume a coefficient of variation (CV) of 7% for body size (average of intrapopulation means of 19 global Holocene hunter-gatherer populations, n = 510, data from JTS) and 3.5% for brain size (from ref. ²⁸ populations, dataset: http://volweb.utk.edu/~auerbach/HOWL.htm; ref. ⁷⁰, see also ref. ³²). Previous research has demonstrated that the range of body size variation in Holocene human populations is larger than any taxonomic unit of earlier hominin and encompasses the range of variation found within earlier hominins⁶ and that sexual dimorphism in size among Mid-Pleistocene hominins is comparable to that of modern humans⁷¹. While there are significant differences in brain shape through recent hominin evolution, the range of size variation within Pleistocene hominin taxa remains comparable to that observed among modern humans⁶⁰. These observations suggest that modeling the intrapopulation variation among hominin taxa upon modern human coefficients of variation provides a reasonable estimate of variation within hominin taxa that are often presented only by much smaller sample sizes. To create a synthetic dataset that has a mean and a variance as close to the fossil dataset, we introduced taxonomic size differences (β₁ in Eq. (2)) that is based on the taxonomic differences in the mean size. This difference was estimated directly from the fossil dataset.

The procedure to generate a single synthetic dataset is as follows. First, we selected a relationship strength, e.g., strong, and calculated the slope and intercepts. For each synthetic data point, we:

1. 1.

   Looked up the age and added a randomly sampled error (±5% or ±10%).

2. 2.

   Looked up the fossil site and selected the climate record from the previously calculated time series for that location and sampled age.

3. 3.

   Added a randomly sampled error (i.e., S.D. of 2 K or 20%) to the climate record. This is now the X value.

4. 4.

   Multiplied X with the slope β₂ and added the intercept β₁ with the respective taxonomic correction. This translates the climate record X into a size estimate Y.

5. 5.

   Added a random term to Y based on the CV, i.e., 3.5% for brain and 7% for body size.

6. 6.

   Repeated steps (1)–(5) for each fossil record and saved all locations, ages, Xs, and Ys to a file. This is a single synthetic dataset in the same format as the original fossil dataset.

We repeated this N times to generate N synthetic datasets and repeated the same procedure for the other relationship strengths, i.e., medium, and weak and for all other climate variables. Panels of exemplary synthetic datasets for body and brain sizes in comparison with the original data are shown in Supplementary Figs. 10 and 11.

We use the same thinning approach as described in the main text (n = 1000) for the synthetic datasets. These are then used for a power analysis to test whether the linear relationship between any climate variable and body or brain size can be detected. We fitted both the LM-TC model, Eq. (2) (in which the slope defining the relationship between climate and size is the same for the three taxonomic groups, which can differ in their intercept), and the LM-T*C model, Eq. (3) (different slopes and intercepts for the three groups). A climate effect was deemed present if the null model had a higher AIC value compared to either of the alternative models, LM-TC or LM-T*C, (ΔAIC > 2), i.e., LM-T, Eq. (1), in which the three groups differ in size but there is no effect of climate. Figure 2 in the main text shows the power to detect a true relationship between size and climate. Individual records are color-coded according to the AIC difference between the LM-T and the alternative models, LM-TC or LM-T*C, ranging from −2 (red) to +15 (blue) with 2 as midpoint (white).

All statistical tests were undertaken in Python version 3.8.5 using the following Python packages: statsmodels 0.12 (for linear models), pandas 1.1.3 (for dataframes, reading/writing CSV/Excel files), netCDF4 1.5.3 (reading NetCDF files), matplotlib 3.3.2 (for plotting), and numpy 1.19.2 (numerics).

Reporting summary

Further information on research design is available in the Nature Research Reporting Summary linked to this article.

Source: Ecology - nature.com