Experimental selection for body size at age modifies early life-history traits and muscle gene expression in adult zebrafish

The domestication process is of fundamental importance to both aquaculture production and conservation efforts that involve restocking habitats with hatchery-reared fish because selection on traits that are beneficial in culture may reduce fitness in the wild (Christie et al., 2012). Domestication is driven by epigenetic and genetic mechanisms occurring in each successive generation of culture. Adaptation to the novel captive environment (unintended or inadvertent selection) occurs because of differences in food availability, and physical and ecological factors relative to the wild, producing profound effects on numerous behavioural, physiological, morphological and life-history traits (Álvarez and Nicieza, 2003; Robison and Rowland, 2005), sometimes in as little as a single generation (Mayer et al., 2011; Christie et al., 2012). For example, satiation feeding, the surface presentation of food and the absence of predators has been linked to a reduction in stress, loss of the startle response, and an increase in aggression and growth rate in cultured relative to wild populations (Lepage et al., 2000; Yamamoto and Reinhardt, 2003). Breeding programmes involving artificial selection of individuals or families are designed to direct the domestication process to produce stains with desirable characteristics for farming, including high growth rate, delayed sexual maturation, increased disease resistance, etc. (Gjedrem, 2000). For Atlantic salmon (Salmo salar), 10–15% improvements in growth rate per generation have been achieved (Gjedrem, 2000). After 16 generations of selection for rapid growth, domesticated coho salmon (Oncorhynchus kisutch) grew faster than unselected strains with satiation feeding, possibly due to a higher feed intake and feed conversion efficiency (Neely et al., 2008). Strong selection for fast growth in Atlantic silverside (Menidia menidia) was shown to increase food intake relative to unselected lines in only four generations (Lankford et al., 2001; Chiba et al., 2007).

Somatic growth is a complex trait involving all the metabolic pathways that control energy acquisition (food ingestion) and energy expenditure (food assimilation, locomotion, reproduction and maintenance). Heritable variations in the coding and regulatory sequence of hormones, growth factors and their downstream signalling pathways are good candidates for explaining differences in growth rate between populations. Growth hormone (GH) produced by the hypophysis has both direct and indirect effects on growth (Reinecke et al., 2005). The indirect effects of GH are mainly produced by the insulin-like growth factor (IGF) pathway, in which igf1 and igf2 are the ligands. Binding of the IGF ligands to igf1r leads to phosphorylation of the phosphatidylinositol 3-phosphate kinase/v-akt murine thymoma viral oncogene/mammalian target of rapamycin (PI3K/AKT/mTOR) pathway, which is responsible for the activation of expression of growth-related genes and repression of transcription of genes that function in protein degradation by the ubiquitin proteasome system (reviewed in Johnston et al., 2011). Binding of IGFs to igf2r results in lysosomal degradation of the ligand and is probably a mechanism for regulating the concentration of circulating IGFs (Lau et al., 1994; Wang et al., 1994). The availability of circulating IGFs and binding to the respective receptors is also regulated by six IGF-binding proteins (IGFBPs) (Wood et al., 2005). There is evidence that changes in the GH–IGF axis might explain some of the effects of domestication observed in selective breeding programmes in fish culture. For example, channel catfish (Ictalurus punctatus) selected for fast growth for two generations had lower plasma cortisol levels and higher muscle igf2 mRNA than fish selected for slow growth, leading to the conclusion that the differences in growth were partially explained by the GH–IGF and stress axis (Peterson et al., 2008). A microarray experiment with coho salmon (O. kisutch) showed similar transcriptional responses in a GH-transgenic strain and fish domesticated for 12 generations, which differed from wild-type fish, leading to the hypothesis that domestication affected similar downstream components to those in GH-transgenic fish (Devlin et al., 2009). The IGF–mechanistic target of rapamycin (mTor) signalling pathway integrates environmental signals to regulate body size and growth rate in diverse animal lineages (Lupu et al., 2001; McCulloch and Gems, 2003). Adaptive modification in IGF–mTor pathway genes, including mTOR and 4e-bp-1 (two proteins that regulate protein synthesis) may underlie the rapid evolution of dwarfism from large-bodied populations in the Arctic charr (Salvelinus alpinus) (Macqueen et al., 2011).

Maternally deposited gene transcripts are involved in the control of zygotic gene activation and drive development prior to the mid-blastula transition (Mathavan et al., 2005). In zebrafish, some maternal transcripts were found to persist until the end of the blastula stage whilst others were degraded more slowly and were still present at the segmentation stage, indicating some functional role after the start of zygotic transcription (Mathavan et al., 2005). Ovulatory ageing was used to manipulate the content of maternal mRNA in rainbow trout oocytes (Aegerter et al., 2004). Interestingly, igf1rb content was negatively correlated with the occurrence of morphological abnormalities and igf-I, igf2 and igf1rb levels were positively correlated with embryonic survival (Aegerter et al., 2004). To our knowledge, no previous studies have investigated whether domestication and/or artificial selection results in changes in the deposition of maternal mRNAs in the yolk.

In the present study, we used a wild-derived stock of zebrafish to produce replicate laboratory-reared lineages over four generations without or with artificial selection for body size at 90 days post-fertilization (d.p.f.). The zebrafish (Danio rerio) is an excellent model for studies of artificial selection and domestication because of its short generation time, our extensive knowledge on its development and physiology, and the ready availability of molecular tools (Robison and Rowland, 2005; Wright et al., 2006). We tested the hypothesis that laboratory culture and artificial selection for body size at age have pervasive effects on transcriptional components of the IGF system in 2-cell embryos and adult fish, providing a potential mechanism for effects on growth at various stages of the life cycle.

This study was conducted on the F3 generation of a wild-caught zebrafish population (D. rerio, Hamilton) from Mymensingh, Bangladesh (27 males and 28 females). Adult fish were reared in 25 l tanks at 27±0.3°C (mean and range) under a 12 h light:12 h dark photoperiodic regime in a filtered freshwater recirculation system and fed bloodworms to satiation twice daily. Replicated unselected (U), small (S) and large (L) lineages were bred from a random selection of 18–19 fish per lineage (sex ratio ~1 female:1 male) derived from the founder population (supplementary material Table S1). All breeding was conducted at ~120 d.p.f. when fish were sexually mature. For each lineage, previously separated male and female fish were introduced into a breeding tank and eggs were collected each morning for 3 days. Eggs were immediately cleaned and maintained in glass tanks (1 l) under the same environmental conditions as the main recirculation system. After 7 d.p.f., the larva were fed ZM-100 (Fish Food, Winchester, Hampshire, UK) and microworms; 50% of the water was changed twice daily until 30 d.p.f., when fish were transferred to the main recirculation system and fed to satiation with ZM-200 (Fish Food) and bloodworms (Ocean Nutrition, Essen, Belgium) twice daily. Three rounds of experimental selection were conducted based on body size at ~90 d.p.f. (supplementary material Table S1). For the S-lineage, fish of standard length (SL, tip of snout to last vertebrae) greater than 75% of the mean SL for the population were removed from the breeding population at each generation. The L-lineages were generated by removing fish with SL less than 125% of the mean SL for the population. A third lineage was produced in which fish were not selected for body size at age (U-lineage), but which were nevertheless undergoing adaptation to the captive environment, i.e. domestication. The percentage of each lineage selected for breeding in the next generation was 46–52% (L-lineage), 63–75% (S-lineage) and 100% (U-lineage) (supplementary material Table S1) and this had little impact on the sex ratio of the populations. The number of fish used to produce the 2nd, 3rd and 4th generations ranged from 24 to 78 per replicated lineage (supplementary material Table S1). As part of routine husbandry procedures SL, fork length (FL), total length (TL), maximum body depth (H) and body mass (M_b) were measured periodically to access the growth and health of all populations in the colony.

Embryos from the 4th generation of the experimental selection were collected in the first hour after fertilization and photographed using an AxioCam CCD camera (Zeiss, Göttingen, Germany) and a Leica Wild M3Z stereo microscope (Leica, Heerbrugg, Switzerland) at 10× magnification. Egg and yolk size were calculated by measuring the ferret diameter of 100 eggs and their yolk from three different spawnings from each lineage, using ImageJ V. 1.42i software (National Institutes of Health, Bethesda, MD, USA). Deformities in the animal pole representing developmental abnormalities were also recorded. The embryos were then reared in glass tanks as described above until 30 d.p.f., and mortality was checked twice per day. The TL of 50 larva from each lineage was measured at 6 d.p.f., corresponding to the time when the yolk was almost completely assimilated and the larvae were free swimming.

Maternal transcript levels were measured by quantitative PCR (qPCR) in the eggs of the 4th generation of experimental selection using 12 replicates of 60 embryos per lineage with SYBR II chemistry (Stratagene, La Jolla, CA, USA). Embryos at the 2- to 4-cell stage were collected at 1 h after fertilization and snap frozen in liquid nitrogen for total RNA extraction using Tri-reagent (Sigma, St Louis, MO, USA) and subsequent first strand cDNA synthesis using a Quantitect Reverse Transcription Kit (Qiagen, Hilden, Germany) with 1 μg of RNA. cDNA at 40-fold dilution was used as the working solution. Primer pairs for the 15 known genes of the IGF system in the zebrafish (Amaral and Johnston, 2011) and myogenic regulatory factors (MRFs) were as described previously (Amaral and Johnston, 2012). Primers were designed to amplify the ‘fecundity genes’ (bmp15 and gdf9) and their receptors (bmpr1aa, bmpr1ab, bmpr1ba, bmpr1bb, bmpr2a, bmpr2b), and GH (gh1) and its receptors (ghra and ghrb) (Table 1). qPCR procedures were compliant with MIQE guidelines (Bustin et al., 2009) and have been described in detail previously (Amaral and Johnston, 2011). In order to establish the best normalization strategy, the expression of 13 reference genes (Amaral and Johnston, 2011) was analysed across lineages using Genorm v3.5 (Vandesompele et al., 2002) with the expression stability measure M set to <1.5. qPCR efficiency (E) was calculated using LinRegPCR V. 12.5 software (Heart Failure Research Center, Amsterdam, The Netherlands) (Table 1). Transcripts levels, expressed in arbitrary units, were calculated using the mean efficiency of 72 reactions with posterior normalization to the level of the two most stable reference genes (tomm20b and ef1α, M=0.129). The specificity of each qPCR assay was validated by direct sequencing of the PCR product.

Fish from the fifth generation (G5) were used to investigate the effect of the experimental selection on the transcriptional response to nutrient levels in skeletal muscle. Two replicate populations from the S- and L-lineages were reared to adult stage in the same conditions as described for the F3 population; fish densities during development were the same for all replicates (N=160 fish per tank). For logistical reasons, fish from the U-lineages were not used in this part of the study. At 9 months of age, 50 male fish were randomly collected from each lineage and transferred to a separate tank, where they were fed bloodworms to satiety twice daily for the 2 week acclimation period. Food was withdrawn for 1 week of fasting and was then resumed. Six fish from each lineage were collected at −170 h (just before the fasting period), 0 h (just before refeeding), and 1, 3, 6, 24 and 48 h after resuming feeding. The feeding schedule after the 0 h time point was the same as that during the acclimation period. Fish were killed by an overdose of ethyl 3-aminobenzoate methanesulphonate salt (MS-222, Fluka, St Louis, MO, USA) and had their TL and M_b measured. Fast skeletal muscle was dissected from the dorsal epaxial myotomes, flash frozen in liquid nitrogen and stored at −80°C prior to total RNA extraction. The digestive tract was dissected and fixed in 4% (m/v) paraformaldehyde for later quantification of intestine content to the nearest mg. The transcription levels of the 15 genes of the IGF system, two ubiquitin ligases, 12 other nutritionally responsive genes, and the MRFs in skeletal muscle were measured by qPCR. Total RNA extraction, first strand cDNA synthesis and dilution, and qPCR procedures were as described previously (Amaral and Johnston, 2011).

All data were tested for equality of variance and normal distribution. Growth patterns were modelled using Growth II software (Pisces Conservation, New Milton, Hampshire, UK). The data on body size at age of the different lineages were fitted to six growth models (von Bertalanffy, exponential, 3 parameters Gompertz, 3 parameters logistic, 4 parameters Gompertz and 4 parameters logistic) and the one with the lowest Akaike information criterion (AIC) and Schwarz criterion was chosen. To facilitate the interpretation of the relationship between SL and M_b in the different populations, these measurements were transformed by raising them to the power 0.33 and log transformation, respectively. A general linear model (GLM) was used to test the effect of selection on the proportionality among body size measurements (SL, FL, TL, H and M_b), using PASW Statistics 18 software (SPSS Inc., Chicago, IL, USA). t-test in PASW Statistics (SPSS Inc.) was used to compare the body size at age between sexes. The chi-square test in R V. 2.10.0 software was used to analyse the sex ratio among the zebrafish populations at sexual maturity. ANOVA followed by Tukey's post hoc test using PASW Statistics 18 software (SPSS Inc.) was employed to analyse normally distributed data. The normally distributed data on maternal transcripts were analysed by principal component analysis (PCA) and differences in the extracted principal components among the lineages were tested by ANOVA followed by Tukey's post hoc test using PASW Statistics. Transcripts with a weighting value higher than 0.75 were considered to be an important contributor to the extracted principal component. Differences in the normally distributed data on transcript level in skeletal muscle of selected adult fish were tested by MANOVA in PASW Statistics. The test set-up was as follows: lineage as factor, gut content and time as co-variates, and transcript level as dependent variable. Heat maps of transcript levels were produced using PermutMatrix software (http://www.lirmm.fr/~caraux/PermutMatrix/EN/index.html), with gene expression normalized for rows and McQuitty's method used for hierarchical clustering.

After four generations of selection the mean SL of the S-lineage was 2% lower than that of the founding population (mean of 0.6% loss in SL per generation), while SL was 10% higher in the L-lineage (gain in SL of 3.3% per generation). The relationship between various measures of body length (SL, FL, TL and H) and M_b was independent of sex and lineage (GLM, P>0.05). However, a significant difference in body size at age between sexes was observed, with females having a larger body size than males in all lineages. The sex ratio of adult fish among the replicate populations from each lineage was close to 1:1 (chi-square test, P=0.6). This allowed for comparison of the body size and growth rate of the lineages without considering sex as a confounding variable. Based on the AIC, which measures the closeness of the experimental points to the model, the best model for growth was a 4-parameter Gompertz equation. In this model, SL follows a sigmoid curve with an inflection point (when growth rate starts to decrease) at 1/3 of age of when the fish reach the predicted maximum body size. There was a statistically significant difference in body size at age among the three different populations observable from 120 d.p.f. (not shown). Body size measurements between the replicates of a given lineage were not statistically different when adult stage was reached and therefore were pooled to facilitate the interpretation of results (supplementary material Fig. S1). Fish from the L-lineages reached a maximum SL of 32.9 mm at 390 d.p.f., which corresponded to 6.8% and 12.3% greater SL than fish from U- and S-lineages (Fig. 1A). These differences of body size at age between the two selected lines from G4 are in accordance with the change in SL per generation (~12%). The mean M_b of the L-lineages (0.52 g) was 16.3% and 41.0% greater than that of the U- and S-lineages at 390 d.p.f. (Fig. 1B). The data on body size at age only allowed the calculation of mean growth rates for each lineage, which were 0.194, 0.181 and 0.173 mm day⁻¹ for the L-, U- and S-lines, respectively, during the exponential phase of growth between 6 and 120 d.p.f.

The production of eggs per breeding batch in the 4th generation was five times greater in the L- than in the S-lineage (Table 2). Eggs from the S-lineage showed a small but significant reduction in diameter (6.4%) compared with eggs from the L-lineage, corresponding to 13% less yolk (Table 2). Three waves of mortality were observed in all lineages peaking at <24 h (embryos), 2–4 days (yolk-sac larvae) and 8–15 days (free-swimming larvae). Embryos from the L-lineage had the highest rates of deformity and mortality in both the first (~12%) and second waves (40%), followed by the U- and S-lineages (Table 2). In the third wave, however, this scenario was reversed, with the S-lineage having the highest mortality (~64%) when compared with the other lineages. At the end of the juvenile stage these differences in mortality seem to have been equalized as the survival rate at 30 d.p.f. was around 30% for all three lineages (Table 2). The mean TL of larvae at 6 d.p.f. just prior to complete yolk absorption was ~3.7, 3.6 and 3.4 mm for the L-, U- and S-lineages, respectively (Table 2).

The PCA identified two main trends in the data. The level of transcripts that contributed to the first trend was higher in the L-lineage than in the S- and U-lineages, whereas the second trend includes transcripts that were at a higher level in the U-lineages than in the other lineages (Fig. 2).

Transcript abundance for some genes of the GH–IGF system in 2- to 4-cell embryos was higher in the L-lineage than in either the U- or S-lineages including igf2a (5.6-fold higher, P<0.01) (Fig. 2B). Three of the five IGF and GH receptors (igf1ar, igf2r and ghra) were similarly higher in the L-lineage than in the S-lineage by an average of 1.4-fold (Fig. 2B) whereas there was no difference in igf1br and ghrb transcripts between lineages (supplementary material Fig. S2). GH mRNA was 1.4- and 1.9-fold higher in eggs from the U-lineage than from the S- and L-lineages, respectively (Fig. 2B). igfbp5a and igfbp3 were 51% and 29% higher whereas igfbp1b was 5-fold lower in the L-lineage than in the U- or S-lineages, respectively (Fig. 2B; supplementary material Fig. S2). igfbp6a mRNA increased in the series U->L->S-lineage (supplementary material Fig. S2). Transcripts of igf1b, igfbp2a, igfbp5b and igfbp6b were not detected in 2- to 4-cell embryos.

No difference in transcript abundance between the lineages was detected for the members of the MRF family (supplementary material Fig. S2). The ligands bmp15 and gdf9 (so-called ‘fecundity genes’) showed no significant differences between lineages (supplementary material Fig. S2). In contrast, bmpr1aa transcripts were 1.4-fold higher in the L- than in the S-lineage. bmpr2b transcript level increased in the series S->L->U-lineage across a 6-fold range of abundance, while no difference was recorded for bmpr1ba and bmpr2a transcripts (Fig. 2B; supplementary material Fig. S2).

The S- and L-lineages had a similar feeding response across the fasting and refeeding experiment, with no significant differences in gut fullness (Fig. 3A). The gut food content decreased from 1.4% M_b at −170 h to 0% after 7 days of fasting. After feeding was resumed, the gut food content increased to 4.0%, 3.3% and 4.0% at 1, 3 and 6 h, respectively. At 24 and 48 h, the gut food content was similar to that found before the fasting period (~1% of body mass), indicating that fasting resulted in a transient increase in food intake following refeeding (Fig. 3A).

The expression of 11 out of 15 genes from the IGF system was very similar in the S- and L-lineages in response to changes in nutritional input. Similarly, no significant difference was observed for nutritionally responsive genes and the two ubiquitin ligases [mafbx and trim63 (also known as murf1)] in the two lineages. The transcript levels of these genes in relation to the gut food content were comparable to those reported previously.

Transcript abundance for the myogenic regulatory factors were unaffected by selection for body size (supplementary material Fig. S3). myf5 and myoD were the only two members of the MRF that showed a clear pattern of expression in response to feeding, which was not affected by the selection experiment. myf5 showed an initial small downregulation followed by a 2.5-fold upregulation with re-feeding to satiation (supplementary material Fig. S3). myoD transcript abundance increased around 3-fold 1 h after re-feeding before returning to fed levels (supplementary material Fig. S3).

Differences in gene expression in skeletal muscle of adult fish between the S- and L-lineages were found only for igf1a, igf2r, igfbp1a and igfbp1b (Fig. 3B). igf1a had constitutively higher expression in the L- than in the S-lineage (P=0.04) (Fig. 3B, Fig. 4A). igf2 receptor paralogues increased with fasting and showed higher constitutive expression in the L- than in the S- lineage (Fig. 4B).

Transcripts for IGF-binding proteins 1a and 1b (the main isoforms expressed in zebrafish fast muscle) (Amaral and Johnston, 2011) were generally higher (up to 2-fold) in the S- than in the L-lineages (Fig. 3B, Fig. 5A,B). igfbp1a showed higher constitutive expression in the S- than in the L-lineage (Fig. 5A), whereas igfbp1b exhibited regulated differences in expression with selection. igfbp1b increased more with fasting and decreased more with refeeding in the L- than in the S-lineage (Fig. 5B).

After three rounds of experimental selection, the L-lineage of zebrafish showed increased growth rate and final body size compared with the U- and S-lineages. The change in body size affected the physiology of the fish at many levels. For example, larger fish produced more and larger eggs containing more yolk, which probably explains the larger body size of larvae at the yolk-absorption stage. Deposition of maternal transcripts was also affected by selective breeding as evidenced by higher levels of igf2a, GH and IGF receptors in fertilized eggs from the L-lineage compared with the S-lineage. Experimental selection also affected the transcription of a limited number of genes from the IGF pathway in skeletal muscle of adult fish in response to changing nutritional input, with higher expression of igf1a and igf2r, and lower expression of igfbp1, in the L- than in the S-lines. These findings point to a differential regulation of transcript deposition in eggs and altered transcriptional regulation in adult muscle with selection for body size. The changes in the IGF pathway observed in fish selected for larger body size are similar to those observed in domesticated and GH-transgenic rainbow trout (Devlin et al., 2009). In coho salmon, selection for body size at age resulted in increased food intake and growth with satiating but not with restricted feeding (Neely et al., 2008). Food intake was not measured in the present study, although initial gut fullness in fasted fish that were fed to satiation was similar for the selected lines, which does not suggest major changes in appetite.

Experimental selection was found to have marked effects on early life-history traits in zebrafish. A considerably higher number of embryos were produced by the L- than by the S-lineage (Table 2), using the same number of fish and the same sex ratio per spawning. The difference in embryo production could simply be due to differences in body mass as a strong positive relationship was previously described between body mass and egg production in the zebrafish (Forbes et al., 2010). However, increased egg production in the L-lineage could also be due to a higher number of females ready to lay eggs, which could be investigated in the future by single-pair breeding experiments. In a previous study, the yolk volume of zebrafish embryos was experimentally manipulated and a clear effect of yolk volume on larval body size was found (Jardine and Litvak, 2003). The larger body size of larvae from the L-lineage (Table 2) can therefore be explained at least in part by the higher amount of energy available for development and growth in the yolk.

During oocyte maturation, maternal transcripts are deposited in the oocyte and remain detectable until ~6 h post-fertilization with varying degrees of degradation (Pelegri, 2003; Mathavan et al., 2005; Tadros and Lipshitz, 2009). Maternally deposited mRNA is necessary for normal embryo development prior to the activation of zygotic transcription, which occurs at the mid-blastula transition (Pelegri, 2003; Dosch et al., 2004; Lubzens et al., 2010). Maternal mRNA is also thought to partially serve as a nutritional reserve (Hunter et al., 2010). Research in this area has tended to focus on observations of developmental abnormalities in zebrafish embryos triggered by chemical mutations produced by exposing the parents to N-ethyl-N-nitrosurea to cause random point mutations. Although important, this approach lacks the power to investigate the changes in deposition of non-lethal transcripts, which could have milder but important and persistent effects on fitness. For example, morpholino knockdown of maternal transcripts for the oestrogen receptor 2a negatively affected embryonic and larval development in zebrafish (Celeghin et al., 2011). Maternally deposited transcripts of radar, a TGF-β signalling molecule, were shown to be essential for the ventral fate of cells during development, with loss-of-function causing lethal dorsalized phenotypes (Sidi et al., 2003). Thus, it is possible that variation in maternal mRNA deposition represents a trade-off in egg quality, with permanent effects on early life-history traits. The levels of most of the transcripts investigated in the present study have recently been reported to belong to gene clusters with very few changes between the 1- and 512-cell developmental stages in the zebrafish embryo as assayed by RNA-seq (Vesterlund et al., 2011). We investigated transcript abundance in 1- to 2-cell embryos of genes involved in fertilization (gdf9, bmp15 and bmp receptors), growth (GH–IGF axis) and myogenesis (MRFs). No difference was observed in maternal deposition of gdf9 and bmp15 transcripts among the three selected lineages. These two genes are closely related members of the TGF-β superfamily with important roles in oocyte maturation in the zebrafish (Liu and Ge, 2007; Peng et al., 2009), and are thought to be involved in ovulation and fecundity. A previous study found no correlation between gdf9 and bmp15 levels in mature oocytes and fecundity in zebrafish populations fed four different diets (Forbes et al., 2010). It is possible, however, that ovulation and fecundity are regulated not solely at the ligands but also at the receptor level. For example, polymorphisms in bmpr1b correlated well with distinct breeds of sheep in which different prolificacy rates were observed (Chu et al., 2011). However, bmp receptors also recognize and bind a number of ligands from the TGF-β family (Koenig et al., 1994; Penton et al., 1994). The present study shows that deposition of some paralogues of the receptors for BMPs was differentially regulated among the three selected lineages, with the L-lineage having significantly more transcripts of bmpr1aa and bmpr2b (Fig. 2B).

Embryos from the L-lineage had more igf2a, igf1ar and ghra transcripts than those from the S-lineage (Fig. 2B). There is no information available on the effect that changes in the levels of these maternal transcripts might have on zebrafish development, but experiments with morpholinos in zebrafish embryos revealed the importance of igf1 receptors for the normal development of the embryo (Schlueter et al., 2007). Information on the effect of knockdown of igf2r is currently lacking. A differential deposition of igfbp1a and igfbp1b in the selected lineages was observed, with igfbp1a being lower in S- than in L-lineages, while the opposite was observed for igfbp1b (Fig. 2B; supplementary material Fig. S2). Sub-function partitioning between these two igfbp1 transcripts in the zebrafish has been reported before, with overexpression of either binding protein causing developmental retardation (Kamei et al., 2008). However, care must be taken when interpreting these results as there is no experimental evidence that the different growth phenotypes obtained here are caused by changes in maternal deposition of these transcripts.

We also tested the hypothesis that IGF pathway gene expression is differentially regulated between the two selected lines in adult stages. Following 7 days of fasting, zebrafish were fed to satiation and transcript abundance was measured over 48 h during the transition from a catabolic to an anabolic state. The results of transcriptional regulation in response to fasting were very similar to those described previously (Amaral and Johnston, 2011). Twenty-three of the genes investigated showed no difference in transcript abundance or expression pattern between the S- and L-lineages. myf5, a member of the MRF family involved in muscle differentiation and proliferation (Chen and Tsai, 2008), showed a clear response to fasting–refeeding with downregulation during fasting and gradual upregulation with feeding (supplementary material Fig. S3). This pattern of myf5 expression was independent of the zebrafish lineage and has been previously described in Atlantic salmon (Bower and Johnston, 2010). Likewise, a transient increase MyoD1b expression with re-feeding has previously been reported in Atlantic salmon (Bower and Johnston, 2010).

In the zebrafish, igf1a is downregulated during fasting with a peak of expression in response to a single meal, while its receptors (igf1ar and igf1br) show the opposite pattern of expression, with upregulation during fasting and subsequent downregulation during refeeding (Amaral and Johnston, 2011). The selection regime affected the level of expression of igf1a, with the L-lineage having a higher basal expression of igf1a across the fasting–refeeding experiment (Fig. 4A). IGF promotes cell growth after binding to its receptors, activating the AKT/PI3K/mTOR pathway (Laplante and Sabatini, 2009). Assuming changes in transcript abundance are transmitted to the protein level, the present study shows that the L-lineage might be more responsive to the same level of energy as a result of a higher expression of igf1a.

Transcripts of igf2r and igf2a were upregulated with fasting and downregulated with refeeding, with the igf2r transcripts showing a very gradual change in transcript level with respect to igf2a (Fig. 4B). This pattern of expression and the correlation between igf2r and igf2 expression were as observed previously in Atlantic salmon (Bower et al., 2008; Bower and Johnston, 2010). The constitutive expression of igf2r was higher in the L- than in the S-lineage and igf2a showed a greater increase between fed and fasted states (Fig. 4B). In mice, the absence of a functional allelic copy of the maternal igf2r leads to perinatal overgrowth and lethality, with elevated levels of igf2 in embryos due to the absence of the degrading igf2 function of the igf2r (Lau et al., 1994; Wang et al., 1994). In the zebrafish, two copies of the igf2 gene exist, with distinct transcriptional regulation during embryonic and adult phases (Zou et al., 2009; Nelson and Van Der Kraak, 2010), and with differential regulation in response fasting and refeeding (present study).

It is known that overexpression of igfbp1a/b in zebrafish embryos under normoxic conditions retards development and growth (Kajimura et al., 2005; Kamei et al., 2008). The growth-inhibiting action of igfbp1 occurs by rendering igf1 less available to tissues. This negative regulation of growth by igfbp1 genes might explain the downregulation of these genes during growth stimulation by satiation feeding (Amaral and Johnston, 2011) and suggests a putative role of igfbp1 in the differential response to growth stimuli in zebrafish lineages selected for divergent body size at age (Fig. 5A,B).

In conclusion, the experimental selection regime reported in this study successfully produced three lineages of zebrafish with divergent body size. Selection for larger body size had positive effects on offspring during embryonic and adult stages in terms of growth (during both phases) and absolute number of fish produced per spawning. A change in expression of a limited number of genes from the IGF pathway was demonstrated in the present study, which points to a better growth opportunity for fish from the L-lineage in response to a growth stimulus (higher basal expression of igf1a and lower transient expression of igfbp1 genes in response to fasting–refeeding). The current study provides a good starting point for an investigation of the potential epigenetic and genetic mechanisms underlying domestication and artificial selection in fish and suggests that zebrafish may be an excellent model species for this purpose.

We thank Dr Vera Vieira-Johnston for help with the fish husbandry.

 
• bmp15

  bone morphogenic factor 15

 
• bmpr1aa

  bone morphogenic factor receptor 1aa

 
• bmpr1ab

  bone morphogenic factor receptor 1ab

 
• bmpr1ba

  bone morphogenic factor receptor 1ba

 
• bmpr1bb

  bone morphogenic factor receptor 1bb

 
• bmpr2a

  bone morphogenic factor receptor 2a

 
• bmpr2b

  bone morphogenic factor receptor 2b

 
• d.p.f.

  days post-fertilization

 
• FL

  fork length

 
• gdf9

  growth differentiation factor 9

 
• GH

  growth hormone

 
• ghra

  growth hormone receptor a

 
• ghrb

  growth hormone receptor b

 
• H

  maximum body depth

 
• igf1a

  insulin-like growth factor 1a

 
• igf1ar

  insulin-like growth factor receptor 1a

 
• igf1b

  insulin-like growth factor 1b

 
• igf1br

  insulin-like growth factor receptor 1b

 
• igf2a

  insulin-like growth factor 2a

 
• igf2r

  insulin-like growth factor receptor 2

 
• igfbp1a

  insulin-like growth factor binding protein 1a

 
• igfbp1b

  insulin-like growth factor binding protein 1b

 
• L-lineage

  zebrafish lineage experimentally bred to have a larger body size at age compared with the founding population

 
• M_b

  body mass

 
• myf5

  myogenic factor 5

 
• myf6

  myogenic factor 6 (also known as mrf4, myogenic regulatory factor 4)

 
• myoD

  myogenic differentiation 1

 
• myog

  myogenin

 
• PI3K/AKT/mTOR

  phosphatidylinositol 3-phosphate kinase/v-akt murine thymoma viral oncogene/mammalian target of rapamycin

 
• SL

  standard length

 
• S-lineage

  zebrafish lineage experimentally bred to have a smaller body size at age compared with the founding population

 
• TGF-β

  transforming growth factor β

 
• TL

  total length

 
• U-lineage

  control lineage from which no fish were excluded