SERT gene polymorphisms are associated with risk-taking behaviour and breeding parameters in wild great tits

A combination of genes, the environment and experience-based learning is known to be responsible for variation in animal (including human) behaviour (Carere and Maestripieri, 2013). However, the molecular mechanisms underlying behavioural differences between individuals in natural animal populations, and to what extent variation in these mechanisms cause fitness differences, are not well understood. Although the gene–behaviour associations are environment and condition dependent (Adriaenssens and Johnsson, 2010), behavioural traits consistently differ between individuals of the same population (for a review, see Carere and Maestripieri, 2013) and have a significant genetic basis (for a review, see Laine and van Oers, 2017). Hence, it is important to conduct experiments in different situations and habitats, and across different years to evaluate the repeatability of behaviours and assess the generality of previously reported genetic effects on behavioural variation. Moreover, in order to understand evolutionary change in phenotypes, a connection between genes, behaviour and individual fitness should be tested in the same study system (see review by van Oers et al., 2005). Growing evidence has suggested that variation in many behavioural traits such as fear, anxiety and novelty-seeking are associated with dopamine- and serotonin-related signalling affected by genes (Korsten et al., 2010; van Oers et al., 2005; Serretti et al., 2006).

The role of serotonin transmission in shaping fear-related behavioural variation is widely studied in humans and laboratory animals (Canli and Lesch, 2007; Gryglewski et al., 2014) and also in free-living populations (Champoux et al., 2002; Serretti et al., 2006). Alterations in the expression of neurotransmitter-related genes could change serotonergic signalling and thus ultimately shape behaviours in which serotonin has been implicated (O'Leary and Cryan, 2010). A candidate gene that is central to the regulation of extracellular and synaptic serotonin concentration is the serotonin transporter gene (SERT). Variation in SERT gene expression has been associated with differences in the SERT protein availability and in modulating active reuptake of serotonin between synapses (Lesch et al., 1996). SERT gene polymorphisms in humans have been repeatedly linked with different psychiatric states such as depression, anxiety and suicidal behaviour (Murphy et al., 2008). In addition, serotonin transport affects parental behaviour through behavioural decisions and by prolactin secretion controlled by neural pathways containing serotonin (e.g. Bakermans-Kranenburg and van Ijzendoorn, 2008). However, the relationships found between behavioural traits and genetic differences in serotonin-related genes in humans have not always been consistent (Savitz and Ramesar, 2004).

In this study we used the candidate gene approach to explore the association of SERT gene polymorphisms with behaviour in a natural population of breeding great tits (Parus major Linnaeus 1758). For birds, SERT gene variation was found to be different in urban compared with rural populations of blackbirds (Turdus merula) (Müller et al., 2013). Also, Riyahi et al. (2015) have provided evidence that these polymorphisms in great tits differ between environments. Finally, the SERT gene in dunnocks (Prunella modularis) has been shown to affect female behaviour, where homozygous females took fewer risks than heterozygotes (Holtmann et al., 2016). We hypothesize that the SERT gene plays an important role in affecting great tit individual novelty-seeking and fear-related behaviours. Moreover, because individual differences in behavioural decisions could affect fitness, we explored the link between genetic polymorphisms and breeding parameters. For this, we determined the polymorphic regions of the SERT gene (exons) and the promoter region. From these, we identified nine single nucleotide polymorphisms (SNPs) for testing the relationships between behavioural phenotypes and genotypes. We used samples originating from two different years in order to examine the generality of the effect. Specifically, we first predicted that the parental latency to feed nestlings at the presence of different novel objects varies between genotypes. Previous studies have shown that the delay in nest visitation and nestling provisioning rates are repeatable characteristics in great tits (Pagani-Núñez and Senar, 2013; Vrublevska et al., 2015). Hence, this trait can be regarded as a potential candidate trait for assessing personality in the wild. Second, we predicted that variation in breeding parameters such as the first egg date, clutch and brood size are associated with these genetic variants, probably via consistent behavioural differences between individuals. This assumption is based on previous studies reporting that reproductive success covaries with serotonergic signalling and maternal care in mice (Lerch-Haner et al., 2008) and in relation to personality variation (Smith and Blumstein, 2008).

The study was conducted in southwest Estonia (58°7′N, 25°5′E) at the Kilingi-Nõmme study site during springs 2012 and 2014 in the great tit breeding season (April to June). The area (approximately 50 km²) consists of deciduous (Betula pendula, Alnus spp. and Salix spp.) and coniferous (Pinus sylvestris and Picea abies) forest patches. Great tits in our study population breed in nest boxes mounted on tree trunks at a height of 1.5–1.8 m.

The nest boxes were checked throughout the breeding period to obtain data on lay date of the first egg, clutch size, hatching date (hatch date=day 0 post-hatch) and number of hatchlings and nestlings (15 days old). Adults (both male and female) were caught using nest-box traps during the second half of the nestling period (between days 7 and 15 post-hatch). Traps were checked every 15 min. Adults were ringed and weighed with a Pesola spring balance (precision 0.1 g). Captured adults were 1–2 years of age. Tarsus and wing length were measured using sliding calipers (with precision to 0.1 mm). Blood samples (∼70 µl) were taken from the brachial vein with a sterile lancet and collected in a heparinized capillary tube. The blood was immediately stored at +4°C, and afterwards centrifuged to separate cells from plasma and then stored at −80°C until DNA analysis. Birds were ringed under Estonian Department of Environment Licence No. 11 and blood samples were taken under Animal Procedures Committee of the Estonian Ministry of Agriculture Licence No. 100.

DNA was extracted from the blood samples (2012) with a commercial kit (High Pure PCR Template Preparation Kit 18, Roche, Basel, Switzerland) according to the manufacturer's instructions. Blood samples collected in 2014 were treated with the Puregene DNA Purification Kit (Qiagen, Hilden, Germany). The concentration of received DNA was approximately 1 μg 1 ml^–1. Primers were designed using the great tit genome (assembly Parus_major1.1.; accession number: SRS1185780; Laine et al., 2016). Primer sequences for SERT gene exonic regions (exons 1–13) were designed de novo (see Table 1).

Primers for promoter region were used from Riyahi et al. (2015). The PCR reaction mixture (25 μl) consisted of 2 μl (consisting of 2 ng genomic DNA), 5 μl recombinant Taq Fire polymerase, 0.5 μl of each primer and water (17 μl). Amplification was performed sequentially as follows: 95°C for 3 min, 38 cycles at 95°C for 30 s, 52°C for 30 s, 72°C for 1 min and finally 72°C for 4 min. All the samples were sequenced at the Estonian Biocentre (Tartu, Estonia), where an Applied Biosystems sequencing platform was used. Received sequences were aligned using ChromasPro (Technelysium Pty Ltd, South Brisbane, Australia) and visually examined for polymorphisms (SNPs). The SNPs were named according to the position in the exonic region and/or were taken from a previous study (Riyahi et al., 2015). Construction of haplotypes followed the DNAsp programme (Librado and Rozas, 2009) and the linkage disequilibrium data were analysed with DNAsp 5.10 and its online calculator (Gaunt et al., 2007).

To test for Hardy–Weinberg equilibrium, chi-square tests were used. In order to take multiple testing into account, we calculated critical P-values using Holm's correction for the number of independent SNPs to test for deviations from Hardy–Weinberg equilibrium (P_critical=0.0062) and associations of genotype with latency behaviour and reproductive parameters (P_critical=0.0125, after omitting four SNPs that deviated from Hardy–Weinberg equilibrium). Sample sizes (see Tables 1 and 2) vary between analyses because not all individuals were successfully trapped and bled, and were therefore not genotyped.

The response of adult birds to a novel object was tested when nestlings were 7–13 days old. A novel object (Eppendorf box) was placed on the rooftop of the nest box at the beginning of the experimental phase (see Timm et al., 2015 for details). The experimental phase was preceded by a control phase, where the normal feeding rate of birds was recorded for 15 min. We measured the feeding delay (latency in seconds to enter the nest box from the moment when the novel object was placed). The experimental phase lasted 30 min. Overall, we tested the adult (both parents) response to a novel object on 113 breeding pairs (2012 and 2014). The experiment (performed randomly in the population) was carried out between early morning (06:00 h) and early afternoon (14:00 h), during the times of highest feeding activity. The behaviour of birds was recorded with a digital video camera placed on a tripod approximately 5–10 m from the nest box.

In 2014, in addition to the response to the novel object, we also recorded the delay of birds to enter the trap (116 individuals). When nestlings were at least 8 days old, a trap was placed in the inner side of the nest-box hole, allowing birds to enter, but not to leave. The latency to enter the nest box started from the moment the bird was in sight of the camera (approximately 1–2 m from the nest box; observed from the video recordings). All traps were checked after every 15 min. The trap was considered to be a more hidden novel object, as it is difficult to see from a distance, but imposing significantly stronger effect at proximity.

Data were analysed using R statistical software version 3.1.2 (https://www.r-project.org/). Linear mixed-effect models (LMMs) with restricted maximum likelihood (REML), type III sum of squares and Satterthwaite approximation for degrees of freedom were conducted using ANOVA from the car package and the lmer function in the lme4 package. F-statistics were obtained using the lmerTest package. The R package compute.es was used to calculate effect sizes (Fisher's z).

Responses to a novel object and a trap during the breeding time are potentially affected by several confounding factors. First, in the analysis using data over two different years, we examined the effect of genotype, parental sex, genotype by sex interaction, sex order (two levels: the first or second parent to enter the nest box), treatment (two levels: control, novel object), habitat (deciduous versus coniferous) and brood size (as a covariate) on the latency to enter the nest box. Year, individual ID (to account for repeated measures) and brood identity (female and male parent form a breeding pair) were also included in this model as random factors. Second, in the analysis on the latency to enter the nest box in response to the nest-box trap, we examined the effect of genotype, parental sex, sex order, habitat (deciduous versus coniferous), brood size (as a covariate) and brood identity (as a random factor) and we included the interaction between genotype and sex. In none of these above models were sex order and the habitat type (two levels: deciduous or coniferous habitat) significant, and thus they were omitted from final models. Third, using data on different breeding parameters (separate models for the start of egg-laying, clutch size, number of hatched young, number of fledged young), we examined the effect of genotype, parental sex, genotype by sex interaction, habitat, brood size as fixed factors and individual ID, brood identity and year as random factors.

Pairwise post hoc comparisons were carried out with Tukey HSD post hoc multiple comparison tests (R package lsmeans; Lenth, 2016). Latency to enter the nest box was log(base) transformed before the analysis to meet normality assumptions of residuals from the models. Breeding parameters were normalised using Box–Cox power transformation (R package MASS) (Venables and Ripley, 2002). Standardised β values were obtained from models with standardised variables. In order to cover allelic and genotypic effects, we used a general genetic model (that retains the three distinct genotype categories), an additive model (assuming that each additional copy of the variant allele increases the response; genotypes coded as a continuous predictor) and an over-dominant model (wherein heterozygous individuals are tested against homozygous individuals; genotype categories are coded as heterozygotes versus homozygotes). In order to take into account multiple testing, we calculated critical P-values using Holm's correction for the number of independent SNPs (P_critical=0.0125).

The SERT gene is polymorphic in our great tit population. Three SNPs were identified in the promoter region (SNP136, SNP290 and SNP478), one in exon 1 (SNP187), two in exon 3 (SNP253 and SNP278) and three in exon 8 (SNP197, SNP407, SNP457). The other exonic regions (exons 2, 4, 5, 6, 7, 9–13) studied were monomorphic. SNP197 in exon 8 is non-synonymous (complex data) and could affect behaviour through changes in amino acids. SNP457 in exon 8 is also non-synonymous and results in an amino acid change from proline (Pro) to leucine (Leu). SNPs found in the promoter region are the same as the SNPs detected by Riyahi et al. (2015), but the number of significant SNPs varies between the populations, as Riyahi and colleagues found more polymorphic nucleotides (11) in this gene region.

Four SNP variants in our population deviated significantly from Hardy–Weinberg equilibrium even after correcting for multiple testing (Table 2), and were removed from further analysis. All the other allele frequencies tended to stay similar over the years in all SNPs. Underrepresented genotypes were present in low frequencies in both 2012 and 2014 (Table 2). We found no complete linkage disequilibrium between SNPs (all R²<0.2±0.2). Confidence intervals (CI) for promoter regions were as follows: Kelly’s Z_nS=0.24; Rozas’ Z_A=0.24; Rozas’ ZZ=0.03. For linkage disequilibrium between exon and promoter, the CIs were: Kelly’s Z_nS=0.00; Rozas’ Z_A=0.00; Rozas’ ZZ=0.00.

In the complete sample (combined data of 2012 and 2014), no effects of sex and treatment on the delay in nest visitation were detected (Table 3). There was a positive correlation between the latency in the control and novel object phase [β=0.36±0.07 (s.e.m.), P<0.001, N=87 nests; corrected for sex, brood size and sex order]. We also detected significant intra-individual repeatability for the latency between the control and the novel object phase (repeatability=0.67, F_1,381=5.13, P=0.024; sex: F_1,381=1.53, P=0.22). SNP253 in exon 3 was significantly associated with the delay in nest visitation (Table 3, Fig. 1). Birds with a CC genotype resumed feeding significantly earlier than those with a CT genotype. There was also a tendency for an association between SNP197 in exon 8 and the delay in nest visitation, but this relationship was not significant after correcting for multiple testing (Table 3). The explained variance in the feeding delay by SNP253 in exon 3 was 11% (R² effect size, see Fig. 2 for effect sizes). Owing to the low frequency of the minor allele at SNP253 (exon 3), additive and over-dominant models were not fitted for this site (Table 3).

In 2014, we also measured the time it took birds (N=50 nests) to enter the nest box in response to the nest box trap. The latency in the control phase was positively related to the latency in the novel object phase (β±s.e.=0.40±0.10) and trapping phase (β=0.48±0.12). A positive association was also found between the latency in the novel object and the trapping phase (β=0.24±0.11, N=50 nests), while correcting for sex, brood size and sex order. The latency to enter the nest box differed significantly between all three treatment phases (F_2,216.0=19.84, P<0.001), being the longest in the trapping (mean±s.e.m.=847±59 s), intermediate in the novel object on the nest box (494±53 s) and the shortest in the control (281±58 s) phase (all post hoc tests P<0.01). In the trapping phase, females delayed nest visitation for a longer time than males (F_2,17,53=13.52, P=0.002) and the latency differed between SNP197 genotypes (exon 8) (F_2,19,92=6.06, P=0.008; post hoc tests: CC versus CT, P=0.029; CC versus TT, P=0.011, CT versus TT, P=0.99; Fig. 3). We also found an association between SNP253 in exon 3 and the delay in nest visitation (F_2,43=3.83, P=0.029) in response to the nest box trap, although this relationship did not remain significant after correcting for multiple testing.

Egg-laying started significantly earlier in deciduous than in coniferous forests (F_1,100.1=18.25, P<0.001). SNP457 in exon 8 was significantly associated with the start of egg-laying in the additive model [Fisher's z (effect size)=0.36; confidence limits (CL): lower=0.14, upper=0.58; Table 4]. Females with the AG (mean±s.e.m.=26.62±0.31, N=61; post hoc P=0.004) and AA genotype (26.64±0.52, N=21; post hoc P=0.024) started egg-laying earlier than those with the GG genotype (28.00±0.50, N=23), whereas no difference was observed between the AG and AA genotypes (P=0.99). Clutch size was on average larger in deciduous than in coniferous forests (F_1,106,85=8.50, P=0.003, N=113 broods), but was not related to any SERT polymorphism (all P>0.21; Table 4).

Brood size at hatching did not depend on habitat (F_1,91.35=3.26, P=0.07, N=94 broods) and was not associated with any SNP in a general or an additive genetic model (Table 4). However, in an over-dominant model, SNP197 in exon 8 was associated with the number of hatched young. This relationship was sex specific (sex: F_1,42.7=1.09, P=0.30; sex×SNP: F_1,44.16=4.28, P=0.044) and significant in females only (F_1,73.0=7.43, P=0.008; Table 4). Post hoc analysis among females showed that CT heterozygotes had larger broods at hatching than CC and TT homozygotes combined (mean±s.e.m.=10.79±0.56 versus 9.86±0.54 offspring). Brood size on day 15 post-hatch (almost fledging) was not related to habitat (F_1,95,45=4.80, P=0.029, P_critical=0.0125, N=98 broods), or to any of the genetic polymorphisms (all P>0.08; Table 4).

In this study, we determined whether SNPs in the coding regions of the SERT gene could predict behavioural differences in a natural population of great tits. We found that SERT indeed plays an important role in individual behavioural response, and that this effect is consistent over sexes and contexts (e.g. contact with novel object, reaction to predator, etc.). Also, the SERT gene could be related to fitness-related traits, as one of the SNPs was significantly associated with breeding parameters. In this study, genetic variability in four polymorphic exonic regions and the promoter area in the SERT gene was associated with variation in the delay in feeding in response to a mild stressor (camera set-up in the control phase), a moderate stressor (novel object) and a strong stressor (nest box trap) and three fitness-related traits (start of egg-laying, clutch size and brood size). Of those five regions, SNPs in exons 3 and 8 were significantly associated with behavioural traits and SNPs in exon 8 with breeding parameters. SNP197 and SNP457 in exon 8 are non-synonymous and have the potential to affect behaviour through changes in amino acids. The polymorphism in exon 3 is synonymous, and thus does not alter the amino acid sequence; therefore, it is less likely that the relationship between the SNP253 in exon 3 and parental reproductive behaviour represents a causal relationship. However, this effect could be mediated via linkage with other non-synonymous polymorphisms, or effects on splicing or transcription (Chamary et al., 2006).

The Hardy–Weinberg equilibrium was violated in some cases, indicating that some genotypes are constantly missing or underrepresented in both years (Table 2). There is a probability that sequencing errors may be the cause of this. However, these differences tended to stay similar over both years, despite good-quality DNA and repeated sequencing. The likely cause of Hardy–Weinberg equilibrium deviations can be complex, especially in a wild, outbred population. Mate choice, mutation rates, selection, genetic drift and gene flow are not stable over time as mating is not random, there is migration in open populations and some of the individuals may be more successful in passing on their genes to the next generation than others (Hartl et al., 1997). When we compare our results with an earlier study in which the SERT promoter region was sequenced (Riyahi et al., 2015), we see that the number of detected polymorphisms in the promoter area differs between studies. We assume that the main differences are caused by the differences in PCR product length, where especially the ends of sequences could vary owing to the variation in sequence quality. Moreover, some SNPs found in a Spanish population (Riyahi et al., 2015) were not polymorphic in our population.

We found that two newly discovered SNPs in exons 3 and 8 potentially play an important role in novelty-seeking and risk-taking behaviours of wild great tits. In our experiment, we manipulated nestling feeding behaviour by exposing parents to different novel stimuli. As a result, we found that the delay in feeding was the longest in the trapping phase, intermediate in the novel object phase and the shortest in the control phase. Given that feeding interruption after exposure to a stressor is a reliable measure of parental fear level (Tilgar et al., 2011; Van der Veen and Sivars, 2000), we have reason to expect parents to be differently stressed when facing a novel or threatening object. We also found that behavioural latencies measured in the presence of different stressors were moderately correlated, suggesting that this trait is repeatable over contexts. In an earlier study, a polymorphism found in the SERT gene promoter region was found to be associated with exploration behaviour in great tits (Riyahi et al., 2015), indicating a more general role for SERT in great tit personality traits.

We suggest that genotype–behaviour associations are relatively consistent over different situations. Despite varying fear levels, genotypic effects on the latency response remained the same in different experimental phases. Heterozygous individuals carrying CT genotype at the SNP253 (exon 3) tended to have significantly longer delays in feeding nestlings in all situations compared with homozygous CC genotype. In the case of SNP197 (exon 8), significant differences between genotypes emerged when the latency response during the most stressful manipulation (trapping phase) was measured. We found that homozygous individuals carrying CC genotype at the SNP197 had longer feeding delays than CT and TT genotypes.

We hypothesized that genetic differences in SERT may affect breeding success via changes in behavioural decision-making. Here, we show for the first time that different SERT genotypes are associated with breeding time as well as offspring number at hatching.

Firstly, starting egg-laying earlier could be beneficial, especially in seasonal environments, where the selection pressure on timing of reproduction exists as the nestling phase falls within the peak of caterpillar abundance (Visser et al., 1998). Because SERT mainly affects behaviour-related traits, we suggest that the impact on breeding time could be linked to serotonin effects on feeding behaviour (Magalhães et al., 2010), which is important for reaching breeding condition. Alternatively, it could be mediated via serotonin effects on hormonal levels, such as prolactin. In this respect, it has been shown that prolactin is secreted by anterior pituitary cells under the stimulatory control of vasoactive intestinal polypeptide, the activity which is controlled by neural pathways containing dopamine, serotonin and opioids (Freeman et al., 2000).

Secondly, the number of hatchlings was higher in CT heterozygotes compared with TT and CC homozygotes combined (over-dominant model). In the great tit, a female incubates eggs alone and the male feeds her during incubation (Gosler and Clement, 2007). Thus, we can assume that heterozygotes can recover more rapidly from disturbing events and resume parental activities, such as incubation and seeking food, after a shorter delay than homozygotes. Genetic polymorphisms in these loci are non-synonymous, potentially affecting SERT enzymatic activity and neuronal information processing in the central nervous system (Freeman et al., 2000). In order to understand these effects in the future, we have to keep in mind that the potential impact of SERT and the polymorphisms are linked with fitness through behavioural decisions.

The impact of genotype on reproduction could indeed occur via individual behavioural decisions (Dingemanse and Wolf, 2010). This may especially be true for reproductive traits that are related to parental sensitivity to environmental stressors through decision-making (e.g. incubation or provisioning behaviour), rather than the traits that directly depend on the parental quality and food availability (e.g. clutch size and egg size). Hence, in a heterogeneous environment, individually consistent behaviours can be adaptive in a particular environment (e.g. reduced fearfulness can be related to increased parental care in stressful environments), but maladaptive in another (for example, shy individuals can be less effective in finding novel food sources; Bell et al., 2007). Thus, it is likely that fluctuating natural selection translates genetically driven behavioural patterns into individual differences in fitness components. This could lead to the maintenance of genetic polymorphisms in natural populations.

We kindly thank Jüri Parik for help with sequencing and analysing data, and Kaarin Koosa for help in the field. We also thank all the reviewers for reading and commenting on the manuscript.