SUMMARY

Behavior and physiology are regulated by both environment and social context. A central goal in the study of the social control of behavior is to determine the underlying physiological, cellular and molecular mechanisms in the brain. The African cichlid fish Astatotilapia burtoni has long been used as a model system to study how social interactions regulate neural and behavioral plasticity. In this species, males are either socially dominant and reproductively active or subordinate and reproductively suppressed. This phenotypic difference is reversible. Using an integrative approach that combines quantitative behavioral measurements, functional genomics and bioinformatic analyses, we examine neural gene expression in dominant and subordinate males as well as in brooding females. We confirm the role of numerous candidate genes that are part of neuroendocrine pathways and show that specific co-regulated gene sets (modules), as well as specific functional gene ontology categories, are significantly associated with either dominance or reproductive state. Finally, even though the dominant and subordinate phenotypes are robustly defined, we find a surprisingly high degree of individual variation in the transcript levels of the very genes that are differentially regulated between these phenotypes. The results of the present study demonstrate the molecular complexity in the brain underlying social behavior, identify novel targets for future studies, validate many candidate genes and exploit individual variation in order to gain biological insights.

INTRODUCTION

Brain and behavior are sculpted by a dynamic interplay between genotype and environment. Changes in social status mediate changes in a wide range of molecular, neuroanatomical, endocrine and behavioral characteristics(Oliveira et al., 2002; Fernald, 2004; Sapolsky, 2005) at different time scales (Hofmann, 2003). Nervous systems respond by modulating neural properties and neuroendocrine responses. Such rapid, short-term changes will often lead to sustained,longer-term changes in the brain and behavior via differential gene expression, leading to subsequent structural and physiological changes. However, few studies have attempted to determine the molecular underpinnings of socially regulated phenotypes within a comprehensive genomic framework,especially in vertebrates. In the present study, we apply a molecular systems approach in a highly social fish species to ask how social context regulates the molecular and physiological substrates that, in turn, sculpt subsequent behavior.

The African cichlid fish, Astatotilapia burtoni (Taxonomic Authority) (formerly Haplochromis burtoni) has become an important model system to study the mechanisms underlying socially mediated behavioral change. In this species, 20–30% of males are dominant (D), slow growing,brightly coloured and actively defend territories for mating. The remaining subordinate (S) males mimic females by schooling and displaying cryptic coloration, while experiencing faster growth(Fernald, 1977; Hofmann et al., 1999a). Subordinate males show little aggression and territoriality and, importantly,have regressed gonads and are thus not reproductive(Fernald and Hirata, 1977a; Fernald and Hirata, 1977b; Francis et al., 1993). These behavioral and physiological characteristics are plastic and influenced by the immediate social environment, such that an individual male switches between the D and S phenotypes several times during its life depending upon its relative ability to obtain and maintain access to a territory through encounters with other males (Hofmann et al., 1999a). Environmental conditions, availability of territorial shelters, relative body size and physiological condition influence the probability of acquiring and maintaining a territory. The phenotypic switch occurs over a timescale of minutes to days to weeks in both the field and the laboratory (White et al.,2002; Hofmann,2003; Burmeister et al.,2005).

In the laboratory, A. burtoni has been the focus of hormonal and molecular studies related to a broad range of phenotypic traits that are affected by the transition between the two male phenotypes(Fig. 1) (for reviews, see Hofmann and Fernald, 2000; Fernald, 2002; Hofmann, 2003; Fernald, 2004). Variation in several different components of the endocrine system reflects the complexity underlying the corresponding phenotypic switch between these dramatically different male phenotypes. Neuroendocrine pathways regulating androgen production (Parikh et al.,2006), growth (Hofmann et al.,1999a; Hofmann et al.,1999b) and stress response(Fox et al., 1997) change in a complex fashion as animals undergo phenotypic change(Fig. 1). Parikh et al.(Parikh et al., 2006)suggested that the higher levels of testosterone (T) and 11–ketotestosterone (11–KT) that they measured in D males might promote aggressive behavior, as has been shown using androgen manipulation in other fish species. Production and release of neuropeptides and neuromodulators such as gonadotropin-releasing hormone (GnRH1)(White et al., 2002) and somatostatin (Hofmann and Fernald,2000; Trainor and Hofmann,2006; Trainor and Hofmann,2007) are higher in dominant males. According to White et al.,GnRH1 mRNA levels and gonadosomatic index (GSI) (a measure of gonadal development) are positively correlated(White et al., 2002). The higher level of GnRH1 in D males coincides with territory acquisition and mating opportunity; as such, it is similar to GnRH changes in seasonally reproducing species (Amano et al.,1995; Dawson et al.,2001; Nelson,2005; Hofmann,2006). Fluctuations at the molecular level also co-vary with the observed behavioral switch of male phenotype. In the case of GnRH(Parhar et al., 2005; Au et al., 2006) and somatostatin (Trainor and Hofmann,2006), specific receptor sub-types are regulated according to social status. Furthermore, the membrane properties of GnRH1-expressing neurons, reduced excitability in S males correspond with the differences in peptide release between D and S males(Greenwood and Fernald,2004).

Fig. 1.

Schematic representation of (A) phenotypic and (B) physiological characteristics associated with dominant (D) and subordinate (S) male phenotypes in A. burtoni. The graphs are based on the following studies: gonad size (Hofmann and Fernald,2000); growth (Hofmann et al.,1999a); GnRH1 (White et al.,2002); testosterone (Francis et al., 1993); somatostatin(Hofmann and Fernald, 2000);growth hormone (Hofmann et al.,1999b); and cortisol (Fox et al., 1997).

Fig. 1.

Schematic representation of (A) phenotypic and (B) physiological characteristics associated with dominant (D) and subordinate (S) male phenotypes in A. burtoni. The graphs are based on the following studies: gonad size (Hofmann and Fernald,2000); growth (Hofmann et al.,1999a); GnRH1 (White et al.,2002); testosterone (Francis et al., 1993); somatostatin(Hofmann and Fernald, 2000);growth hormone (Hofmann et al.,1999b); and cortisol (Fox et al., 1997).

Manipulation of the social environment allows experimental control of the phenotypic switch (Francis et al.,1993; White et al.,2002; Burmeister et al.,2005). The ease of experimental manipulation is paired with a wealth of ecological and evolutionary information available for haplochromine cichlids in general (for a review, see Kocher, 2004; Salzburger et al., 2005). For example, there is thought to be a trade-off between reproduction and survival,such that territorial males, the sole reproducers, are also subject to higher mortality rates through predation, probably due to their conspicuous coloration (Fernald and Hirata,1977a; Fernald and Hirata,1977b; Maan et al.,2008). Furthermore, A. burtoni is currently undergoing whole-genome sequencing, along with three other African cichlid species (see www.broad.mit.edu/models/tilapia/),and thus offers an unrivaled laboratory-based model system for the genomic analysis of complex and ecologically relevant phenotypes.

cDNA microarrays for transcript profiling, have become a powerful tool when applied to species of behavioral, ecological or evolutionary interest [e.g. alternative life histories (Aubin-Horth et al., 2005), cooperative breeding(Aubin-Horth et al., 2007),social behavior (Grozinger and Robinson,2002; Whitfield et al.,2003; Robinson et al.,2005), response to heat stress(Buckley et al., 2006),response to environmental estrogens(Martyniuk et al., 2006) and physiology of drug addiction (Rhodes and Crabbe, 2005)]. In the present study, we employ a microarray platform that contains many known candidate genes as well as ∼4000 brain-derived cichlid cDNAs (Renn et al.,2004). Such a combined candidate gene and genomic strategy allows hypothesis-driven and discovery-based experiments on a single platform. The obviously complex nature of behavioral traits – such as this socially regulated, reversible switch between D and S phenotypes – requires a discovery-based approach in order to identify the many genes involved.

In the present study, we analyze gene transcript patterns for reproductively active D males and reproductively suppressed S males as well as for brooding females. We then compare the expression patterns for each of the three phenotypes in order to identify gene sets associated with reproduction or dominance behavior, providing insight into the molecular modularity underlying these phenotypes. Next, we annotate the array features according to gene ontology (GO), with the goal of identifying gene regulation within molecular categories free of a priori expectations or experimenter bias. Finally, we examine variation in gene expression patterns between individual animals within a social phenotype and ask whether any of these variable genes are also those that are differentially expressed between the social phenotypes.

MATERIALS AND METHODS

Animals

Astatotilapia burtoni (Günter) from a wild-caught stock population were kept in aquaria under conditions mimicking their natural environment (Fernald and Hirata,1977b): pH 8.0, 28°C water temperature and 12 h:12 h light:dark cycle. Gravel substrate and terracotta shelters allowed the establishment and maintenance of territories necessary for reproduction(Fernald and Hirata, 1977a). Fish were fed every morning with cichlid pellets and flakes. All work was carried out in compliance with Animal Care and Use Guidelines.

Behavioral experiments

Males were marked using colored bead combinations attached near to the dorsal fin. Nine groups of 2–3 males with 2–3 gravid females were established in half of a 100 liters aquarium. Each group was visually isolated from neighboring fish. Ten-minute behavioral observations were made approximately twice every week for five weeks by directly observing chasing,threat, display, border threat, courtship, flee, schooling and carousel (as described by Fernald, 1977). Behavioral measures were used to calculate a dominance index (DI), i.e. the sum of all aggressive behaviors observed (threat, chasing, border threat and carousel) minus the number of submissive behaviors observed (flee) and a reproductive index (RI), i.e. the total of all reproductive behaviors (court,dig, spawn). On the final day of the experiment, males of each status were taken for gene expression profiling only if they had continuously expressed either the D or S phenotype for the past 28 days. Six mouth-brooding females were also selected for expression profiling. The 18 fish used for gene expression analysis were obtained from a total of seven different groups on the basis of stable phenotype. For the purpose of hybridization design the animals were treated as equivalent within phenotype, though comparisons between individuals can be made indirectly (see below). Close inspection of the behavioral and microarray data did not reveal any covariation between animals derived from the same social group. Standard length, body mass and gonad mass were measured for each fish in order to calculate the GSI as gonad mass/body mass. Condition factor (CF) was calculated based on the residuals from the regression of body mass on standard length for before and after the experiment (r2=0.95; P=2.79×10–6), and growth rate (GR) was calculated as the relative change in standard length over the course of the experiment. Whole brains were dissected and stored in RNAlater (Ambion,Austin, TX, USA) within 5 min of initial tank disturbance.

The array used in the present study has numerous redundant features, i.e. two or more features represent the same gene. We exploited this property for quality control purposes and to assess the sensitivity of the approach. For example, the array includes four independent features that represent the neuropeptide GnRH1 (two previously cloned cDNAs and two obtained independently from the cDNA library), known to be expressed in only ∼300 neurons and upregulated in D males (Davis and Fernald,1990; White et al.,2002). All four features were similarly significantly differentially expressed, being upregulated with dominance, and thus serve as both a biological and a technical control demonstrating the sensitivity even with whole-brain RNA.

Microarray analysis

Brains were homogenized (Tissue tearor, Biospec products, Bartlesville, OK,USA) and total RNA was extracted according to standard Trizol (Invitrogen,Carlsbad, CA, USA) and phase-lock gel (Eppendorf; Westbury, NY, USA)protocols. RNA integrity was determined on the Bioanalyzer (Agilent; Santa Clara, CA, USA) and spectrophotometer (Agilent) prior to indirect RNA labeling protocol, starting with 2 μg of total RNA according to Renn et al.(Renn et al., 2004). Briefly,each RNA sample was labeled twice, once with Cy3 and once with Cy5. After purification from unincorporated label, each sample was divided in two,combined with the appropriate samples and every individual was compared with two individuals from each of the other two phenotypes in a balanced loop design incorporating a dye-reversal (see below for how to access annotation of submitted data at NCBI's GEOdatabase for further details). Targets were hybridized to the brain-specific cDNA array from A. burtoni(Renn et al., 2004) at 65°C for 12–16 h. Arrays were scanned on an Axon 4000B arrays scanner (Genepix 4.0; Molecular Devices, Sunnyvale, CA, USA). The loop design allows for the direct comparison of samples of interest, thus offering greater statistical power with fewer replicates(Churchill, 2002).

After filtering for bad feature morphology, hybridization artifacts and low intensity (<2 s.d. above local background), raw data were imported into R software [v.1.0 R-Development team, 2004, Vienna, Austria(R Development Core Team 2006)] and normalized using the Linear Models for Microarray Data package [LIMMA v.1.6.6 (Smyth et al.,2004)]. Background-subtracted mean intensities were normalized using a within-array printtip-lowess normalization and used to calculate ratios for a Bayesian analysis of gene expression levels [BAGEL v.3.6(Townsend and Hartl, 2002)]. BAGEL takes advantage of additional information obtained from transitive comparisons of individuals in loop designs experiments(Townsend and Hartl, 2002; Churchill, 2002). Genes represented by more than one feature on the array were only counted as significant if at least one representative passed the significance threshold and the full complement survived Fisher's test of combining probabilities from multiple tests of significance [p.794 in Sokal and Rohlf(Sokal and Rohlf, 1995]. All raw and processed data are available at GEOdatabase(www.ncbi.nlm.nih.gov/projects/geo/),sample numbers GSM267785–819 of series GSE10624.

Functional annotation

The clone templates for PCR amplification were end sequenced(Salzburger et al., 2008),resulting in 4258 expressed sequence tags (ESTs), 3670 of which were deemed to be of high quality and have been submitted to GenBank (accession number CN472211–CN468542) and are maintained at the Dana-Farber Cancer Institute as a GeneIndex(http://compbio.dfci.harvard.edu/tgi/cgi-bin/tgi/gimain.pl?gudb=a_burtoni)such that 1280 clones combine into 399 tentative contigs (TC), leaving 2381 singleton sequences for a total of 2780 unique sequences. Half (49.5%) of the sequences (1408 out of 2842) could be annotated with a `best hit' at a threshold of e–12 or better according to BLAST alignment to UniProt (v.5.2). GO terms were applied to 869 unique cichlid genes by transitive annotation, meaning that the GO annotations for a cichlid gene's best hit were collected and used in the present study for further analysis. In order to avoid species bias, we collected GO terms from all genes with the same name as the best hit annotation. GO annotations also include confidence codes. The less reliable annotations derive from `inferred sequence similarity(ISS)' or `inferred electronic annotation (IEA)'(The Gene Ontology Consortium,2000), therefore we excluded those annotations when transitively applying GO terms to the ESTs represented on the A. burtonimicroarray. The resulting GO graphs [referred to as directed acyclic graphs(DAGs)] were then `slimmed' to 183 terms and a total of 4102 total annotations. For the `slimming' process, the leaf-most nodes that were selected to contain a minimum of 10 annotated cichlid sequences, and parent nodes were retained only if an additional 10 cichlid sequences were annotated at that level.

Over- and under-representation of GO terms for a regulated set of genes was determined in Cytoscape (Shannon et al.,2003) using the Biological Network Gene Ontology tool, BiNGO(Maere et al., 2005), which relies upon hypergeometric statistical significance. As GO categories are highly non-independent, the statistical treatment of these terms is still under discussion (Ge et al.,2003). Also, owing to the small number of genes for each ontology term and the relatively small number of genes that are regulated, there is less statistical power to identify significantly under-represented GO terms. For these reasons, we use GO analysis as a hypothesis-generating tool and report only uncorrected hypergeometric P-values.

Clustering

Prior to clustering, features representing replicate ESTs were collapsed by combining probabilities from multiple tests of significance [p. 794 in Sokal and Rohlf (Sokal and Rohlf,1995)], and the mean expression level was determined for each set of features. A hierarchical clustering analysis was applied to the list of genes that were significantly regulated according to D, S and brooding female phenotypes. The estimated gene expression levels were used to obtain the dissimilarity matrix by applying Euclidian distance measure, which integrates the effects of amplitude of ratios as well as direction (correlation) in patterns. Clustering analysis of gene expression patterns of each individual was performed using the hclust function in R software v.2.0.1. Clustering was based on dissimilarity measures obtained using the dist functions in the stats package. The consensus tree and bootstrap confidence values for each tree node were obtained with the consensus function in the maanova package(Wu et al., 2002). The consensus tree and confidence values were calculated as the proportion of trees obtained with bootstrapped datasets that agreed with the original tree. Each bootstrapped tree was based on the Euclidian distance matrix calculated for each of 1000 permutated gene expression profile datasets obtained by resampling with replacement. Alternate clustering methods and different measures of distance are available and are similarly appropriate for gene expression analysis. Hierarchical clustering, based upon all features on the full array rather than on regulated genes only, provided a similar tree with reduced confidence values at each node (not shown). The heatmapfunction and colour options in the package gplots were used to visualize clusters of gene expression, the z-transformed expression ratios were grouped by k-means function in the stats package of R and ordered as such, while the samples were ordered according to the consensus hierarchical cluster. GO terms provide a means to address the possible functional relationship of a cluster of genes that are coordinately regulated. However, no statistically significant over- or under-representation of GO terms was seen for any of the gene groups identified by the k-means clusters according to a hypergeometric test (not shown).

RESULTS AND DISCUSSION

Male dominance behavior and gene expression differences

Individually tagged males were placed into community tanks with females such that naturalistic dominance hierarchies were established among the males according to standard protocol (Hofmann et al., 1999a). Regular observations over the course of five weeks identified males of definitive dominant (D) and subordinate (S) phenotypes. The territory-holding males (N=6) behaved more aggressively than S males (N=6), as reflected by their significantly higher DI (D males,15.7±1.7; S males, –3.5±0.8; t-test, P<0.001) and RI (D males, 2.12±0.57; S males,0.0±0.0; t-test, P=0.004). The GSI varied continuously across the two groups, showing the expected (yet, in this study,non-significant) trend for higher GSI in dominant males (D males,0.85±0.23; S males, 0.56±0.04; t-test, P=0.25). There was no significant difference in growth rate (GR) over the course of the experiment (D males, 1.25±0.215; S males,1.39±0.285; t-test, P=0.69) or CF (residuals) at the time of sacrifice (D males, 0.01±0.053; S males,–0.01±0.061; t-test, P=0.58). Note, however,that CF differed marginally when the social communities were established at the beginning of the experiment (D males, 0.03±0.051; S males,–0.03±0.068; t-test, P=0.07), which could suggest that increased condition may be advantageous when ascending to dominance in the first place. Specific behaviors also varied significantly and predictably between the D and S males (Fig. 2): increased presence of an eye-bar in D (P<0.001);increased schooling behavior of S (P=0.0005); more fleeing behavior by S (P=0.004); and more frequent chasing behavior(P=0.0003), border threat behavior (P=0.08) and sexual behavior (P=0.003) by D. These results are representative of D and S males in A. burtoni as previously described (e.g. Francis et al., 1993; Hofmann and Fernald, 2000; White et al., 2002).

For the first investigation of our microarray data, we conducted Bayesian analysis of gene expression (Townsend and Hartl, 2002) by treating the six individuals of each phenotype as biological replicates and calculated gene expression differences at the phenotype level (D, S or brooding female). When comparing the different phenotypes, array features representing 3647 cichlid unique sequences passed quality filters in a sufficient number of hybridizations to be considered in the results. D and S phenotypes showed significant differences in gene expression for 171 genes [Bayesian posterior probability, BPP>0.99) (87 upregulated; 84 downregulated with dominance (supplementary material Tables S1 and S2)]. Permutation analysis of the data shows that, at this significance threshold of BPP>0.99, only ∼17% of these 171 genes could be expected to show significant differential expression by chance, a reasonable rate of false positives (Table 1). Therefore, almost five percent of the genes studied were differentially regulated in the brain according to male dominance phenotype. This percentage is considerably smaller than that found in typical honey bee colonies (39%),where nurse and forager phenotypes are not only distinguished by their social role but also differ in their daily behavioral routines and the environment in which they move about (Whitfield et al.,2003). These results also show phenotype-specific regulation of a considerably smaller portion of the genome compared with the 15% observed difference between alternative mating tactics in a study of Atlantic salmon males, which similarly compared sexually mature and immature males(Aubin-Horth et al., 2005).

Table 1.

False discovery rate (FDR) estimated by randomization of hybridization ratios both between genes and between arrays using permutation with replacement and resampling for analyses conducted at the level of phenotype and sex

Bayesian probability threshold
Phenotypes comparedP⩾0.99P⩾0.995P⩾0.999P⩾0.9995
Male dominant>Male subordinate 89 61 41 37 
Male dominant>Female 273 218 148 127 
Male subordinate>Male dominant 84 61 26 23 
Male subordinate>Female 277 220 136 116 
Female>Dominant male 205 154 96 78 
Female>Subordinate male 219 164 95 83 
Males>Female 316 269 174 147 
Female>Males 253 211 119 103 
FDR % 17.08 10.76 2.07 
standard error 0.61 0.57 0.39 0.39 
Bayesian probability threshold
Phenotypes comparedP⩾0.99P⩾0.995P⩾0.999P⩾0.9995
Male dominant>Male subordinate 89 61 41 37 
Male dominant>Female 273 218 148 127 
Male subordinate>Male dominant 84 61 26 23 
Male subordinate>Female 277 220 136 116 
Female>Dominant male 205 154 96 78 
Female>Subordinate male 219 164 95 83 
Males>Female 316 269 174 147 
Female>Males 253 211 119 103 
FDR % 17.08 10.76 2.07 
standard error 0.61 0.57 0.39 0.39 

The number of features identified by analysis of the permutated dataset as significant for each threshold is reported as the average (N=5)percentage of the number of features identified at that threshold by analysis of the original dataset

Candidate genes

The microarray was designed to include candidate genes previously studied in the context of social dominance and other behavioral contexts in A. burtoni (Hofmann, 2003; Fernald, 2004). The inclusion of known candidate genes allowed us to test multiple hypotheses and also offered validation of the microarray results by comparison with previous studies for some of these genes (summarized in Fig. 1). For instance,peptidergic neurons in the pre-optic area (POA) and other brain regions express several neurohormones (e.g. arginine vasotocin, GnRH1, galanin) as well as neurohormone receptors and steroid receptors, which have previously been shown in separate studies to play a role in the regulation of social behavior. Below, for the first time, we provide a combined analysis of these neuroendocrine pathways in A. burtoni.

Fig. 2.

Social behavior varied significantly and predictably between the dominant(D) and subordinate (S) males. y-axis represents (A,C) the percent time or (B,D,E,F) the mean number of observed events per 10 min focal observation. The numbers identify individual males (filled circles, dominants;open squares, subordinates).

Fig. 2.

Social behavior varied significantly and predictably between the dominant(D) and subordinate (S) males. y-axis represents (A,C) the percent time or (B,D,E,F) the mean number of observed events per 10 min focal observation. The numbers identify individual males (filled circles, dominants;open squares, subordinates).

As predicted from previous studies using ribonuclease protection assays and in situ hybridization, among the three GnRH neuropeptide genes that are expressed in the brain of fish, only GnRH1, the form expressed in the POA(White et al., 1994; White et al., 2002), showed highly significant regulation in the microarray results (BPP=0.9998), with D males having higher levels. Given the small number of cells expressing this neuropeptide (∼300) (Soma et al.,1996; Munz, 1999),confirmation of GnRH1 regulation by our microarray analysis provides an important cross-validation and confirms the sensitivity of the array even when using whole-brain RNA. As predicted from previous studies(White et al., 1994; White et al., 2002), our results also confirmed that the other two forms of GnRH, GnRH2, (midbrain) and GnRH3 (terminal nerve), are not regulated according to male social phenotype(BPP=0.366 and 0.700, respectively). None of the GnRH receptor sub-types on the array were significantly regulated, although studies have demonstrated their regulation in the pituitary relative to sexual maturity and social status (Parhar et al., 2005; Au et al., 2006).

Galanin, a neuropeptide that links metabolic activity and reproduction through regulation of GnRH release (reviewed by Kageyama et al., 2005; Tortorella et al., 2007), was marginally upregulated in D males (BPP=0.9592). There is considerable evidence from mammals that galanin reduces nociception(Wiesenfeldhallin et al.,1992), increases food intake(Schick et al., 1993) and stress reactivity (Holmes et al.,2002), plays a role in regulation of sexual behavior, and is itself regulated by GnRH and estrogen(Gabriel et al., 1993). Specifically, in fish, galanin is thought to play a role in regulation of food intake and is widely distributed in the brain, being localized to the olfactory bulb, telencephalon, hypothalamus, midbrain and posterior brain(reviewed by Volkoff et al.,2005), as well as the pituitary(Jadhao and Pinelli, 2001) and peripheral tissues (Johnsson et al.,2001). The results of the present study suggest the intriguing possibility that galanin might be upregulated in D males as a response to reduced food intake and constant challenges to their social status by other individuals. Future studies will test these novel hypotheses.

Arginine vasotocin (AVT; represented by multiple clones on the array), the non-mammalian homolog of arginine vasopressin (AVP), was among the most strongly regulated genes in this study, being upregulated in the brains of D males (BPP>0.9999). AVP/AVT has been implicated in the regulation of social behavior across vertebrates, including aggression and social affiliation(Goodson, 1998; Goodson and Adkins-Regan,1999; Winslow et al.,1993). In teleost fish, AVT is known to play a role in male mating tactics [peacock blenny (Grober et al.,2002; Carneiro et al.,2003); midshipman (Goodson and Bass, 2001)], as well as in the behavioral regulation of sex change and the associated territory acquisition [bluehead wrasse(Semsar and Godwin, 2003; Semsar and Godwin, 2004)]. AVT is also associated with territorial aggression [damselfish(Santangelo and Bass, 2006)]as well as dominant and territorial behavior in both the male and female of a breeding pair when compared with their subordinate helpers in another cichlid species, the cooperative breeding Neolamprologus pulcher(Aubin-Horth et al., 2007). Our data from A. burtoni suggest a role for AVT in regulation of dominance and are consistent with an in situ hybridization-based study (Greenwood et al.,2008). The AVT V1a receptor, which plays a fundamental role in affiliative behaviors in voles (Lim et al., 2004), was not represented on the array.

The enzyme aromatase, which converts testosterone to estrogen, is important in sex determination (Nakamura and Kobayashi, 2005), sex change in fish(Black et al., 2005; Marsh et al., 2006) and regulation of social behavior (Hallgren et al., 2006). There are two isoforms of aromatase, one localized to the brain and the other to the gonads, and both are represented on the microarray. In the present study, D males showed increased neural expression of the brain form (two features on the array, BPP=0.9914; 0.9997) but not the gonad form of aromatase (BPP=0.3192). This result suggests that the elevated testosterone levels found in D males(Parikh et al., 2006) may affect aggression, courtship or dominance through aromatization and subsequent action via estrogen receptors in the brain. In birds, aromatase activity increases during the territorial period and correlates with aggression (e.g. Soma et al.,2003; Silverin et al.,2004). Blocking brain aromatase reduces male courtship in guppies(Hallgren et al., 2006),further suggesting estrogen-mediated neuroendocrine regulation of reproductive behavior for some species. However, studies on gonadal sex change in fish(Black et al., 2005; Marsh et al., 2006) suggest the opposite relationship between brain aromatase and male aggression and thus a more complex mechanism, possibly involving differences in receptor expression, binding proteins or anatomical localization. Estrogen receptors did not show differential regulation on the array, although they may have been expected to, according to Burmeister et al.(Burmeister et al., 2007). The inability to reliably detect differences in receptor gene expression is probably due to small localized effects that are masked by whole-brain gene expression levels.

Novel genes

We bioinformatically annotated the ESTs obtained from the cichlid cDNA library features represented on the microarray (see Materials and methods). Several of the genes thus identified fall into categories that represent candidates likely to play a role in the social regulation of a complex phenotype (Table 2). In the present study, these genes are considered to be `novel genes' rather than`candidate genes' because the annotation process does not involve rigorous manual curation of genes a priori, which was employed for the candidate genes discussed above. In addition to many genes involved in cellular metabolism that are differentially regulated between the two social phenotypes, we found genes encoding structural proteins, cell-cycle regulators, specific transcription factors, a plethora of neuropeptides,components of the neurosecretory machinery and neurotransmitter receptors.

Table 2.

Best hits for dominant (D) and subordinate (S) regulated genes

Dominant (D) greater than subordinate (S) (D>S)Subordinate (S) greater than dominant (D) (S>D)
Peptides, neurotransmitters and growth hormones Prolactin 1 Orphan nuclear hormone receptor 
 Prolactin 2 Kainate receptor 
 Somatotropin Cholecystokinin 
 Somatolactin C-type natriuretic peptide 
 Proopiomelanocortin  
 Arginine vasotocin  
 GABA (receptor β subunit)  
 Tilapia growth hormone  
 Glycoprotein α  
 GnRH1  
 Brain aromatase  
Granins Chromogranin  
 Neurogranin  
 Secretogranin II  
Cell structure β-actin Clasp2 (regulation of microtubule dynamics) 
 β-tubulin  
 Septin 7  
 Elongation factor 1a  
Synaptic vesicle Similar to synaptophysin  
 Q5RCZ2_PONPY (secretory vesicle electron transport)  
Transcription factor ETS family efr1 bHLH 
  HLH 
Axonal growth Neuromodulin  
 Neuroserpin  
Cell-cycle regulation Schip-1 Rbm5 protein putative tumor suppressor 
 Pdcd4  
Cellular metabolic enzymes Enolase 1 Many ATPases 
 Enolase 2 NADH dehydrogenase and oxidoreductase 
 Uridine kinase Glutathione S-transferase 
 Lactate dehydrogenase Glutamine synthetase 
  Betaine aldehyde dehydrogenase 
  Leucine carboxyl methyltransferase 
Dominant (D) greater than subordinate (S) (D>S)Subordinate (S) greater than dominant (D) (S>D)
Peptides, neurotransmitters and growth hormones Prolactin 1 Orphan nuclear hormone receptor 
 Prolactin 2 Kainate receptor 
 Somatotropin Cholecystokinin 
 Somatolactin C-type natriuretic peptide 
 Proopiomelanocortin  
 Arginine vasotocin  
 GABA (receptor β subunit)  
 Tilapia growth hormone  
 Glycoprotein α  
 GnRH1  
 Brain aromatase  
Granins Chromogranin  
 Neurogranin  
 Secretogranin II  
Cell structure β-actin Clasp2 (regulation of microtubule dynamics) 
 β-tubulin  
 Septin 7  
 Elongation factor 1a  
Synaptic vesicle Similar to synaptophysin  
 Q5RCZ2_PONPY (secretory vesicle electron transport)  
Transcription factor ETS family efr1 bHLH 
  HLH 
Axonal growth Neuromodulin  
 Neuroserpin  
Cell-cycle regulation Schip-1 Rbm5 protein putative tumor suppressor 
 Pdcd4  
Cellular metabolic enzymes Enolase 1 Many ATPases 
 Enolase 2 NADH dehydrogenase and oxidoreductase 
 Uridine kinase Glutathione S-transferase 
 Lactate dehydrogenase Glutamine synthetase 
  Betaine aldehyde dehydrogenase 
  Leucine carboxyl methyltransferase 

Dominant (D) and subordinate (S) regulated genes were grouped according to presumed functional category (see Materials and methods for details). For array feature identification and BPP for differential regulation see Tables S1 and S2 in supplementary material

Genes coding for structural proteins, such as tubulin and actin, and proteins that bind scaffold elements, such as septin 7 and ELF-1a, were more highly expressed in D males reflecting the observed differences in soma size between D and S for pre-optic neurons expressing GnRH1(Francis et al., 1993) and somatostatin (Hofmann and Fernald,2000). Furthermore, genes involved in axonal growth, neuromodulin[also known to play a role in modeling of sex-specific brain regions(Simerly, 2002)] and neuroserpin (Miranda and Lomas,2006), were also upregulated in D males. Taken together, the regulation of this gene set strongly suggests increased neuronal re-wiring in D males not previously reported and possibly similar in scale to the massive remodeling of neural circuits seen in seasonal accession to territoriality and mating accompanied by increased testosterone levels in song birds(Devoogd and Nottebohm, 1981;reviewed in Arnold, 1992). It is particularly intriguing that neuroserpins may play a role in anxiety and sexual behaviors. Specifically, neuroserpin-deficient rats showed decreases in exploratory behavior together with increases in anxiety and neophobia(Madani et al., 2003). In swordtail fish, neuroserpin expression increased in the brain of females exposed to an attractive male compared with females exposed to a non-attractive male (Cummings et al.,2008). This association of neuroserpin with social behavior is intriguing in that it may enable dominant males to approach and interact with novel stimuli such as competitors and potential mates.

Similarly, several cell-cycle regulators(Table 2) were significantly regulated in D and S phenotypes, suggesting that the extent of neurogenesis and subsequent cell death may also differ between these phenotypes, a hypothesis consistent with the finding that cell proliferation in the brain is correlated with high social status in rainbow trout(Sorensen et al., 2007). While there is currently no other evidence for plasticity of this kind in neuroanatomical structures outside the POA in A. burtoni, gross neuroanatomical differences that correspond to species' typical reproductive strategies have been identified in other cichlid species(Pollen et al., 2007).

Genes encoding neuropeptides and protein hormones that have not been previously studied in this system were perhaps the most striking, although not unexpected, class of genes regulated according to social status. In addition to the neuropeptides GnRH1, AVT and galanin discussed above, we found somatotropin, prolactin and somatolactin [all members of the growth hormone(GH) family of genes], as well as proopiomelanocortin (POMC), to be upregulated in D males. Interestingly, a similar pattern of endocrine gene regulation (GH, prolactin, somatolactin, POMC) is observed in Atlantic salmon,such that the expression profile for the early maturing `sneaker' male compared with immature males matches the profile observed in the present study for reproductive D males, suggesting conserved function of these pathways(Aubin-Horth et al., 2005). We found cholecystokinin (CCK) and natriuretic peptide to be upregulated in S males.

Since the activity of pituitary somatotrophes is associated with testis maturity and is stimulated by high levels of GnRH in several fish species(reviewed by Legac et al.,1993; Yu and Peter,1991), we suggest that the increased expression of growth-related genes in D males is probably related to gonad maturation. Somatolactin, which thus far has been found only in teleost fish, is involved in both growth(Forsyth and Wallis, 2002) and color change (Fukamachi and Sugimoto,2004), two plastic traits associated with social dominance in A. burtoni. Importantly, the observed increase in expression of GH is consistent with the previous finding that circulating GH levels are higher in D males (Hofmann et al.,1999b). Additionally, the growth hormone-releasing hormone(GHRH)/GH axis facilitates territorial behavior in A. burtoni(Hofmann et al., 1999b; Trainor and Hofmann, 2006)(note that GHRH was not represented on the array). Finally, antagonists of the neuropeptide somatostatin (which inhibits GH production and release in the pituitary; this gene was not represented on the array) inhibit aggressive behavior in A. burtoni males without affecting sexual behavior(Trainor and Hofmann,2006).

In addition to neuropeptide genes, we found many genes involved in production, maturation, release and reception of neuropeptides and neurotransmitters to be differentially expressed between social phenotypes. For example, secretory granule proteins, such as a homolog of synaptophysin as well as members of the granin family of acidic proteins, were upregulated in D males. These are notably found in a wide variety of endocrine and neuro-endocrine cells (for reviews, see Gerst, 1999; Helle, 2004), and the regulation pattern found in the present study may simply be a consequence of increased neuroendocrine activity in D males.

Fig. 3.

Analysis, at all Gene Ontology (GO) levels, revealed 22 categories that were statistically over- (red) or under- (blue) represented. The three separate GO vocabularies (molecular function, biological process and cellular location) provide overlapping information, as many of the genes are annotated according to each. P-values represent uncorrected results for the hypergeometric test. GO number and name are according to 200605 releases.

Fig. 3.

Analysis, at all Gene Ontology (GO) levels, revealed 22 categories that were statistically over- (red) or under- (blue) represented. The three separate GO vocabularies (molecular function, biological process and cellular location) provide overlapping information, as many of the genes are annotated according to each. P-values represent uncorrected results for the hypergeometric test. GO number and name are according to 200605 releases.

Although few neurotransmitter receptors were present on the array, there is one very intriguing result to report: a GABA-(A) receptor was upregulated in D males whereas at least two subunits of the kainate-type glutamate receptor were upregulated in S males. Clearly, these differences were not predicted nor do we have, at this point, a clear understanding of their specific roles in regulating social dominance. However, the fact that these receptors are affected by (or affect) social phenotypes probably reflects the interconnected nature of neurotransmitter systems and suggests that regulation at all levels of neural circuits may underlie the transition from one social phenotype to another. Interestingly, an antagonistic relationship between GABA and kainate signaling has been suggested in the regulation of reproductive physiology in other systems (e.g. Sagrillo et al.,1996; Chu and Moenter,2005; Clarkson and Herbison,2006; Liu et al.,2006) (see Molecular modules underlying dominance and reproduction).

Molecular functions, biological processes and cellular locations

The GO annotation scheme applied to the cichlid microarray allows rigorous statistical analysis for over- and under-representation of particular molecular functions, biological processes, and cellular locations in genes that are differentially expressed in each male phenotype. Despite the biased nature of the GO terms due to their origin and application in model organisms and directed research, this tool offers a mechanism for statistical analysis of microarray results according to function(Shaw et al., 1999). These terms, unlike specific gene names, avoid experimenter bias and cross-referencing between experiments and even between species and relate experimental results between organisms and platforms. Of the 171 features regulated according to male social phenotype, GO terms could be applied to 85. Analysis at all GO levels revealed 22 categories that are statistically over-or under-represented among the genes that are regulated by social phenotype(P<0.05; compared with their representation among all genes above threshold) (Fig. 3). Using permutation analysis, we determined that only five GO terms were expected to show significant over- or under-representation by chance alone (hypergeometric test P<0.05, only one GO term at P<0.01). Importantly,this unbiased statistical approach confirms our observation discussed earlier that cytoskeleton/structural molecules as well as hormone signaling are upregulated in D males (Table 3). Furthermore, cation/potassium transport pathways appear to be important building blocks for each male phenotype. In D males, the biological processes of GTP binding, iron ion binding and motor activity were significantly enriched, whereas in S males potassium ion transport, regulation of cellular cation transport and ligand-gated ion channel function were activated. Although difficult to interpret directly, this bioinformatic approach results in a considerable data reduction, facilitates comparisons across species and platforms, and provides a framework of hypotheses for future studies on the molecular underpinnings of socially regulated brain function.

Table 3.

Gene Ontology (GO) terms that were significantly over- or under-represented among the genes regulated according to male phenotype

  A. Biological process  
GO: 0006118 Electron transport D>S over-represented P=1.01×10–3 (eight of 55)    
TC400 0.9991 Q5RCZ2_PONPY Similar to secretory vesicle-specific electron transport protein 
HH_NB_Harvardcol_000005749 0.9914 Cloned fragment Brain aromatase 
HH_AB_BRAIN2000_000003171 0.9999 Q803D7_BRARE Sb:cb825 protein (fragment) 
HH_AB_BRAIN2000_000001582 0.9997 Q9DEZ3_ORENI Brain aromatase (fragment) 
HH_AB_BRAIN2000_000001970 0.9925 Q6GLW8_XENLA Hypothetical protein 
TC59 0.9988 Q94SU5_9TELE Cytochrome c oxidase polypeptide II 
TC2 >0.9999 COX1_GADMO NADH dehydrogenase 1 
TC9 >0.9999 CYB_TROMR Cytochrome b family 
GO:0006813 Potassium ion transport S>D over-represented P=4.83×10–4 (five of 25)    
TC15 >0.9999 Q90279_CARAU Kainate receptor alpha subunit 
TC34 0.9972 Q7ZU25_BRARE Atp1a1a.1 protein 
TC16 0.998 Q90278_CARAU Kainate receptor β subunit 
TC17 >0.9999 Q90278_CARAU Kainate receptor β subunit 
TC139 0.9909 Q9DGL2_BRARE Na+/K+ ATPase β subunit 
GO:0000079 Regulation of cell-cycle S>D over-represented P=4.76×10–2 (three of 30)    
TC109 0.9971 SERF2_HUMAN Small EDRK-rich factor 2 
HH_AB_BRAIN2000_000001151 0.9919 RAS_CARAU Ras-like protein 
HH_AB_BRAIN2000_000005191 0.9965 RBM5_HUMAN RNA-binding protein 5 
GO:0006464 Protein modification D>S under-represented P=1.7×10–2 (none of 148)    
GO:0043283 Biopolymer metabolism D>S under-represented P=3.47×10–2 (two of 161)    
TC42 0.9999 PAB1_HUMAN Polyadenylate-binding protein 1 
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
GO:0009056 Catabolism D>S under-represented P=3.38×10–2 (none of 81)    
  B. Cellular component  
GO:0005737 Cytoplasm D>S over-represented P=1.27×10–2 (23 of 383)    
TC40 0.9997 Q8AVI3_XENLA Arbp-prov protein 
TC236 >0.9999 LDHB_FUNHE l-lactate dehydrogenase B chain 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
TC2 >0.9999 COX1_GADMO NADH dehydrogenase 1 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
TC32 0.9998 Q6PC12_BRARE Enolase 1 (α) 
HH_AB_BRAIN2000_000001970 0.9925 Q6GLW8_XENLA Hypothetical protein 
TC348 >0.9999 RB6A_HUMAN Ras-related protein Rab-6A 
TC118 0.9998 Q9DG82_BRARE Prosaposin 
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
HH_AB_BRAIN2000_000002804 0.9995 Q6P4P6_BRARE Zgc:73142 protein (HMP19 protein) 
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC42 0.9999 PAB1_HUMAN Polyadenylate-binding protein 1 
TC399 0.9999 Q802G7_BRARE Selenoprotein M 
HH_AB_BRAIN2000_000003171 0.9999 Q803D7_BRARE Sb:cb825 protein (fragment) 
TC37 0.9944 Q6GQM9_BRARE Enolase 2 
HH_AB_BRAIN2000_000003205 0.9947 Q6GNE7_XENLA MGC82852 protein ribosome receptor activity 
HH_AB_BRAIN2000_000003279 0.9905 Q9YHD1_ANGAN Vacuolar-type H+ transporting ATPase B2 subunit 
TC59 0.9988 Q94SU5_9TELE Cytochrome c oxidase polypeptide II 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
TC9 >0.9999 CYB_TROMR Cytochrome b family 
HH_AB_BRAIN2000_000003055 0.9989 Q5U4S4_XENLA LOC495448 protein 
GO:0005856 Cytoskeleton D>S over-represented P=2.85×10–2 (seven of 77)    
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
HH_AB_BRAIN2000_000001469 0.9967 SEPT7_MOUSE Septin 7 
GO:0005856 Cytoskeleton D>S over-represented P=2.85×10–2 (seven of 77)    
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
TC12 >0.9999 ACT1_FUGRU cytoplasmic 1 (β-actin 1) 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
GO:0005576 Extracellular region D>S over-represented P=5.25×10–6 (15 of 106)    
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
TC49 >0.9999 SOMA_ORENI Somatotropin precursor 
HH_AB_BRAIN2000_000003361 0.99 O70534_RAT ZOG protein 
TC124 0.9998 Q6ZM74_BRARE Similar to synaptophysin 
HH_AB_BRAIN2000_000005345 >0.9999 Q9W7F7_OREMO Proopiomelanocortin 
HH_AB_Stanfordcol_000005677 0.9998 Cloned fragment Arginine vasotocin 
TC290 0.9986 GBRA5_HUMAN GABA receptor α-5 subunit 
TC118 0.9998 Q9DG82_BRARE Prosaposin 
TC194 >0.9999 NEUS_CHICK Neuroserpin precursor 
HH_AB_BRAIN2000_000001923 >0.9999 PRL2_OREMO Prolactin-2 precursor 
TC325 >0.9999 PRL1_OREMO Prolactin-1 precursor 
TC331 >0.9999 SOML_PAROL Somatolactin precursor 
TC310 0.9987 Q8JHD7_HAPBU Arginine vasotocin 
HH_AB_Stanfordcol_000005722 >0.9999 Cloned fragment Arginine vasotocin 
TC134 >0.9999 Q9DEI5_OREMO Glycoprotein α subunit 
GO:0005743 Mitochondrial inner membrane S>D over-represented P=3.07×10–2 (four of 44)    
TC36 0.9929 Q8JHI0_BRARE Solute carrier family 25 member 5 protein 
TC13 >0.9999 Q8M522_HAPBU NADH dehydrogenase 
TC21 0.9998 Q8HLY8_9SMEG NADH-ubiquinone oxidoreductase 
TC4 >0.9999 Q7YDM6_9LABR ATP synthase 6 
GO:0005634 Nucleus D>S under-represented P=2.3×10–2 (five of 275)    
TC12 >0.9999 ACT1_FUGRU Cytoplasmic 1 (β-actin 1) 
TC32 0.9998 Q6PC12_BRARE Enolase 1 (α) 
HH_AB_BRAIN2000_000001906 0.9910 O42448_ONCMY Id2 protein 
HH_AB_BRAIN2000_000003955 0.9919 Q5MBG6_BRARE ETS family transcription factor efr1 
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
  C. Molecular function  
GO:0005525 GTP binding D>S over-represented P=5.23×10–3 (seven of 56)    
HH_AB_BRAIN2000_000001469 0.9967 SEPT7_MOUSE Septin 7 
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
HH_AB_BRAIN2000_000005599 0.9959 Q6PBW0_BRARE Hypothetical protein zgc:73218 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
TC348 >0.9999 RB6A_HUMAN Ras-related protein Rab-6A 
GO:0005179 Hormone activity D>S over-represented P=1.08×10–9 (nine of 16)    
TC49 >0.9999 SOMA_ORENI Somatotropin precursor 
HH_AB_BRAIN2000_000001923 >0.9999 PRL2_OREMO Prolactin-2 precursor 
TC325 >0.9999 PRL1_OREMO Prolactin-1 precursor 
HH_AB_BRAIN2000_000005345 >0.9999 Q9W7F7_OREMO Proopiomelanocortin 
TC331 >0.9999 SOML_PAROL Somatolactin precursor 
TC310 0.9987 Q8JHD7_HAPBU Arginine vasotocin 
TC134 >0.9999 Q9DEI5_OREMO Glycoprotein α subunit 
HH_AB_Stanfordcol_000005722 >0.9999 Cloned fragment Arginine vasotocin 
HH_AB_Stanfordcol_000005677 0.9998 Cloned fragment Arginine vasotocin 
GO:0005506 Iron ion binding D>S over-represented P=1.05×10–2 (five of 35)    
TC400 0.9991 Q5RCZ2_PONPY Similar to secretory vesicle-specific electron transport protein 
GO:0005506 Iron ion binding D>S over-represented P=1.05×10–2 (five of 35)    
HH_NB_Harvardcol_000005749 0.9914 cloned fragment Brain aromatase 
HH_AB_BRAIN2000_000001582 0.9997 Q9DEZ3_ORENI Brain aromatase (fragment) 
TC2 >0.9999 COX1_GADMO NADH dehydrogenase 1 
TC9 >0.9999 CYB_TROMR Cytochrome b family 
GO:0003374 Motor activity D>S over-represented P=1.91×10–2 (three of 15)    
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
TC12 >0.9999 ACT1_FUGRU Cytoplasmic 1 (β-actin 1) 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
GO:001711 Nucleoside-triphosphatase activity D>S over-represented P=2.50×10–2 (seven of 75)    
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
HH_AB_BRAIN2000_000005599 0.9959 Q6PBW0_BRARE Hypothetical protein zgc:73218 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
HH_AB_BRAIN2000_000003279 0.9905 Q9YHD1_ANGAN Vacuolar-type H+ transporting ATPase B2 subunit 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
TC348 >0.9999 RB6A_HUMAN Ras-related protein Rab-6A 
GO:0005198 Structural molecule activity D>S over-represented P=6.53×10–4 (11 of 93)    
TC40 0.9997 Q8AVI3_XENLA Arbp-prov protein 
TC236 >0.9999 LDHB_FUNHE l-lactate dehydrogenase B chain 
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
HH_AB_BRAIN2000_000001469 0.9967 SEPT7_MOUSE Septin 7 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
TC12 >0.9999 ACT1_FUGRU Cytoplasmic 1 (β-actin 1) 
HH_AB_BRAIN2000_000005599 0.9959 Q6PBW0_BRARE Hypothetical protein zgc:73218 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
HH_AB_BRAIN2000_000003205 0.9947 Q6GNE7_XENLA MGC82852 protein ribosome receptor activity 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
GO:0008324 Cation transporter activity S>D over-represented P=2.84×10–4 (nine of 84)    
TC13 >0.9999 Q8M522_HAPBU NADH dehydrogenase 
TC15 >0.9999 Q90279_CARAU Kainate receptor alpha subunit 
TC1 >0.9999 NU1M_BRARE NADH-ubiquinone oxidoreductase chain 1 
TC21 0.9998 Q8HLY8_9SMEG NADH-ubiquinone oxidoreductase 
TC4 >0.9999 Q7YDM6_9LABR ATP synthase 6 
TC34 0.9972 Q7ZU25_BRARE Atp1a1a.1 protein 
TC16 0.998 Q90278_CARAU Kainate receptor β subunit 
TC17 >0.9999 Q90278_CARAU Kainate receptor β subunit 
TC139 0.9909 Q9DGL2_BRARE Na+/K+-ATPase β subunit 
GO:0005230 Extracellular ligand-gated ion channel activity S>D over-represented P=2.78×10–3 (three of 11)    
TC15 >0.9999 Q90279_CARAU Kainate receptor α subunit 
TC16 0.998 Q90278_CARAU Kainate receptor β subunit 
TC17 >0.9999 Q90278_CARAU Kainate receptor β subunit 
GO:0003676 Nucleic acid binding D>S under-represented P=3.47×10–2 (six of 280)    
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC32 0.9998 Q6PC12_BRARE Enolase 1 (α) 
TC40 0.9997 Q8AVI3_XENLA Arbp-prov protein 
TC42 0.9999 PAB1_HUMAN Polyadenylate-binding protein 1 
HH_AB_BRAIN2000_000003955 0.9919 Q5MBG6_BRARE ETS family transcription factor efr1 
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
GO:0016740 Transferase activity D>S under-represented P=1.25×10–2 (one of 150)    
HH_AB_BRAIN2000_000004040 0.9925 UCKL1_HUMAN Uridine/cytidine kinase-like 1 
GO:0008270 Zinc ion binding D>S under-represented P=1.13×10–2 (none of 106)    
  A. Biological process  
GO: 0006118 Electron transport D>S over-represented P=1.01×10–3 (eight of 55)    
TC400 0.9991 Q5RCZ2_PONPY Similar to secretory vesicle-specific electron transport protein 
HH_NB_Harvardcol_000005749 0.9914 Cloned fragment Brain aromatase 
HH_AB_BRAIN2000_000003171 0.9999 Q803D7_BRARE Sb:cb825 protein (fragment) 
HH_AB_BRAIN2000_000001582 0.9997 Q9DEZ3_ORENI Brain aromatase (fragment) 
HH_AB_BRAIN2000_000001970 0.9925 Q6GLW8_XENLA Hypothetical protein 
TC59 0.9988 Q94SU5_9TELE Cytochrome c oxidase polypeptide II 
TC2 >0.9999 COX1_GADMO NADH dehydrogenase 1 
TC9 >0.9999 CYB_TROMR Cytochrome b family 
GO:0006813 Potassium ion transport S>D over-represented P=4.83×10–4 (five of 25)    
TC15 >0.9999 Q90279_CARAU Kainate receptor alpha subunit 
TC34 0.9972 Q7ZU25_BRARE Atp1a1a.1 protein 
TC16 0.998 Q90278_CARAU Kainate receptor β subunit 
TC17 >0.9999 Q90278_CARAU Kainate receptor β subunit 
TC139 0.9909 Q9DGL2_BRARE Na+/K+ ATPase β subunit 
GO:0000079 Regulation of cell-cycle S>D over-represented P=4.76×10–2 (three of 30)    
TC109 0.9971 SERF2_HUMAN Small EDRK-rich factor 2 
HH_AB_BRAIN2000_000001151 0.9919 RAS_CARAU Ras-like protein 
HH_AB_BRAIN2000_000005191 0.9965 RBM5_HUMAN RNA-binding protein 5 
GO:0006464 Protein modification D>S under-represented P=1.7×10–2 (none of 148)    
GO:0043283 Biopolymer metabolism D>S under-represented P=3.47×10–2 (two of 161)    
TC42 0.9999 PAB1_HUMAN Polyadenylate-binding protein 1 
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
GO:0009056 Catabolism D>S under-represented P=3.38×10–2 (none of 81)    
  B. Cellular component  
GO:0005737 Cytoplasm D>S over-represented P=1.27×10–2 (23 of 383)    
TC40 0.9997 Q8AVI3_XENLA Arbp-prov protein 
TC236 >0.9999 LDHB_FUNHE l-lactate dehydrogenase B chain 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
TC2 >0.9999 COX1_GADMO NADH dehydrogenase 1 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
TC32 0.9998 Q6PC12_BRARE Enolase 1 (α) 
HH_AB_BRAIN2000_000001970 0.9925 Q6GLW8_XENLA Hypothetical protein 
TC348 >0.9999 RB6A_HUMAN Ras-related protein Rab-6A 
TC118 0.9998 Q9DG82_BRARE Prosaposin 
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
HH_AB_BRAIN2000_000002804 0.9995 Q6P4P6_BRARE Zgc:73142 protein (HMP19 protein) 
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC42 0.9999 PAB1_HUMAN Polyadenylate-binding protein 1 
TC399 0.9999 Q802G7_BRARE Selenoprotein M 
HH_AB_BRAIN2000_000003171 0.9999 Q803D7_BRARE Sb:cb825 protein (fragment) 
TC37 0.9944 Q6GQM9_BRARE Enolase 2 
HH_AB_BRAIN2000_000003205 0.9947 Q6GNE7_XENLA MGC82852 protein ribosome receptor activity 
HH_AB_BRAIN2000_000003279 0.9905 Q9YHD1_ANGAN Vacuolar-type H+ transporting ATPase B2 subunit 
TC59 0.9988 Q94SU5_9TELE Cytochrome c oxidase polypeptide II 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
TC9 >0.9999 CYB_TROMR Cytochrome b family 
HH_AB_BRAIN2000_000003055 0.9989 Q5U4S4_XENLA LOC495448 protein 
GO:0005856 Cytoskeleton D>S over-represented P=2.85×10–2 (seven of 77)    
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
HH_AB_BRAIN2000_000001469 0.9967 SEPT7_MOUSE Septin 7 
GO:0005856 Cytoskeleton D>S over-represented P=2.85×10–2 (seven of 77)    
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
TC12 >0.9999 ACT1_FUGRU cytoplasmic 1 (β-actin 1) 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
GO:0005576 Extracellular region D>S over-represented P=5.25×10–6 (15 of 106)    
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
TC49 >0.9999 SOMA_ORENI Somatotropin precursor 
HH_AB_BRAIN2000_000003361 0.99 O70534_RAT ZOG protein 
TC124 0.9998 Q6ZM74_BRARE Similar to synaptophysin 
HH_AB_BRAIN2000_000005345 >0.9999 Q9W7F7_OREMO Proopiomelanocortin 
HH_AB_Stanfordcol_000005677 0.9998 Cloned fragment Arginine vasotocin 
TC290 0.9986 GBRA5_HUMAN GABA receptor α-5 subunit 
TC118 0.9998 Q9DG82_BRARE Prosaposin 
TC194 >0.9999 NEUS_CHICK Neuroserpin precursor 
HH_AB_BRAIN2000_000001923 >0.9999 PRL2_OREMO Prolactin-2 precursor 
TC325 >0.9999 PRL1_OREMO Prolactin-1 precursor 
TC331 >0.9999 SOML_PAROL Somatolactin precursor 
TC310 0.9987 Q8JHD7_HAPBU Arginine vasotocin 
HH_AB_Stanfordcol_000005722 >0.9999 Cloned fragment Arginine vasotocin 
TC134 >0.9999 Q9DEI5_OREMO Glycoprotein α subunit 
GO:0005743 Mitochondrial inner membrane S>D over-represented P=3.07×10–2 (four of 44)    
TC36 0.9929 Q8JHI0_BRARE Solute carrier family 25 member 5 protein 
TC13 >0.9999 Q8M522_HAPBU NADH dehydrogenase 
TC21 0.9998 Q8HLY8_9SMEG NADH-ubiquinone oxidoreductase 
TC4 >0.9999 Q7YDM6_9LABR ATP synthase 6 
GO:0005634 Nucleus D>S under-represented P=2.3×10–2 (five of 275)    
TC12 >0.9999 ACT1_FUGRU Cytoplasmic 1 (β-actin 1) 
TC32 0.9998 Q6PC12_BRARE Enolase 1 (α) 
HH_AB_BRAIN2000_000001906 0.9910 O42448_ONCMY Id2 protein 
HH_AB_BRAIN2000_000003955 0.9919 Q5MBG6_BRARE ETS family transcription factor efr1 
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
  C. Molecular function  
GO:0005525 GTP binding D>S over-represented P=5.23×10–3 (seven of 56)    
HH_AB_BRAIN2000_000001469 0.9967 SEPT7_MOUSE Septin 7 
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
HH_AB_BRAIN2000_000005599 0.9959 Q6PBW0_BRARE Hypothetical protein zgc:73218 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
TC348 >0.9999 RB6A_HUMAN Ras-related protein Rab-6A 
GO:0005179 Hormone activity D>S over-represented P=1.08×10–9 (nine of 16)    
TC49 >0.9999 SOMA_ORENI Somatotropin precursor 
HH_AB_BRAIN2000_000001923 >0.9999 PRL2_OREMO Prolactin-2 precursor 
TC325 >0.9999 PRL1_OREMO Prolactin-1 precursor 
HH_AB_BRAIN2000_000005345 >0.9999 Q9W7F7_OREMO Proopiomelanocortin 
TC331 >0.9999 SOML_PAROL Somatolactin precursor 
TC310 0.9987 Q8JHD7_HAPBU Arginine vasotocin 
TC134 >0.9999 Q9DEI5_OREMO Glycoprotein α subunit 
HH_AB_Stanfordcol_000005722 >0.9999 Cloned fragment Arginine vasotocin 
HH_AB_Stanfordcol_000005677 0.9998 Cloned fragment Arginine vasotocin 
GO:0005506 Iron ion binding D>S over-represented P=1.05×10–2 (five of 35)    
TC400 0.9991 Q5RCZ2_PONPY Similar to secretory vesicle-specific electron transport protein 
GO:0005506 Iron ion binding D>S over-represented P=1.05×10–2 (five of 35)    
HH_NB_Harvardcol_000005749 0.9914 cloned fragment Brain aromatase 
HH_AB_BRAIN2000_000001582 0.9997 Q9DEZ3_ORENI Brain aromatase (fragment) 
TC2 >0.9999 COX1_GADMO NADH dehydrogenase 1 
TC9 >0.9999 CYB_TROMR Cytochrome b family 
GO:0003374 Motor activity D>S over-represented P=1.91×10–2 (three of 15)    
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
TC12 >0.9999 ACT1_FUGRU Cytoplasmic 1 (β-actin 1) 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
GO:001711 Nucleoside-triphosphatase activity D>S over-represented P=2.50×10–2 (seven of 75)    
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
HH_AB_BRAIN2000_000005599 0.9959 Q6PBW0_BRARE Hypothetical protein zgc:73218 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
HH_AB_BRAIN2000_000003279 0.9905 Q9YHD1_ANGAN Vacuolar-type H+ transporting ATPase B2 subunit 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
TC348 >0.9999 RB6A_HUMAN Ras-related protein Rab-6A 
GO:0005198 Structural molecule activity D>S over-represented P=6.53×10–4 (11 of 93)    
TC40 0.9997 Q8AVI3_XENLA Arbp-prov protein 
TC236 >0.9999 LDHB_FUNHE l-lactate dehydrogenase B chain 
TC291 0.9987 ACT2_FUGRU Actin, cytoplasmic 2 (β-actin 2) 
HH_AB_BRAIN2000_000001469 0.9967 SEPT7_MOUSE Septin 7 
TC19 >0.9999 TBA6_HUMAN α-tubulin 6 
TC12 >0.9999 ACT1_FUGRU Cytoplasmic 1 (β-actin 1) 
HH_AB_BRAIN2000_000005599 0.9959 Q6PBW0_BRARE Hypothetical protein zgc:73218 
TC28 >0.9999 Q9DFT6_9PERC β-tubulin 
HH_AB_BRAIN2000_000001815 0.993 Q9I921_CARAU β-actin 
HH_AB_BRAIN2000_000003205 0.9947 Q6GNE7_XENLA MGC82852 protein ribosome receptor activity 
TC48 0.9934 Q9DFT6_9PERC β-tubulin 
GO:0008324 Cation transporter activity S>D over-represented P=2.84×10–4 (nine of 84)    
TC13 >0.9999 Q8M522_HAPBU NADH dehydrogenase 
TC15 >0.9999 Q90279_CARAU Kainate receptor alpha subunit 
TC1 >0.9999 NU1M_BRARE NADH-ubiquinone oxidoreductase chain 1 
TC21 0.9998 Q8HLY8_9SMEG NADH-ubiquinone oxidoreductase 
TC4 >0.9999 Q7YDM6_9LABR ATP synthase 6 
TC34 0.9972 Q7ZU25_BRARE Atp1a1a.1 protein 
TC16 0.998 Q90278_CARAU Kainate receptor β subunit 
TC17 >0.9999 Q90278_CARAU Kainate receptor β subunit 
TC139 0.9909 Q9DGL2_BRARE Na+/K+-ATPase β subunit 
GO:0005230 Extracellular ligand-gated ion channel activity S>D over-represented P=2.78×10–3 (three of 11)    
TC15 >0.9999 Q90279_CARAU Kainate receptor α subunit 
TC16 0.998 Q90278_CARAU Kainate receptor β subunit 
TC17 >0.9999 Q90278_CARAU Kainate receptor β subunit 
GO:0003676 Nucleic acid binding D>S under-represented P=3.47×10–2 (six of 280)    
TC18 >0.9999 Q8UW60_ORENI Elongation factor 1a 
TC32 0.9998 Q6PC12_BRARE Enolase 1 (α) 
TC40 0.9997 Q8AVI3_XENLA Arbp-prov protein 
TC42 0.9999 PAB1_HUMAN Polyadenylate-binding protein 1 
HH_AB_BRAIN2000_000003955 0.9919 Q5MBG6_BRARE ETS family transcription factor efr1 
HH_AB_BRAIN2000_000004166 0.9992 Q6NX86_BRARE High mobility group box 1 
GO:0016740 Transferase activity D>S under-represented P=1.25×10–2 (one of 150)    
HH_AB_BRAIN2000_000004040 0.9925 UCKL1_HUMAN Uridine/cytidine kinase-like 1 
GO:0008270 Zinc ion binding D>S under-represented P=1.13×10–2 (none of 106)    

Each of the three independent ontologies is given. Significance is determined to be P<0.05 for the uncorrected hypergeometric P-value. TC number refers to `tentative contig' number assigned in release 1.0 A. burtoni cichlid gene index(compbio.dfci.harvard.edu)

Molecular modules underlying dominance and reproduction

The notion that biological entities (e.g. cognitive tasks, developmental programs, neural circuits, metabolic pathways) operate as functional and discrete (i.e. largely non-overlapping) units, or modules, is not new(Fodor, 1983; Redies and Puelles, 2001; Schlosser and Wagner, 2004; Op de Beeck et al., 2008). In molecular systems biology, a module can simply be defined as a set of co-regulated genes or proteins (Segal et al., 2004). Many such modules may serve as building blocks for the assembly of more complex processes. To date, most studies in this area have primarily been concerned with molecular and cellular networks and pathways in simple unicellular systems (Hartwell et al., 1999; Wolf and Arkin,2003). However, the ultimate challenge in the biology of complex systems is the integration across many levels of biological organization, from molecules to whole organisms, in an ever-changing environment. In the following, we make an initial attempt at such an integration of genomic data with physiological and behavioral phenotypes to provide a comprehensive conceptual framework for understanding phenotypically plastic traits. Specifically, we examine male dominance phenotypes in terms of molecular modules of socially controlled traits such as aggression, territoriality,reproduction and growth.

As we also obtained neural expression profiles of brooding females, we examined variation in transcript levels in relation to sex. 569 genes were regulated according to sex (316 male upregulated, 253 female upregulated)(Table 1), a number that is far greater than that observed for dominance phenotypes within males (171 genes or∼5%). This suggests that while the switch between dominance phenotypes is multifaceted (see Fig. 1), the difference in brain gene expression profiles between males and females is even more dramatic, affecting ∼16% of the genes on the array. Interestingly,this proportion is comparable to the 15% observed difference between alternative mating tactics in males of Atlantic salmon(Aubin-Horth et al., 2005). A considerable proportion of the 171 genes associated with male social phenotypes was also regulated according to sex (supplementary material Tables S1 and S2). However, among the social status-regulated genes, there were similar proportions of female-enriched and male-enriched genes(Fig. 4). Among the 87 dominant upregulated genes, were 11 female-enriched and 20 male-enriched genes(hypergeometric test; P=0.09), and among the 84 subordinate upregulated genes there were 16 female-enriched genes and 21 male-enriched genes (hypergeometric test; P=0.13). In other words, S males do not appear to be molecularly feminized nor are D males simply `super-males'. This result makes sense in light of the reproductive state of these animals. Although female behavior is, in many ways, similar to that of S males, both brooding females and D males are reproductively active. Furthermore, just as the metabolic demands of maintaining a territory are associated with reduced growth in D males (Hofmann et al.,1999a), mouth-brooding females starve while incubating their offspring and exhibit a marked reduction in body mass(Mrowk, 1984). By contrast, S males do not reproduce, and metabolic energy is directed toward growth(Hofmann et al., 1999a). Thus,the 11 genes (including synaptophysin, neuroserpin and GABA-receptor)upregulated in both D males and brooding females may be part of a module facilitating reproduction and/or reducing growth and may not necessarily be involved in dominance behavior per se.

Eighteen of the 316 male-enriched genes were even more highly expressed in D males compared with S males. We suggest that these genes, which include several structural proteins (actin, tubulin, ELF-1) and hormones (GnRH1, AVT,somatolactin, POMC) (supplementary material Table S1) are part of a`super-male' module, which probably plays a role in the aggressive and/or sexual behaviors typical of the D male(Fig. 4). Even though we know from previous studies using other approaches that the absolute levels of GnRH1(Francis et al., 1993; White et al., 2002) and GH(Hofmann et al., 1999b) are higher in D males, it is important to note that cDNA microarrays measure only relative transcript levels. Future studies will have to determine whether these genes must be activated in order to produce a dominant phenotype or whether their expression needs to be suppressed for the subordinate phenotype to occur. Following a similar rationale, we suggest that the 21 male-enriched genes (including neuroD and kainate receptor) that were further upregulated in S males relative to D males represent a masculinizing module in some sense. In this case, reduced expression would be critical for the manifestation of behaviors linked to aspects of social dominance and reproduction. Consequently, the 21 genes in this masculinizing module could be considered to be complementary to the super-male module. Finally, as brooding females and S males all display submissive behaviors and very little aggression, it is tempting to suggest that the 16 genes (including natriuretic peptide)expressed in a similar pattern in these two phenotypes are part of a`subordination module' (Fig. 4).

The opposing pattern of regulation for receptor expression in two classic neurotransmitter systems that we observed between male social phenotypes (see above) is also maintained in relation to sex. The gene encoding a GABA-(A)receptor was upregulated in D males and in females whereas the kainate-type glutamate receptor was upregulated in males in general and particularly in S males. There is a wealth of research that ties the GABA-(A) receptor to the regulation of the hypothalamic–pituitary–gonadal axis viaintegration of steroid feedback to GnRH neurons (for a review, see Sagrillo et al., 1996). In mammals, GABA has mixed inhibitory and excitatory effects on the release of GnRH, due in part to a developmental switch from GABA-(A) receptor depolarization to hyperpolarization(Clarkson and Herbison, 2006). In fish, depending on the species, GABA has either excitatory or inhibitory effects on GnRH release (Trudeau et al.,2000). Interestingly, glutamate-controlled GABA release has been implicated in GnRH regulation (Chu and Moenter, 2005; Clarkson and Herbison, 2006). The kainate system has also been proposed to underlie observed sex differences in the mechanisms of the neural glucocorticoid/stress response. Female mice show less atrophy of hippocampal neurons in response to elevated glucocorticoid levels, possibly due to the increased expression of NMDA, AMPA and kainate glutamate receptor sub-types(Liu et al., 2006). In A. burtoni, the increased kainate receptor expression seen in S males could similarly provide a neuroprotective effect against the elevated cortisol levels seen in S males during specific social situations(Fox et al., 1997). Future pharmacological and neurohistochemical experiments will elucidate the mechanistic interactions of these neurotransmitter systems in relation to sex and social behavior.

Fig. 4.

Venn diagram depicting the relationship of sexually regulated and socially regulated genes. These relationships subdivide the gene classes to indicate modules of gene expression that potentially underlie reproduction (orange),submissive behavior (lavender), and super-male dominance (green) and opposing super-male dominance (teal). Numbers indicate total unique sequences and unsequenced array features. Gene names given represent best-hit blast annotation for available sequences. The Venn diagram indicates regulation at a BPP of >0.99, (the specific BPP, down to 0.80, for regulation is indicated in supplementary material Tables S1 and S2)

Fig. 4.

Venn diagram depicting the relationship of sexually regulated and socially regulated genes. These relationships subdivide the gene classes to indicate modules of gene expression that potentially underlie reproduction (orange),submissive behavior (lavender), and super-male dominance (green) and opposing super-male dominance (teal). Numbers indicate total unique sequences and unsequenced array features. Gene names given represent best-hit blast annotation for available sequences. The Venn diagram indicates regulation at a BPP of >0.99, (the specific BPP, down to 0.80, for regulation is indicated in supplementary material Tables S1 and S2)

Taken together, the molecular systems analysis in the present study supports the notion that transcript patterns may indeed be organized in a modular fashion and can be strongly associated with behavioral and/or physiological traits associated with social phenotypes or sex in either concordant or contrasting ways. Additionally, we can exploit expression variation between phenotypes for tentatively annotating gene function and predicting functional roles of these genes.

Individual variation in gene expression

To appreciate the importance of variation in gene expression for phenotypic plasticity, we need to evaluate the expression differences between individuals. To determine the extent to which individuals of the same phenotype differ in their expression profiles, we estimated transcript levels for each individual separately. Individual profiles were then clustered for similarity according to the estimated transcript levels using the gene list that had been identified as significantly regulated between any two phenotypes(BPP>0.99) (Euclidian distance matrix based on resampling for bootstrap confidence levels) (Fig. 5).

Similarities in expression profiles of individuals

Sex was a strong factor in the clustering of individual expression profiles, as five of the six females clustered separately from all males. While males of the same social phenotype tended to have similar brain expression profiles, two D individuals (males 4 and 5) clustered with the S males. The analysis of behavioral data of these males did not suggest any obvious difference from other D males (Fig. 2). It is unlikely (since we controlled for age and size and based on our observations) that these two D males were experiencing dominance for the first time. Also, at the behavior level, these two D males were not similar to each other. While male 5 did show the fewest chases and the most border threats (Fig. 2) of all D males, the behaviors of male 4 were not in any way different from those of the other D males. While it is possible that there is more than one molecular substrate for constructing and maintaining a dominant male phenotype, the expression profiles of these two D males are no more similar to each other than to S males.

When individual profiles were clustered according to estimated transcript levels using all genes on the array that passed filtering for every fish,rather than only those genes that showed significant regulation, similar results were obtained (data not shown). While the bootstrapped confidence values were lower, the same four D males formed a cluster, suggesting that gene regulation according to sex and social status account for the greatest variation among individuals in this study. Principal component analysis (PCA)corroborated this conclusion (supplementary material Fig. S1).

Variation among individuals within and between male phenotypes

While, overall, the clustering resulted in a robust separation of the three phenotypes according to sex and social status, the individual variation,already apparent at this level, prompted further inquiry into expression variation between individuals of the same phenotype. In voles, expression patterns of oxytocin receptor for females and vasopressin receptor for males is correlated with individual variation in social and anxiety-related behaviors (Olazabal and Young,2006). Similarly, male mice show individual variation in estrogen receptor distribution that correlates with aggressive behavior(Trainor et al., 2006). Previous studies in A. burtoni and other cichlid species have also reported strong covariation patterns within a social phenotype between the expression of candidate genes and specific phenotype measures [somatostatin correlated with aggression (Trainor and Hofmann, 2006); AVT correlated with hormone titers(Aubin-Horth et al., 2007)]. We therefore asked whether significant differences in gene expression existed between individuals of the same phenotype. Furthermore, we tested whether significant differences in gene expression between individuals of the same phenotype could be found among those genes that are differentially regulated between phenotypes.

Fig. 5.

Hierarchical clustering of phenotypes (A) and individuals (B) (based on expression profiles for the genes regulated at level of phenotype). The genes have been ordered according to k-means clustering. The heatmaps (red,upregulated, green, downregulated) show estimated gene expression levels. Heat values are relative only within, not across, genes. The numbers identify individual males consistent with Fig. 2. Confidence values at the nodes were obtained by bootstrap analysis (1000 permutations with resampling). The colors indicate S,subordinate males in teal; D, dominant males in gold; and F, brooding females in pink.

Fig. 5.

Hierarchical clustering of phenotypes (A) and individuals (B) (based on expression profiles for the genes regulated at level of phenotype). The genes have been ordered according to k-means clustering. The heatmaps (red,upregulated, green, downregulated) show estimated gene expression levels. Heat values are relative only within, not across, genes. The numbers identify individual males consistent with Fig. 2. Confidence values at the nodes were obtained by bootstrap analysis (1000 permutations with resampling). The colors indicate S,subordinate males in teal; D, dominant males in gold; and F, brooding females in pink.

In order to investigate the degree of individual variation in gene expression, we determined the number of genes that varied significantly in expression between individuals of the same phenotype and compared it with the number that varied between individuals of different phenotypes. Although statistical power was lower because we had only four technical replicates per individual as opposed to six biological replicates in the analyses above(Clark and Townsend, 2007), we were able to measure the average number of genes significantly regulated(BPP>0.99) in each possible pairwise comparison of two individuals within and between phenotypes. For intra-phenotype variations among males, a mean(±s.e.m.) of 82.8±4.5 (2.3% of all genes analyzed) genes varied between any two D males, while a mean of 92.3±4.3 (2.6%) genes varied between any two S males. The identity of these genes was substantially different for each pairwise comparison, such that 38% of all the array features varied in at least one intra-phenotype comparison. Interestingly, for inter-phenotype variation between individual males, a mean of 132±8.4(3.6%) genes varied between any two males of differing social phenotype, which was not significantly different from the intra-phenotype variation(t-test, P=0.13). Although it is difficult to set an equivalent threshold for significant variation between individuals and between phenotypes, this high degree of individual variation is consistent with other studies that have examined genome-scale individual differences in gene expression in order to study the molecular basis of natural variation. Whitehead and Crawford found that 69% of the metabolic pathway genes showed significant variation between individuals within a population, while only 12%were significantly regulated between populations adapted to different temperatures (Whitehead and Crawford,2006). Similarly, in yeast, up to 50% of the expressed genes show significantly different levels of expression among individual strains(Brem et al., 2002). Not only do we find a similar number of genes to be regulated between individuals of the same or different phenotypes, we also find no significant difference in coefficient of variation of expression level for sets of these individually regulated genes. This result indicates that absolute gene expression levels vary between individuals of the same phenotype as much as between phenotypes.

While low variation within a phenotype for those genes that define that phenotype can be expected, an alternative hypothesis posits that those genes that define the social phenotype vary across individuals displaying that phenotype in a manner associated with variation in physiological and behavioral traits (e.g. Trainor and Hofmann, 2006; Aubin-Horth et al., 2007; Cummings et al.,2008). In support of this notion, we found a statistically significant over-representation of intra-phenotype regulated genes among those that were regulated by social status (Table 4). Approximately 72% of the 87 genes upregulated in the dominant phenotype were also significantly regulated among individuals within a phenotype (41 among D and 39 among S, 24 of which are shared). Similarly, 64%of the 84 genes that were upregulated in the S phenotype were also significantly regulated among individuals (47 among D and 52 among S, 38 of which are shared). This variation cannot be explained by technical variation in array hybridizations: for a given animal, an array feature must show consistent results across four dye-reversed hybridizations before it can be identified as regulated across individuals according to our statistical analysis. Rather, the results of the present study show that even considerable and potentially important within-phenotype variation in gene expression can give rise to reliable and readily identifiable between-phenotype differences. Future integrative studies will help determine whether the observed variation between individuals is caused by, or causes, subtle phenotypic differences or represents a dramatic, alternative molecular mechanism for constructing the same phenotype.

Table 4.

Enrichment of intra-phenotype gene expression variation is significant for genes that are both upregulated and downregulated for that phenotype in comparison with the other phenotype

Upregulated by phenotype
DominantSubordinate
Regulated by individual  87 84 
   w/in dominant 824 41 (P=2.9×10–747 (P=2.76×10–11
   w/in subordinate 925 39 (P=4.3×10–552 (P=1.2×10–12
Upregulated by phenotype
DominantSubordinate
Regulated by individual  87 84 
   w/in dominant 824 41 (P=2.9×10–747 (P=2.76×10–11
   w/in subordinate 925 39 (P=4.3×10–552 (P=1.2×10–12

Hypergeometric test (3598 total genes analyzed)

CONCLUSIONS

By including candidate genes on the microarray, we have validated the discovery-based approach that identifies the gene expression patterns for reproductively active D males and reproductively suppressed S males as well as for brooding females. While the regulation of neuroendocrine genes was predicted from previous research, the unexpected and novel discovery of opposing roles for two classic neurotransmitter systems opens exciting new avenues for future research. Furthermore, the increased neuronal remodeling activity in D males suggests that the known neuroanatomical changes observed in pre-optic GnRH and somatostatin neurons may extend to other cell types and/or additional brain regions. Using the GO framework to interrogate the data for statistical significance reinforces the importance of hormonal regulation and highlights the hitherto under-appreciated roles of cytoskeletal components in addition to neurotransmitter pathways. Interestingly, expression profiles vary among individuals within a male phenotype for roughly the same number of genes and with similar magnitudes of transcript abundance as seen between male phenotypes. Specifically, there is a surprisingly high level of variation among those genes that molecularly define the very phenotypes in the first place. Taken together, our genome-scale analysis of molecular systems in the brain has identified complex patterns of gene expression that are associated with a socially regulated switch in behavioral phenotype.

LIST OF ABBREVIATIONS

     
  • AMPA

    alpha-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid

  •  
  • AVP

    arginine vasopressin

  •  
  • AVT

    arginine vasotocin

  •  
  • BPP

    Bayesian posterior probability

  •  
  • CF

    condition factor

  •  
  • CCK

    cholecystokinin

  •  
  • DAG

    directed acyclic graphs

  •  
  • DI

    dominance index

  •  
  • EST

    expressed sequence tags

  •  
  • FDR

    false discovery rate

  •  
  • GABA

    gamma-aminobutyric acid

  •  
  • GH

    growth hormone

  •  
  • GHRH

    growth hormone-releasing hormone

  •  
  • GnRH

    gonadotrophin-releasing hormone

  •  
  • GO

    gene ontology

  •  
  • GR

    growth rate

  •  
  • GSI

    gonadosomatic index

  •  
  • IEA

    inferred electronic annotation

  •  
  • ISS

    inferred sequence similiarity

  •  
  • KT

    ketotestosterone

  •  
  • NMDA

    N-methyl-D-aspartic acid

  •  
  • POA

    pre-optic area

  •  
  • POMC

    proopiomelanocortin

  •  
  • RI

    reproductive index

  •  
  • T

    testosterone

  •  
  • TC

    tentative contigs

Acknowledgements

All authors contributed to the design, implementation, analysis and presentation of this experimental data set. We are grateful to Josiah Altschuler and Melinda Snitow for animal care, Sarah Annis for assistance with the behavioral experiments, Christian Landry for programming and photography,Amir Karger for bioinformatics advice, and Victoria Zero, Brian Dias, Lin Huffman and Kim Hoke for comments on earlier versions of this manuscript. We thank the H ofmann laboratory and the members of the Bauer Center for Genomics Research for stimulating discussions. This research was supported by an NRSA post-doctoral fellowship to S.C.P.R., a postdoctoral fellowship from FQRNT(Fonds Québécois de la Recherche sur la Nature et les Technologies) and a postdoctoral fellowship from the Natural Science and Engineering Research Council of Canada (NSERC) to N.A.H., and by National Institutes of Health grant NIGMS GM068763, the Bauer Center for Genomics Research (H.A.H.).

References

Amano, M., Hyodo, S., Kitamura, S., Ikuta, K., Suzuki, Y.,Urano, A. and Aida, K. (
1995
). Short photoperiod accelerates preoptic and ventral telencephalic salmon GnRH synthesis and precocious maturation in underyearling male masu salmon.
Gen. Comp. Endocrinol.
99
,
22
-27.
Arnold, A. P. (
1992
). Hormonally-induced alterations in synaptic organization in the adult nervous-system.
Exp. Gerontol.
27
,
99
-110.
Au, T. M., Anna, K. G. and Fernald, R. D.(
2006
). Differential social regulation of two pituitary gonadotropin-releasing hormone receptors.
Behav. Brain Res.
170
,
342
-346.
Aubin-Horth, N., Landry, C. R., Letcher, B. H. and Hofmann, H. A. (
2005
). Alternative life histories shape brain gene expression profiles in males of the same population.
Proc. R. Soc. Lond. B
272
,
1655
-1662.
Aubin-Horth, N., Desjardins, J. K., Martei, Y. M., Balshine, S. and Hofmann, H. A. (
2007
). Masculinized dominant females in a cooperatively breeding species.
Mol. Ecol.
16
,
1349
-1358.
Black, M. P., Balthazart, J., Baillien, M. and Grober, M. S.(
2005
). Socially induced and rapid increases in aggression are inversely related to brain aromatase activity in a sex-changing fish, Lythrypnus dalli.
Proc. Biol. Sci.
272
,
2435
-2440.
Brem, R. B., Yvert, G., Clinton, R. and Kruglyak, L.(
2002
). Genetic dissection of transcriptional regulation in budding yeast.
Science
296
,
752
-755.
Buckley, B. A., Gracey, A. Y. and Somero, G. N.(
2006
). The cellular response to heat stress in the goby Gillichthys mirabilis: a cDNA microarray and protein-level analysis.
J. Exp. Biol.
209
,
2660
-2677.
Burmeister, S. S., Jarvis, E. D. and Fernald, R. D.(
2005
). Rapid behavioral and genomic responses to social opportunity.
PLoS Biol.
3
,
1996
-2004.
Burmeister, S. S., Kailasanath, V. and Fernald, R. D.(
2007
). Social dominance regulates androgen and estrogen receptor gene expression.
Horm. Behav.
51
,
164
-170.
Carneiro, L. A., Oliveira, R. F., Canario, A. V. M. and Grober,M. S. (
2003
). The effect of arginine vasotocin on courtship behaviour in a blenniid fish with alternative reproductive tactics.
Fish Physiol. Biochem.
28
,
241
-243.
Churchill, G. A. (
2002
). Fundamentals of experimental design for cDNA microarrays.
Nat. Genet.
32
,
490
-495.
Chu, Z. G. and Moenter, S. M. (
2005
). Endogenous activation of metabotropic glutamate receptors modulates GABAergic transmission to gonadotropin-releasing hormone neurons and alters their firing rate: a possible local feedback circuit.
J. Neurosci.
25
,
5740
-5749.
Clark, T. A. and Townsend, J. P. (
2007
). Quantifying variation in gene expression.
Mol. Ecol.
16
,
2613
-2616.
Clarkson, J. and Herbison, A. E. (
2006
). Development of GABA and glutamate signaling at the GnRH neuron in relation to puberty.
Mol. Cell. Endocrinol.
254
,
32
-38.
Cummings, M. E., Larkins-Ford, J., Reilly, C. R. L., Wong, R. Y., Ramsey, M. and Hofmann, H. A. (
2008
). Sexual and social stimuli elicit rapid and contrasting genomic responses.
Proc. Biol. Sci.
275
,
393
-402.
Davis, M. R. and Fernald, R. D. (
1990
). Social control of neuronal soma size.
J. Neurobiol.
21
,
1180
-1188.
Dawson, A., King. V. M., Bentley, G. E. and Ball, G. F.(
2001
). Photoperiodic control of seasonality in birds.
J. Biol. Rhythms
16
,
365
-380.
Devoogd, T. and Nottebohm, F. (
1981
). Gonadal-hormones induce dendritic growth in the adult avian brain.
Science
214
,
202
-204.
Fernald, R. D. (
1977
). Quantitative behavioural observations of Haplochromis burtoni under semi-natural conditions.
Anim. Behav.
25
,
643
-653.
Fernald, R. D. (
2002
). Social regulation of the brain: sex, size and status. In
Genetics and Biology of Sex Determination
no. 244 (ed. Novartis Foundation), pp.
169
-186. Chichester: John Wiley and Sons.
Fernald, R. D. (
2004
). Social influences on the brain.
Horm. Behav.
46
,
129
-130.
Fernald, R. D. and Hirata, N. R. (
1977a
). Field study of Haplochromis burtoni: quantitative behavioral observations.
Anim. Behav.
25
,
964
-975.
Fernald, R. D. and Hirata, N. R. (
1977b
). Field study of Haplochromis burtoni: habitats and co-habitants.
Environ. Biol. Fishes
2
,
299
-308.
Fodor, J. A. (
1983
).
The Modularity of Mind
. Cambridge, MA: The MIT Press.
Forsyth, I. A. and Wallis, M. (
2002
). Growth hormone and prolactin-molecular and functional evolution.
J. Mammary Gland Biol. Neoplasia
7
,
291
-312.
Fox, H. E., White, S. A., Kao, M. H. F. and Fernald, R. D.(
1997
). Stress and dominance in a social fish.
J. Neurosci.
17
,
6463
-6469.
Francis, R. C., Soma, K. and Fernald, R. D.(
1993
). Social regulation of the brain pituitary-gonadal axis.
Proc. Natl. Acad. Sci. USA
90
,
7794
-7798.
Fukamachi, S. and Sugimoto, M. (
2004
). Identification of the gene for morphological body-color change: Somatolactin.
Zool. Sci.
21
,
1224
.
Gabriel, S. M., Koenig, J. I. and Washton, D. L.(
1993
). Estrogen stimulation of galanin gene-expression and galanin-like immunoreactivity in the rat and its blockade by the estrogen antagonist keoxifene.
Regul. Pept.
45
,
407
-419.
Ge, Y. C., Dudoit, S. and Speed, T. P. (
2003
). Resampling-based multiple testing for microarray data analysis.
Test
12
,
1
-77.
Gerst, J. E. (
1999
). SNAREs and SNARE regulators in membrane fusion and exocytosis.
Cell. Mol. Life Sci.
55
,
707
-734.
Goodson, J. L. (
1998
). Territorial aggression and dawn song are modulated by septal vasotocin and vasoactive intestinal polypeptide in male field sparrows (Spizella pusilla).
Horm. Behav.
34
,
67
-77.
Goodson, J. L. and Adkins-Regan, E. (
1999
). Effect of intraseptal vasotocin and vasoactive intestinal polypeptide infusions on courtship song and aggression in the male zebra finch(Taeniopygia guttata).
J. Neuroendocrinol.
11
,
19
-25.
Goodson, J. L. and Bass, A. H. (
2001
). Social behavior functions and related anatomical characteristics of vasotocin/vasopressin systems in vertebrates.
Brain Res. Rev.
35
,
246
-265.
Greenwood, A. K. and Fernald, R. D. (
2004
). Social regulation of the electrical properties of gonadotropin-releasing hormone neurons in a Cichlid fish (Astatotilapia burtoni).
Biol. Reprod.
71
,
909
-918.
Greenwood, A. K., Wark, A. R., Fernald, R. D. and Hofmann, H. A. (
2008
). Arginine vasotocin mediates both dominant and subordinate behavior via two distinct preoptic pathways.
Proc. Biol. Sci.
doi:
Grober, M. S., George, A. A., Watkins, K. K., Carneiro, L. A. and Oliveira, R. F. (
2002
). Forebrain AVT and courtship in a fish with male alternative reproductive tactics.
Brain Res. Bull.
57
,
423
-425.
Grozinger, C. M. and Robinson, G. E. (
2002
). Microarray analysis of pheromone-mediated gene expression in the honey bee brain.
Integr. Comp. Biol.
42
,
1237
-1237.
Hallgren, S. L. E., Linderoth, M. and Olsen, K. H.(
2006
). Inhibition of cytochrome p450 brain aromatase reduces two male specific sexual behaviours in the male Endler guppy (Poecilia reticulata).
Gen. Comp. Endocrinol.
147
,
323
-328.
Hartwell, L. H., Hopfield, J. J., Leibler, S. and Murray, A. W. (
1999
). From molecular to modular cell biology.
Nature
402
,
C47
-C52.
Helle, K. B. (
2004
). The granin family of uniquely acidic proteins of the diffuse neuroendocrine system: comparative and functional aspects.
Biol. Rev.
79
,
769
-794.
Hofmann, H. A. (
2003
). Functional Genomics of Neural and Behavioral Plasticity.
J. Neurobiol.
54
,
272
-282.
Hofmann, H. A. (
2006
). GnRH signaling in behavioral plasticity.
Curr. Opin. Neurobiol.
16
,
343
-350.
Hofmann, H. A. and Fernald, R. D. (
2000
). Social status controls somatostatin neuron size and growth.
J. Neurosci.
20
,
4740
-4744.
Hofmann, H. A., Benson, M. E. and Fernald, R. D.(
1999a
). Social status regulates growth rate: Consequences for life-history strategies.
Proc. Natl. Acad. Sci. USA
96
,
14171
-14176.
Hofmann, H. A., Le Bail, P. Y. and Fernald, R. D.(
1999b
). Social control of growth and growth hormone levels in an African cichlid fish.
Abstr. – Soc. Neurosci.
29
,
348
.16.
Holmes, A., Yang, R. J. and Crawley, J. N.(
2002
). Evaluation of an anxiety-related phenotype in galanin overexpressing transgenic mice.
J. Mol. Neurosci.
18
,
151
-165.
Jadhao, A. and Pinelli, C. (
2001
). Galanin-like immunoreactivity in the brain and pituitary of the `four-eyed' fish, Anableps anableps.
Cell Tissue Res.
306
,
309
-318.
Johnsson, M., Axelsson, M. and Holmgren, S.(
2001
). Large veins in the Atlantic cod (Gadus morhua)and the rainbow trout (Oncorhynchus mykiss) are innervated by neuropeptide-containing nerves.
Anat. Embryol.
204
,
109
-115.
Kageyama, H., Takenoya, F., Kita, T., Hori, T., Guan, H. L. and Shioda, S. (
2005
). Galanin-like peptide in the brain: effects on feeding energy metabolism and reproduction.
Regul. Pept.
126
,
21
-26.
Kocher, T. D. (
2004
). Adaptive evolution and explosive speciation: The cichlid fish model.
Nat. Rev. Genet.
5
,
288
-298.
Legac, F., Blaise, O., Fostier, A., Lebail, P. Y., Loir, M.,Mourot, B. and Weil, C. (
1993
). Growth-Hormone (GH) and Reproduction – a review.
Fish Physiol. Biochem.
11
,
219
-232.
Lim, M. M., Hammock, E. A. and Young, L. J.(
2004
). The role of vasopressin in the genetic and neural regulation of monogamy.
J. Neuroendocrinol.
16
,
325
-332.
Liu, H. H., Payne, H. R., Wang, B. and Brady, S. T.(
2006
). Gender differences in response of hippocampus to chronic glucocorticoid stress: Role of glutamate receptors.
J. Neurosci. Res.
83
,
775
-786.
Maan, M. E., Eshuis, B., Haesler, M. P., Schneider, M. V., Van Alphen, J. J. M. and Seehausen, O. (
2008
). Color polymorphism and predation in a Lake Victoria cichlid fish.
Copeia.
2008 (3),
621
-629.
Madani, R., Kosloz, S., Akhmedov, A., Cinelli, P., Kinter, J.,Lipp, H. P., Sonderegger, P. and Wolfer, D. P. (
2003
). Impaired explorative behavior in neophobia in genetically modified mice lacking or overexpressing the extracellular serine protease inhibitor neuroserpin.
Mol. Cell. Neurosci.
23
,
473
-494.
Maere, S., Heymans, K. and Kuiper, M. (
2005
). BiNGO: a Cytoscape plugin to assess overrepresentation of Gene Ontology categories in Biological Networks.
Bioinformatics
21
,
3448
-3449.
Marsh, K. E., Creutz, L. M., Hawkins, M. B. and Godwin, J.(
2006
). Aromatase immunoreactivity in the bluehead wrasse brain, Thalassoma bifasciatum: Immunolocalization and co-regionalization with arginine vasotocin and tyrosine hydroxylase.
Brain Res.
1126
,
91
-101.
Martyniuk, C. J., Xiong, H. L., Crump, K., Chiu, S., Sardana,R., Nadler, A., Gerrie, E. R., Xia, X. H. and Trudeau, V. L.(
2006
). Gene expression profiling in the neuroendocrine brain of male goldfish (Carassius auratus) exposed to 17 alpha-ethinylestradiol.
Physiol. Genomics
27
,
328
-336.
Miranda, E. and Lomas, D. A. (
2006
). Neuroserpin: a serpin to think about.
Cell. Mol. Life Sci.
63
,
709
-722.
Mrowk, A. (
1984
). Brood care motivation and hunger in the mouth brooding cichlid Pseudochrenilabrus multicolor.
Behav. Processes
9
,
181
-190.
Munz, H. (
1999
). GnRH-systems in the forebrain of cichlid fish.
Eur. J. Morphol.
37
,
100
-102.
Nakamura, M. and Kobayashi, Y. (
2005
). Sex change in coral reef fish.
Fish Physiol. Biochem.
31
,
117
-122.
Nelson, R. J. (
2005
).
An Introduction to Behavioral Endocrinology. 3rd edn.
Sunderland,MA: Sinauer Associates.
Olazabal, D. E. and Young, L. J. (
2006
). Oxytocin receptors in the nucleus accumbens facilitate `spontaneous' maternal behavior in adult female prairie voles.
Neuroscience
141
,
559
-568.
Oliveira, R. F., Hirschenhauser, K., Carneiro, L. A. and Canario, A. V. M. (
2002
). Social modulation of androgen levels in male teleost fish.
Comp. Biochem. Physiol. B, Biochem. Mol. Biol.
132
,
203
-215.
Op de Beeck, H. P., Haushofer, J. and Kanwisher, N. G.(
2008
). Interpreting fMRI data: maps, modules and dimensions.
Nat. Rev. Neurosci.
9
,
123
-135.
Parhar, I. S., Ogawa, S. and Sakuma, Y. (
2005
). Three GnRH receptor types in laser-captured single cells of the cichlid pituitary display cellular and functional heterogeneity.
Proc. Natl. Acad. Sci. USA
102
,
2204
-2209.
Parikh, V. N., Clement, T. S. and Fernald, R. D.(
2006
). Androgen level and male social status in the African cichlid, Astatotilapia burtoni.
Behav. Brain Res.
166
,
291
-295.
Pollen, A. A., Dobberfuhl, A. P., Scace, J., Igulu, M. M., Renn,S. C. P., Shumway, C. A. and Hofmann, H. A. (
2007
). Environmental complexity and social organization sculpt the brain in Lake Tanganyikan cichlid fish.
Brain Behav. Evol.
70
,
21
-39.
R Development Core Team (
2006
).
R:A language and environment for statistical computing
. R. Foundation for Statistical Computing, Vienna, Austria.〈 http://www.R-project.org
Redies, C. and Puelles, L. (
2001
). Modularity in vertebrate brain development and evolution.
BioEssays
23
,
1100
-1111.
Renn, S. C. P., Aubin-Horth, N. and Hofmann, H. A.(
2004
). Biologically meaningful expression profiling across species using heterologous hybridization to a cDNA microarray.
BMC Genomics
5
,
42
.
Rhodes, J. S. and Crabbe, J. C. (
2005
). Gene expression induced by drugs of abuse.
Curr. Opin. Pharmacol.
5
,
26
-33.
Robinson, G. E., Grozinger, C. M. and Whitfield, C. W.(
2005
). Sociogenomics: Social life in molecular terms.
Nat. Rev. Genet.
6
,
257
-270.
Sagrillo, C. A., Grattan, D. R., McCarthy, M. M. and Selmanoff,M. (
1996
). Hormonal and neurotransmitter regulation of GnRH gene expression and related reproductive behaviors.
Behav. Genet.
26
,
241
-277.
Salzburger, W., Mack, T., Verheyen, E. and Meyer, A.(
2005
). Out of Tanganyika: genesis, explosive speciation,key-innovations and phylogeography of the haplochromine cichlid fishes.
BMC Evol. Biol.
5
,
15
.
Salzburger, W., Renn, S. C. P., Steinke, D., Hofmann, H. A.,Braasch, I. and Meyer, A. (
2008
). Annotation of expressed sequence tags for the East African cichlid fish species Astatotilapia burtoni and evolutionary analyses of cichlid ORFs.
BMC Genomics
9
,
96
-110.
Santangelo, N. and Bass, A. H. (
2006
). New insights into neuropeptide modulation of aggression: field studies of arginine vasotocin in a territorial tropical damselfish.
Proc. Biol. Sci.
273
,
3085
-3092.
Sapolsky, R. M. (
2005
). The influence of social hierarchy on primate health.
Science
308
,
648
-652.
Schick, R. R., Samsai, S., Zimmermann, J. P., Eberl, T., Endres,C., Schusdziarra, V. and Classen, M. (
1993
). Effect of galanin on food-intake in rats-involvement of lateral and ventromedial hypothalamic sites.
Am. J. Physiol.
264
,
R355
-R361.
Schlosser, G. and Wagner, G. P. (ed.) (
2004
).
Modularity in Development and Evolution
. Chicago:University Of Chicago Press.
Segal, E., Friedman Koller, D. and Regev, A.(
2004
). A module map showing conditional activity of expression modules in cancer.
Nat. Genet.
36
,
1090
-1098.
Semsar, K. and Godwin, J. (
2003
). Social influences on the arginine vasotocin system are independent of gonads in a sex-changing fish.
J. Neurosci.
23
,
4386
-4393.
Semsar, K. and Godwin, J. (
2004
). Multiple mechanisms of phenotype development in the bluehead wrasse.
Hormones and Behavior
45
,
345
-353.
Shannon, P., Markiel, A., Ozier, O., Baliga, N. S., Wang, J. T.,Ramage, D., Amin, N., Schwikowski, B. and Ideker, T. (
2003
). Cytoscape: a software environment for integrated models of biomolecular interaction networks.
Genome Res.
13
,
2498
-2504.
Shaw, D. R., Ashbumer, M., Blake, J. A., Baldarelli, R. M.,Botstein, D., Davis, A. P., Cherry, J. M., Lewis, S., Lutz, C. M., Richardson,J. E. et al. (
1999
). Gene ontology: a controlled vocabulary to describe the function, biological process and cellular location of gene products in genome databases.
Am. J. Hum. Genet.
65
,
A419
.
Silverin, B., Baillien, M. and Balthazart, J.(
2004
). Territorial aggression, circulating levels of testosterone, and brain aromatase activity in free-living pied flycatchers.
Horm. Behav.
45
,
225
-234.
Simerly, R. B. (
2002
). Wired for reproduction:organization and development of sexually dimorphic circuits in the mammalian forebrain.
Annu. Rev. Neurosci.
25
,
507
-536.
Smyth, G. K., Thorne, N. P. and Wettenhall, J.(
2004
).
LIMMA: Linear Models for Microarray Data, Version 1.6.6, User's Guide
[http://bioinf.wehi.edu.au/limma/].
Sokal, R. R. and Rohlf, F. J. (
1995
).
Biometry: The Principles and Practice of Statistics in Biological Research. 3rd edn.
New York: W.H. Freeman.
Soma, K. K., Francis, R. C., Wingfield, J. C. and Fernald, R. D. (
1996
). Androgen regulation of hypothalamic neurons containing gonadotropin-releasing hormone in a cichlid fish: Integration with social cues.
Horm. Behav.
30
,
216
-226.
Soma, K. K., Schlinger, B. A., Wingfield, J. C. and Saldanha, C. J. (
2003
). Brain aromatase, 5 alpha-reductase, and 5 beta-reductase change seasonally in wild male song sparrows: relationship to aggressive and sexual behavior.
J. Neurobiol.
56
,
209
-221.
Sorensen, C., Overli, O., Summers, C. H. and Nilsson, G. E.(
2007
). Social regulation of neurogenesis in teleosts.
Brain Behav. Evol.
70
,
239
-246.
The Gene Ontology Consortium (
2000
). Gene ontology: tool for the unification of biology.
Nat. Genet.
25
,
25
-29.
Tortorella, C., Neri, G. and Nussdorfer, G. G.(
2007
). Galanin in the regulation of the hypothalamic-pituitary-adrenal axis (Review).
Int. J. Mol. Med.
19
,
639
-647.
Townsend, J. P. and Hartl, D. L. (
2002
). Bayesian analysis of gene expression levels: statistical quantification of relative mRNA level across multiple strains or treatments.
Genome Biol.
3
,
research0071.1
-0071.16.
Trainor, B. C. and Hofmann, H. A. (
2006
). Somatostatin regulates aggressive behavior in an African cichlid fish.
Endocrinology
147
,
5119
-5125.
Trainor, B. C. and Hofmann, H. A. (
2007
). Somatostatin and somatostatin receptor gene expression in dominant and subordinate males of an African cichlid fish.
Behav. Brain Res.
179
,
314
-320.
Trainor, B. C., Kyomen, H. H. and Marler, C. A.(
2006
). Estrogenic encounters: How interactions between aromatase and the environment modulate aggression.
Front. Neuroendocrinol.
27
,
170
-179.
Trudeau, V. L., Spanswick, D., Fraser, E. J., Lariviere, K.,Crump, D., Chiu, S., MacMillan, M. and Schulz, R. W. (
2000
). The role of amino acid neurotransmitters in the regulation of pituitary gonadotropin release in fish.
Biochem. Cell Biol.
78
,
241
-259.
Volkoff, H., Canosa, L. F., Unniappan, S., Cerda-Reverter, J. M., Bernier, N. J., Kelly, S. P. and Peter, R. E. (
2005
). Neuropeptides and the control of food intake in fish.
Gen. Comp. Endocrinol.
142
,
3
-19.
White, S. A., Bond, C. T., Francis, R. C., Kasten, T. L.,Fernald, R. D. and Adelman, J. P. (
1994
). A 2nd gene for gonadotropin-releasing-hormone-cdna and expression pattern in the brain.
Proc. Natl. Acad. Sci. USA
91
,
1423
-1427.
White, S. A., Nguyen, T. and Fernald, R. D.(
2002
). Social regulation of gonadotropin-releasing hormone.
J. Exp. Biol.
205
,
2567
-2581.
Whitehead, A. and Crawford, D. L. (
2006
). Variation within and among species in gene expression: raw material for evolution.
Mol. Ecol.
15
,
1197
-1211.
Whitfield, C. W., Cziko, A. M. and Robinson, G. E.(
2003
). Gene expression profiles in the brain predict behavior in individual honey bees.
Science
302
,
296
-299.
Wiesenfeldhallin, Z., Xu, X. J., Langel, U., Bedecs, K.,Hokfelt, T. and Bartfai, T. (
1992
). Galanin-mediated control of pain-enhanced role after nerve injury.
Proc. Natl. Acad. Sci. USA
89
,
3334
-3337.
Winslow, J. T., Hastings, N., Carter, C. S., Harbaugh, C. R. and Insel, T. R. (
1993
). A role for central vasopressin in pair bonding in monogamous prairie voles.
Nature
365
,
545
-548.
Wolf, D. M. and Arkin, A. P. (
2003
). Motifs,modules and games in bacteria.
Curr. Opin. Microbiol.
6
,
125
-134.
Wu, H., Kerr, M. K., Cui, X. and Churchill, G. A.(
2002
). MAANOVA: a software package for the analysis of spotted cDNA microarray experiments. In
The Analysis of Gene Expression Data: Methods and Software
(ed. G. Parmigiani, E. S. Garett, R. A. Irizarry and S. L. Zeger), pp.
313
-341. New York:Springer-Verlag.
Yu, K. L. and Peter, R. E. (
1991
). Changes in brain levels of gonadotropin-releasing-hormone and serum levels of gonadotropin and growth-hormone in goldfish during spawning.
Can. J. Zool.
69
,
182
-188.