The genetic landscape of animal behavior

Although most animal behaviors are associated with some form of heritable genetic variation we do not yet understand how genes sculpt behavior across evolution, either directly or indirectly. To address this, I here compile a dataset comprised of over 1,000 genomic loci representing a spectrum of behavioral variation across animal taxa. Comparative analyses reveal that courtship and feeding behaviors are associated with genomic regions of significantly greater effect than other traits, on average three fold greater than other behaviors. Investigations of whole-genome sequencing and phenotypic data for 87 behavioral traits from the Drosophila Genetics Reference Panel indicate that courtship and feeding behaviors have significantly greater genetic contributions and that, in general, behavioral traits overlap little in individual base pairs but increasingly interact at the levels of genes and traits. These results provide evidence that different types of behavior are associated with variable genetic bases and suggest that, across animal evolution, the genetic landscape of behavior is more rugged, yet predictable, than previously thought.


Introduction
Nearly all behaviors are associated with some form of heritable genetic variation (Kendler and 48 Greenspan 2006). This interplay between genetic and other forces that shape behavior is complex and 49 disentangling it occupies an array of research endeavors, spanning disciplines from evolutionary biology 50 to psychiatry. Accordingly, recent years have seen reasonable progress toward understanding the genetic 51 architecture of certain behavioral traits using model systems (Reaume and Sokolowski 2011). The general 52 conclusion from this research in mice, flies, worms, and humans is that the genetic architectures of 53 behaviors generally fit an exponential distribution, with a small number of loci of moderate to large effect 54 and a larger number of loci with small effects (Robertson 1967;Flint and Mackay 2009). However, owing 55 to limits in data and methods, the extent to which genetic architectures vary across a full spectrum of 56 behaviors and animal taxa has remained largely unexplored. 57 Behaviors can exhibit considerable variation in genetic influence. Comparative analyses reveal 58 that behaviors vary substantially in heritability estimates, most often ranging between 10% and 50% 59 (Kendler and Greenspan 2006; Mousseau and Roff 1987;Meffert et al. 2002). Analyses of individual 60 behaviors reveal even greater diversity. For example, a single retro-element is responsible variation in a 61 courtship song between Drosophila species (Ding et al. 2016) while other traits, such as deer mouse 62 burrowing, have modular genetic architectures comprised of multiple interacting loci (Weber et al. 2013). 63 Furthermore, the structure and effect of genetic architectures may vary with behavioral traits, as suggested 64 by the preponderance of large effect loci found for insect courtship traits across multiple species 65 (Arbuthnott 2009). Despite these observations the extent to which behavioral traits may systematically 66 vary across species and behaviors remains unknown. Understanding this could provide insights into how 67 behaviors respond to evolutionary processes, the prospects for finding general principles in the genetic 68 evolution of behavior, and even potentially why there has been such variable success in the mapping of 69 human neuropsychiatric traits. 70 Here, using reports associating behavioral variation with the genes for specific traits across 71 diverse species, I assemble a comparative behavior genetics resource composed of 1,007 significant 72 4 genomic loci from 114 QTL studies conducted in 30 species across 5 taxonomic classes. These data 73 exploit advances in sequencing and genetic marker design that have accelerated reports using quantitative 74 trait locus (QTL) mapping to identify genomic regions that are associated with behavioral variation 75 (Lander and Botstein 1989;Flint and Mackay 2009). With the compiled dataset I compare the genetic 76 architecture of behavioral types across animal taxa. I then corroborate these observations and assay 77 genetic processes involved in the early stages of behavioral differentiation in a natural population using 78 whole genome data from the Drosophila Genetic Resource Panel (DGRP). These analyses provide insight 79 into the genetic architecture of behavior across animals and the interplay between specific behavioral 80 traits and their genetic influence through evolutionary history.

82
I performed a comprehensive analysis of results aggregated from 114 QTL studies conducted in 83 30 species across 5 taxonomic classes to assemble a comparative behavior genetics resource composed of 84 1,007 significant genomic loci (Database S1). The species examined represent over 500 million years of 85 evolutionary divergence and over a broad spectrum of phylogenetic data (Fig 1a). For each locus I 86 annotated the trait measured and its associated effect size (percent phenotypic variation explained), the 87 reported measure of significance (e.g., LOD score), genomic locus, and study sample size. I focused the 88 analyses on the reported effect sizes to allow comparison of the genomic architecture of traits across 89 studies similar to previous meta-analyses of behavioral QTL in mice and flies (Flint 2003;Flint and 90 Mackay 2009). 91 I found that the distribution of effect sizes in the dataset is similar to that found in these previous 92 studies (Fig 1b). In the majority of loci (89.51%) the effect sizes are less than 20% with a mean effect size 93 of 9.54%, suggesting that the genetic bases of most behaviors assayed are complex and composed of 94 many loci of moderate effect.

95
Though these results support a model of many loci with small effects for behavior overall, I then 96 asked whether genetic architecture might vary across types of behavior. I identified ten behavioral 97 categories for which traits had been measured in at least two species (See supplementary methods). My 98 5 null hypothesis was that individual categories would likely reflect the overall distribution seen across the 99 dataset, consistent with previous observations that QTL have relatively similar effect sizes across mouse 100 and fly phenotypes (Flint 2003;Flint and Mackay 2009). Surprisingly, I found instead that behaviors 101 differed significantly in their effect sizes. Specifically, loci associated with courtship (n=124) explained 102 significantly more phenotypic variance than all other behaviors combined (Kruskal-Wallis p = 6.7 x 10 -29 ) 103 and had a mean effect size three times larger than found in all other categories (Fig. 1c). Loci associated 104 with feeding behaviors (n=11) also explained significantly more phenotypic variance than all other 105 behaviors combined (p = 6.8x10 -13 ) while emotion and social behaviors explained significantly less (p = 106 8.6 x 10 -33 ; p = 2.5 x 10 -21 , respectively). These data suggest that, across species, courtship and feeding 107 behaviors possess genetic architectures different from those of other traits. To assess whether these observations arose from differences in the behavioral traits, as 111 byproducts or experimental artifacts I controlled for factors that might have contributed bias. I first 112 considered the effect of intraspecific (within species) compared to interspecific (between species) crosses 113 used for the QTL mapping, a known source of influence in QTL studies (Broman 2001). I indeed found   After eliminating sources of potential biases inherent to individual datasets, I next considered the 129 possibility that the detection of courtship and feeding behaviors as outliers was a trivial outcome of our 130 own classification method for grouping single behaviors into ten categories. Minimally assuming that the 131 categorizations of courtship and feeding traits were correct, it is possible that the binning of traits into the 132 other eight categories may have masked a real signal from some biologically relevant categorization. To test this possibility, I compared the distribution of effect sizes for the courtship and feeding 134 categories to the distribution for all other behaviors combined (Supplementary methods). I found that 135 courtship behaviors explained significantly more variation (p <0.05) than 89% of non-courtship behaviors 136 while feeding behaviors explained more variation than 46% of non-feeding behaviors (Fig 2b; Fig S4b). I 137 complemented this test with a bootstrap analysis that created a null distribution from 10,000 permutations 138 of the non-courtship/feeding trait effect sizes. The observed mean adjusted effect size for both courtship 139 and feeding fell significantly outside the bootstrap null distribution created for each comparison (p < 5 x 140 10 -200 ) (Fig 2c; Fig S4c). These findings reject the notion that there may be another categorization of non-141 courtship and feeding behaviors missed by our schema that explains substantially more variation of effect.

142
My results suggest that courtship behaviors, and to a lesser extent feeding, may respond to 143 evolutionary pressures differently than other behavioral traits. Consistent with this notion, previous 144 analyses of the QTL behavior literature in insects found that a majority of courtship traits are associated 145 with few loci of particularly strong effect that play a potential role in rapid speciation through prezygotic 146 isolation (Arbuthnott 2009). In addition, theoretical work has suggested that traits controlling local 147 adaptation during speciation, such as courtship and feeding, evolve more rapidly if they are associated 148 with a smaller number of loci (Gavrilets et al. 2007). Given the importance of behavior's role in the early 149 stages of speciation it may be possible that for the organisms and traits analyzed here, courtship and 150 feeding traits with simpler genetic components of large effect were selected for during the evolution of 151 these lineages. These observations led me to hypothesize that, in a naturally interbreeding population,

164
After running GCTA 20 behavioral traits passed a p-value threshold of 0.05, indicating that autosomal 165 SNPs could explain more trait variation than by chance in these cases ( Fig. 3a; Supplementary methods).

166
The majority of these traits were enriched for involvement in courtship and feeding: 30% (6/20) were 167 associated with courtship and 50% (10/20) were either involved in olfactory behavior or feeding. Notably, 168 for a number of these traits the vast majority of phenotypic variation could be explained by genome-wide In addition to an increase in genomic heritability, my QTL analyses also showed that the genomic 174 architectures of courtship and feeding traits may be simpler and of higher effect. To test this I performed 175 a separate GWA experiment for each trait across all lines with available phenotypic data and filtered for 176 SNPs with a nominal p-value of 5 x 10 -6 (Supplementary methods). At this threshold I found 25,919 SNPs 177 ( Fig. 3b; Table S1).

178
I re-ran GCTA for each trait using only SNPs identified at p < 5 x 10 -6 from the GWAS 179 (supplementary methods). This test is more conservative compared to genome-wide GCTA since it uses  (Fig. 3d). Many of these SNPs fell within the same genomic regions. I found 72 genes had at 204 least 2 SNPs associated with multiple traits, several of which contained a multitude of variants (Fig. S6a).  Table S3). In addition, I found 81 intergenic SNPs that 207 each occurred within 20kb of their nearest gene -26 genes in total -suggesting potential regulatory roles 208 for these SNPs (Fig. S6b).   The connection between male courtship behaviors and body size has long been recognized in laboratory 241 strains of Drosophila though with little evidence of a molecular basis for this effect (Ewing 1961). In 242 general I find extensive evidence of both directional (G→P 1 →P 2 ) and general (P 1 ←G→P 2 ) pleiotropic 243 effects between traits in the DGRP, supporting the notion that the early stages of behavioral 244 diversification involve the role of genes that can effect multiple types of traits. Furthermore, I observe that 245 while variation in behavior across trait categories is associated with non-overlapping variants these may 246 occur in common genes and molecular pathways with pleiotropic effects, reflecting suggestions of the 247 existence of phenotypic "hotspots" that are recurrently used by evolution to sculpt phenotypes (Stern & 248 Orgogozo 2008).

249
Taken together these results suggest that behavioral traits may respond to evolutionary processes 250 with greater variation than previously appreciated. For example, researchers may now anticipate that 251 assaying a courtship ritual will likely yield a higher genetic effect than, say, variation in a personality   Table S1.

319
Phylogeny 320 I used the phylogenetic relationships reported in Ponting 2008 as a template for our phylogeny of species 321 examined (Fig. 1a). I added unrepresented species and adjusted dates of evolutionary divergence using the  For the analyses plotted in Fig. 2a-2c and Figs. S1-4 I summed the effect sizes of all loci associated with a 344 specific behavioral measure for each study. This was done to allow for a comparison of the maximum 345 amount of phenotypic variance explained for each trait in order to allow for conservative test between comparison between all categories (Fig. S1) and between courtship/non-courtship (Fig. 2a) and 362 feeding/non-feeding (Fig. S4a). The bootstrap comparisons in Figs. 2c and S3c were done using the custom R function 379 bootstrap.2independent which is available on the Fernald lab website. For these tests I 380 permuted the non-courtship/non-feeding residual effect sizes 10,000 times (with replacement) to create a 381 null distribution against I which I tested the observed median residual effect size for each trait. A p-value 382 for each test was calculated by dividing the sum of instances in which the permuted medians were greater 383 than the observed by 10,000. All plots were produced using base graphics in R and adjusted for design in 384 Adobe Illustrator. Associations were then filtered for a p-value < 5 x 10 -6 . SNPs associated with multiple traits were  (Table S4). This matrix could then be directly queried for comparison 446 of the effect sizes associated with a certain set of SNPs across traits of interest. In order to assess the 447 overall structure of this data set I used Spearman rank correlations to test the associations between all 448 possible trait pairs. The results of this test were visualized using the clustering functionality of 449 heatmap2 in R (Fig. S7).

450
Tests for trait pair directionality 451 Directionality in the relationships between trait pairs was tested by first obtaining pairwise rank 452 correlations for each trait pair in which both traits were associated with >3 significant SNPs (60 traits).

453
For traits x and y, s 1 is the vector of SNPs significantly associated with trait x and s 2 is the vector of SNPs The resulting p-values were adjusted using Bonferroni correction.