Evolutionary analysis of human prosociality: an approach based on the study of HPA axis related genes in informative ancient populations

1. Introduction

From the early 20th century, research has been conducted on human evolution and comparative neuroscience[1],[2]. One of the propositions that have attained the greatest consensus is the idea that primates developed large brains to face their unusually complex social life, which is the basis for the social brain hypothesis[3],[4]. Although there are different interpretations of this hypothesis, all of them state that the increase in selective pressures has favored the development of a brain with greater social and cognitive skills in H. sapiens[5],[6]. One of the arguments supporting this interpretation is the increase of human population throughout evolution, which forced individuals to coexist with a larger number of people, thereby increasing the complexity of social interactions and favoring the development of a more advanced brain in social and cognitive terms[7]-[9].

Similar propositions about the behavioral evolution of humans claim that one of the selective pressures that shaped H. sapiens was the decrease of certain aggressive behaviors[10], which may have led our species to the development of more prosocial behaviors (or vice versa), with better sociocognitive skills[11]-[14]. A relevant explanation for this decrease of aggressiveness is the hypothesis of human self-domestication (HSD), which postulates that the selection against aggression operated in human evolution, producing secondary effects (morphological, physiological and behavioral) similar to those characteristic of domestic animals[11]. According to this hypothesis, the more social the individual, the greater the reproductive advantage in complex social environments, thus favoring the selection of prosocial traits.

This HSD hypothesis, although attractive and supported by several lines of evidence[11], continues to be under debate[15],[16], with authors underlining the importance of other factors in the behavioral evolution of H. sapiens in addition to prosociality, such as climate, ecology, food, and other factors associated with competition and survival[17]-[19]. However, regardless of the involvement or not of these or other factors, the selective pressure on genes that regulate stressful and/or aggressive behavioral phenotypes throughout the evolution of H. sapiens is undeniable[20]. The extensive search of these genes has revealed that many of them belong to the hypothalamic-pituitary-adrenal axis (HPA)[20]-[23]. This is not surprising, as this axis regulates human behavior in the face of different environmental stressful stimuli, not only modulating the release of cortisol, the main physiological stress hormone (see Figure 1), but also interacting with oxytocin, hormone related to social behaviors and attachment[24]-[28]. The involvement of both hormones, cortisol and oxytocin, in this axis enables, among other pathways (see some of them in Figure 3), the regulation of human behavior, emotions and cognitive function[23],[29]-[36].

Figure 1. Functioning of the hypothalamic-pituitary-adrenal axis (HPA)

Reception of stressful stimuli from the environment activates, in various ways, the paraventricular nucleus (PVN) of the hypothalamus, which secretes CRH (Corticotropin Releasing Hormone). CRH is taken by the pituitary gland (or hypophysis), where it binds its receptors and triggers the production of ACTH (Adrenocorticotropic hormone). ACTH travels to the adrenal glands, where it binds its receptors, activating the secretion of glucocorticoids (GC), such as cortisol. These GC produce the physiological effect of stress. In turn, to ensure that the axis is not over-activated, GC induce a negative feedback on the pituitary gland and hypothalamus. Source: [29].

1.1. Aims of the present study

In this context, our hypothesis is that, during the Upper Palaeolithic, humans acquired complex sociocognitive skills, which may be in part due to the selection of genes involved in the regulation of the HPA axis, enabling the development of relevant prosocial behaviors, which enhanced the biological and social adaptation of H. sapiens. These Upper Palaeolithic hunter-gatherers underwent ecological, technological and cultural changes during the Last Glacial Maximum (LGM) (between 20,000 and 22,000 years ago[37]), reflected in the notorious raise in portable and parietal art in the Magdalenian the last period of the Upper Palaeolithic[38],[39]. Adverse climatic conditions may also affected our species relationship with other animals as wolves, whose domestication occurred as early as the Late Magdalenian, with the identification of one of the oldest dogs known to date in the Erralla site (Gipuzkoa, Basque Country), which dates back ~17,000 years[40].

1.2. Sample and gene selection

In this study, the possible genetic basis of prosocial behaviors (and greater cognition), which may be present in H. sapiens since the beginning of the Upper Paleolithic, was explored by analysing the nuclear genomes available in s corresponding to European and Asian individuals from ~47,000-1,000 BP. These individuals were grouped into two chronological periods: i) those from the Upper Paleolithic (including the oldest individuals analysed in this study, from the Magdalenian and earlier periods), who led the domestication of the wolf; and ii) those from the Epipaleolithic and/or later chronologies. The bioinformatic analysis was focused on a few genes in which single-nucleotide polymorphism mutations (SNPs) have been described with phenotypes associated with prosocial and cognitive behaviors. Among the 7 studied genes, the following are worth highlighting: OXTR (Oxytocin Receptor), AVPR1B (Arginine Vasopressin Receptor 1B) and GABRA6 (Gamma-Aminobutyric Acid type A Receptor subunit Alpha 6), which play a direct and dual role in the regulation of the HPA axis, regulating the release of the main hormones (CRH and ACTH, see Figure 1) that activate the axis and, therefore, the behavioral responses to stress; and OPRM1 (Opioid Receptor µ 1) and COMT (Catechol-O-Methyltransferase), which indirectly and dually modulate the activation of the HPA axis, and thus the presence of catecholamines in the system, favoring the release of the main hormones that activate the HPA axis (see Figure 3 and Discussion for further details).

2. Materials and methods

2.1 Materials: Online s and bioinformatic tools

For this study, we used nuclear genomes of Homo sapiens corresponding to individuals of Europe and Asia from ~47,000-1,000 BP, which can be found in the online European National Archive [ENA, ENA Browser (ebi.ac.uk)]. Selected SNPs were analysed after a search in s of scientific articles, such as PubMed (Advanced Search Results - PubMed (nih.gov)), Connected Papers (Connected Papers | Find and explore academic papers), dbSNP of NCBI (Home - SNP - NCBI (nih.gov)) and ENSEMBL (Ensembl genome browser 111), for the determination of allele frequencies and their associated phenotypes.

For the treatment of the genomic data, we used the servers of the Human Evolutionary Biology research group (HP Proliant DL580 GEN8 and Dell Precision Tower 7920) at the University of the Basque Country/Euskal Herriko Unibertsitatea (UPV/EHU).

2.2. Methods

2.2.1. Bibliographic search

We carried out a bibliographic search of SNPs related to cognition, prosociality and empathy, building a . The resources used are, in addition to the s mentioned above, the selected publications in Medline and specialized on line libraries, where keywords, such as: prosociality, prosocial behavior, oxytocin, HPA, domestication, etc. pronted ~250 manuscripts, which were filtered considering the abstract content. Finally, 80 papers were reviewed, eliminating those focused on clinical aspects. The most informative papers selected for the objectives of this study are cited in the reference section.

This research, led to the identification of other behavioral traits that might be related to those previously mentioned, such as the behaviors included within Autism Spectrum Disorder (ASD), since ASD behaviors have also been linked to some of the selected SNPs related to prosocial skills depending on the allele that is present in those positions (see Table 5). The selected genes and their respective SNPs are shown in Table 1 and 5. The latter describes the phenotypes associated with the different allelic variants of the studied SNPs.

Table 1. Genes and SNPs selected for the present study

Moreover, a selection of nuclear genomes of Homo sapiens was performed (from the Mousterian to late Antiquity), generating another with all the information associated with them. The 92 ancient individuals selected are listed in Table 2.

Table 2. Ancient Homo sapiens selected for the present study

2.2.2. Bioinformatic analysis

Through the ENA , the genomic data of the 92 individuals selected (Table 2) were downloaded and subjected to bioinformatic analysis. To this end, we followed a process of: i) quality control, ii) alignment and editing, iii) search of selected variants (SNPs) (Table 1), and iv) recording of the genetic information. The bioinformatic work flow is detailed in Table 3. The reference human genome used was “Genome assembly GRCh38.p14” [Homo sapiens genome assembly GRCh38.p14 - NCBI - NLM (nih.gov)[52]].

Table 3. Followed work flow for the bioinformatic analysis

The genetic data of the SNPs of interest were recorded in another , indicating the total number of reads and the number of reads of each allele (nucleotides) per individual and position.

Table 4. Individuals selected after applying the quality criteria

With the aim of increasing the reliability of the study and preventing the overestimation of the presence/absence of the SNPs of interest, only those positions with 5 or more reads were considered. This selection reduced the number of individuals selected down to 20 (Table 4). A total of 14 SNPs were analysed, corresponding to 7 genes (Table 5).

Table 5. SNPs analysed in the current study and their associated phenotype

In the second column, in addition to the identifier of the SNP, the following is indicated in brackets: “chromosome:nucleotide position, reference allele > alternative allele”.

2.2.3. Limitation of the analysis

Due to the scarce preservation of ancient DNA (which is one of its main characteristics, especially at the nuclear genome level; de la Rúa, C. & Hervella, M, 2013), and the screening of the recorded genetic information (raw data) based on the established quality criteria, the number of individuals included in this study was very small. Similarly, the set of individuals that met said quality criteria present a large amount of missing data. This scarcity of data implies a limitation when taking into consideration possible demographic processes that affect the continuity or not of a specific human genetic lineage, as with so few individuals and so many missing values they could go unnoticed. This condition prevented the calculation of the allele frequencies in most of the analysed SNPs; therefore, it was necessary to work with qualitative data of presence and/or absence of these allele variants. All this hindered the implementation of statistical analyses, such as the imputation of missing values and the use of multivariate statistics that allow managing the missing values. Likewise, we do not know with certainty the proportion of the variation that genetics can explain in the behavioral traits analysed in this study. Thus, it would be interesting to carry out a future study that includes not only a larger sample size, but also the calculation of the weight that genetics implies in the whole variation and expression of this trait (that also includes environment) through the f statistic (central tool in the analysis of variance, ANOVA).

2.2.4. Heatmap

Given the limitations of the data managed in this study, it was decided to create a heatmap that allow visualising and understanding the data of the gathered SNPs in a simple and qualitative manner. To this end, the genetic data were coded as follows: i) missing values = 0 (represented with white colour in the heatmap); ii) presence of the alternative allele in a percentage lower than 26% of the total number of reads per individual and position –defined as the absence of the alternative allele– = 1 (represented with blue colour in the heatmap); and iii) presence of the alternative allele in a percentage equal to or higher than 26% of the total number of reads per individual and position –defined as the presence of the alternative allele– = 2 (represented in red colour in the heatmap). To establish the genotypes of the analysed individuals (Tables 7-9), the following were considered: i) homozygous for the reference allele is present at ≤26%; ii) homozygous when the alternative allele is present at ≥76%; and iii) heterozygous when the alternative allele is present at 27-75%. For the elaboration of the heatmap, statistical software R was used, through the R-Studio interface, employing the basic statistical packages R (“base”) and “gplots”. Once the heatmap was completed, the individuals were divided into two different chronological groups: the individuals of the Upper Paleolithic (the oldest individuals analysed in this study, from the Magdalenian and earlier periods) are represented in green, and the individuals of the Epipaleolithic and/or later chronologies are represented in blue.

2.2.5. Fisher’s exact test

After applying the quality criteria, when a SNP presented 10 or more individuals who did show reads in its nucleotide position (SNP3, SNP8 and SNP9 of the OXTR gene, SNP1 of the OPRM1 gene, and SNP1 of the COMT gene), Fisher’s exact test was performed. To execute this test, R was used, through the R-Studio interface (“base” statistical packages). This statistic allows determining the existence of significant differences between the frequency of the alternative allele and that of the reference allele in the two different chronologies: individuals of the Upper Paleolithic (the oldest individuals analysed in this work) and those of the Epipaleolithic and later chronologies. To this end, the number of individuals carrying the alternative allele at a percentage of ≤26% (assigned as genotype 0, absence of the alternative allele) on the one hand, and at a percentage of >26% (assigned as genotype 1, presence of alternative allele in homozygosis and heterozygosis) in the other hand, was counted in both chronological groups for each SNP. Next, we present the hypotheses from which the p-value and the Odds ratio were calculated under a significance level of 5% (or 95% confidence interval):

• Alternative hypothesis or H₁: there are significant differences in the frequency of the alternative allele between the two groups of individuals.

Fisher’s test is particularly useful for small sample sizes, as it does not rely on specific distributional assumptions or minimum sample requirements. Like other statistical tests, it assumes independence of observations within each category. While kinship structure and shared evolutionary history can challenge this assumption in genetic studies, our dataset—comprising individuals from diverse geographic regions and distant chronological periods—minimizes the risk of strong genetic correlations. This broad distribution of samples ensures that Fisher’s test remains a suitable and robust tool given the limitations of our study.

3. Results and discussion

3.1. Results and discussion of the heatmap

As it can be seen in Figure 2, there are different scenarios. Firstly, the alternative alleles of SNP1 of AVPR1B gene, SNP1 of OXTR gene, and SNP1 of OPRM1 gene are not present in prehistoric humans from the Upper Palaeolithic onwards considered in this analysis, which could mean that prosocial skills related to these SNPs may not have been present during the studied time lapse. Secondly, the alternative alleles of SNP3 and SNP6 of OXTR gene, SNP1 of GABRA6 and SNP1 of COMT, that have been proposed to be related to prosocial behaviors, are present in most of the Upper Paleolithic hunter-gatherers considered in this study, and it seems they could be preserved, with slight variations, throughout late prehistoric times. Thirdly, there are other SNPs analysed in this work that seem to have different frequencies in the individuals from the Upper Paleolithic compared to those from the Epipaleolithic and later prehistoric periods, which could be the case of SNP8 and SNP9 of OXTR gene. Finally, SNP1 of CD38 gene and SNP1 of FKBP5 gene are only represented in a few of the individuals analysed in this study, which makes it difficult to interpret the meaning of the presence/absence of these SNPs in prehistoric humans. These results are described and discussed in detail below, for each of the genes analysed.

Figure 2. Heatmap of the presence/absence of the SNPs selected in this study

Presence of the alternative allele (red): when the alternative allele is present at ≥ 26%. Absence of the alternative allele (blue): when the alternative allele is present at <26%. The rest of the colours correspond to the absence of data or a number of reads below 5 for each SNP. The ancient humans analysed are shown in chronological order, from greater (Konstenki14) to lower (I10866) antiquity, grouped into two groups: green >~15,000 BP, and blue <~11,000 BP.

3.1.1. AVPR1B gene

Firstly, regarding the AVPR1B gene, it was observed that the alternative allele C of the SNP (rs28373064) analysed in this gene –which has been related to a possible role in emotional empathy and in prosocial behaviors24 (see Table 5 and Figure 3)– is not present in any of the analysed individuals that did have reads at this position (Figure 2). In fact, the thorough analysis of the reads of each of these individuals (Table 6) show that the alternative allele C was not detected in any of the reads of the analysed humans. These data suggest a low frequency or the absence of this allele in the prehistoric chronologies, at least in the ancient individuals analysed in this study. Although the frequency of the alternative allele in the current global population is also low (information obtained from dbSNP), it is important to highlight that, from the analysed sample, only 5 individuals did have reads at this position in their nuclear genome. It is necessary to expand the sample size in order to obtain a broader and more accurate view of the presence or absence of the alternative allele of this SNP.

Table 6. Information of the homozygosity and heterozygosity of the analysed AVPR1B gene’s SNP1

The number of reads presented by the analysed individuals for this position is shown in parenthesis, next, the number of times that each of the alleles appeared, indicating in red alternative allele and in blue the reference allele. Homozygosity of the reference and alternative allele are shown in blue and red boxes, respectively. Heterozygosity is indicated in yellow boxes.

3.1.2. OXTR gene

Secondly, there is the OXTR gene, which encodes the receptor of oxytocin, a hormone that plays a key role in the development of selective social bonds[64],[65] (see Figure 3). With regard to the analysed SNPs of the OXTR gene, different scenarios were detected (Figure 2). The alternative allele A of SNP1 (rs59190448), possibly associated with prosocial behaviors[25], is absent in all the analysed individuals that did have reads at this position, being detected in none of their reads (Table 7), although only 4 individuals presented reads at this position. Given that this alternative allele is also found in low frequency in the current global population (information obtained from dbSNP), this result suggests the same comment as in the case of the SNP rs28373064 of the AVPR1B gene. A similar scenario was presented by the haplotype constituted by SNP3 (rs237887), SNP8 (rs1042778) and SNP9 (rs2268490) from the same OXTR gene, which present reads at these positions in many of the analysed individuals (Figure 2) (SNP3, n= 10; SNP8, n= 10; SNP9, n=14). The combination of the reference alleles of these three SNPs (GGC, respectively) has been associated with possible tendencies toward altruistic behaviors[26]. The heatmap shows that the alternative allele is the most represented in SNPs 3 and 8 in all the analysed chronologies (except for SNP8 in Chan_Meso and SNP3 in I11614, which presented the reference allele in homozygosis), being present in all individuals from over ~14,000 BP, in both homozygosis (especially in the case of SNP8) and heterozygosis (especially in the case of SNP3) (Figure 2, Table 7). In the case of SNP9, the alternative allele, as in the case of SNP3 and SNP8, is present in most of the ancient individuals from the Upper Paleolithic (Figure 2), especially in heterozygosis (Table 7). On the other hand, its reference allele is present in homozygosis in all individuals from <~15,000 BP indicated in the blue chronological group in the heatmap (Figure 2, Table 7). Therefore, none of the analysed individuals in this study presented the haplotype that has been related to altruistic behaviors. Similarly, since the alternative allele of SNP3 has been significantly related to Autism Spectrum Disorder[53] (Table 7), in regard with SNP3, the results of this study suggest that part of the analysed prehistoric population may have presented behaviors characteristic of Autism Spectrum Disorder.

Table 7. Information of the homozygosity and heterozygosity of the analysed SNPs in the OXTR gene

The number of reads presented by the analysed individuals for this position is shown in parenthesis, next the number of times that each of the alleles appeared, indicating in red alternative allele and in blue the reference allele. Homozygosity of the reference and alternative alleles is shown in blue and red boxes, respectively. Heterozygosity is indicated in yellow boxes.

Continuing with the analysed SNPs in the OXTR gene, we have SNPs 2 (rs7632287), 4 (rs2268491) and 5 (rs2254298). All these SNPs present their alternative alleles in many of the oldest individuals analysed in the current study (Kostenki14, Sunghir3, Sunghir2 and Sunghir 4, among others). However, these alleles are practically absent in the rest of individuals from <~35,000 BP who did have reads at these positions. The alternative alleles of these SNPs have been associated with different phenotypes that are contrary to prosocial, empathetic and altruistic behaviors, with the reference alleles showing a greater tendency toward behaviors related to greater social skills[55],[56] (see Table 5). Consequently, these results, along with the low frequency of these alternative alleles in the current global population [SNP2: G_reference  = 0.74252, A_alternative = 0.2575; SNP4: C_reference = 0.8716, T_alternative = 0.1284; SNP5: G_reference = 0.8602, A_alternative = 0.1398; information obtained from dbSNP], suggest that the evolutionary pattern of these SNPs points to the decrease of alternative alleles along time and, therefore, a greater probability of developing prosocial and empathetic behaviors and better communication skills (see phenotype of the reference alleles of SNPs 2, 4 and 5 in Table 5)[24],[26],[55],[56]. However, it is important to highlight that the number of individuals who presented reads at these SNPs in this study was very small (SNP2, n=5; SNP4, n=7; SNP5, n=7), thus it is necessary to expand the sample size in order to safely confirm this possible tendency.

Lastly, with regard to SNP6 (rs53576) in the OXTR gene, it was observed that its alternative allele G is present in all the individuals who did have reads at this position, although these were only four individuals (Figure 2). The great functional variability (or variability of phenotypic consequences) represented by this alternative allele –related to better sociocognitive skills (see phenotypes in Table 5)–, along with the high frequency of the alternative allele present in the current global population (A_reference = 0.3265, G_alternative = 0.6735; information obtained from dbSNP), could be indicative of the great relevance of said SNP in human evolution. Moreover, although this position presents reads in few of the analysed individuals, all of them (from both chronological groups indicated in the heatmap) presented the alternative allele in homozygosity (Table 7). Thus, it could be inferred that, regarding this SNP, prehistoric H. sapiens may have developed prosocial and empathetic behaviors, with lower reactivity to stress, since at least 35,000 BP. Nevertheless, as in the previous cases, it is necessary to increase the number of analysed individuals to confirm this possibility.

Figure 3. Hypothalamic-Pituitary-Adrenal axis (HPA) and its relationship with the genes analysed in this study

The hypothalamus comprises, among other nuclei, the paraventricular nucleus (PVN). The neurons of the PVN, which secrete CRH (Corticotropin Releasing Hormone) (responsible for the activation of the HPA axis), present the following genes: OXTR (oxytocin receptors), CD38 (Cluster of Differentiation 38), OPRM1 (Opioid Receptor µ 1), GABRA6 (Gamma-Aminobutyric Acid type A Receptor) and FKBP5 (FK506 Binding Protein 51). All these genes regulate the release of CRH, but enzyme FKBP5 affects the sensitivity of the reception of GC (Glucocorticoids) and, therefore, modulates the negative feedback of the axis exerted by the GC (see Figure 1). CRH, when released, travels to the hypophysis, where the GABRA6, FKBP5 and AVPR1B (Arginine Vasopressin Receptor 1B) genes are also expressed, regulating the release of ACTH (Adrenocorticotropic Hormone). ACTH then binds the receptors of the adrenal glands, inducing the release of GC, such as cortisol, into the blood stream, generating the physiological response to stress. The cells of the adrenal gland present, among others, enzyme COMT (Catechol-O-MethylTransferase), which metabolises catecholamines, that are adrenal hormones that exert a positive feedback on the activation of the axis by travelling to the hypothalamus through the blood stream: the greater the presence of catecholamines, the greater the activation of the HPA axis. Sources: [23], [29]-[36].

3.1.3. CD38 gene

Thirdly, there is gene CD38. This gene encodes a cell-surface glycoprotein present in some immune cells, and, among other functions, it regulates the secretion of oxytocin (Figure 3). The reference allele G of the analysed SNP (rs3796863) in this gene has been associated with a lower parental contact or involvement in upbringing[30],[60]. In the population analysed in this study, the alternative allele T is in all the individuals from the Upper Paleolithic (~23,000-40,000 BP) that did have reads at this position (Figure 2), in both heterozygosis and homozygosis (Table 8). In turn, with regard to the 3 epipalaeolithic and earlier individuals who presented reads at this position, from 7,000-11,000 BP, the alternative allele T is absent in Chan Meso and Canes1 Meso, and present in La Braña1. Although these results could indicate a possible decrease in the frequency of the alternative allele from the beginning of the Upper Paleolithic and, consequently, a possible human inclination toward greater parental or upbringing involvement and/or attachment[30], the number of individuals (n=7) was very small to generate such approximations. Given the possible tendency of this SNP toward possible prosocial and more sensitive behaviors, it is necessary to increase the sample size in order to characterise, more accurately, the possible involvement of this SNP in the development of prosocial behaviors in the evolution of H. sapiens.

3.1.4. GABRA6 gene

Fourthly, there is the SNP (rs3219151) analysed in gene GABRA6, which is one of the main regulators of the central nervous system (Figure 3). Under the Trier Social Stress Test (TSST) –a psychosocial test used in psychological and medical research that induces acute stress in the participants to evaluate the effects of stress during different activities[66]–, the alternative allele T of this SNP has been associated with greater responses to ACTH (adrenocorticotropic hormone) and cortisol in both homozygosis and heterozygosis[31]. Therefore, this phenotype could be associated with a possible overactivation of the HPA axis and stress under certain circumstances[23],[29]. Although this nucleotide position only presents reads in few of the analysed individuals (n=5), the alternative allele T was identified in all of them, in both homozygosis and heterozygosis (Table 8). These data thus suggest a possible greater sensitivity to ACTH and cortisol of the analysed H. sapiens under stressful circumstances. However, due to the small sample size, it is not possible to establish greater approximations regarding this SNP.

Table 8. Information of the homozygosity and heterozygosity of the analysed SNPs in CD38 gene (a), GABRA6 gene (b) and FKBP5 gene (c)

The number of reads presented by the analysed individuals for this position is shown in parenthesis, next the number of times that each of the alleles appeared, indicating in red alternative allele and in blue the reference allele. Homozygosity of the reference and alternative alleles is shown in blue and red boxes, respectively. Heterozygosity is indicated in yellow boxes.

3.1.5. FKBP5 gene

Fifthly, we have the FKBP5 gene. The SNP analysed in this gene (rs1360780) presents the reference allele T associated with poorly impulsive behaviors in the making of decisions that involve the choice between short-term and long-term rewards[61], although also a greater probability of suffering disorders related to anxiety and greater reactivity to fear[61],[62]. This position presents reads in very few of the prehistoric humans analysed in the current study (Figure 2). Kostenki14 and Sunghir3, i.e., the oldest individuals (~38,000 BP and ~35,000 BP), present the reference allele T in heterozygosis, whereas La Braña1 (~7.000 BP) presents it in homozygosis (Table 8). The allele frequency of this SNP in the current global population indicates a lower frequency of allele T with respect to the alternative allele C (T_reference = 0.3140, C_alternative = 0.6859; information obtained from dbSNP), which could suggest that the reference allele T granted a lower evolutionary advantage than allele C, thereby reflecting a possible human evolution toward a lower reactivity to fear and stress. However, although it is obvious that the reference allele T was present in the human population since, at least, the early Upper Paleolithic, the scarcity of individuals analysed in this study who did have reads at this position (n=3) does not allow establishing any hypothesis based on this SNP.

3.1.6. OPRM1 gene

Sixthly, there is the OPRM1 gene. This gene encodes the opioid receptor µ 1, which is expressed in some neurons. When this receptor binds β-endorphins, the secretion of CRH is inhibited by these neurons, and, since this hormone activates the HPA axis (thus generating stress), the role of this receptor is to reduce the behaviors associated with stress[31]. In vitro studies have reported that the reference allele A is related to a change in the morphology of the receptor, which then binds β-endorphins with lower affinity, resulting in a lower inhibition and a possible activation of the HPA axis[31]. Regarding the presence of the alternative and reference alleles of the SNP analysed in this gene (rs1799971) in the human population included in this study (Figure1), this position presented reads in some of the analysed individuals (n=11). The alternative allele G is present in only one of the individuals who presented reads at the position of this SNP in the nuclear genome (ELT002, ~5,700 BP), appearing in heterozygosis (Table 9). The rest of the 9 older individuals (~38,000-5,800 BP) showed the presence of the reference allele A in homozygosis, as well as I10866, i.e., the latest individual included in the study. These results, along with the low allele frequency of the alternative allele G of this SNP in the current global population (A_reference = 0.8612, G_alternative = 0.1388; information obtained from dbSNP), could indicate a greater evolutionary benefit of allele A along time in general, possibly due to some other reason that has not been considered in this work, although it is necessary to increase the sample size in order to corroborate the low frequency of allele G in the prehistorical individuals.

Table 9. Information of the homozygosis and heterozygosis of the SNPs analysed in OPRM1 gene (a) and COMT gene (b)

The number of reads presented by the analysed individuals for this position is shown in parenthesis, next the number of times that each of the alleles appeared, indicating in red alternative allele and in blue the reference allele. Homozygosity of the reference and alternative alleles is shown in blue and red boxes, respectively. Heterozygosity is indicated in yellow boxes.

3.1.7. COMT gene

Lastly, there is the SNP (rs4680) analysed in the COMT gene. This position presented reads in many of the ancient humans analysed in this study (Figure 2) (n=11), from all the chronological groups established; only two (Sunghir2 and Chan_Meso, from very different chronologies) of the 11 individuals who did have reads at this SNP presented the reference allele G in homozygosis (see Table 9). Taking into account that COMT metabolises catecholamines (e.g., noradrenaline and adrenaline), which activate the HPA axis by stimulating the release of CRH and ACTH[31],[63], it has been described that the alternative allele A causes a decrease in the degradation rate of these chemical substances compared to the variant G[31],[63]. As a result, the levels of catecholamine may be relatively high in individuals with the allele variant A with respect to those with allele G, leading to an increase of HPA activation and, therefore, stress generation. The fact that only two individuals presented the reference allele G in homozygosis indicates that allele A was widely distributed throughout prehistory, or at least in the individuals analysed in this work. Thus, these results suggest the possibility that the reduced enzymatic activity of COMT under certain circumstances was common in the analysed prehistoric individuals, with a greater probability of HPA overactivity and its consequences.

3.2. Results and discussion of Fisher’s exact test

Fisher’s exact statistic was calculated in 5 of the 14 SNPs studied in this work, which correspond to 3 genes: i) SNPs 3, 8 and 9 of the OXTR gene; ii) SNP1 of the OPRM1 gene; and iii) SNP1 of the COMT gene.

Firstly, in the case of SNPs 3 and 8 of the OXTR gene, SNP1 of the OPRM1 gene and SNP1 of the COMT gene, the results indicate (Table 10) that there is insufficient evidence to reject the null hypothesis, which states that the presence of the alternative allele is similar in the group of individuals from the Upper Paleolithic and in the group of individuals from the Epipaleolithic and later periods. Likewise, the estimation of the Odds ratio suggests that there is no substantial difference in the presence of the alternative allele between the two groups. However, it is necessary to expand the number of ancient humans analysed in order to determine the allele presence, frequencies and their phenotypic consequences more accurately.

Table 10. Comparison between individuals presenting the alternative allele in ≤26% and >26% of the reads for the selected SNPs that did present reads at those nucleotide positions in 10 or more individuals

The table indicates the number of individuals presenting the alternative allele in a percentage >26% in those SNPs that conserved 5 or more reads for the SNPs indicated. Comparison is made between two chronological groups (Upper Palaeolithic and Epipalaeolithic and later periods). The p-values are also shown, which were calculated from Fisher’s exact test (5% significance level), as well as the Odds ratio (95% confidence interval).

Secondly, for SNP9 of the OXTR gene, the results (Table 10) indicate that there are significant differences in the presence of the alternative allele between the individuals of the Upper Paleolithic and those of the Epipaleolithic and later periods, since the p-value is far below the significance level (5%) established (Table 10). Similarly, the estimation of the Odds ratio shows that the presence of the alternative alleles is greater in the group of individuals from the Upper Paleolithic compared to the group of individuals from the Epipaleolithic and later periods. Therefore, these data suggest a possible evolutionary advantage of the reference allele of SNP9, i.e., one of the three alleles (along with the reference alleles of SNP3 and SNP8 of OXTR gene) that constitute the haplotype associated with more altruistic behaviors (Table 5). This would strengthen the probability of the GGC haplotypic combination with altruistic behaviors in the H. sapiens from the Epipaleolithic and later periods. Nevertheless, it is important to point out that only 3 individuals provided information for this haplotype for Epipalaeolithic and later humans in the case of SNP8. Thus, as was previously mentioned, it is necessary to increase the sample size in order to corroborate these results.

4. Conclusions

This study analysed human prosociality from an evolutionary perspective. Although the genomic data of prehistoric individuals are scarce, making it difficult to understand human prosocial behaviors, we explored possible gene variants underlying some of the behaviors that could favored the biological and social adaptation of our species in the last ~40,000 years.

The analysis of 14 SNPs (single nucleotide polymorphisms of the DNA sequence) –selected based on their relationship with behaviors regulated by the HPA axis– in a sample of 20 individuals of Homo sapiens gathered from different genomic s, provided some data of interest. These results on the variation of these SNPs in relation to the prosocial behaviors of our species in individuals from different periods may be synthesized into three different scenarios:

1.- There are SNPs whose alternative alleles are not present in the groups of prehistoric humans considered in this analysis from the Upper Paleolithic onwards, such as the alternative alleles of SNP1 of the AVPR1B gene, SNP1 of the OXTR gene, and SNP1 of the OPRM1 gene. The alternative alleles of these SNPs have been related to prosocial behaviors and/or lower HPA-axis activity, and thus also lower physiological stress response. Nevertheless, their allele frequency in the current global population is low. Therefore, these alleles may not only be associated with social behaviors but could also have other phenotypical consequences, explaining this low present day frequency. However, in order to interpret the meaning of the variations of these SNPs, some of the limitations mentioned above would have to be overcome, such as: i) to increase the sample size; and ii) to study in depth the possibility that there were other phenotypic consequences of these SNPs that might have caused them to be unfavorable for the lifestyle of our species, which would explain why its frequency would not have increased over the course of human evolution.

2.- There are SNPs whose alternative alleles are present in the Upper Paleolithic hunter-gatherers considered in this study, and it seems that they are preserved, with slight variants, throughout late prehistoric times. This is the case of some SNPs of the OXTR gene, such as SNP3 and SNP6. However, the phenotypic consequence of both mentioned SNPs is the opposite, since the alternative allele of SNP3 is associated with Autism Spectrum Disorder behaviors, whereas the alternative allele of SNP6 is associated with greater social skills. In the case of the COMT and GABRA6 genes, there are SNPs whose alternative alleles are present in H. sapiens since the Upper Paleolithic throughout late prehistoric times; but these alleles are associated with greater reactivity to stress under certain circumstances.

Therefore, there are allele variants present in the first representatives of our species in Eurasia from the early Upper Paleolithic that are associated with diverse behaviors, and even opposite behaviors, in terms of the development and promotion of social interactions. These results show, on the one hand, the complexity of the phenotypes of “prosocial behaviors” and, on the other hand, the importance of other external factors, such as the gene-environment interaction (which could not be assessed in the present study), in the genesis of these phenotypes.

3.- There are other SNPs analysed in this work that seem to have different frequencies in the individuals from the Upper Paleolithic and in those from the Epipaleolithic and later prehistoric periods. We detected a statistically significant decrease in the alternative allele presence in favor of the reference allele (SNP 9 of OXTR gene) when these two chronological groups were compared. In the case of SNPs 2, 4 and 5 of the OXTR gene, chronological differences could not be statistically corroborated.

Altogether, in this third scenario, there is evidence of a temporal variation in the frequencies of the SNPs of some genes. Considering that the reference allele of SNP9 of OXTR gene constitute, along with the reference allele of SNP3 and SNP8 (OXTR gene), a haplotype associated with the possible development of altruistic behaviors, this results could support the increase of behaviors of greater social skills (e.g., altruistic behaviors) in post-Paleolithic periods. Therefore, these allele variants, whose frequency has increased since the end of the Upper Paleolithic, may have been beneficial, given the adaptive advantages of their phenotypes for the human group, thereby increasing the development of behaviors that would favor social cohesion.

Lastly, SNP1 of the CD38 gene and SNP1 of the FKBP5 gene were only represented in a few of the individuals analysed in this study, which made it impossible to classify them in any of the three mentioned scenarios.

To sum up, this study highlights the complexity of the evolution of human prosociality, not only due to the limitation posed by the scant preservation of ancient DNA, but also as a result of the possible influence of other factors (e.g., genes-environment) on behavioral phenotypes, which demonstrates the need to contextualise the variants in the allele frequencies with the environment that the ancient H. sapiens were exposed to. However, although no categorical assertions can be established, this preliminary analysis shows temporal changes in the allele frequencies of some genes, which provides support for our methodological proposal to analyze the evolution of sociocognitive skills of Homo sapiens. Moreover, further research on this topic could help to understand the evolution of one of the most important human characteristic as a social species. A greater comprehension on this topic would allow knowing, more accurately, the influence that the evolution of human prosociality could have had on the interaction with other animal species, favoring their domestication, as is the case with the origin of the dog.

5. Acknowledgments

This research was supported by a Basque Government Grant for Research Groups in the Basque University Education System (IT 1693-22) (IP: C.R), and the grant to A.A.Z. from the Basque Government pre-doctoral fellowship program. We would like to thank two referees for their valuable contributions to the review of this manuscript, as well as the editors of this special volume of RIEV (N. Ibarretxe y R.Jimeno) for their support and encouragement to complete this issue.

6. References