Transcriptional analysis of abdominal fat in chickens divergently selected on bodyweight at two ages reveals novel mechanisms controlling adiposity: validating visceral adipose tissue as a dynamic endocrine and metabolic organ

Background Decades of intensive genetic selection in the domestic chicken (Gallus gallus domesticus) have enabled the remarkable rapid growth of today’s broiler (meat-type) chickens. However, this enhanced growth rate was accompanied by several unfavorable traits (i.e., increased visceral fatness, leg weakness, and disorders of metabolism and reproduction). The present descriptive analysis of the abdominal fat transcriptome aimed to identify functional genes and biological pathways that likely contribute to an extreme difference in visceral fatness of divergently selected broiler chickens. Methods We used the Del-Mar 14 K Chicken Integrated Systems microarray to take time-course snapshots of global gene transcription in abdominal fat of juvenile [1-11 weeks of age (wk)] chickens divergently selected on bodyweight at two ages (8 and 36 wk). Further, a RNA sequencing analysis was completed on the same abdominal fat samples taken from high-growth (HG) and low-growth (LG) cockerels at 7 wk, the age with the greatest divergence in body weight (3.2-fold) and visceral fatness (19.6-fold). Results Time-course microarray analysis revealed 312 differentially expressed genes (FDR ≤ 0.05) as the main effect of genotype (HG versus LG), 718 genes in the interaction of age and genotype, and 2918 genes as the main effect of age. The RNA sequencing analysis identified 2410 differentially expressed genes in abdominal fat of HG versus LG chickens at 7 wk. The HG chickens are fatter and over-express numerous genes that support higher rates of visceral adipogenesis and lipogenesis. In abdominal fat of LG chickens, we found higher expression of many genes involved in hemostasis, energy catabolism and endocrine signaling, which likely contribute to their leaner phenotype and slower growth. Many transcription factors and their direct target genes identified in HG and LG chickens could be involved in their divergence in adiposity and growth rate. Conclusions The present analyses of the visceral fat transcriptome in chickens divergently selected for a large difference in growth rate and abdominal fatness clearly demonstrate that abdominal fat is a very dynamic metabolic and endocrine organ in the chicken. The HG chickens overexpress many transcription factors and their direct target genes, which should enhance in situ lipogenesis and ultimately adiposity. Our observation of enhanced expression of hemostasis and endocrine-signaling genes in diminished abdominal fat of LG cockerels provides insight into genetic mechanisms involved in divergence of abdominal fatness and somatic growth in avian and perhaps mammalian species, including humans. Electronic supplementary material The online version of this article (doi:10.1186/s12864-017-4035-5) contains supplementary material, which is available to authorized users.


(Continued from previous page)
Conclusions: The present analyses of the visceral fat transcriptome in chickens divergently selected for a large difference in growth rate and abdominal fatness clearly demonstrate that abdominal fat is a very dynamic metabolic and endocrine organ in the chicken. The HG chickens overexpress many transcription factors and their direct target genes, which should enhance in situ lipogenesis and ultimately adiposity. Our observation of enhanced expression of hemostasis and endocrine-signaling genes in diminished abdominal fat of LG cockerels provides insight into genetic mechanisms involved in divergence of abdominal fatness and somatic growth in avian and perhaps mammalian species, including humans.

Background
The domestic chicken (Gallus gallus domesticus) is a widely used biomedical model and serves as a major source of high-quality dietary protein for humans. Decades of intensive genetic selection have led to the remarkable growth rate and feed efficiency of commercial broiler chickens. However, this rigorous genetic selection for growth rate has led to increased adiposity, skeletal abnormalities, and disorders of metabolism and reproduction [1][2][3][4]. Few studies have attempted to identify biological pathways and gene networks that promote abdominal fatness and related metabolic disorders in the chicken. Previous studies of global gene transcription in different models of chicken growth have concentrated on either skeletal muscle [5] or the hypothalamus [6] for identification of candidate genes responsible for differences in growth rate. Another study compared intra-muscular adipose tissue between two lines of chickens selected for fast growth or slow growth [7]. This microarray analysis showed that expression of several differentially expressed (DE) genes was correlated with increased or decreased growth of breast muscle and intramuscular fat, which suggests that adipose tissue per se could regulate the rate of muscle growth, albeit no mechanisms related to abdominal fatness were uncovered.
A large volume of research has been published on growth and metabolic characteristics of the Virginia Tech (VT) population [8,9] of chickens, which were divergently selected on body weight at 56 days of age (8 wk) [10][11][12][13][14][15][16][17][18][19]. The VT low-weight strain (LWS) chickens exhibit diminished growth, anorexia [15], impaired endocrine signaling [13,15,19,20], and higher expression of genes involved in catabolism of lipid [17]. A genomic scan of the VT high-weight (HWS) and low-weight (LWS) selected chickens revealed several quantitative trait loci (QTL) and their associated polymorphic genes related to divergent selection for body weight at 8 wk [9]. Among the candidate genes identified by combined analyses of genome sequence variation in the VT HWS and LWS chickens were glucagon (GCG), insulin-like growth factor binding protein-2 (IGFBP2) and endothelin-1 (EDN1). Furthermore, we found that these three genes (GCG, IGFBP2 and END1) were differentially expressed in abdominal fat of genetically fat (FL) and lean (LL) chickens [21,22]. Our transcriptional study of abdominal fat in the FL and LL chickens also revealed over-expression of numerous lipogenic genes in the FL, which suggests that adipose tissue has a more significant involvement in the synthesis and metabolism of lipids than previously thought [21]. More importantly, we discovered over-expression of an unusually large number of genes involved in hemostasis, endocrine signaling and lipid catabolism [22] in the diminished abdominal fat of LL chickens.
The present study focuses on the transcriptional analysis of abdominal fat in a random-bred population of Bresse-Pile (meat-type) chickens, which were divergently selected by Ricard [23] for either high growth (HG) or low growth (LG) body weight at two developmental ages [juvenile (8 wk) and adult (36 wk)]. The diverse growth curves [24] and muscle growth patterns [25][26][27][28][29] of the HG and LG chickens have been described in detail. These unique growth models were chosen for our original functional genomics project aimed at identification of genomic regions and gene networks, which control growth and metabolism of the broiler chicken [30][31][32]. Despite genetic selection for a large difference in bodyweight (3.2-fold at 7 wk), these divergent lines exhibit an even greater difference in abdominal fat weight (19.6-fold at 7 wk). Furthermore, we have published several papers from quantitative trait loci (QTL) analyses using the F2 resource population created from an intercross of HG and LG chickens [33][34][35][36]. These genetic analyses have revealed numerous genomic loci containing positional candidate genes associated with several metabolic and growth traits. Another study of this unique F2 resource population has identified a cis expression QTL (eQTL) controlling β-carotene 15, 15′-monooxygenase (BCMO1) expression, which consequently determines the extent of yellow coloring in chicken breast muscle, an important meat quality trait [37]. The LG chickens, which carry the inactive BCMO1 gene, exhibit a slower growth rate, greatly-reduced abdominal (visceral) fatness, lower plasma insulin, and hyperglycemia when compared to the HG chickens [34].
Thus, the present study aimed to identify differential gene expression in abdominal fat of HG and LG chickens, which show a greater divergence in abdominal fatness that is incidental to their selection on body weight at two ages (juvenile and adult). Our analyses reveal overexpression of genes in abdominal fat of the HG birds that are involved in increased adiposity, including transcriptional regulators and metabolic (lipogenic) enzymes, throughout juvenile development (1-11 wk). Conversely, LG chickens shown up-regulation of several energy generating processes (i.e., peroxisomal β-oxidation, mitochondrial β-oxidation, ketogenesis and oxidative phosphorylation) early in juvenile development which are likely responsible for their extreme leanness. The RNA-Seq analysis of abdominal fat at 7 wk. revealed up-regulation of several hemostatic factors in the LG cockerels that could contribute to their extreme leanness. Furthermore, this transcriptional study of visceral adiposity in HG and LG meat-type (Bresse-Pile) chickens also serves as cross-validation of abdominal fat as a dynamic endocrine and metabolic organ as indicated in our genetically fat (FL) and lean (LL) chicken lines [21,22].

Animal management and tissue preparation
The chickens used in this study were divergently selected from a population of Bresse-Pile (meat-type) chickens by Ricard [23,38] for an extreme difference in body weight at two ages: 8 (juvenile) and 36 (adult) weeks of age (wk). Chickens were bred and raised at INRA UE1295 Pôle d'Expérimentation Avicole de Tours, F-37380, Nouzilly, France. At hatching, HG and LG cockerels were wing-banded and vaccinated against Marek's disease virus. Birds were provided with ad libitum access to water and fed a conventional starter ration (22% crude protein and 3050 kcal ME/kg) from hatching to 3 wk and then with a grower pelleted ration from 3 to 11 wk. (20% crude protein, and 3100 kcal). The HG birds were separated from LG birds for the first 3 weeks (at which time LG chickens were provided crushed feed pellets) to increase early survival of the LG; afterwards, both lines were placed together and raised in floor pens (4.4 m × 3.9 m). Continuous incandescent light was provided for the first two days followed by a maintenance of a 14 h light /10 h dark cycle (14 L:10D). Infrared gas heaters provided supplemental heat and the ambient temperature was decreased progressively from 32°C at hatching, until 22°C was reached at 22 days. At 1, 3 5, 7, 9 and 11 wk, eight fed cockerels from each genetic line (HG and LG) were randomly selected, weighed and bled into heparinized syringes prior to cervical dislocation, and the excision and weighing of abdominal fat mass. Abdominal adipose tissue samples were immediately snap frozen in liquid nitrogen and stored at −75°C until further processing for RNA analysis.

RNA extraction
Abdominal fat aliquots from forty-eight individuals (4 HG and 4 LG per age at 1, 3, 5, 7, 9 and 11 wk) were homogenized and total cellular RNA extracted using guanidine thiocyanate and CsCl gradient purification [39], followed by DNase I treatment. The quality of RNA was determined with an RNA 6000 Nano Assay kit and the Model 2100 Bioanalyzer (Agilent Technologies; Palo Alto, CA). All samples used for RNA analyses had an RNA integrity number (RIN) greater than 9.0.

Microarray analysis and statistical analysis of microarray data
The Del-Mar 14 K Chicken Integrated Systems Microarrays (Geo Platform # GPL1731) described earlier [31] were used for transcriptional profiling of four abdominal fat samples from each genotype (HG and LG), across 11 weeks of juvenile development (48 total individuals). Methods used for microarray preparation including labeling, hybridization, and image acquisition were described earlier [21]. Briefly, twenty-four Del-Mar 14 K Chicken Integrated Systems Microarrays were hybridized with 48 labeled samples using a balanced block design [40], where half of the birds from each genotype and age were labeled with Alexa Flour® 647 (red dye) and the other half with Alexa Flour® 555 (green dye); see Additional file 1 for experimental design. These hybridized microarrays were scanned with a GenePix 4000B scanner using GenePix Pro 4.1 software (Molecular Devices, Union City, CA) at wavelengths of 635 nm (Alexa® 647-labeling) and 532 nm (Alexa® 555-labeling) producing a combined TIFF image file for each slide. Laser power was set at 100% with the photomultiplier tube (PMT) setting adjusted for each scan producing a PMT count near unity. All slides were checked manually for quality, and all spots with inadequacies in signal, background or morphology were eliminated. The image analysis results were merged with Excel files (in GPR format) containing clone identification, spot location on slide, and most current gene name/function (based on BLAST score). The gene list were first annotated by clone ID, GenBank ID as determined by BLASTN or BLASTX analysis using the GeneBase function on our laboratory website [41].
The GPR files were used to determine differential expression in abdominal fat of HG and LG chickens. Log2 transformed median intensity values (for each dye) were normalized using a global LOWESS transformation (without background subtraction) to remove dye bias within microarray [42]. A two-way ANOVA was used to determine main effects of age (A) and genotype (G), and interactions between age and genotype (A x G); differences between genotypes at each age were also determined. The Benjamini-Hochberg procedure [43] was used to control experiment-wise false discovery rate (FDR) associated with multiple testing. The differentlyexpressed (DE) gene lists from the time course (1-11 wk) microarray analysis were functionally annotated using the GeneBase function on our laboratory website [41] and finally the Ingenuity Knowledge Base [44]. Expression values of 25 DE genes at 7 wk were retrieved from the microarray analysis for comparison across the three methods. For these 25 genes, a Student's T-test was used to determine significant differences between genotypes.

RNA-sequencing and statistical analyses
The same eight RNA samples (4 HG and 4 LG at 7 wk), originally used for the microarray analysis, were also used for construction of indexed (bar-coded) sequencing libraries. Libraries were made from 2 μg of total adipose RNA with the Illumina TruSeq® Stranded mRNA library preparation kit following standard Illumina protocols. All eight barcoded libraries were pooled and paired-end sequenced (101-bp reads) in duplicate lanes on an Illumina HiSeq 2000 Sequencing System (Illumina, Inc., San Diego, CA) at the Delaware Biotechnology Institute, University of Delaware (Newark, DE).
Sequences were trimmed for quality using a combination of custom Perl scripts and Btrim64 software [45]. Boxplot graphing of pre-and-post trimming reads confirmed the absence of outlier samples based on read count. After trimming, reads were mapped to the chicken genome assembly Galgal4.0 (down-loaded from Ensembl) using Tophat (version 1.3.3), followed by assembly and quantitation using Cufflinks software (v1.3.0). The fragments per kilobase of exon per million fragments mapped (FPKM) threshold for detection of a gene was set at FPKM > 0.5. The resulting gtf files were merged with cuffmerge, and differential expression was assessed using Cuffdiff. The two-sided P-value was corrected using the false discovery rate (FDR) which accounts for multiple testing procedures [43]. Genes with a FDR adjusted P-value (P ≤ 0.05) and fold change ≥ 1.2 were considered to be differentially expressed (DE) transcripts. This detection threshold was based on our extensive transcriptional analyses of multiple tissues in the chicken.

Ingenuity Pathway Analysis of differentially expressed genes from microarray and RNA-Seq analyses
Differentially-expressed (DE) gene datasets from the time-course microarray analyses were submitted to Ingenuity® Pathway Analysis (IPA) [44] for functional annotation as "Analysis Ready" (AR) genes according to annotated mammalian genes and proteins accrued in the Ingenuity® Knowledge Base. The unique and commonly shared (intersect) gene sets were then used for Ingenuity Up-Stream Regulator analysis to identify transcription factors and their direct target genes.

Quantitative RT-PCR analysis
Candidate DE genes were selected for verification of expression by quantitative RT-PCR (qRT-PCR) analysis from both the time-course microarray analysis (1-11 wk) and RNA-Seq analysis (7 wk). Superscript III reverse transcriptase (Invitrogen) and an oligo (dT) primer were used to prepare cDNA from 1 μg of RNA. Primers were designed using Primer Express v2.0 software (Applied Biosystems, Foster City, CA). Detailed information for each primer pair including gene name, gene symbol, primer sequences (forward and reverse), GenBank accession number and amplicon size is provided in Additional file 2.
The qRT-PCR assays were performed in an ABI Prism Sequence Detection System 7900HT using 50 ng of cDNA, Power SYBR green PCR master mix (Applied Biosystems, Foster City, CA) and 400 nM of each primer (forward and reverse; Sigma-Aldrich, St. Louis, MO) in duplicate wells. Disassociation curves were analyzed to confirm specific amplification and to verify absence of primer dimers. Resulting PCR products were subjected to agarose gel electrophoresis to compare PCR product size to expected amplicon size. To verify gene expression from the microarray and RNA-Seq analyses, the cycle time (Ct) for each sample was normalized to the corresponding sample geometric mean of two housekeeping genes: cytochrome c oxidase subunit VIIa polypeptide 2 like (COX7A2L) and ribosomal protein L14 (RPL14). The housekeeping genes were selected based on invariability in the microarray and RNA-Seq analyses. Their stable expression in qRT-PCR analysis was determined by Biogazelle qbase+ software [46]. The 2 -(ΔΔCt) formula was used to calculate relative transcript abundance [47]. The statistical analysis was performed using a general linear model procedure in SAS v9.3 and differences between genotypes at each age were determined using Tukey's multiple comparisons test. For genes only analyzed at 7 wk., a Student's T-test was used to detrmine differential expression. The significance level for statistical analysis was set at P ≤ 0.05.

Independent bioinformatic analysis of RNA-Seq analyses of abdominal fat across four divergent genotypes
An independent analysis was conducted on two deposited RNA-Seq datasets of abdominal fat (7 wk) in four distinct genotypes/phenotypes, which were created by divergent genetic selection on juvenile and adult BW (HG vs. LG; NCBI GEO Series Accession # GSE49121) or abdominal fatness at the same BW [fat line (FL) vs. lean line (LL); NCBI GEO Series Accession # GSE42980)]. The systems biology analysis of these two RNA-Seq datasets was performed by the Animal Systems Biology Analysis and Modeling Center (ASBAMC), University of California at Davis, CA. The components of the custom bioinformatics pipeline used for systems biology analysis of domestic animals are described in detail on the project website [48]. The main purpose of this independent metaanalysis of abdominal fat transcriptomes across four distinct genotypes (HG-LG; FL-LL) at the same age (7 wk) was to identify a set of common and unique DE gene across the four genotypes for further functional and pathway analyses using IPA.

Microarray and RNA-Seq analyses of abdominal fat gene expression
Differentially-expressed (DE) genes were defined as those having a significant adjusted P-value and false discovery rate (FDR≤0.05). Statistical analysis of the time-course microarray study provided significant gene sets from three contrasts: the main effect of genotype (312) DE genes) or age (2918 DE genes), and the interaction between genotype and age (718 DE genes). These annotated DE gene sets from the time-course microarray study are presented as HG/ LG expression ratios in Additional file 3. These DE gene sets were used as input files for IPA and functionally annotated as "Analysis Ready" (AR) DE genes according to the Ingenuity Knowledge Base, which is largely based on annotations accrued from the human and murine biomedical literature.
RNA sequencing of abdominal fat at 7 wk. yielded 65.8 million (M) paired-end (101 bp) reads in HG cockerels and 66.6 M paired-end reads in the LG ( Table 1). The percentage of reads mapped to chicken transcripts was 80.7% for the HG and 82.6% for the LG birds. A power analysis (Additional file 4) was conducted on the RNA-Seq dataset, using the webbased software program "Scotty" [49,50], to demonstrate sufficient biological sample size and sequencing depth for detection of DE genes. Power was calculated using the average of 50 M reads per sample at three levels of fold-change (≥1.5-, 2-, or 3-fold change) between HG and LG chickens at a significance of P≤0.05. The "Scotty" program also provided a hierarchical cluster analysis using the Spearman correlation as the distance metric to demonstrate relatedness among the eight individual (4 HG and 4 LG) birds used for RNA-Seq analysis of abdominal fat at 7 wk. The 4 HG and 4 LG abdominal fat samples are tightly clustered according to their genotype. Statistical analysis of the RNA-Seq dataset identified 2410 DE (FDR ≤ 0.05) genes in abdominal fat of HG versus LG (HG/LG log2 expression ratio) chickens at 7 wk. (Additional file 5).
From the time-course microarray analysis, IPA annotated 329 "Analysis Ready" (AR)-DE genes from the main effect of genotype, 559 AR-DE genes from the interaction of genotype x age, and 647 AR-DE genes in a non-redundant dataset combined from the main effect of genotype and genotype x age interaction ( Fig. 2-a). The second Venn diagram ( Fig. 2-b) presents unique and commonly shared AR-DE genes from the time-course microarray study (genotype and age x genotype interaction) and 2026 AR-DE genes from RNA-Seq analysis of abdominal fat in HG and LG chickens at 7 wk. The number of DE genes considered by IPA as "Analysis Ready" (AR-DE) is shown in parenthesis.The Venn diagram ( Fig. 2-a) presents unique and commonly shared DE gene sets found in HG and LG abdominal fat from the time-course microarray study. Only "Analysis Ready" (AR) gene sets, annotated with the Ingenuity Knowledge Base, are represented in the Venn diagram. An additional non-redundant gene set (combined-unique, 647 DE genes) was assembled by combining the genotype (329 DE genes) and age x genotype (559 DE genes) data sets and removing duplicated genes (cDNAs) printed on the microarray. A fourth dataset of 492 commonly-shared DE genes, found in these intersects, was also used for Ingenuity® Pathway Analysis (IPA). This commonly-shared DE gene dataset seems highly enriched with genes controlling the divergence in growth and abdominal fatness traits in HG and LG cockerels.
The RNA-Seq analysis ( Fig. 2-b) of abdominal fat from the same 4 HG and 4 LG birds at 7 wk. provided 2026 DE genes that were annotated by IPA as "Analysis Ready" (AR). The intersection of the genotype and age x genotype datasets from the microarray study shows 126 commonly shared AR-DE genes. The RNA-Seq AR-DE gene set shared 86 common genes with the main effect of genotype AR-DE genes, and another 130 genes in common with the age x genotype AR-DE gene set. Only 47 AR-DE genes were shared among all three DE gene sets, derived from microarray and RNA-Seq analyses.

IPA of gene interaction networks and functional pathways
Time-course microarray data analysis of abdominal fat in HG and LG cockerels (1-11 wk) Significant gene (FDR≤0.05) lists from the microarray analysis were first annotated using the GeneBase tool on our laboratory website [41], which provides protein IDs (from GenBank or Swiss-Prot) derived from BLASTX/ BLASTN analysis of 18,240 cDNA probes printed on the Del-Mar 14 K chicken array. Annotated microarray data files ( Fig. 2-a), containing the Ensembl protein ID and  LG (red triangles) cockerels at six ages (1-11 wk). Each symbol represents the mean ± SE of 8 individual chickens from each genotype. Average abdominal fat content (b) (% bodyweight, %BW) of 8 individual birds from the HG and LG chickens; four chickens per age and genotype were randomly selected for transcriptional analyses. Significant differences between genotypes at each age were determined using a one-way analysis of variance (ANOVA) and Tukey's multiple comparisons procedure at a significance level of P≤0.01 (**) or P≤0.001 (***). Fold differences between genotypes (HG/LG) in body weight (kg) and abdominal fat (%BW) from 1 to 11 wk. (c). The maximum divergence in body weight and abdominal fatness occurs at 7 wk., which was the age selected for deeper RNA-Seq analysis log2 fold-difference for each gene, were submitted for Ingenuity® Pathway Analysis [44]. IPA provides functional annotation from the Ingenuity Knowledge Base [indicated as "Analysis Ready" (AR)], mapping to canonical metabolic/regulatory pathways, gene interaction networks and Ingenuity® Upstream Regulator Analysis, which provides transcription factor interaction networks and predicted interactions between transcription factors and their target genes. A summary of the IPA of the time-course (1-11 wk) microarray dataset of 647 unique AR-DE genes condensed from the main effect of genotype and interaction of genotype x age datasets is presented in Table 2. The top five canonical pathways identified by microarray analysis were "Oxidative Phosphorylation" (19 AR-DE genes; Additional file 6), "Mitochondrial Dysfunction" (22 AR-DE genes), "Eukaryotic translation initiation factor (EIF2) Signaling" (21 AR-DE genes), "NRF2-mediated Oxidative Stress Response" (19 AR-DE genes), and "Phagosome Maturation" (13 AR-DE genes). The top five "Upstream Regulators" were TP53 (105 direct target Table 2. Ingenuity Pathway Analysis (IPA) was used for functional annotation and mapping of 647 DE genes from the time-course microarray analysis of abdominal fat in HG and LG cockerels (1-11 wk) that were "Analysis Ready" (AR-DE). This unique (non-redundant) gene set was compiled from AR-DE genes found in the main effect of genotype or the genotype x age interaction other hand, PPARD and many of its direct target genes (10 AR-DE genes) are up-regulated in the slower-growing and leaner LG genotype, with the exception of two genes (SCD and CD36), which were highly expressed in the HG. The two opposing transcription factors (PPARG and PPARD) also share several down-regulated target genes in the HG, including PDK4, ALDH9A1, ACSL, SLC27A1 and or chaperonin (HSPD1) and hydroxysteroid (11-beta) dehydrogenase 1 (HSD11B1), which degrades glucocorticoid, interact with direct targets of PPARD (i.e., VLDLR, THBS1 and CD36), whereas, THBS2 and collagen type VI were up-regulated in visceral fat of the HG chickens. Interestingly, an equal number of unique ATPases were associated with either PPARG in the HG or PPARD in the LG chickens.
Numerous direct targets of PPARG were identified by Ingenuity Upstream Regulator Analysis from DE genes identified by the time-course (1-11 wk) microarray study of abdominal fat of HG and LG chickens ( Fig. 3-b). Twenty-seven direct DE genes of PPARG were expressed at higher levels in abdominal fat of the HG chickens, whereas 13 direct gene targets were more abundant in the LG chickens. This detailed view of PPARG target genes predicts (based on literature accrued in the Ingenuity Knowledge Base and observed expression values) that PPARG up-regulates or activates 31 direct gene targets (orange arrows) in the HG chickens. Among these up-regulated genes, most control lipogenesis, adiposity and energy metabolism (SCD, FASN, FN1, DGAT2, CEBPA, LPL, FABP4, INSIG2, IRS1, PCK1, SCD, SOD1, SERPINA1, and ACSL1). Whereas, PPARD expression was higher in the LG chickens and predicted to be inhibited (blue color arrows) and downregulated in the HG chickens (or higher in the LG), which agrees with our observed expression values ( Fig. 3-b). Eight down-regulated genes (PDK4, PLP1, SIRT5, SLC22A3, THBS1, VLDLR, ALDH9A1, ANGPTL3, and CSAD) were predicted to be inhibited by PPARD (indicated by the blue arrows) in visceral fat of the HG; these genes were actually expressed higher in the slow-growing LG chickens as indicated by the green-colored gene symbols.
Ingenuity Upstream Regulator Analysis identified 18 AR-DE genes as direct targets of CEPBA ( Fig. 4-b), 13 lipogenic AR-DE genes were up-regulated in the HG genotype (FASN, DAGT2, LPL, CEBPA, FABP4 and IRS2). Whereas, only five AR-DE genes were downregulated in the HG birds (ALB, KLF5, HSD11B1, HLA-B, and GBP1) in the subset of 252 commonly-shared genes from microarray contrasts (see Fig. 2-a). The Ingenuity Upstream Regulator Analysis predicts that CEBPA was activated in the HG chickens and led to activation of eight genes (orange arrows) or inhibition (blue blunted line) of a single gene, major histocompatibility complex, class I, B (HLA-B). Based on our observed expression values and expected responses from the mammalianbased Ingenuity Knowledge Base, the Upstream Analysis predicts that KLF2 would be inhibited (blue gene symbol), while clearly its expression was more abundant in HG abdominal fat. KLF2 was predicted to inhibit three major genes (CEBPA, FABP4, and MYH10), highly expressed in the HG, and would further inhibit the embryonic form of hemoglobin, epsilon 1 (HBE1). The yellow arrows point to a contradiction between information accrued in the Ingenuity Knowledge Base and the expected state (inhibited, blue) of two target DE genes [thrombospondin receptor (CD36) and ASS1], although their actual state was upregulated (activated, red symbol). There is predicted uncertainty (grey arrows) about KLF2 inhibition of ID3 (blunted-yellow line), which is down regulated in the HG. Ingenuity correctly predicts FBXW7 to be inhibited (blue symbol), since it was expressed lower in HG, and that this ubiquitin protein ligase would directly block up-regulation of three ligogenic genes (CD36, DGAT2 and FABP4) (blunt orange line), along with uncertainty about the ability of FBXW7 to block (blunt yellow line) the down-regulated KLF5 (green symbol). It is clear from these observations that CEBPA is major lipogenic transcription factor and a major contributor to the higher growth and greater fat accretion found in the HG cockerels.   The other hemostatic gene network revealed by RNA-Seq analysis (Fig. 5-b) was composed of several coagulation factors (F2, F3, F5, F8, F10, SERPINB2, and FGB), which were all over-expressed in abdominal fat of LG cockerels. The complement component (C1S), calreticulin (CALR), the anion exchanger SLC4A1, kinesin family member 7 (KIF7), carbonic anhydrase II (CA2) and regulator of G-protein signaling 19 (RGS19) were also expressed higher in the LG. Several additional genes in this direct interaction network were up-regulated in the HG, including mitogenactivated protein kinase 1 (MAPK1), parvin alpha (PARVA), the monocarboxylate transporter SCL16A7, coiled-coil domain containing 2 (GCC2), CDC42 binding protein kinase alpha (CDC42BPA), LDL receptor related protein 6 (LRP6), frizzled-related protein (FRZB), wingless-type MMTV integration site family member 16 (WNT16), protein kinase C and casein kinase substrate in neurons 3 (PACSIN3), fibrinogen-like 2 (FGL2) and mitochondrial fission regulator 1 (MTFR1). This gene interaction network was functionally annotated by IPA as "Hematological Disease".
Sixteen of the 17 DE genes were over expressed in the IPA canonical coagulation system, which includes the extrinsic and intrinsic prothrombin activation pathways ( Fig. 6; Additional file 7). The up-regulated hemostasis genes found in the LG chickens include six coagulation factors (F2, F3, F5, F8, F10 and F13B), three fibrinogen subunits (FGA, FGB and FGG), two serpin peptidase inhibitors (SERPINC1 and SERPINF2), plasminogen (PLG), kininogen 1 (KNG1), bradykinin (BDK), protein C (PROC), and von Willebrand factor, C and EGF domains (VWCE). In contrast, plasminogen activator, urokinase (PLAU) was the only coagulation-related gene that was expressed higher in HG abdominal fat. Clearly, the diminished abdominal fat mass of slower-growing LG cockerels reflects a highly prothrombotic state; although the consequence on endocrine signaling and/or lipid metabolism remains unknown. However, microarray analysis of liver in the HG and LG cockerels does show higher expression of four coagulation genes (FGA, FGB, SERPINF2, and KNG1) in the LG birds (Cogburn, LA unpublished data). Presently, we do not know if local activation of the coagulation system in adipose tissue (or liver) affects the systemic hemostatic mechanism, which normally would rely on availability of coagulation precursors synthetized mainly in the liver.
One metabolic gene network, annotated by IPA as "Lipid Metabolism, Small Molecule Biochemistry and Energy Production", was centered on interactions of two transcription factors (SREBF1 and SREBF2) and their numerous activated target genes, most of which were highly expressed in abdominal fat of the faster-growing and fatter HG cockerels at 7 wk ( Fig. 7-a). The two highest-expressed genes found in this lipogeneic network were SCD (10.7-fold higher) and 24-dehydrocholesterol reductase (DHCR24) [5. Ingenuity® Upstream Regulator Analysis identified 36 AR-DE genes that are direct targets of SREBF1 (22 genes were expressed higher in the HG, which is consistent with activated SREBF1) and 16 AR-DE genes that are direct targets of SREBF2 (12 genes were expressed higher in the HG, which is consistent with activation of SREBF2) (Fig. 7-b). As direct targets of SREBF1, the four most abundant genes in abdominal fat of the HG were myosin light chain 1 (MYL1), SCD, troponin I type 2, skeletal, fast (TNNI2) and FADS2, whereas MAT1A, APOA5, TF and ADH1C were the highest-expressed LG cockerels. Panel a shows an interaction network functionally annotated by IPA as "Cell-to-Cell Signaling, Hematological System Development and Function". The second panel (b) shows a gene network, annotated by IPA as "Hematological Disease", which involves interaction of a large cluster of clotting factors with multiple cell signaling components target genes of SREBF1 found in the LG at 7 wk. The Ingenuity Upstream Analysis predicts that both SREBF1 and SREBF2 are activated (orange gene symbols) and that 19 DE target genes would be activated (orange arrows). Twelve of these highly-expressed genes are activated by both transcription factors (ACSL1, CYB5A, DHCR7, FADS2, FASN, INSIG1, LSS, MSMO1, PCYT1A, SC5D, SCD and SQLE). SREBF1 appears to inhibit only three genes [alcohol dehydrogenase 1C (class I), gamma polypeptide (ADH1C), apolipoprotein A-V (APOA5) and androgen receptor (AR)], which were over-expressed in the LG chickens. Only 13 direct target genes of SREBF1 and/or SREBF2 were expressed higher in abdominal fat of the LG chickens. Perilipin 2 (PLIN2) was the only direct target of SREBF2 that was expressed higher in visceral fat of the LG cockerels.
RNA-Seq analysis revealed interactions between the vitamin D receptor (VDR) and several direct targets of SREBF1 including FASN, ACACA and LPIN1, which were upregulated in the HG birds ( Fig. 8-a) Ingenuity Upstream Regulator Analysis predicts that the expression of VDR (Fig. 8-b) would be inhibited (blue gene symbol), which would lead to downregulation of four direct target genes [serpin peptidase inhibitor, clade C (antithrombin), member 1 (SER-PINC1), transmembrane protease, serine 2 (TMPRSS2), ATP binding cassette subfamily C member 3 (ABCC3), and cadherin 1 (CDH1). Of the 17 DE direct targets of the VDR, three genes [angiotensinogen converting enzyme (ACE), myosin, heavy chain 3, skeletal muscle, embryonic (MYH3) and phosphate regulating endopeptidase homolog, X-linked (PHEX)] were expressed at higher levels in abdominal fat of the HG chickens. The up-regulated expression of six genes in the LG birds (AGT, AGTR1, F3, FGG, PNPLA2 and SELE) was inconsistent with the expected state of these genes (blunted yellow edges) from the Ingenuity Knowledge Base. However, the effect of the VDR was not predicted for three AR-DE genes [the inflammatory chemokine, C-X-C motif chemokine ligand 8 (CXCL8), phosphate regulating endopeptidase homolog, X-linked (PHEX) and perilipin 2 (PLIN2).

Higher Expression in HG AbFAT
Another gene network functionally annotated by IPA as controlling lipid metabolism was centered on the androgen receptor (AR), which was expressed higher in abdominal fat of the LG cockerels, and the interaction of its direct target genes with three additional transcription factors [SREBF2, PPARG coactivator 1 beta (PPARGC1B) and MKL1/myocardin like 2 (MKL2)] ( Fig. 9-a). Eight genes are direct targets of the lipogenic transcription factor SREBF2, where CYB5A, FADS2, DHCR7, MSMO1, SQLE, LSS and SCD were up-regulated in HG visceral fat, and only solute carrier family 4 member 4 (SLC4A4) was expressed higher in the LG. Two direct targets of the AR [methylsterol monooxygenase 1 (MSMO1) and 7dehydrocholesterol reductase (DHCR7)] are also direct targets of SREBF2, while DHCR24 is a direct target of the AR and PPARGC1B, both of these transcription regulators were up-regulated in LG chickens. PPARGC1B also shares three genes that are direct targets of SREBF2 LG cockerels at 7 wk. (Fig. 9-b). Thirty-five of the known direct targets of the AR were upregulated in HG cockerels, while 50 DE genes were over expressed in visceral fat of the LG birds. The predicted inhibition of the AR and 26 of its target genes is based on the observed down-regulation of the AR (blue symbol and arrows) in the HG abdominal fat (or up-regulation in LG birds as indicated by green gene symbols) from the RNA-Seq LG cockerels at 7 wk. shows a gene interaction network centered on interactions between SREBP2 and the AR (a). SREBP2 and its direct targets are expressed higher in the HG cockerels, while the AR and 8 target genes were up-regulated in the LG. Ingenuity Upstream Analysis predicts that the AR should be activated in the HG (orange gene symbol), based on the up-regulated condition (orange arrows and red gene symbols) of 27 target genes (b). An additional 22 DE genes, known direct targets of the AR, were up-regulated (green gene symbols) in abdominal fat of the LG birds. Ingenuity predicted activation is indicated by orange arrows), while predicted inhibition is shown by blunted blue lines

Verification of differential gene expression by quantitative RT-PCR (qRT-PCR)
Based on biological function, candidate DE genes were selected from the microarray (1-11 wk) and RNA-Seq (7 wk) analyses for verification of expression by qRT-PCR assay. Three genes (ME1, DIO3, and scGH), not identified by either microarray or RNA-Seq analysis, were also examined by qRT-PCR analysis due to special interest. Expression patterns of four transcriptional regulators which directly regulate adipogenesis and/or lipogenesis are shown in Fig. 10-a. Three of these transcriptional regulators [PPARG, CEBPA and thyroid hormone responsive spot 14, alpha (THRSPA)] show very similar patterns of expression, being up-regulated in HG chickens. Similarly, SREBF1 was 4-fold higher in HG chickens at 1 wk. (P ≤ 0.01), although not different at 3 and 11 wk. Four targets of these transcriptional regulators, controlling lipogenesis, are shown in Fig. 10-b. Interestingly, SCD was the second highest-expressed gene identified in abdominal fat of HG chickens by time-course (1-11 wk) microarray analysis (Table 2), being over~110-and~90-fold higher (P≤0.001) at 1 and 7 wk., respectively. Malic enzyme 1, NADP(+)-dependent, cytosolic (ME1) was up-regulated in the HG at 1, 5 and 7 wk., with a peak difference at 5 wk. (~2.5-fold increase in HG chickens). Both FASN and DGAT2 were also significantly up-regulated in HG chickens (P≤0.001) at 1, 7 and 9 wk., with the greatest fold-change observed at 7 wk. (2.8-and 9.7-fold, respectively).
Transcript abundance was also examined by qRT-PCR for seven genes which appear to be associated with leanness (Fig. 11). The nuclear-hormone receptor Data points represent the mean ± SE of 4 birds/genotype and age. Significant differences between genotypes at each age were determined using a one-way ANOVA and Tukey's multiple comparisons procedure at a significance level of P≤0.05 (*), P≤0.01 (**) and P≤0.001 (***) PPARD (see Fig. 2) was up-regulated in abdominal fat of LG chickens between 1 and 5 wk., with an 8.6-fold difference at 1 wk. Carnitine palmitoyltransferase I (CPT1A), which mediates the transport of long chain fatty acids across the outer mitochondrial membrane, was also higher in the LG during early development (1-7 wk) with an average of a 2.4-fold increase across these ages. Two key enzymes in ketogenesis [HMG-CoA synthase 1, soluble (HMGCS1) and 3-hydroxy-3methylglutaryl-CoA lyase (HMGCL)] also exhibited early up-regulation in LG chickens. Interestingly, the endogenous avian leukosis virus envelope protein (envPr57) was highly expressed in abdominal fat of LG chickens across juvenile development (1-11 wk). The truncated or short chicken growth hormone (scGH) transcript was expressed higher (P≤0.01) in abdominal fat of HG chickens at 1 and 7 wk. Differences in the thyroid hormone-activating enzyme, deiodinase, iodothyronine, type I (DIO1), were seen at 3, 5 and 9 wk. (upregulated in LG chickens at all 3 ages) with the largest difference seen at 5 wk. (7.7-fold higher in LG chickens). Conversely, the gene for the thyroid hormone-deactivating enzyme [deiodinase, iodothyronine, type III (DIO3)] was expressed higher in LG chickens at 1 and 3 wk., whereas DIO3 abundance was greater in visceral fat of the HG at 9 and 11 wk.
The differential expression of 46 genes identified by RNA-Seq analyses at 7 wk. was also examined by qRT-PCR at 7 wk. (Table 4). A subset of 25 DE genes, identified  Fig. 11 Verification of differential expression of genes associated with leanness by qRT-PCR analysis. The abundance of seven genes, expressed higher in LG birds and associated with leanness, was also verified by qRT-PCR analysis. An additional gene, short isoform of chicken GH (scGH), was included in this figure, although its expression was higher in the HG at 7 wk. Data points represent the mean ± SE of 4 birds/genotype and age. Significant differences between genotypes at each age were determined using one-way ANOVA and Tukey's multiple comparisons procedure at significance levels of P≤0.05 (*) and P≤0.01 (**). The ‡ symbol denotes a data point that approaches significance (P≤0. 10) by microarray analysis, was included in the three-way comparison. Eleven of these genes were significantly different (P ≤ 0.05) between HG and LG chickens across all three methods (ALB, ALDOB, A2M, EX-FABP, FADS2, FGA, HSD17B7, PDK4, SCD, THBS2, and TTR). Twelve genes were significantly different (P ≤ 0.05) between HG and LG chickens by qRT-PCR, although they did not reach significance level by either RNA-Seq or microarray analyses (AGTR1, ANXA1, CEBPA, F5, FASN, HPGDS, LDHA, LPIN1, LPL, ND6, PPARG, and PYGL). Two genes were not significantly different by qRT-PCR or RNA-Seq analyses, but reached significance in the microarray analysis (PLG and PGRMC1). Expression ratios of these 25 genes for qRT-PCR analysis versus microarray analysis and RNA-Seq analysis versus microarray analysis comparisons had significant Spearman's rank coefficients [rho (P) = 0.824615 and 0.854395, respectively]. The additional twenty-two genes analyzed by qRT-PCR at 7 wk. were similar by RNA-Seq analysis in magnitude, direction of fold change and significance level, except for four genes (GPD1, LCN15, PDE1C and SELP), which were not significant (P≤0.05) by RNA-Seq analysis but reached significance by qRT-PCR analysis. Expression fold-change ratios of forty-seven genes produced a significant Spearman's rank coefficient [rho (P) = 0.928789] across qRT-PCR and RNA-Seq analyses.
Independent analysis of two RNA-Seq datasets of abdominal fat at 7 wk. in chickens divergently selected on abdominal fatness (FL vs. LL) or growth rate (HG vs. LG) The two datasets from RNA-Seq analysis of abdominal fat in the high-growth (HG) versus low-growth (LG) chickens (NCBI GEO Series Accession # GSE49121) and the fat line (FL) versus lean line (LL) broiler chickens (# GSE42980) were independently analyzed by the USDA-funded Animal Systems Analysis and Modeling Center (ASBAMC) [48]. This metaanalysis of abdominal fat transcriptomes across four distinct genotypes revealed 1500 DE genes [adjusted P-value (≤0.05) and FDR (≤0.05)] in the HG-LG cockerels and 653 DE genes in the FL-LL at 7 wk. These two datasets were used as input files for "Comparison Analysis" in IPA, which identified 97 commonly shared DE adipose genes. The top overrepresented canonical pathways found by IPA were "Acute Phase Response Signaling, Coagulation System, and FXR-RXR Activation". The top gene interaction network represented by the 97 commonly shared DE genes in HG vs. LG (Fig. 12-a) and FL vs. LL cockerels (Fig. 12-b) was functionally annotated by IPA as "Cardiovascular/Hematological Disease and Developmental Disorder". This network was composed of a core of fibrinogenic (FGA, FGB, FGG, HRG, and F10), fibrinolytic (PLG) and transporter (ALB, APOA4, and AHSG) genes which were highly expressed in leaner chickens from either divergent line (i.e., the LG and LL). However, the expression pattern of several other genes in this network differed according to genetic background. These divergent genes include G-protein coupled receptors (CHRM2, TRHR, TAS1R3), estrogen biosynthesis (CYP19A1), gluconeogenesis (PCK1), retinol transport (RBP4), and an upstream regulator of triglyceride metabolism (CREB3L3).

Discussion
Our avian models of growth and metabolism were originally derived from a population of Bresse-Pile broiler-type chickens, which were divergently selected by Ricard [23] for either high (HG) or low (LG) body weight (BW) at 8 and 36 wk. Divergent genetic selection on juvenile and adult BW for more than 30 generations has resulted in a 2.7-fold difference in BW of HG and LG cockerels between 1 and 11 wk. (see Fig. 1-a). Perhaps more remarkable was achievement of an even greater difference in visceral fatness, where HG cockerels are 8-fold fatter than the LG birds during juvenile development ( Fig. 1-b and c). Thus, the divergent HG and LG chickens serve as unique avian models to unravel the genetic basis of extremes in fatness and leanness, which are incidental to their primary divergence in BW [30,32]. Divergent genetic selection on abdominal fat alone in the FL and LL chickens has resulted in 2.5-fold difference in visceral fatness at the same body weight [51]. From our current analyses, we found DE genes that were unique to abdominal fat transcriptomes of the HG and LG cockerels, while a common set of DE genes, related to the divergence in abdominal fatness of the HG and LG chickens, was shared with our previous models of divergently selected FL and LL chickens.
For example, our recent transcriptional analyses of abdominal fat in FL and LL chickens revealed overexpression of numerous hemostasis genes in abdominal fat of the leaner LL chickens [21,22]. The present study has extended this novel finding to slower-growing and leaner LG chickens, which also show over-expression of coagulation and fibrinogenic genes in their diminished visceral fat. Hemostatic processes were amongst the top canonical pathways represented by DE genes identified by RNA-Seq and IPA analyses of abdominal fat in the LG cockerels. We found marked up-regulation of several hemostatic genes in visceral fat of LG chickens that belong to both the intrinsic and extrinsic coagulation pathways (see Fig. 6). Other overexpressed components of the coagulation system and fibrinolysis include serine protease inhibitors that control coagulation (SERPINA1, SERPINB2 and SERPINF1) and major transport proteins found in the bloodstream (A2M, ALB, GC, and TTR).   These plasma transporters could be important in maintaining leanness. For example, A2M is known to inhibit proteases (including coagulation factors) and transport cytokines [52]. Likewise, ALB, GC, and TTR are part of the same family of proteins that transport metabolites associated with leanness including: steroids, fatty acids, vitamin D [53], thyroid hormone [54] and retinoid [55]. Numerous studies of obese humans and murine models have demonstrated direct links between hemostasis and development of metabolic disorders (reviewed in [56][57][58][59][60][61][62]). Obese mammals present with a prothrombotic and hypercoagulable state, which is marked by elevated levels of thrombin in circulation. Enhanced production of inflammatory cytokines and adipokines in abdominal fat of obese mammals leads to insulin resistance and impaired fibrinolysis [60,62]. In contrast, we have discovered over-expression of numerous coagulation genes and "ectopic" endocrine factors in reduced abdominal fat of LL chickens, but not the FL [21,22]. We have proposed that these highly expressed coagulation factors (mainly proteases and protease inhibitors) could serve a novel role within visceral fat of genetically lean chickens by controlling the proteolytic activation or deactivation of adipokines or other endocrine factors. The present transcriptional analysis of abdominal fat in divergently selected HG and LG chickens clearly shows up-regulation of several hemostatic genes in the diminished visceral fat of the leaner and slower-growing LG chickens, which validates our idea that the enhanced expression of prothrombic and fibrogenic genes in genetically leaner (LL or LG) chickens contributes to their extreme leanness. Furthermore, we have no indication that enhanced expression of coagulation factors in visceral fat of genetically lean chickens would adversely affect systemic hemostasis, which depends upon hepatic synthesis of most clotting precursors. A recent study of three breeds of chickens with distinct differences in abdominal fatness indicates that increased lipid catabolism and reduced lipid synthesis within adipose tissue largely determine leanness in the chicken [63]. The independent meta-analysis of two RNA-Seq datasets from abdominal fat of the HG and LG lines versus the FL and LL lines at 7 wk. has revealed a common gene interaction network of 97 shared genes that are involved in blood coagulation, endocrine signaling, and energy metabolism (Fig. 12). The leaner (LG and LL) chickens over-express a core set of prothrombotic (FGA, FGB, FGG, F10 and HRG) genes, the fibrinolytic gene (PLG), acute-phase proteins (AHSG and ORM1), and several transporter (ALB, APOA4, HDL and LDL) genes in their diminished abdominal fat. This common "core" network also contains a clusters of divergently-expressed genes related to G-protein-coupled signaling (CHRM2, TRHR, MLNR and TAS1R3). The sweet taste receptor (TAS1R3) was expressed higher in abdominal fat of HG cockerels; whereas, the TAS1R3 was up-regulated in abdominal fat of the LL chickens at 7 wk. Previously, we found by microarray analysis that the sweet taste receptor (TAS1R1) was expressed higher in the hypothalamus of FL chickens at 1 wk. [64], although in visceral fat the TAS1R1 was over-expressed in the LL chickens [21]. Furthermore, TAS1R3 knockout mice fed an adipogenic diet present with decreased adiposity and increased bone mass [65]. However, the mechanism by which the sweet taste receptors contribute to diet induced obesity in mice was not elucidated. Clearly, the role of sweet taste receptors (TAS1R1 and TAS1R3) in regulation of energy utilization and adiposity of the chicken and other avian species will require further investigation.
The Del-Mar 14 K Chicken Integrated Systems microarray was used for time-course gene expression profiling in abdominal fat of HG and LG cockerels (1-11 wk). In particular, this custom chicken microarray [30,31] was well suited for functional gene discovery in abdominal fat, since one-quarter (26%) of cDNAs printed on this 18 K-feature array were derived from abdominal fat cDNA clones. Similar to our transcriptional analyses of fat accretion in the divergent FL and LL cockerels [22], RNA-Seq analysis of abdominal fat was also completed on four HG and four LG cockerels at 7 wk., which was the age with the largest divergence in body weight and visceral fatness. Our time-course microarray analysis clearly shows engagement of several pairs of DE transcription factors in genetic divergence of growth and abdominal fatness traits in HG and LG cockerels. For example, PPARG has a larger number of direct gene targets that support increased synthesis and storage of fatty acids in visceral fat of HG cockerels, while PPARD exerts an opposing action by inhibiting lipid synthesis (SCD and FASN), while promoting lipolysis and fat catabolism in the LG chickens. RNA-Seq analysis revealed two additional opposing pairs of up-stream regulators [SREBF1 vs. VDR and SREBF2 vs. AR], where the AR and VDR seems to promote expression of catabolic genes and activation of the acute phase response in the LG cockerels. In contrast, other pairs of transcription factors [SREBF1 vs. SREBF2; and CEBPA vs. CEBPB] seem to be complementary or even synergistic in promoting lipogenesis and angiogenesis in abdominal fat of HG cockerels. Interestingly, the PGR appears to act similar to ESR1 in promoting growth and fat accretion in the HG cockerels, whereas PPARD, AR and VDR seem to support slower growth and diminished visceral fat in the LG cockerels.
The differential expression of numerous transcription factors and their direct target genes in abdominal fat of juvenile HG and LG cockerels could favor either the fast-growing fatter phenotype (HG) or the slower-growing leaner LG phenotype (see Additional files 6 and 7). The major DE upstream regulators that contribute to increased bodyweight and abdominal fatness phenotype of HG chickens include PPARG, SREBF1, SREBF2, CEBPA, E2F4, FOXO1, MSC, KLF9, KLF13 and THRSPA. On the other hand, the up-regulated transcription regulators that support reduced bodyweight and diminished visceral fatness of the LG phenotype were PPARD, CEBPB, CALR, PITX2, CREM, STAT5B, AR, PGR, FBXW7, KLF5, MYC and NKX2-1.
Divergent selection for slow growth in meat-type chickens (i.e., LG) substantially inhibits fat deposition, while enhancing catabolism of abdominal fat (see Fig. 1). The most remarkable changes found in expression of abdominal fat genes appear early in post-hatching development of the LG chickens (1-5 wk.; see Fig. 11). The long-chain fatty acid transporter (SLC27A1) was upregulated in abdominal fat of LG chickens throughout juvenile development. This gene is highly expressed in human adipose tissue and muscle and much lower in tissues with lower metabolic activity, and is barely detectable in liver [66]. The SLC27 proteins are responsible for the cellular up-take of long chain fatty acids for storage [67] or, more likely, for β-oxidation and enhanced energy production in LG chickens. Correspondingly, the rate-limiting step of β-oxidation is controlled by CPT1A, which is highly expressed in abdominal fat of LG chickens. The over-expression of CPT1A is driven by increased energy expenditure in chickens, where hepatic expression is up-regulated by fasting [68] or down-regulated after the embryo-hatch transition [32]. Furthermore, several single nucleotide polymorphisms were identified in CPT1A of Yup'ik Eskimos, which exhibit fasting-lipid and obesity phenotypes [69]. After transport of long chain fatty acids to the mitochondria, several steps of β-oxidation are catalyzed by the tri-functional protein hydroxyacyl-CoA dehydrogenase (HADH), which produces acetyl-CoA for entry into the TCA cycle. Both subunits of the tri-functional protein (HADHA and HADHB) were over expressed in LG chickens from 1 to 5 wk., which corresponds to increased expression of mitochondrial malic enzymes (ME2 and ME3) and several subunits of enzyme complexes in the electron transport chain. Although most long-chain fatty acids undergo mitochondrial β-oxidation, several substrates including very-long-chain fatty acids, polyunsaturated fatty acids and several others must be broken down by peroxisomal β-oxidation. The transcript for a major enzyme involved in peroxisomal β-oxidation, HSD17B4, was also upregulated in abdominal fat of LG chickens. Deficiency of this enzyme in humans causes a severe developmental syndrome due to cellular build-up of long chain fatty acids and leads to death soon after birth; whereas, knockout of the HSD17B4 gene in mice blocks the peroxisomal βoxidation of long chain fatty acids, without causing neonatal death [70].
Ultimately, the TCA cycle appears to be unable to handle the overload of acetyl-CoA produced by β-oxidation of fatty acids in abdominal fat of LG chickens, causing the over-expression of HMGCS2 and HMGCL (see Fig. 11), the two main enzymes in mitochondrial ketone body synthesis. The expression of HMGCS2 transcripts is greatly induced by fasting in the liver of suckling piglets [71] and 4-week old chickens [68]. In humans, both HMGCS2 and HMGCL exhibit alternative splicing which could regulate their tissue-specific expression [72]. Whether or not this complex regulation exists among different tissues in chickens has not been determined. Proprotein convertase subtilisin/kexin type 1 (PCSK1) is another candidate gene, whose expression was 2.2-fold higher in visceral fat of the LG cockerels (Additional file 5). Inhibitors of another members of this gene family, PCSK9, show great promise as a new therapeutic intervention to control hypercholesterolemia and promote cardiovascular health in humans [73]. Taken together, these findings support the idea that LG chickens increase cellular fatty acid up-take, peroxisomal β-oxidation, and transport to/into the mitochondria for β-oxidation and usage in the TCA cycle during early (1-5 wk) development, driving their divergence in abdominal fatness.
The LG chickens appear to down-regulate the mitochondrial breakdown of pyruvate through up-regulation of PDK4, which inhibits the pyruvate dehydrogenase (PDH) complex. PDH activity increases with increasing skeletal muscle stimulation (increasing contraction intensity) to a much greater extent in PDK4 knockout mice than in controls, which suggests that PDK4 is essential for preventing over-activation of the PDH complex [74]. Along with exercise, PDK4 is over-expressed in rodents subjected to high fat diet or fasting [75], which is similar to what is seen in chickens where both fasting and acute insulin immunoneutralization cause up-regulation when compared to the fed controls [76]. Furthermore, PDK4 appears to be crucial pre-hatch for the acquisition of stored energy, since its expression is quite high between embryonic day e16 and e20, then drops dramatically at hatch and remains low from day 1 through 9 post hatch [32].
Gene expression of type 1 iodothyronine deiodinase (DIO1), the major enzyme responsible for the conversion of the pro-hormone T 4 to metabolically-active T 3 , and DIO3, which degrades active T 3 and converts T 4 to inactive reverse T 3 (rT 3 ) [79] were examined by qRT-PCR analysis (see Fig. 11). Interestingly, DIO1 appears to be over-expressed in LG chickens during early juvenile development (3-7 wk), while expression of the degrading enzyme DIO3 is higher in HG chickens at later ages (9 and 11 wk). This could provide higher levels of T 3 in abdominal fat of LG chickens throughout juvenile development, which would contribute to their leanness [54]. The nuclear hormone receptor PPARD (see Fig. 2) has similar activation and inhibition targets and directly up-regulates several additional genes associated with leanness (i.e., ANGPTL3, PLP1 and PDK4). Selective knock-in of PPARD in adipose tissue of Lep db/db mice reduces lipid accumulation and prevents obesity in mice [80]. Early (1 to 5 wk) up-regulation of PPARD in abdominal fat of LG chickens and its direct action on target genes suggest that PPARD is a major regulator of lipolysis in these leaner and slower growing animals.
Presently, we show that HG cockerels over-express many genes controlling adipogenesis and lipogenesis in visceral fat (see Figs. 2 and 4; Table 4), which seems against the convention that abdominal fat only makes a minimal contribution to lipogenesis in the chicken [81]. This finding is supported by our earlier studies that show enhanced expression of genes controlling these processes in abdominal fat of fat line (FL) chickens [21]. Upregulation of adipogenesis and lipogenesis in abdominal fat of HG chickens is likely supported by increased activity and interaction of several transcription factors, including PPARG [82], CEPBA [83], SREBF1 [84] and THRSPA [85,86]. These transcriptional regulators are over-expressed throughout juvenile development and control common lipogenic targets, including DGAT2, FADS2, FASN and SCD (discussed below). Furthermore, IPA predicts the up-regulation of an additional 13 transcriptional regulators associated with fatness of HG chickens, which also directly increase the expression of lipogenic genes and other genes associated with fatness (see Additional files 6 and 7). Microarray analysis of liver in transgenic mice over-expressing SREBF1 or SREBF2 show that they have many common targets and both regulators play a crucial role in synthesis and metabolism of fatty acids and cholesterol [87]. Currently, our microarray and RNA-Seq analyses of abdominal fat in the HG and LG chickens show that SREBF1 and SREBF2 support increased growth rate and abdominal fatness of HG cockerels. These ligand-activated transcription factors enhance lipogenesis as indicated by their over-expression in HG chickens and direct redundant activation of numerous target genes controlling lipid metabolism (CYB5A,  DHCR7, FADS2, FASN, INSIG1, LSS, MSMO1, PCYT1A, SC5D, SCD and SQLE) (see Fig. 7-b).
INSIG1, also higher in HG chickens, may be a major regulator of fatty acid and cholesterol synthesis through the direct regulation of the SREBF chaperone (SCAP), and HMG-CoA reductase [88]. Two additional genes involved in cholesterol synthesis (DHCR24 and HSD17B7) are commonly up-regulated by divergent selection for high body weight or high abdominal fatness [21]. The expression of DHCR24 in whole blood was highly correlated and decreased with weight loss after bariatric surgery in humans [89]. Two other genes which are over-expressed in adipose tissue of HG chickens, ACACA and ACSS2, were identified as predictors of weight loss in another group of patients who underwent bariatric surgery, where decreased expression in adipose tissue was associated with a decrease in weight and hip circumference [90]. Acetyl CoA carboxylase (ACACA or ACC1) is critically important for the generation of malonyl CoA and synthesis of long chain fatty acids. This was demonstrated in human hepatoma HepG2 cells where inhibition of ACC1 by soraphen A reduces de novo lipogenesis by attenuating the formation of malonyl CoA and long chain fatty acids [91]. Overexpression of fatty acid elongases, ELOVL5 and ELOV6, or desaturases, FADS1 and FADS2 were not successful in reversing the effect of soraphen A on production of long chain fatty acids, which supports a crucial role of ACC1.
Perhaps one of the most significant findings of abdominal fat transcriptomes in the HG and LG chickens was a very large difference in expression of SCD, the rate-limiting enzyme responsible for conversion of saturated fatty acids into monounsaturated fatty acids. At 7 wk., SCD was the highest-expressed gene found in abdominal fat of HG chickens (see Fig. 10-b), where its expression was approximately 140-fold greater than that of the LG. This large difference was present throughout juvenile development (1-11 wk), where SCD transcripts were nearly undetectable in abdominal fat of LG chickens. Mutation of the SCD gene in mice produces a lean and hyper-metabolic phenotype, while ob/ob mice that also carrying mutations in SCD exhibit reduced fatness and elevated metabolism [92]. Another desaturase expressed very highly in abdominal fat of HG chicken is FADS2, which is down-regulated in liver in response to fasting in both pigs [93] and chickens [68], as well as in the pre-hatch chick embryo [32]. The major multi-enzyme protein responsible for the synthesis of fatty acids, FASN, is also very highly expressed in HG chickens at several ages, with the peak difference observed at 7 wk., the age of maximal distinction in abdominal fat between the two genotypes. The essential role of FASN in the fatty acid biosynthesis is highlighted in mouse embryos where FASN knockout mice die before birth and heterozygous knockout mice die at various stages of development [94]. An incredibly similar expression pattern (as FASN) was observed for DGAT2, which was almost 10-fold higher in HG chickens at 7 wk. DGAT2 catalyzes the final committed step in triacylglycerol synthesis and is critical for formation of adipose tissue [95,96]. Overexpression of DGAT2 in mammalian HEK293 cells significantly increases triglyceride (TG) synthesis while inhibition of DGAT2 by compound 122 decreases TG synthesis in a dose dependent manner [97]. Further, DGAT2 is up-regulated in goose hepatocytes exposed to a mixture of long chain fatty acids [98]. Correspondingly, ACSL1, which generates fatty acyl-CoA for entry into the triacylglycerol synthesis pathway, and AGPAT9, which catalyzes the conversion of glycerol-3phosphate to lysophosphatidic acid, the initial step of triacylglycerol synthesis [99] were both overexpressed in HG chickens. Abdominal fat of HG chickens also over expresses GPD1, which serves as a link between carbohydrate metabolism and lipid metabolism by catalyzing the conversion of dihydroxyacetone phosphate to glycerol-3-phosphate, the major component of glycerophospholipids. In fact, HG chickens up-regulate several genes involved in glycolysis including hexokinase 1 (HK1), which phosphorylates glucose to produce glucose-6-phosphate, the first step in the glycolytic pathway, and ALDOB, which converts fructose 1,6-bisphosphate into glyceraldehyde 3-phosphate and dihydroxyacetone phosphate (substrate for GPD1). Taken together, HG chickens up-regulate genes required for production of acetyl-CoA used in biosynthesis of fatty acids, cholesterol and triglycerides, which are highly upregulated processes in the abdominal fat of these animals.
Another important observation from the present transcriptional study was the predominate overexpression of avian leukosis virus (ALVE or envPr57 polyprotein) transcripts in the LG chickens across all ages (1-11 wk) and tissues (pituitary, abdominal fat, liver and breast muscle) examined [32]. The abundant overexpression of endogenous envPR57 transcripts in abdominal fat of LG cockerels during juvenile development (1-11 wk) was verified by qRT-PCR analysis (Fig. 11). It is of particular interest that ALVE transcripts are also highly-expressed in the hypothalamus [10], brain and peripheral tissues of the lowweight strain (LWS) chickens developed at Virginia Tech, particularly in immature pullets [6,11], where a greater number of ALVE integration sites are found in LWS pullets when compared to high-weight (HWS) females. These authors concluded that ALVE sequences might directly affect growth or could be linked to loci regulating growth. The relationship between the multiple ALVE loci and genetic selection of production traits was examined in several meat-type chickens [100,101]. Although our HG and LG lines were divergently selected from the Bresse-Pile breed in France [23] and the Virginia HWS and LWS lines were independently selected from the White Plymouth Rock breed. It is not known whether the higher abundance of ALVE in the low-growth (and low abdominal fatness) lines could reflect inheritance from a common ancestor or a direct response of selection for low growth. It is possible that the high frequency of ALVE insertion sites in the genome of meat-type chickens could disrupt key functional genes controlling growth or metabolism traits.
Lastly, we found evidence for a divergence in growth hormone (GH) signaling in the divergently selected HG and LG cockerels. Plasma levels of pituitary GH are elevated 2.5-fold in LG chickens, whereas plasma IGF-1 levels are 3-fold higher in the faster-growing and fatter HG birds (Duclos MJ, Simon J and Cogburn LA; unpublished observations). The LG cockerels also exhibit diminished GH binding to hepatic GH receptors (GHR), which seems to reflect the subsequent lack of negative feedback control on GH secretion from the pituitary gland. Thus, the reduced growth rate of LG chickens reflects a disruption in GHR signaling as evidenced by a reduction in hepatic synthesis of IGF-I and low plasma IGF-I levels [29]. Presently, we found a greater abundance of the small chicken GH (scGH) in abdominal fat of HG chickens at 7 wk. when compared to LG birds (see Fig. 11). The scGH was originally identified in the eye of chick embryos [102] and apparently exerts a local action in ocular tissue via an autocrine or paracrine action, since the signal sequence is absent [103][104][105]. Earlier, we identified a similar truncated cGH transcript in adipose tissue of the chicken from an exhaustive expressed-sequence tag (EST) sequencing project [106]. Insulin-like growth factor (IGF) signaling could be up-regulated in visceral fat of HG chickens, since IGFBP2 is over-expressed in abdominal fat of HG birds. Interestingly, IGFBP2 is primarily secreted from adipocytes and diminishes with childhood obesity [107]. Clearly, further studies on GH and IGF signaling in chicken adipose tissue are needed to fully understand how GH and IGF-1 contribute to lipogenesis and adiposity of chickens. Our earlier studies clearly show that exogenous GH increases abdominal fatness in juvenile chickens without affecting growth rate or plasma IGF-1 levels [108][109][110][111][112], which is unlike the typical mammalian response to exogenous GH [113][114][115].

Conclusions
Divergent selection for bodyweight at two ages (juvenile and adult) in meat-type chickens has profound effects on abdominal fatness and body weight traits (phenotypes). Throughout juvenile development, visceral fat of HG chickens over-express several transcription factors, which promote lipogenesis and adipogenesis. These regulators of transcription appear to be responsible for the up-regulation of several processes (i.e., biosynthesis of fatty acids, cholesterol and triglycerides) that increase abdominal fatness in HG cockerels. These lipogenic genes are not only differentially expressed, but are also among the most abundant transcripts found in abdominal fat of the HG. During early juvenile (1-5 wk) development, the slower growing and leaner LG chickens increase expression of several genes that promote energy expenditure, which likely contributes to their extreme leanness (i.e., oxidative phosphorylation, peroxisomal β-oxidation, mitochondrial β-oxidation and ketogenesis). Also, the expression of functional energy-consuming genes could be altered by the ALVE protein itself or random viral ALVE insertions, found throughout the chicken genome, could disrupt gene expression. Up-regulation of hemostatic factors (proteases and protease inhibitors) appears critical to development of extremes in leanness found in LG and LL chickens derived from different genetic backgrounds. This lends support to the idea that genetic selection for a divergence in growth or fatness traits also provides insight into mechanisms controlling fatness and leanness phenotypes. Metaanalysis of two RNA-Seq datasets from abdominal fat in divergently selected chickens from two distinct genetic backgrounds provides novel evidence that adipose genes encoding coagulation proteins are associated with the lean phenotype in meat-type chickens. While the exact mechanisms by which hemostatic genes contribute to leanness in chickens have not been clearly defined, the present transcriptional study of abdominal fat in HG and LG cockerels validates our previous reports that describe involvement of the hemostasis system in limiting adiposity of the leaner (LL) chickens [21,22]. Furthermore, the present study demonstrates that abdominal (visceral) fat of the chicken is a dynamic metabolic and endocrine tissue, which has an unappreciated capacity for in situ lipogenesis.