pathogen disgust, contamination avoidance, and hygiene behaviours that contribute to survival in the face of potential pollutants.
This simulated trait was allowed to affect up to 0.54% of yearly mortality, which corresponds to 10% of all deaths due to gastrointestinal illness in a hunter-gatherer population 25 . Mate-finding, reproduction, and the remaining 99.46% of mortality were governed by models and estimates from studies on hunter-gatherers. Inheritance could occur genetically, as the product of sexual reproduction, either with or without germ cell mutations. In addition, cultural inheritance was modelled as a simple Markov process, resulting in an undirected random walk (geometric Brownian motion) across generations 26, 27 . This model of a socially transmittable behavioural immune system is plausible, because disgust is similar between parents and their offspring 28, 29 , and there is empirical evidence of parental moderation of disgust in their children 30 .
Between the simulations (natural selection for the behavioural immune system) and their controls (mortality unrelated to the trait), cultural transmission ( Figure 1A ) was significantly different from year 3900 [Welch's t=2.00 -24.17, p=0.048 -p<1e -16] (Figure 2A , yellow line). By contrast, the polygenic phenotype ( Figure 1B ) had only evolved to be different from control from year 38400 [Welch's t=2. 25 -3.37 , p=0.027 -0.001]. Not only was the genetic trait slower to develop than the cultural, it also developed to a lesser magnitude ( Figure 1C) . Furthermore, the p values associated with the natural selection and control simulations were around significance threshold. One alternative approach was comparing each type of simulation against the trait's starting value of 0.75 using one-sample t-tests. This produced more reliable results: culture was significantly different from year 2100 [t(39)=-2.22 --41.92, p=0.032 -p<1e- 16] , and phenotype from year 25600 [t(39)=-2.22 --4.34, p=0.032 -9.8e-5]; whereas only three such differences occurred by chance in the control condition [phenotype in years 600-800, p=0.022 -0.042]. In sum, these results showed that cultural evolution is not only a plausible alternative to biological evolution, but that it can operate more quickly when a trait is beneficial for survival.
The inter-generational drift in cultural traits due to copying errors allowed it a wider variance than polygenic phenotypes. We thus also conducted simulations in which the cultural transmission was artificially constrained to match polygenic phenotype variance. In these simulations, phenotype was consistently different between the real and control conditions from year 7900 [Welch's t=2. 13 -4.99 , p=0.036 -3.72e-6], and culture from year 19500 [Welch's t=2.01 -4.21, p=0.048 -7.02e-5] (Figure 2A ). These results indicate that cultural evolution of the behavioural immune system is not due to its wider variance, as it also occurs when trait variance is artificially reduced to an unrealistically low value. 
While polygenic inheritance of the behavioural immune system is perhaps the most plausible based on other psychological traits (e.g. anxiety 31 ), its actual genetics are unknown. We thus simulated biological evolution for several genetic models in isolation, without cultural evolution occurring alongside. Populations started with equal numbers of high and low alleles, where the lower allele was more beneficial for survival. When inheritance was polygenic, the proportion of low alleles was significantly different from 0.5 from year 2500 [t (19) =2.43 -6.44, p=0.025 -3.55e-6] ( Figure 2B ). When inheritance was monogenic with dominance for the high allele, this proportion was significantly different from 0. 5 Figure 2B ). This suggests that either polygenic and monogenic inheritance could have shaped the behavioural immune system, regardless of whether the allele most beneficial to survival is dominant. However, monogenic adaptation comes at the cost of reduced genetic variance, with the high allele being practically eliminated after 60-80 thousand years in the no-mutation simulations ( Figure 1D ), and a clearer bias in the simulations that did incorporate mutations towards alleles that aid survival ( Figure 2C ).
Our agent-based simulation approach is not limited in applicability to pathogen avoidance, and can serve more widely as an example of how cultural evolution can outpace genetics under natural selection. The results presented here support the hypothesis that the uniquely high-strung qualities of the human behavioural immune system evolved as a "cognitive gadget" 7 , a mental mechanism passed down generations through cultural learning.
A total of 16 types of simulations were run to test a variety of hypotheses. Each simulation lasted 100 000 years, and was repeated 20 times. We manipulated whether natural selection occurred for the behavioural immune system (simulations 1-4 were controls in which it did not), culture and polygenic inheritance in tandem, only genetic inheritance, only cultural inheritance, different types of genetic inheritance, and restricted cultural inheritance so that its variance was equal to that of a polygenic phenotype. A full overview is provided in Table S1 . 
Each simulation started with 1000 individuals at year 0, with age distribution 25 and sex ratio (51.8% male) 32 based on research in hunter-gather communities. How cultural and genetic traits were initialised is outlined below. After initialisation, every simulated year progressed as follows:
• A proportion of fertile pairs produced offspring. The birth probability per pair was initialised at 0.292, estimated from a reported average of 7.0 children born to Savannah Pumé women aged 40 or over (see 32 , Table 1 , "completed fertility"). Each child inherited cultural and genetic traits from their parents via the processes described in further detail below. Sex ratio at birth was kept at 1.076 (slightly skewed towards male), based on adult sex ratios among Savannah Pumé (see 32 , Figure 1c ). This is slightly higher than in industrialised nations, for which birth sex ratios range from 1.031 to 1.067 33 .
• A proportion of the population died, based on an established mortality model 34 , using parameters estimates from hunter-gatherer societies Hadza, Ache (forest communities), Hiwi, !Kung, and Agta (see 25 , Table 2 , "average hunter-gatherer"). Mortality depended on age, and was highest for newborns and elderly (for details, see under "Mortality" below). Deceased individuals were removed from the population.
• Every individual's age was incremented by 1 year.
• The probability of individuals to reach fertility was computed and applied. The age of fertility onset was sampled from a normal distribution with mean 16.0 and standard deviation 2.5 for female, or mean 19.5 and standard deviation 3.4 for male individuals (based on Savannah Pumé, 32 , Table 1 , "age at first birth"). While the referenced population has their menarche and marriage earlier than 16.0 or 19.5 years, "fertility" here is very narrowly defined as "producing offspring at the current birth rate".
• The probability of fertile individuals to become infertile was computed and applied. Age of fertility offset was sampled from a normal distribution with mean 45.0 (similar to 35 , Figure  1 ; and a rough average of values reported by 36 and 37 , Table 1 , "OLR") and standard deviation 2.5 (matching that of fertility onset). For male individuals, fertility offset was sampled from a normal distribution with mean 55.0 (matches !Kung and Hadza, 37 , Table 1 , "OLR") and standard deviation 3.4 (again mimicking onset). The 10-year sex difference in mean fertility offset resembled the average difference in oldest recorded age at last reproduction among hunter-gatherers !Kung, Hadza, and Agta (see 37 , Table 1 , "OLR"). Infertility as a consequence of anything other than age was not modelled.
• Pairs were made from the pool of fertile and unpaired individuals, which included both both newly fertile individuals and those whose partner had died. Pairing was random, and stopped when no further pairs could be made. This resulted in a low number of unpaired individuals, reflecting high marriage rates among e.g. the !Kung 38 . Polygeny was not modelled, because of its low prevalence in many hunter-gatherer societies [37] [38] [39] .
• Newborns produced in the first step were now added to the population, starting at age 0.
• If the population size decreased this year, the birth probability per pair per year was increased by 0.002. If the population decreased or was over 5000, the birth probability decreased by 0.002. This step was introduced to prevent the population from exponentially growing or fatally collapsing under an unchanging birth rate.
At every age, but particularly in infancy and late life, humans face a likelihood of dying. Crucially, our interest was in how the behavioural immune system affects mortality. Examples of situations where this could occur include a foraging individual deciding to forego feeding on food that inspired disgust, to avoid close contact with bodily waste or a corpse due to implicit pathogen disgust, and more complex behaviours such as washing or cooking victuals. This is only a small proportion of total mortality, here defined as γ (Equation 1).
Gastrointestinal illness was the reported cause of 5.5% of deaths among Ache (forest communities) hunter-gatherers (see 25 , Table 5 , "All ages, gastrointestinal"). The behavioural immune system was allowed to affect 10% of these deaths, ultimately resulting in a γ parameter of 0.0054. (Note: that this is 0.0054 instead of 0.0055 is the result of a manual copying error in the simulation source code; but it was consistent across all runs.)
Where m(x) is the function that describes total mortality at age x; γ is the proportion of total mortality affected by the behavioural immune system, described by function f(x); and h(x) describes the probability of death at age x according to the Siler hazard function 34 :
Where α 1 is the initial infant mortality rate, β 1 the rate of infant mortality decline, α 2 the general mortality (unaffected by age), α 3 the initial adult mortality rate, and β 3 the increase in adult mortality rate. We employed parameters values that were estimated for hunter-gatherers Hadza, Ache (forest communities), Hiwi, !Kung, and Agta (see 25 , Table 2 , "Average hunter-gatherer"): α 1 =0.422, β 1 =1.131, α 2 =0.013, α 3 =1.47e-4, and β 3 =0.086. In practice, this resulted in a high infant mortality, a relatively low mortality throughout early adulthood, and an increase in mortality between four and eight decades of life.
The trait behavioural immune system is computed as the unweighted average of culture and phenotype if both are simulated, or either culture or phenotype if only one is simulated. It is made age-dependent, with a normally distributed onset (Equation 3). The assumption behind this gradual increase was that children are still minded during early childhood, meaning their own potentially preventable mortality is initially determined by their parents and not by their own actions.
(
Where τ is the trait behavioural immune system, and Φ(x) is a cumulative density function for a normal distribution with mean 5 and standard deviation 3 years , N(5,3) . These estimates were based on the notions that children in industrialised nations develop disgust between ages 3 and 5 years 4 , that Hadza children engage in food acquisition and processing from as early as 3.5 years 40 , and that children in hunter-gatherer societies function as "helpers at the nest" 41 .
The behavioural immune system was computed as the average of the cultural trait and the (genetically determined) phenotype, or fully determined by either in runs where only cultural or only genetic transmission was simulated. This combined trait (τ) was represented as a numerical value that weighted only the trait-related mortality (Equation 3), which was 0.54% of total mortality (γ, Equation 1). A value of 0 would eliminate all trait-related mortality, a value between 0 and 1 would reduce it, and values over 1 would increase it.
In control simulations, inheritance of cultural and/or genetic traits occurred, but was not allowed to affect mortality. This was achieved by setting the γ parameter (Equation 1) to 0, thereby ensuring that cultural and genetic inheritance occurred, but that they did not contribute to survival (and thus were not subject to natural selection).
It is likely that the behavioural immune system is at least in part a socially transmittable trait, because disgust is similar between parents and their offspring 28, 29 , and given empirical evidence of parental moderation of disgust in their children 30 . It would thus be subject to cultural inheritance, which occurs when children copy traits from the generation before them. Because they cannot copy perfectly, this inheritance is affected by a certain error rate. This process can be modelled as a simple Markov process, originally described as Equation 2 26 :
Where Y is the trait of interest, t is time (in generations), c is half the error rate, and N(0,1) is a normal distribution with mean 0 and standard deviation 1. It can be rewritten in the form of Equation 5 27 :
Where s is the culturally transmitted trait of interest in the population, t is time (in generations), and ε(t) is the normally distributed copying error rate with mean 0 and standard deviation σ, or N(0,σ) .
Note that children do not always copy their parents, but could instead conform to a social norm in the population. Such a bias can be incorporated in the model, resulting in Equation 6:
Where λ is the conformist bias, and s (t) is the population average for trait s in generation t. In practical terms, this means a proportion (λ) of children in the next generation will not inherit trait s from their parents, but will instead conform to the population average. When doing so, they will still be subject to copying error ε. This model has been validated for various culturally inherited traits, including the size of projectile points 26, 27 . Crucially, the error rate has been empirically determined to be 5% 27, 42 , and the conformist bias λ as 0.38 27 . These estimates come from the social transmission of technical knowledge, which are more likely to be copied more precisely and from a respected authority (e.g. an arrowhead maker instead of a parent) than the traits that constitute the behavioural immune system. Hence, they are likely conservative estimates for our purposes.
At initialisation, each individual's cultural behavioural immune system trait value was sampled from a normal distribution with mean 0.75, and standard deviation of 0.025 (half the error rate). In simulations with artificially constrained cultural evolution, the standard deviation was 0.0025 (half the restricted error rate), and the conformist bias to 0.10. This led to an artificial equalisation of variance between polygenic and cultural traits (phenotype at year 0: M=0.750, SD=3.32e-3; culture at year 0: M=0.750, SD=3.47e-3), and was used to test whether cultural evolution was dependent on its wider variance. At birth, an individual's cultural trait was randomly selected to be inherited from their mother or father (each with a probability of half 1-λ), or the population average (with probability λ).
Genetic inheritance of the behavioural immune system was simulated as the result of sexual reproduction. Each individual carried two alleles per locus, one passed down from the father and one from the mother. In reproduction, offspring received one randomly selected allele per locus from each parent. At year 0, the population was initialised with two alleles: 0.73 and 0.77. Each occurred with equal probability, and one was randomly selected for each locus in each individual.
An individual's phenotype was modelled in three different ways. It could emerge as a monogenic trait, with either the highest or the lowest allele being dominant. Alternatively, it could emerge as a polygenic trait, with phenotype computed as the average of all alleles on 10 loci.
Polygenic inheritance is perhaps the most likely option for the behavioural immune system, as it would resemble other psychological traits such as anxiety 31 .
Half the genetic simulations did not allow for any mutations to occur, whereas the other half allowed for germ cell mutations to occur when offspring was generated. The probability of such a mutation was set to 1.1e-8 per position per haploid genome 43 , which corresponds to roughly 70 de novo mutations in an individual's whole genome. Using the smallest gene (AURKB, ~6000 bases) implicated in polygenic anxiety 31 as a reference, the chance of at least 1 mutation occurring in an allele was set at (1-(1-1.1e-8 ) 6000 =) 6.6e-05. This is likely a high estimate, because de novo mutations are not equiprobable across the genome 44 . If a germ cell mutation occurred in an allele, its new value was randomly sampled from 11 variants (0.70, 0.71, …, 0.79, 0.80). Note that this resulted in an unchanged allele in 1 out of every 11 mutations (by randomly selecting the value that an allele already had before mutation).
Control conditions were incorporated to observe cultural and genetic inheritance without the pressure of natural selection. Side-by-side comparisons show that the cultural component of the behavioural immune system ( Figure S1 ) systematically drifted away from its starting value only when it was relevant to survival. The same was true for the polygenic component ( Figure S2) , particularly in the conditions where cultural transmission was artificially constrained. Table S2 shows the averages of a number of population characteristics across all runs for each simulation type. The only differences seem to have existed between control simulations and simulations in which the behavioural immune system did contribute to survival. In the control runs, the birthrate and number of children born to women aged 40 years and over were slightly lower, and the average age of those aged 5 years and over was slightly higher. No obvious differences in basic population characteristic existed between the non-control simulations. 
Inferential statistics in population simulations are somewhat arbitrary: because the whole population is known, its properties do not have to be inferred from a sample. However, one could consider a single simulation as a random sample from a wider population of stochastically varying simulations. Using this approach, we conducted two types of tests: Welch's t-tests for differences between the simulations with natural selection and their controls without, and one-sample t-tests for differences between simulated trait means and their starting value.
These tests were conducted every 100 years of the simulation. We define "consistent differences" as a consecutive period of 500 years in which statistically significant differences occur, using a threshold of α=0.05. The main text reports t and p values in the first year of a consistent difference, and at its maximum value. Furthermore, p values are reported in Figure 2A -B, allowing them to be compared against more stringent significance thresholds. 
We modelled the behavioural immune system as a dampening factor on mortality through gastrointestinal illness, but remain agnostic about the exact roles of specific emotions and behaviours. For example, one could model the contributions of pathogen disgust and hygiene independently. This would be a fruitful approach if these traits' exact contributions to survival were known, or could be realistically estimated. In the absence of such data from early homo sapiens, we think that our less granular approach with estimates from current hunter-gatherers was a good approximation. In addition, ours was a more conservative approach than introducing more degrees of freedom and uncertainty by guessing at specific processes.
We did not model potential sex differences in the behavioural immune system. While some report small sex differences in self-reported pathogen disgust sensitivity 45 , behavioural observations yield mixed results. Women are more likely than men to use soap in hand-washing in a public bathroom with an observer, and are more sensitive to external prompting 46 ; but no sex differences appear in oculomotor disgust avoidance 47 . It thus remains unclear whether observed sex differences are genuine, or reflect differential demand effects. Furthermore, if real underlying differences exist, it is unknown whether they emerge through genetics or socialisation. Hence, for the sake of simplicity, no sex differences were incorporated in this model of early human evolution.
In addition, our simulations did not include intersex individuals, non-heterosexual relationships, or offspring produced outside of an established pair. These were not included to retain relative simplicity, but also because it is unclear how they would impact or be impacted by the behavioural immune system.
To many researchers, the concept of disgust comprises more than just a behavioural extension to the immune system, and incorporates avoiding reminders of mortality, sexual selection, and even morality 11 . It is beyond the scope of this study to assess whether these functions actually build on core disgust, or simply share a superficial overlap grounded in common parlance ("disgusting" is a popular description for unsavoury individuals, actions, and ideas). However, the notion of disgust affecting sexual selection is relevant in an evolutionary context. Unfortunately, existing theories merely assert that disgust could be invoked by mates who display "low genetic compatibility and quality" 5 . From a computational point of view, disgust against incompatible individuals does not add explanatory value over sexual attraction for compatible individuals 48 . Furthermore, mate selection among hunter-gatherers requires economic viability 38 and is often arranged 39 , thereby rendering "disgust" for genetically incompatible mates potentially irrelevant. Due to the lack of computationally viable theories, and for simplicity's sake, a potential role for disgust in sexual selection was not modelled here.
We elected not to include cultural mutations. These could be conceived of as novel ideas and innovations that act as a cultural analogue to genetic mutations 49 . However, whether such inventions are truly new or simply recombinations of existing concepts is a matter of debate 49, 50 . Furthermore, the employed models of cultural evolution are parsimonious and well-validated 26, 27 , rendering their use without modification the most conservative approach.
Finally, our estimations for demographic variables come from contemporary hunter-gathers. While these populations are of a different culture to industrialised nations, they are modern humans, and not a perfect model for early homo sapiens.
Simulation code and data is available from https://github.com/esdalmaijer/2020_disgust_evolution