Transcriptional profiling of liver during the critical embryo-to-hatchling transition period in the chicken (Gallus gallus)
BMC Genomics volume 19, Article number: 695 (2018)
Although hatching is perhaps the most abrupt and profound metabolic challenge that a chicken must undergo; there have been no attempts to functionally map the metabolic pathways induced in liver during the embryo-to-hatchling transition. Furthermore, we know very little about the metabolic and regulatory factors that regulate lipid metabolism in late embryos or newly-hatched chicks. In the present study, we examined hepatic transcriptomes of 12 embryos and 12 hatchling chicks during the peri-hatch period—or the metabolic switch from chorioallantoic to pulmonary respiration.
Initial hierarchical clustering revealed two distinct, albeit opposing, patterns of hepatic gene expression. Cluster A genes are largely lipolytic and highly expressed in embryos. While, Cluster B genes are lipogenic/thermogenic and mainly controlled by the lipogenic transcription factor THRSPA. Using pairwise comparisons of embryo and hatchling ages, we found 1272 genes that were differentially expressed between embryos and hatchling chicks, including 24 transcription factors and 284 genes that regulate lipid metabolism. The three most differentially-expressed transcripts found in liver of embryos were MOGAT1, DIO3 and PDK4, whereas THRSPA, FASN and DIO2 were highest in hatchlings. An unusual finding was the “ectopic” and extremely high differentially expression of seven feather keratin transcripts in liver of 16 day embryos, which coincides with engorgement of liver with yolk lipids. Gene interaction networks show several transcription factors, transcriptional co-activators/co-inhibitors and their downstream genes that exert a ‘ying-yang’ action on lipid metabolism during the embryo-to-hatching transition. These upstream regulators include ligand-activated transcription factors, sirtuins and Kruppel-like factors.
Our genome-wide transcriptional analysis has greatly expanded the hepatic repertoire of regulatory and metabolic genes involved in the embryo-to-hatchling transition. New knowledge was gained on interactive transcriptional networks and metabolic pathways that enable the abrupt switch from ectothermy (embryo) to endothermy (hatchling) in the chicken. Several transcription factors and their coactivators/co-inhibitors appear to exert opposing actions on lipid metabolism, leading to the predominance of lipolysis in embryos and lipogenesis in hatchlings. Our analysis of hepatic transcriptomes has enabled discovery of opposing, interconnected and interdependent transcriptional regulators that provide precise ying-yang or homeorhetic regulation of lipid metabolism during the critical embryo-to-hatchling transition.
The developing mammalian embryo is completely dependent upon its mother and her placenta to supply nutrients, exchange respiratory gases, and to remove nitrogenous waste products. Quite the opposite is true among amniotes (i.e., birds and reptiles), which must undergo independent embryonic development in a pre-formed and completely enclosed system—the cleidoic egg. In chickens, the principal source of energy for late embryonic growth is derived from yolk lipids stored in the egg by the hen. The embryo exhibits exponential growth between days 13 to 18 of embryonic development (E13 to E18), when the total energy required for growth comes from β-oxidation of fatty acids derived from yolk lipids [1, 2]. The final phase of embryonic development is marked by a dramatic accumulation of cholesterol esters in the liver (> 30% of its dry mass) . Just prior to hatching, the yolk sac and its remaining lipids are retracted into the embryo’s visceral cavity to fuel the ensuing metabolic switch to endothermy in hatchling chicks. Yolk lipids are utilized by the newly hatched chicken from direct transfer into circulation via transport from the yolk stalk into the small intestine . Although fatty acids are readily absorbed from the lipid drenched intestine, the absorption of carbohydrates and amino acids from ingested high-energy feed is delayed until expression of intestinal enzymes and cotransporters (glucose/sodium/amino acid) reach adequate levels [4,5,6,7]. A recent transcriptional analysis, using serial analysis of gene expression (SAGE), of the chick embryo yolk sac between E13 and E21  provides the first view of genome-wide transcriptional changes that occur during transport of lipid from the yolk sac into the liver prior to hatching. The high abundance of yolk lipids and their transporters (i.e., lipoproteins) appear to interfere with absorption and utilization of nutrients from the intestine of the immediate post-hatch chick [4,5,6]. Thus, the newly hatched chick must undergo a dramatic metabolic shift from ectothermy, with an exclusive dependence on yolk lipids, to endothermy which depends on utilization of ingested feed—mainly carbohydrates and protein . Early hepatic expression of several lipogenic enzymes [9,10,11] provides the chick with the ability to convert dietary carbohydrate into fat stores . This sudden switch from a dependence on fat stored in yolk (embryos) to utilization of nutrients in feed (hatchling) requires major shifts in nutrient transport and metabolism that are controlled by yet uncharted genetic pathways. Although newly hatched chicks can survive for several days on residual yolk, delayed feeding prevents full recovery and normal growth, despite a brief compensatory growth spurt . Prolonged starvation during the immediate post-hatch period retards body and muscle growth, which is irreversible due to cellular changes in skeletal muscle . Although hatching is perhaps the most abrupt and profound metabolic challenge that a chicken must experience, there have been no attempts to functionally map the metabolic pathways induced by nutrients in the chick’s first meal. Furthermore, we know very little about transcriptional control of multiple metabolic and regulatory factors that control metabolism during the critical immediate post-hatch period.
The chicken embryo develops as an aquatic ectotherm with respiratory gas exchange via the chorioallantoic membrane. On the last day of embryonation (E21), the fully-developed embryo must suddenly convert to pulmonary respiration after pipping through the eggshell. During the peri-hatch period, the precocious chick becomes an endotherm with a very high and self-sustaining metabolic rate maintained by a fully-functional thermoregulatory system. The thyroid axis plays a major role in regulating metabolism and energy expenditure. Hence, the activity and importance of the thyroid axis has been extensively studied during the peri-hatch period of the chicken [14,15,16,17,18]. The thyroid axis prepares the embryo for this metabolic jump with decreased activity of the inactivating deiodinase (DIO3) at E20, whereas the activity of the activating deiodinases (DIO1 and DIO2) increases in E20 embryos during hatching and afterwards in the hatchling. The pro-thyroid hormone (T4) is readily deactivated to reverse rT3 by DIO3 in liver of E16 and E18 embryos; then with onset of pulmonary respiration, the activity of (DIO1) and (DIO2) increases in liver to meet the ensuing energy demands of hatching and endothermy. The embryo must switch from chorio-allantoic to pulmonary respiration and instantly meet the high metabolic demands of an endotherm—the hatchling chick. Therefore, it is essential to understand the transcriptional networks that integrate and control function of the endocrine system to meet the dynamic metabolic demands of growth imposed upon the liver.
Our first glimpse of the liver transcriptome during the chicken’s peri-hatch period came from an earlier low-resolution transcriptional scan using the original chicken liver (3.2 K) microarray . This initial transcriptional study using our prototype microarray was completed with 24 liver samples taken during the peri-hatch period [embryonic day 16 (E16) to post-hatch day 9 (D9)]. Self-organizing maps (SOMs) clustering of the 3.2 K microarray data identified two distinct gene expression profiles in liver. Cluster A contained differentially-expressed genes (DEGs) with very high expression in hatching liver, while Cluster B DEGs showed higher expression in embryo liver (see Fig. 2 in ). Hepatic expression of several lipolytic genes was higher in embryos and then abruptly decreased after hatching. In contrast, many up-regulated genes in liver of hatchling chicks control lipogenesis and energy metabolism (THRSPΑ, SCD, FASN, ME1, HMGCS, and EVOL6).
The present high-resolution genome-wide transcriptional analysis of the liver during the peri-hatch period has greatly expanded our hepatic repertoire of regulatory and metabolic genes involved in the embryo-to-hatchling transition. Several transcription factors were identified that appear to exert opposing actions, which lead to the predominance of lipolysis in embryos or lipogenesis in hatchlings. Furthermore, we have gained new knowledge on interactive transcriptional networks and metabolic pathways that enable the abrupt switch from ectothermy (embryo) to endothermy (hatchling) in the chicken.
Animal care and use
Fertile eggs from Ross x Ross broiler (meat-type) chickens were purchased from a commercial hatchery (Moyer’s, Quakertown, PA). The eggs were incubated in an automated incubator (Jamesway Incubator Co., Ontario CAN) at 39 °C and 95% relative humidity. At E20, eggs containing viable embryos were transferred to hatching trays to determine the approximate time of hatch (i.e., wetness of down). Chicks hatched within a 4 h period on the 21st day of incubation were wing-banded and held in the incubator for an additional 4 h. The remaining one day-old (D1) chicks were vaccinated against Marek’s disease virus and transferred to a heated battery brooder (Petersime Incubator Co., Gettysburg OH) maintained at 33 °C. The chicks received a commercial starter ration (22% crude protein and 3100 kcal ME/kg) and water ad libitum throughout the experiment. After continuous light for the first two days, the birds were maintained on a 20 h light: 4 h dark (20 L:4D) photoperiod. All protocols used in this study were in accordance with the United States Department of Agriculture guidelines on the use of agricultural animals in research and approved by the University of Delaware Animal Care and Use Committee.
Experimental design of hepatic transcriptional profiling study
The experimental design used for the present transcriptional profiling with Affymetrix Chicken Genome arrays was similar to the design utilized in our original 3.2 K microarray analysis of liver in peri-hatch chickens . Four male embryos were killed for collection of liver samples on embryonic day 16 (E16), E18 and E20 and four male hatchlings taken at 1, 3 and 9 day (D) of age. The D1 liver samples were taken from hatched chicks held in the hatching incubator (38 °C) for 4 h. After visual verification of male sex (presence of testes), each liver was quickly removed, placed into a 15 ml Nalgene cryogenic vial and snap frozen in liquid nitrogen. The remaining chicks were provide with feed and water ad libitum and housed in a heated battery brooder maintained at under a 20 L:4D photoperiod. The frozen liver samples were transferred to an ultra-low freezer (− 85 °C) and stored until isolation of total RNA.
Extraction and purification of total RNA
Total RNA was purified from each liver sample using an RNeasy Mini kit with on column DNase I treatment, as recommended by the manufacturer’s (Qiagen, Valencia, CA). Quantitation of total RNA was done using a NanoDrop spectrophotometer (NanoDrop, Wilmington DE) and the RNA was diluted to a standard concentration (250 ng/μL). RNA quality was assessed using an Agilent 2100 Bioanalyzer and an Agilent RNA 6000 Nano Kit (Agilent Genomics, Santa Clara, CA). Liver RNA samples (1 μL each) were analyzed on a RNA 6000 Nano Chip. The ribosomal RNA ratio (28S/18S) was assessed for RNA integrity, where with a RNA Integrity Number (RIN) of 9 to 10 indicates RNA of high quality for use in amplification of RNA for the microarray analysis and for verification by qRT-PCR analysis. Each liver RNA sample was divided into two identical aliquots and stored frozen at − 85 °C until use in the transcriptional analyses below.
RNA amplification, biotin labeling and hybridization of chicken genome arrays
Microarray-based transcriptome profiling was carried out using 24 Affymetrix GeneChip™ Chicken Genome Arrays and one-cycle target labeling protocol according to the manufacturer’s (Affymetrix) protocols. The experimental design for hybridization of the 24 chicken genome arrays was based on four replicate hybridization days, where all six ages were represented by one of four biological replicates (individual embryo or hatchling) per age group (Additional file 1). Briefly, 5 μg of total RNA was reverse-transcribed into first-strand cDNA using a T7-oligo(dT) promoter primers, followed by RNase H-mediated second-strand cDNA synthesis. and subsequently in vitro transcribed (IVT) into biotinylated complementary RNA (cRNA) in the present of T7 RNA Polymerase and biotinylated nucleotide analog/ribonucleotide mix using the Affymetrix IVT Labeling Kit (Affymetrix, Santa Clara, CA). The biotin-labeled cRNAs were cleaned up, fragmented and quality checked on an Agilent 2100 Bioanalyzer and hybridized to Chicken GeneChip™ Arrays for 16 h. The probe arrays were washed and stained on the Affymetrix GeneChip™ Fluidics Station and scanned for raw intensity data on the Agilent GeneChip™ Scanner using GeneChip™ Operating Software (GCOS) on the Affymetrix GeneChip™ Workstation. The raw cel files were then used for pre-processing and normalization as describe below.
Processing, normalization and statistical analysis of genome array data
An initial visual analysis (including hierarchical and K-mean clustering) of the microarray data was completed using GeneSpring GX software (Agilent, Santa Clara, CA). Raw cell files from the Affymetrix GeneChip Chicken Genome Arrays were imported into GeneSpring software and unfiltered differentially-expressed genes (DEGs; FDR ≤ 0.01) identified for a preliminary hierarchical and k-means cluster analysis (Fig. 1 only).
Next, open-source R (version 3.3.3)  software packages, available from Bioconducter , were used for data processing and statistical analysis of the single-channel microarray data (raw cell files) obtained in this study. The software package affylmGUI (version 1.48.0) was used to perform normalization and to determine differential expression . The robust multi-array average (RMA) method  was applied for background correction, quantile normalization and output of gene expression values on the log2 scale. Internally, affylmGUI uses the Linear Model Fit of LIMMA package (version 3.30.13) to assess differential expression of genes . Statistical significance was adjusted for multiple testing by controlling the Benjamini and Hochberg’s false discovery rate , which was set at FDR ≤ 0.05. The Affymetrix chicken annotation data package  was then used to annotate the differentially-expressed genes (DEGs). LIMMA was used to make pairwise comparisons (contrasts) across embryo and hatchling ages and to filter DEGs by the FDR adjusted p-value (P ≤ 0.05) and a log2 ratio ± 0.75 cutoff. This gene expression level threshold of log2 ratio ± 0.75 (or a 1.68-fold difference) has proven useful for statistical analysis of microarray data in model plants .
Ingenuity pathways analyses
An annotated list of DEGs for each pairwise contrast was then used as the input file for Ingenuity® Pathway Analysis (IPA), which provides functional annotation, mapping to canonical pathways and biological processes, and identification of gene interaction networks . The DEGs accepted by IPA are considered as “Analysis Ready” (AR)-DEGs, if the gene is curated and annotated in the Ingenuity® Knowledge Base. However, the Ingenuity Knowledge Base is mainly based on mammalian biomedical literature and rather void of avian-specific genes. We used a combination of Human Genome Organisation (HUGO) gene symbol as the primary gene ID and the RefSeq Protein ID as the secondary gene ID to maximize the number of chicken DEGs accepted by IPA. The AR-DEGs were mapped onto canonical pathways, biological processes and gene interaction networks. IPA uses the Fisher’s Exact Test to indicate significance (P ≤ 0.05) of over-representation of AR-DE genes in canonical pathways and biological processes. The Ingenuity® Upstream Regulator Analysis was used to identify major transcription factors, their interaction with other upstream regulators, and their interactions with direct target genes to predict either activation or inhibition of upstream regulators. Lists of annotated AR-DEGs were generated by IPA and provided in Additional files.
Verification of differential gene expression by quantitative RT-PCR analysis
A second aliquot of each total RNA sample from the 24 livers (above) was used for quantitative real-time PCR (qRT-PCR) analysis. The expression of 15 DEGs identified by microarray analysis was verified by qRT-PCR analysis (see Additional file 2). Selection of these candidate genes was based on their known function in pathways of interest or involvement in lipid metabolism. In additional, a panel of three invariant genes (COX7A2L, PCF11 and PRL14) was included in the qRT-PCR analysis for normalization of qRT-PCR expression levels. First-strand cDNA was synthesized using 1 μg of total RNA, oligo-dT, random hexamer primers and SuperScript III reverse transcriptase kit (Invitrogen Life Technologies, Carlsbad CA). Primers used for qPCR were designed using Primer Express (Applied Biosystems, Foster City, CA). The qRT-PCR assay was performed for each sample in duplicate wells using Power SYBR green PCR master mix (Applied Biosystems, Foster City, CA) and a gene-specific primer-pair (Sigma-Genosys, Woodlands, TX) on an ABI Prism Sequence Detection System 7900HT. Biogazelle qbase+ software  was used for pre-processing raw cycle threshold (Ct) data and the geNorm procedure for normalization of relative gene expression levels. The expression stability of candidate genes and a panel of invariant “housekeeping” genes was based on the geNorm (M) and coefficient of variation (CV) values across duplicate measurements of the 24 liver RNA samples. Ribosomal protein L14 (RPL14) and cytochrome c oxidase subunit VIIa polypeptide 2-like (COX7A2L) were selected as optimal co-reference genes with a mean of M = 0.251 and CV = 0.087. The PROC GLM in Statistical Analysis System (SAS, v.9.4; Cary NC) was used to analyze all pairwise differences among ontogenic stages and control for multiple comparisons using the option ADJUST = TUKEY. The Pearson’s correlation procedure in Excel was used to compare the average relative expression levels of each candidate DEG obtained from microarray analysis with that obtained by qRT-PCR analysis. Relative expression levels obtain by these two independent methods was verified by a significant (P ≤ 0.01) Pearson’s correlation coefficient (r).
Initial visualization of gene expression patterns using GeneSpring® analysis
The preliminary GeneSpring analysis of the 24 Affymetrix arrays identified 4566 differentially-expressed genes (DEGs; FDR ≤ 0.01) in the contrast of 12 embryos vs.12 hatchling chicks. The unsupervised hierarchical clustering of these unfiltered DEGs by GeneSpring provided a visual image (heat map) of the most dramatic changes in hepatic gene expression during the peri-hatch period (Fig. 1a). This heat map shows three distinct patterns of hepatic gene expression during the embryo-to-hatchling transition. “Cluster A” DEGs are expressed at higher levels in embryos (E16, E18 and E20) and at lower levels after hatching, especially in fully-fed D3 and D9 hatchings. In contrast, the DEGs in “Cluster B” have a low abundance in embryos and sharply elevated expression in hatchlings, while “Cluster C” DEGs exhibit a mixed pattern of gene expression. The second panel (Fig. 1b), derived from k-means clustering, provides a higher resolution of gene expression patterns in Cluster A (1166 DEGs) and Cluster B (801 DEGs), with examples of two predominant genes per cluster. Deiodinase 3 (DIO3) and monoacylglycerol O-acyltransferase 1 (MOGAT1) are highly expressed in liver of late embryos, whereas thyroid hormone-responsive Spot 14 protein alpha (THRSPA) and deiodinase 2 (DIO2) were the most abundant DEGs found in liver of hatchling chicks. K-means clustering provides identification of individual genes and a high-resolution view of opposing expression patterns of the deiodinases (DIO2, DIO3), which clearly shows the pivotal role that thyroid hormone metabolism plays in the metabolic jump from ectotherm (embryos) to endotherm (hatchlings) during the peri-hatch period.
Ingenuity pathway analysis (IPA): Identification of analysis ready DEGS
The DEGs within each contrast were considered as “Analysis Ready” (AR-DEGs) if accrued and curated in the Ingenuity® Knowledge Base. The numbers of up-regulated and down-regulated AR-DEGs are presented for the14 pairwise contrasts in a split-bar graph (Fig. 2). An annotated list of AR-DEGs in each pairwise contrast is provided in Additional file 3. A primary contrast was made between all embryos (E) vs. all hatchlings (H), which identified 1272 AR-DEGs. The greatest number of AR-DEGs (2440 hepatic genes) was found in the E18 vs. D3 contrast, where 846 AR-DEGS were higher in E18 embryos (lipid-laden-ectotherms) and twice that number (1594 AR-DEGs) was found in liver of fully-fed D3 hatchlings. When E16 embryos were compared against either E18 (231AR-DEGs) or E20 (814 AR-DEGs) embryos, there were fewer genes compared to E16 embryos versus D1 (1084 AR-DEGs) or D3 (1784 AR-DEGs) hatchlings. The lowest number of AR-DEGs was found at ends of the ontogenic spectrum (Fig. 2), namely E16 vs. E18 (231 AR-DEGs) and D3 vs. D9 hatchlings (403 AR-DEGs).
Venn diagrams provide unique and commonly-shared sets of hepatic AR-DEGs among meaningful contrasts of embryo and hatchling ages (Fig. 3). The Venn diagram in Fig. 3a compared three contrasts: E vs. H, E16 vs. E20 embryos and D1 vs. D9 hatchlings. A group of 604 AR-DEGs were unique to the E vs. H contrast, which shared 173 AR-DEGs with the E16 vs. E20 contrast and 313 AR-DEGs with the D1 vs. D9 contrast. The Venn diagram in Fig. 3b shows larger numbers of AR-DEGs in contrasts made between each embryonic age (E16, E18 and E20) vs. fully-fed D3 hatchlings. Also, this Venn diagram contains the largest number of commonly-shared AR-DEGs (1004 genes) among contrasts between E16, E18 and E29 embryos compared with the D3 hatchlings. In Fig. 3c, E18 embryos were compared against E20 embryos, D1 or D3 hatchlings. Comparisons of E18 vs. E20, D1 or D3 (Fig. 3c) revealed 372 common genes and a large group of 1096 AR-DEGs that are unique to the E18 vs. D3 contrast and an additional 914 AR-DEGs that are shared with the E18 vs. D1 contrast. The fourth Venn diagram (Fig. 3d) represents exclusive contrasts of hatchlings. The D1 vs. D3 contrast (643 AR-DEGs) was made between newly-hatched, fasted chicks (D1) and fed D3 hatchlings. The D1 vs. D9 hatchling contrast yielded 1036 AR-DEGs, 416 genes were shared with the D1 vs. D3 contrast, and 228 AR-DEGs shared with the D3 vs. D9 contrast of 403 AR-DEGs. And as expected, this Venn diagram of hatchling contrasts had the lowest number (42 genes) of commonly-shared AR-DEGs.
Ingenuity pathway analysis (IPA): Functional analysis of each pairwise contrast
Embryo (E) vs. hatchling (H) contrast
A description of the IPA of each pairwise contrast follows the same order in which the pair-wise contrasts are presented in Fig. 2 and Additional file 3. A brief summary of the IPA “Core Analysis” of 1272 AR-DEGs identified in the primary contrast of all embryos (E) vs. all hatchlings (H) is provided in Table 1. The top five canonical pathways over-represented by AR-DEGs from the E vs. H contrast were related to cholesterol biosynthesis, triacylglycerol biosynthesis, or activation of the ligand-activated nuclear receptor, retinoid X receptor (RXR) and its heterodimer partners [liver activated receptor (LXR) and farnesoid activated receptor (FXR)]. The top upstream regulator found in liver of embryos was peroxisome proliferator-activated receptor-alpha (PPARA), which is considered as the master regulator of lipid catabolism, namely lipolysis of yolk lipids. It should be noted that IPA estimates the top up-stream regulators according to the number of target AR-DEGs in the data set rather than whether or not the transcription factor (TF) itself was differentially expressed. PPARA was an AR-DEG that had 94 direct targets among the 1272 AR-DEGs. Sterol-response element binding factors 1 and 2 [SREBF1 (51 AR-DEG targets) and SREBF2 (28 AR-DEG targets) are ligand-activated transcription factors, which regulate sterol biosynthesis. Although not an AR-DEG itself, the tumor protein p53 (TP53) had the highest number of direct gene targets observed in the E vs. H contrast. Similarly, fork-head box O1 (FOXO1), a member of the fork-head family of transcription factors was not an AR-DEG in the E vs. H contrast, although IPA considered FOXO1 a top up-stream regulator due to its large number of direct target genes in this contrast. Actually, 23 transcription factors (TFs) were AR-DEGs in the E vs. H contrast, where 11 TFs were expressed higher in liver of embryos and 12 TFs were expressed higher in hatchlings (see Additional file 4). Of particular interest, a recently discovered gene SERTA domain containing 2 (SERTAD2 or TRIP-Br2) was up-regulated in liver of embryos. SERTAD2 is a transcription co-regulator of lipolysis, thermogenesis and oxidative metabolism in mammals.
The top “Molecular and Cellular Functions”, found in the E vs. H contrast (Table 1), were “Small Molecule Biochemistry” (358 AR-DEGs), “Lipid Metabolism” (284 AR-DEGs), and “Molecular Transport” (231 AR-DEGs). Fifty AR-DEGs were functionally annotated by IPA as involved in “Amino Acid Metabolism”, where 28 genes were overexpressed in liver of embryos compared to 22 genes in this category expressed higher in hatchlings (Additional file 4). “Oxidation of Lipid” was another IPA “Biofunction” category that was over-represented by 40 metabolic AR-DEGs in the E (19 genes) vs. H (21 genes) contrast. “LPS/IL1 Mediated Inhibition of RXR Function” was another canonical pathway occupied by AR-DEGs from the embryo vs. hatchling contrast; 13 genes were more abundant in embryos and 19 genes were higher in liver of hatchlings. Hepatic genes involved in “Cholesterol Biosynthesis” are more numerous in liver of hatchlings (16/17 AR-DEGs), where only HADHB was expressed greater in embryos. “Triacylglycerol Synthesis” was also over- represented with 14 AR-DEGs, 6 genes (mainly acyltransferases) were higher in liver of embryos and 8 genes (3 acyltransferases and 3 fatty acid elongases) were expressed a greater levels in hatchlings. Four additional over-represent canonical pathways populated by AR-DEGs from the E vs. hatchling contrast included farnesoid X-activated receptor (“FXR)-RXR Activation” (25 AR-DEGs), liver activated receptor (“LXR-RXR Activation” (24 AR-DEGs), thyroid hormone receptor (“TR-RXR Activation” (16 genes), and the “Coagulation System” (5 genes expressed higher in embryos and one gene higher in hatchlings).
The importance of hepatic lipid metabolism during the embryo-to-hatchling transition becomes clear by examining the subcellular distribution of 149 AR-DEGS that control “Lipid Synthesis” (Fig. 4.). This graphic overview shows higher expression of 54 AR-DEGs (red symbols) in the liver of embryos, while twice that number (95 AR-DEGs) were expressed higher in hatchlings (green symbols). In the nucleus, seven TFs were expressed higher in liver of embryos PROX1, SIRT1, PPARGC1A, TFC2L1, RGN, PPARA and NR1H4). Likewise, seven TFs were expressed higher in liver of hatchlings (PPARG, SREBF2, BRCA1, KLF11, THRSPA, LPIN1 and CASP8). Thus, Ingenuity® Up-stream Regulator Analysis predicts that lipid synthesis was inhibited in liver of hatchlings when compared to the embryos (red symbols) as indicated by blue-dashed lines. However, the over-abundance of genes with green symbols (62 AR-DEGs) found in the cytoplasm are metabolic enzymes, transporters, kinases and phosphatases, which are actually expressed at higher levels in liver of hatchlings than in embryos (only 22 up-regulated genes). Among the most abundant cytosolic genes found in the embryo were several lipolytic AR-DEGs: MOGAT1, PDK4), FDXR, LPL, HMOX2, PTGDS and DAB1. Notable and highly-expressed lipogenic DEGs were found in the liver of hatchlings, including SCD, FASN, ME1, AGPAT2, LSS, ACACA, ELOVL5 and ELOVL6. In the plasma membrane, 17 DEGs were expressed at high levels in hatchlings (i.e., FADS1, FADS2, DDP4, PLIN2, ENPP2 and LDLR), while the GHR, TNFRSF1B, ENDRA and EGFR were upregulated in embryos. Extracellular space was occupied by 12 gene products expressed higher in the embryo (i.e., APOB, ADIPOQ, IGFBP2, TNFSF10, IPC, F2, HGF, LCAT, CETP, CXCL12 and BMP5) when compared to only 5 hepatic gene products that were higher in extracellular space of the hatchlings (APOA4, EREG, NPC2, ANG and FGF23). This subcellular distribution of 149 AR-DEGs from the E vs. H contrast provides a comprehensive view of how 14 differentially-expressed TFs control lipid synthesis, while revealing major lipolytic genes in embryos and an even larger number of lipogenic genes in liver of hatchlings.
Another interaction network the E vs. H contrast was functionally annotated by IPA as “Lipid Metabolism, Molecular Transport” (Fig. 5). This network shows interactions among three lipogenic TFs (THRSPA, PPARG and SREBF2) and their direct target genes. The ligand-activated nuclear hormone receptor PPARG has a direct action on five target genes that were highly expressed in embryos (IDH1, FABP5, MANBA, PDK4 and CPT1A). In addition, PPARG interacts with several genes that are more abundant in liver of hatchling chicks (SCD, THRSPA, EVOLV6, HMGCR, FABP7, GARNl3, ADIPOR2, VNN1, LPIN1 and FDPS). LPIN1 interacts with LPIN2, PLPP3 and FDPS, which has a direct action on FNTB and GGPS1. Six genes in this network have direct interactions with both PPARG and SREBF2 (FDPS, SCD, THRSPA, ELOVL6, HMGCR and IDH1). Numerous additional direct targets of SREBF2 from the E vs. H contrast were identified by the Ingenuity® Upstream Regulator Analysis (Fig. 5b), which predicts that the majority of these 28 target genes would be inhibited (or down-regulated, as indicated by blue arrows) by SREBF2. Ingenuity only predicts that the LDL receptor related protein 1 (LRP1), which was up-regulated in embryos, would be inhibited by SREBP2 as indicated by the blunt orange line. However, two additional genes that were also highly abundant in the embryo, IDH1 and FABP5, appeared to be inconsistent with IPA’s anticipated state. Interestingly, SREBF2 and its direct target genes are involved in cholesterol metabolism and lipogenesis in particular, which are both enhanced in liver of hatchling chicks. As direct targets of SREBF2, only LRP1, IDH1 and FABP5 were up-regulated in liver of embryos. Furthermore, 26 of the most abundant AR-DEGs found in liver of hatchling chicks (Additional file 4) are involved in either cholesterol biosynthesis (i.e., MSMO1, MVD, CYP8B1 and DHCR7) or lipogenesis (i.e., FADS2, FASN, HMGCS1, LPIN2 and ACACA).
The gene interaction network in Fig. 6a, found in the E vs. H contrast, was functionally annotated by IPA as “Hematological System Development and Function”. This network was centered on interaction of several blood clotting factors (TFPI, F2, F11, PROC and PROZ) with two lipolytic TFs (PPARGC1A and SIRT1) and with components of the glutathione transferase pathway (GSTA1, GSTA2, GSTO1, GSTZ1, MGST3, GPX3 and GPX8). Ingenuity® Up-stream Regulator Analysis also identified 24 target genes of SIRT1 (Fig. 6b), which would either activate several genes in embryos (CPT1A, HGF, NR1H4, PDK4, PPARGC1A, TFP1, ADIPOQ and CDH2) or inhibit other genes in (CDKN2B, FASN, H2AFZ and HMGCR). These genes are highly expressed in liver of embryos and directed at lipolysis and gluconeogenesis, while the AR-DEGs with green symbols are mainly lipogenic genes.
E16 vs. E18 contrast
The E16 vs E18 contrast yielded the smallest number (231) of AR-DEGs (Fig. 2; Additional file 3), which prevented generation of highly-populated gene interaction networks by IPA. Interestingly, the most abundant DEGs expressed in the liver of E16 embryos when compared to E18 embryos were seven transcripts that correspond to feather keratin (FKER). Unfortunately, these avian-specific transcripts were rejected by Ingenuity software and not included in the IPA functional analysis of the E16 vs. E18 contrast.
E16 vs. E20 contrast
The contrast of E16 vs. E20 embryos identified 814 AR-DEGs, which included 514 genes up-regulated in E16 liver (Additional file 3). A summary of the IPA of the liver transcriptome of E16 vs. E20 embryos is presented in Table 2. The top five canonical pathways over-populated by 814 AR-DEGs in this contrast were “LPS/IL-1 Mediated Inhibition of RXR Function, Mitotic Roles of Polo-Like Kinase, Complement System, Tyrosine Degradation I, and Cell Cycle Control of Chromosomal Replication”. The top upstream regulators identified by IPA from this contrast were E2F4, TCF3, PPARA, E2F1 and TCF4. The top “Molecular and Cellular Functions” recognized in the E16 vs. E20 contrast were “Cell Cycle”, “Cellular Assembly and Organization”, “DNA Replication, Recombination, and Repair”, “Cell Death and Survival”, and “Cell Morphology”. The top “Physiological Function” category included the “Reproductive” and “Digestive” systems, “Organ Morphology”, “Organismal Development” and “Behavior”. Under “Toxicity Functions”, IPA assigned several AR-DEGs to “Inhibition of RXR Function”, fatty acid metabolism, renal failure, and transmembrane potentials. Interestingly, the hepatic deiodinases were the top 10 opposing DEGs discovered in the liver of E16 (DIO3; + 3.6 log2 ratio) vs. E20 (DIO2; − 7.2 log2 ratio) embryos. Another gene highly expressed in liver of E16 embryos was angiopoietin-like protein 3 (ANGPTL3), a regulator of lipoprotein lipase. Serpin family A member 3 (SERPINA3), HPX, FABP1 and KLF9 were DEGs that were up-regulated in liver of E20 embryos. Of particular interest, the E18 vs. E20 contrast revealed higher expression of three Kruppel-like factors (KLF9, log2–3.52; KLF11, log2–0.85; KLF13, log2–2.31) in liver of E20 embryos (Additional file 3).
E16 vs. D1 contrast
The contrast of E16 embryos vs. D1 hatchlings recognized 1084 AR-DEGs, where 593 AR-DEGs were more abundant in liver of the E16 embryos (Fig. 2; Additional file 3). The highest DEGs found in embryo liver from the E16 vs. D1 contrast included DIO3 and MOGAT1, while DIO2 and FABP1 were among the highest AR-DEGs found in fasted D1 hatchlings.
E16 vs. D3 contrast
The E16 vs. D3 contrast identified 1784 AR-DEGs with 852 transcripts up-regulated in liver of E16 embryos compared with 932 AR-DEGs in D3 hatchlings. DIO3 and MOGAT1 were also among the highest DEGs in E16 embryos from the E16 vs. D3 contrast whereas DIO2 and THRSPA were overexpressed DEGs found in liver of D3 hatchlings.
E18 vs. E20 contrast
A total of 624 AR-DEGs were identified by the E18 vs. E20 contrast, where 258 AR-DEGs were expressed at higher levels in E18 embryos compared to 366 AR-DEGs overexpressed in the E20 embryos (Table 3). The major canonical pathways highly populated by AR-DEGs from this contrast were “LPS/IL-1 Mediated Inhibition of RXR Function, Mitotic Roles of Polo-Like Kinase, Protein Ubiquitination Pathway, Unfolded Protein Response”, and “Tyrosine Degradation I”. The top five upstream regulators found in the E18 vs. E20 contrast included PPARA, TCF3, ATF6, FOXO3 and RARA. “Small Molecule Biochemistry” and “Cell Cycle” were the highly represented molecular and cellular functions according to IPA. Likewise, “Organ Morphology” and “Organismal Development” were highly populated physiological functions identified by IPA. Among the top toxicology categories were “Fatty Acid Metabolism”, inhibition/activation of RXR and “Cell Cycle”. Similar to the E16 vs. E18 contrast, the thyroid hormone deiodinases were the most abundant DEGs in E18 (DIO3; log2 + 3.8) vs. E20 (DIO2; log2–5.6) contrast. Several additional genes, highly expressed in the E20 liver, are involved in lipid metabolism (SERPINA3, HMGCS1 and KLF9) or fatty acid transport (FABP1). Lists of functionally-annotated hepatic AR-DEGs belonging to these canonical pathways and biological functions are provided in worksheets of Additional file 5. “Lipid Synthesis” was a subcategory under “Lipid Metabolism” that was over-populated by 60 AR-DEGs; 18 DEGs were up-regulated in E18 embryos, whereas 42 hepatic genes were expressed at higher levels in liver of E20 embryos. Under “Amino Acid Metabolism”, 21 AR-DEGs were more abundant in E18 embryos, while 12 genes were expressed higher in E20 embryos. Under the “Insulin Resistance” category, IPA identified 12 AR-DEGs that were more abundant in E18 embryo liver, while 16 genes were more abundant in E20 than in E18 embryos. Nineteen AR-DEGs were assigned by IPA to the “Sirtuin Signaling” pathway, where 14 AR-DEGs were higher in E20 embryos. The “FXR-RXR Activation” pathway was occupied by 11 AR-DEGs, where 7 genes were more abundant in E20 embryos.
This gene interaction network (Fig. 7) from the E18 vs. E20 contrast was functionally annotated by IPA as “Lipid Metabolism”. In liver of late chicken embryos, the interactions among four nuclear hormone receptors (NROB2, ESRRG, PPARA and PPARGC1A) control energy balance and lipid metabolism. Ingenuity predicts PPARGC1A itself is inhibited (blue symbol) and blue arrows show inhibition of target genes in the E18 vs. E20 embryos (Fig. 7b). Fifteen direct targets of PPARGC1A (green symbols) are highly expressed in liver of E20 embryos, whereas red-colored gene symbols indicate up-regulated expression in liver of E18 embryos. These genes (i.e., INSIG1, NR1H4 and PPARA) are involved in fat catabolism, while genes with green-colored symbols indicate higher expression of lipogenic genes in E20 embryos. In addition, Ingenuity predicts activation of the orphan receptor NROB2 (nuclear receptor subfamily 0 group B member 2), which leads to inhibition (blunt blue lines) of seven lipogenic genes (PPARGC1A, ANOX1, CPT1A, CTP8B1, EGR1, ESRRG and HMGCR). Another direct target of NROB2 (an inhibitor of estrogen, thyroid and retinoid nuclear hormone receptors) was the ligand-activated transcription factor NR1H4, which is activated by bile acids and thereby controls expression of genes involved in bile acid synthesis and transport. NR1H4 is also a direct target of PPARGC1A.
E18 vs. D1 contrast
The functional analysis of 1591 AR-DEGs found in the E18 vs. D1 contrast is presented in Table 4. The major canonical pathways identified by IPA in this contrast were “NRF2-mediated Oxidative Stress Response” and “Hereditary Breast Cancer Signaling”. The largest number of direct target genes was associated with two TFs, TP53 and ESR1. The top “Molecular and Cellular Functions” recognized by IPA of the E18 vs. D1 contrast were “Cell Cycle”, “Cellular Assembly and Organization”, “DNA Replication, Recombination, and Repair”, “Cell Death and Survival”, and “Lipid Metabolism”. Under the category “Physiological System Development and Function”, the largest number of DEGs was assigned to “Embryonic Development”. The top “Toxicity Functions” were related to the oxidative stress response, liver necrosis/proliferation and RXR activation. The highest DEGs in liver of E18 embryos were DIO3 and MOGAT1, while in D1 hatchling chicks the most abundant transcripts were involved in lipid metabolism (DIO2, ELOVL5, FABP1 and HMGCS1). Furthermore, four members of the Kruppel-like transcription factors (KLF9, KLF11, KLF13 and KLF15) were abundant DEGs found in liver of fasted D1 hatchlings (Additional file 3).
E18 vs. D3 contrast
The largest number of AR-DEGs was found in the E18 vs. D3 contrast (2440 AR-DEGs), which overpopulated the cholesterol and triacylglycerol biosynthetic pathways, and the “Intrinsic Prothrombin Activation Pathway” where all 12 AR-DEGs were up-regulated in liver of E18 embryos (Table 5). The top “Up-stream Regulators” in the E18 vs. D3 contrast included TP53, PPARA, and FOXO3. Among the top “Molecular and Cellular Functions” identified by IPA were the “Cell Cycle”, “DNA Replication, Recombination, and Repair”, “Small Molecule Biochemistry” (535 AR-DEGs; including 215 AR-DEGs that belong to the “Lipid Synthesis” subcategory; see Additional file 6), and “Cellular Assembly and Organization”. “Organismal Survival” was the most populated physiological function found with this contrast. The top “Toxicity Functions” assigned to AR-DEGs from the E18 vs. D3 include “Cholesterol Biosynthesis”, NRF2-Mediated Oxidative Stress Response, “LXR/RXR Activation”, and “TR/RXR Activation”. The top 10 “Up-regulated Genes” in E18 embryos were DIO3, MOGAT1 and PDK4, whereas lipogenic genes (THRSPA, SQLE, ELOVL5, DIO2, FASN and SCD) were highly expressed in liver of fully-fed D3 hatchlings. The E18 vs. D3 contrast revealed seven members of the Kruppel-like transcription factor family as AR-DEGs (Additional files 3 and 6).
Two gene interaction networks depicted in Fig. 8 were derived from IPA of a subset of 1004 commonly-shared genes identified by the Venn diagram shown in Fig. 3b, which compared contrasts of E16, E18 or E20 embryos vs. D3 hatchlings. This subset of commonly-shared AR-DEGs, enriched with multiple TFs and their direct target genes, was functionally annotated by IPA as ‘Lipid Metabolism, Molecular Transport, Small Molecule Biochemistry’ (see Additional file 3; last worksheet). The top panel (Fig. 8a) shows the interaction of three TFs (PPARA, PPARG and PPARD) with several common target genes (SIRT5, GPAM, PDK4, ALDH9A, ACSL, SLC27A4, ELOVL6, ACOT8, CPT1A and PER3). The ligand-activated TF PPARA and its five direct targets (CPT1A, SIRT5, CFI and LECT) were expressed at higher levels in liver of E18 embryos. In contrast, PPARG and PPARD and 11 target genes are expressed higher in liver of D3 hatchlings, to support lipogenesis and energy production of these fully-fed hatchlings. PPARG has four unique direct target genes (ADIPOR2, MANBA, GARNL3 and the collagen alpha 1 complex), whereas PPARD has only four unique targets (MFSD2A, LPIN2 and ST3GAL5) in this gene network. The second gene network (Fig. 8b) shows the interaction of two major lipogenic TFs (THRSPA and SREBF2; both expressed higher in D3 hatchlings) with the transcriptional coactivator PPARGC1A, which was expressed higher in liver of E18 embryos. Only five genes, directly or indirectly responsive to PPARGC1A were expressed at higher in liver of E18 embryos, whereas SREBF2 and its 11 target genes, including THRSPA, were more abundant in liver of D3 hatchlings. Most of the target genes of SREBF2 are involved in metabolism of cholesterol, sterol and ketones.
The contrast of E18 embryos vs. D3 hatchlings revealed another gene network of six transcription factors (RORA, SMARCD3, ESRRG, NR1H4, PPARGC1A and SREBF2) and their interactions with direct gene targets that control “Lipid Metabolism” (Fig. 9a). Only SREBF2 was expressed higher in D3 hatchlings, where most of its target genes are involved in lipid metabolism. The RAR-related orphan receptor A (RORA; Fig. 9b) was predicted by IPA to be inhibited (blue symbol) and three target genes with blue arrows would be inhibited (down-regulated; green symbols) when compared to E18 embryos, which possessed only five up-regulated target genes (IGFBP1, LPIN2, NTRK2, PPARGC1A and AKR1D1). Ingenuity predicts that NR1H4 would be slightly activated mainly due to presence of seven DEGs (CETP, FABP5, LCAT, LIPC, NOS2, PPARA and PPARGC1A) that were highly expressed in liver of the E18 embryos. Ingenuity also predicts that upregulated NR1H4 would inhibit lipid synthesis (i.e., ACACA, CYP7A1 and FASN).
Another gene interaction network from the E18 vs, D3 contrast (Fig. 10a) was functionally annotated by IPA as involved in “Cell Signaling” and under control of eight interacting transcription factors (THRSP, NFYC, KAT2B, KLF13, EAF1, EAF2, SOX7 and SOX8). The lipogenic transcription factor THRSPA has direct interactions with two other transcription regulators (NFYC and KAT2B). In all, 12 genes in this network were expressed higher in E18 embryos, whereas 18 AR-DEGs were over-expressed in liver of D3 hatchlings. The second gene network (Fig. 10b) involves interactions among six clotting factors (F5, F8, F9, F11, PROS1 and VWF), the growth factor gene MST1, collagen genes (COL1A1, COL3A1and COL18A1) and components of the extracellular matrix (FBLN1–2, LAMB1, LAMB3, LOX, SPARC and SOX17). This connection of the extracellular matrix with F9 and other coagulation factors is mediated by the direct interaction between two enzymes (HELLS and HSPD1). Interestingly, most of the AR-DEGs in this network were also expressed higher in liver of E18 embryos than D3 hatchlings. Furthermore, 12 coagulation factors (AR-DEGs) belonging to the canonical “Intrinsic Prothrombin Activation Pathway” were up-regulated in liver of E18 embryos (see Additional file 6).
E20 vs. D3 contrast
The E20 vs. D3 contrast identified 1775 AR-DEGs, whose functional analysis is summarized in Table 6. The top canonical pathways occupied by these AR-DEGs were related to “Cholesterol Biosynthesis”, the “Complement System”, “LXR/RXR Activation” and “Intrinsic Prothrombin Activation Pathway”. The top transcription regulators in the E20 vs. D3 contrast were TP53, PPARA, E2F4, E2F1 and FOXO3. Within the “Molecular and Cellular Functions” category, IPA identified “Cell Cycle’, “Cellular Assembly and Organization”, DNA Replication, Recombination, and Repair”, “Lipid Metabolism”, and “Small Molecule Biochemistry” as over-represented sub-categories. Over-populated physiological functions subcategories included “Organismal Survival”, “Connective Tissue Development” and “Tissue Morphology”. The top “Toxicity Functions” identified by IPA included “LXR/RXR Activation, Cholesterol Biosynthesis, Fatty Acid Metabolism” and “FXR/RXR Activation”. Among the top 10 “Up-regulated Genes” in the E20 vs. D3 contrast were MOGAT1, PDK4, CPT1A (a major regulator of fatty acid β-oxidation) and IGFBP2. The highest DEGs found in the liver of D3 hatchlings are involved in lipogenesis [THRSPA, SQLE, FASN, SCD, ELOVL5 and FADS2).
D1 vs D9 contrast
Finally, the contrast of D1 vs. D9 hatchlings yielded 1036 AR-DEGs (Table 7; see Additional file 7). The top canonical pathways identified by IPA were related to cholesterol and stearate biosynthesis, the “Complement System” and “LXR/RXR Activation”. The top upstream regulators and their direct target genes found by IPA in the D1 vs. D9 contrast were PPARA, TP53, E2F4, sterol response element binding factor-2 (SREBF2) and TCF3. Under the “Molecular and Cellular Functions” category, the largest numbers of AR-DEGs were assigned to the “Cell Cycle”, “Lipid Metabolism” and “Small Molecule Biochemistry” subcategories. Under the “Physiological Systems” category, the largest number of AR-DEGs in the D1 vs. D9 contrast were found under “Connective Tissue Development and Function” and Tissue Morphology. The top five “Toxicity Functions” were “Fatty Acid Metabolism”, “Cholesterol Biosynthesis”, “LXR/RXR Activation”, “Cell Cycle: G2/M DNA Damage Checkpoint Regulation” and “Liver Necrosis/Cell Death”. The “Lipid Synthesis” subcategory contained 134 AR-DEGs, 53 were up-regulated in fasted D1 hatchlings and 81 genes were expressed higher in liver of fed D9 hatchlings (see Additional file 7). Three over-represented canonical pathways were related to AR-DEGs controlling “RXR Function”; these were “LPS/IL1 Mediated Inhibition of RXR Function”, “LXR-RXR Activation”, “FXR-RXR Activation”, and “TR-RXR Activation”. Among the highest AR-DEGs found in the liver of fasted D1 hatchlings were IGFBP2, PDK4, PTGDS, PPARGC1A, CPT1A, LPL and DIO3 (Table 7). In contrast, the fully-fed D9 hatchling chicks showed higher expression of CDKN2B and several lipogenic genes including SQLE, THRSPA, FADS2, SCD, LSS and PNPLA3.
Two gene interaction networks from the D1 vs. D9 contrast were functionally annotated by IPA as involved in “Lipid Metabolism” (Fig. 11). As shown in Fig. 11a, several lipogenic transcription factors [SREBF2, THRSPA, MID1IP1 (a THRSP paralog), PLAGL, JUN and CREM] control expression of multiple metabolic enzymes (ACACA, HMGCR, FDPS, SQLE, AACS, LSS, CYP8B1, IDI1, MVD, COMT, FBXL12 and ACAT2), kinases (MASTL, SIK1), a phosphatase (LPIN1), and transporters (STARD4 and SLC6A6). The gene network in Fig. 11b shows interactions of five transcription factors (PPARG, PPARGC1A, NR0B1, NCOR1and ZBTB20) with key lipogenic metabolic enzymes (FASN, ME1, CPT1A), growth factors (IGF1), binding proteins (IGFBP1, IGFBP2, IGFALS) and the adiponectin receptor (ADIPOR2). Other direct targets of PPARG shown in this gene network were NDRG1, PLIN2, PER3 and VNN1.
The final gene interaction network (Fig. 12a) was also identified in the D1 vs. D9 contrast and functionally annotated by IPA as belonging to “Energy Production” and “Lipid Metabolism”. This network is centered on interactions among four lipogenic TFs (THRSPA, PPARG, PPARD and KLF15), where PPARG and PPARD share seven target genes (PER3, DUT, CPT1A, THRSPA, LPIN1, ACOX2 and PDK4). The TF PPARG and several of its direct targets are expressed at higher levels in liver of D9 hatchlings, including ADIPOR2, FBXO9, ELOVL6, ELOVL2, ELOVL1, HSD17B12, PEX11A, HECTD2, VNN1 and FDPS. The fasting-induced transcriptional factor KLF15 and four direct gene targets (HADHA, HADHB, EHHADH and PDK4) are expressed higher in liver of fasted D1 hatchlings. The Ingenuity Upstream Regulator Analysis (Fig. 12b) identified 30 AR-DEGs in the D1 vs. D9 contrast as direct targets of PPARG, where only 9 genes were up-regulated in the fasted D1 hatchlings. Consequentially, IPA predicted that PPARG would be inhibited (blue gene symbol and blue arrows), since 21 of its target genes are down-regulated (green symbols) AR-DEGs in this contrast. Based on known relationships in the Ingenuity Knowledge Base, Ingenuity predicts that inhibited PPARG would block the expected activation (blunt orange lines) of CDK1, CPT1A and LPL. Ingenuity predicts that KLF15 would be activated, which would lead to activation (up-regulation) of seven direct targets (ACSL1, EHHADH, FABP5, HADHA, HADHB and PDK4) as indicated by the orange arrows. Interestingly, KLF15, a key transcriptional regulator of gluconeogenesis, was highly expressed in liver of fasted D1 hatchling chicks, as revealed by three pairwise contrasts (E18 vs. D1; D1 vs. D3; and D1 vs. D9).
qRT-PCR analysis and verification of differential expression of 15 candidate DEGs
Hepatic expression patterns of seven DEGs that were expressed higher in embryos and an invariant gene (COX7A2L) during the embryo-to-hatchling transition were examined by qRT-PCR analysis (Fig. 13). Six genes (MOGAT1, PDK4, FZD2, LDHB, DIO3 and ADIPOQ) show higher log-based expression in embryos and a sharp decline after hatching. Surprisingly, hepatic expression of DIO1 was relatively stable, except in E18 embryos where DIO1 was slightly elevated. The invariant gene (COX7A2L) was used by geNorm software for normalization of the qRT-PCR expression levels. The expression patterns of eight lipogenic genes during the embryo-to-hatching transition are depicted in Fig. 14. The hepatic expression of THRSPA and SCD was log-scale and increased exponentially, reaching a plateau in liver of fed hatchlings between D3 and D9. Similar log-scale expression patterns were found for fiv other metabolic genes (ME1, SCD, FASN, ATPCL and ELOVL6), while HMGCL and HMGCS2 only showed a slight increase in transcript abundance in D9 hatchlings.
Figure 15 provides a side-by-side comparison of genome array analysis vs. qRT-PCR analysis for 13 “candidate” genes across 6 developmental ages (E16-D9). In general, transcript abundance determined by genome array analysis was of a lower magnitude than the dynamic log-scale range provided by qRT-PCR analysis, although patterns of gene expression were similar between thetwo analytical methods. Three lipolytic DEGs identified by genome array analysis (Fig. 15a, top) exhibited higher expression in liver of embryos than in hatchings [i.e., DIO3, MOGAT1 and pyruvate dehydrogenase kinase 4 (PDK4)]. The deactivating deiodinase DIO2 and the adipokine ADIPOQ showed a progressive increase in hepatic expression from E16 to a peak in E20 embryos. Although DIO2 levels were lower in D1 and D3 hatchlings, a peak in abundance was reached in D9 hatchlings. The qRT-PCR analysis (Fig. 15a, bottom) revealed high magnitude (log2) expression of DIO3, MOGAT1 and PDK4, which were highest in E16 and E18 embryos and lowest in liver of fully-fed D3 and D9 hatchlings. Likewise, the hepatic abundance of ADIPOQ, FZD2 and LDHB transcripts were also elevated in embryos compared to hatchlings, albeit with lower amplitudes. Seven DEGs identified by genome array analysis (Fig. 15b, top) displayed relatively low expression in liver of embryos with a sharp rise in expression in fully-fed D3 and D9 hatchlings. The developmental pattern and large-amplitude expression of THRSPA, ME1 and ELOVL6 transcripts were almost identical and reached a plateau in liver of fed D3 and D9 hatchlings. HMGCL levels were similar across embryonic ages and only slightly greater in D3 and D9 chicks. Also, the abundance of SCD was low in embryos and fasted D1 hatchlings, with a sharp increase in fed D3 and D9 hatchlings. In contrast, qRT-PCR analysis (Fig. 15b, bottom panel) shows an extremely high abundance of THRSPA and SCD in liver of D3 and D9 hatchlings that was several orders of magnitude greater than in liver of embryos (E16-E20). Five hepatic genes (FASN, ME1, ATPCL, ELOVL6 and HMGCL) had higher expression levels, albeit with a lower amplitude, in D3 and D9 hatchlings. We were surprised to find that hepatic expression of DIO1 was rather stable across embryonic (E16-E20) and hatchling (D1-D9) ages. Furthermore, DIO1 was not detected as a DEG by genome array analysis or qRT-PCR analysis; therefore, DIO1 was considered an invariant gene and not included in the Pearson’s Correlation Analysis.
The Pearson’s correlation coefficient (r = 0.929; 11 degrees of freedom) indicates a highly significant (P ≤ 0.01) correlation between gene expression levels obtained from both microarray and qRT-PCR analyses (Additional file 8). As expected, there was close agreement in hepatic transcript abundance between expression platforms (microarray and qRT-PCR analyses); the major difference was the large-magnitude log-scale expression levels achieved by qRT-PCR analysis, albeit the developmental expression patterns were quite similar.
The present study of liver during the embryo-to-hatchling transition (or peri-hatch period) provides a detailed view of the innate choreography of major transcriptional responses involved in the abrupt switch from a lipid-laden ectotherm (embryo) to a free-living endotherm (hatchling chick). We have identified several transcription factors (PPARA, PPARGC1A, KLF15, SIRT1, SERTAD2 and NR1H4) which interact with each other and their respective target genes in late embryos as they utilize yolk lipids and prepare for hatching on E21. Upon feeding, hatchling chicks abruptly increase expression of several ligand-activated transcription factors (THRSPA, SREBF2, PPARG, PPARD and KLF11) which in turn control numerous down-stream metabolic enzymes, transporters and kinases/phosphatases to support greater rates of energy expenditure for lipogenesis, adipogenesis and somatic growth. Hierarchical clustering revealed two opposing patterns of gene expression in liver during the peri-hatch period, where Cluster A had higher expression in embryos, which sharply declined after hatching (see Fig. 1). In contrast, Cluster B genes exhibited lower expression in embryos, which was followed by an exponential increase in abundance of numerous lipogenic genes in fully-fed D3 and D9 hatchlings. These two distinct patterns represent clusters of functionally-related genes that are similar to those found in our original transcriptional profiling of liver in embryos and hatchlings during the peri-hatch period using a prototype 3.2 K chicken liver microarray [19, 30, 31]. Self-organizing maps (SOMs) analysis of the first transcriptional scans using our 3.2 K liver cDNA array (NCBI GEO Platform No. GPL1742) revealed similar distinct patterns of transcript abundance. Cluster 11 contained 46 genes with higher expression in embryos, whereas Cluster 3 (21 genes) had lower expression in embryos with a sharp elevation in hatchlings . Higher resolution of functionally related gene clusters, unique to embryos and to hatchling chicks, was achieved with our spanning-tree clustering method (see Fig. 3 in ). Our original transcriptional profiling analysis with a low density 3.2 K microarray provided us with a coarse view of the topography of hepatic gene expression during the embryo-to-hatchling transition. The present transcriptional analysis using chicken genome arrays has greatly expanded the repertoire of responsive hepatic genes involved in the embryo-to-hatchling transition and considerable insight into the regulatory and metabolic pathways that control the abrupt switch from embryonic ectotherm to hatchling endotherm.
Transcriptional regulation of metabolism during the peri-hatch period
The present study utilized multiple pairwise contrasts between embryos and hatchlings to gain even greater resolution by populating canonical pathways and biological processes with DEGs, which provide homeorhetic control over metabolism during the peri-hatch period. The primary contrast of all embryos against all hatchlings revealed a large number of hepatic DEGs (284 genes) that are involved in lipid metabolism (Table 1). And of these, 149 DEGs belong to the lipid synthesis pathway, where only one-third (54) of the DEGs are more abundant in liver of embryos. We discovered 14 up-stream regulators that control transcription of a large number of metabolic enzymes, transporters, phosphatases and kinases in liver of embryos and hatchlings (see Fig. 4). Membrane-bound and nuclear ligand-activated receptors and G-protein coupled receptors (GCPRs) fine tune and amplify the signals mediated by a few critical transcription factors. These enhanced metabolic activities are enabled by a robust thyroid axis and deiodinases with opposing actions. In the embryo liver, DIO3 converts the prohormone T4 into inactive rT3 until the moment the embryo pips through its shell and begins pulmonary respiration . Furthermore, the abrupt rise in regulated body temperature (Tb = 42–44 °C) of precocial hatchling chicks is supported by increased hepatic DIO2 activity, which converts T4 into metabolic T3, the major hormone of thermoregulation and energy expenditure. DIO3 was one of the highest DEGs in liver of embryos, whereas DIO2 was the highly-expressed and opposing deiodinase found in liver of hatchlings by multiple contrasts. Initially, DIO1 was considered to be the major hepatic deiodinase providing active thyroid hormone (T3) for rapid growth and metabolism in hatchling chicks [15, 18]. The present study clearly shows that DIO2 expression increases in E20 liver to provide metabolically-active thyroid hormone (T3) for hatching and reaches maximum expression in fully-fed rapidly-growing D9 chicks. The importance of DIO3 in deactivating T4 was demonstrated by the D1 vs. D9 contrast, where DIO3 was highly expressed in liver of fasted D1 hatchlings and sharply depressed in fed D3 and D9 hatchings. The opposing actions of DIO3 and DIO2 on lipid catabolism or lipogenesis (respectively) were evident across multiple contrasts of embryo and hatchling liver (see Additional file 3).
As expected from our previous transcriptional scans [19, 30,31,32,33,34,35], major lipogenic transcriptional factors (THRSPA, SREBP2, PPARG, PPARD and KLF11) interact with each other to control numerous target genes encoding metabolic enzymes, transporters and kinases/phosphatases that support lipogenesis and the thermogenic drive, which sharply increases after hatching and consuming high-energy high-protein feed. Of particular interest, was the over-representation of AR-DEGs (27 genes) observed in the “THR-RXR Activation” pathway (see Additional file 4), where the most abundant DEGs expressed in embryo liver was DIO3, which deactivates both prohormone (T4) and metabolically-active (T3) thyroid hormone, whereas DIO2 activates the thyroid axis by converting T4 to T3, the active thyroid hormone and major regulator of energy metabolism and lipogenesis. Furthermore, T3 bound to its ligand-activated receptor (THRB) sequentially activates triple T3-THRB binding sites on the THRSP promoter to provoke an exponential increase in transcription of the THRSP gene . In the present study, THRSPA followed a distinct sigmoidal expression pattern during the embryo-to-hatchling transition (see Fig. 15b), as demonstrated by qRT-PCR analysis. The increased generation of T3 by DIO2 in liver of hatchlings provides ligand for binding to and activating the THRB, an essential step in the progressive activation of THRSPA. Furthermore, activated nuclear receptor (T3-THRB) also serves as a binding partner (heterodimer) for several other nuclear receptors (FXR, LXR, RXR and VDR). For example, the THRSP promoter in the rat has three T3-THR response elements which act synergistically to provide maximal transcription of the THRSP gene . Transcription of THRSP is regulated by several other ligand-activated nuclear receptors: SREBP1c via sterol response elements (SRE) , retinoic acid-RXR , NR1I3, pregnane X receptor (PXR), CAR-RXR complex , estrogen via estrogen response elements (ERE) , and several metabolites: [glucose bound to carbohydrate regulatory element binding protein (ChREBP) , polyunsaturated fatty acid  and prostaglandin . In the present study, we identified a lipogenic gene network in the D1 vs. D9 contrast (see Fig. 11a) where THRSP and its paralog MID1 interacting protein 1 (MID1IP1) [also known as THRSPL/S14R] were up-regulated in liver of D9 hatchlings with direct interactions with SREBF2 and ACACA, which encodes the rate limiting enzyme in fatty acid synthesis. Apparently, both THRSP and MID1IP1 are glucose-responsive genes, which have overlapping roles in regulating expression of several lipogenic genes . These observations strongly support activation of the thyroid axis in hatchlings to increase energy expenditure to meet demands of endothermy and rapid somatic growth. The present study affirms and extends our original idea that THRSPA controls transcription of several metabolic enzymes and transporters, which contribute to lipogenesis and the exponential growth exhibited by newly-hatched (meat-type) chicks. Previously, we have shown abundant expression of THRSPA in abdominal fat and liver of the chicken and demonstrated its involvement in adipogenesis and lipogenesis [19, 30, 31, 33, 34, 45, 46]. Furthermore, we discovered a 9 bp insertion/deletion polymorphism in the DNA binding domain, which is associated with the abdominal fatness trait in multiple populations of broiler chickens. Others have shown that mutations in the THRSPA gene are associated with fatness phenotypes in chickens [47, 48], ducks  and geese . Taken together, THRSPA is obviously a major regulatory of lipogenic genes expressed in liver and visceral fat of the chicken and other poultry, like ducks and geese. A recent study in mice provides direct evidence that enhanced THRSP expression is controlled via a LXR-mediated, SREBP-dependent mechanism . These extensive findings strongly support the full recognition of THRSPA as a principal transcriptional regulator of lipogenesis and adipogenesis in the chicken and other domestic birds.
The transcriptional choreography of the metabolic switch from embryo to hatchling includes interaction of several ligand-activated nuclear receptors  namely, farnesoid X receptor (FXR or NR1H4; activated by bile acids and cholesterol derivatives), liver X receptor (LXR; activated by oxysterols and other cholesterol derivatives), the orphan nuclear receptor small heterodimer partner (SHP or NROB2), which regulates metabolism of cholesterol, bile acids, fatty acids and glucose , and PPARA. The overlapping roles and interplay of the ligand-activated nuclear receptors LXR and RXR has been elegantly described as a “ying and yang relationship”, which precisely regulates bile acid, cholesterol and triglyceride metabolism . Presently, most pairwise contrasts of embryo and hatchling liver transcriptomes have identified activation of LXR-RXR, and FXR-RXR heterodimers as over-represented metabolic pathways. Furthermore, we discovered other examples of opposing (or ying-yang) transcriptional regulators (i.e., PPARA vs. PPARG; SIRT1 vs. PPARG; SIRT1 vs. SREBF2; PPARG vs. PPARGC1A; SREBF2 vs. PPARA; NR0B2 vs. NR1H4; THRSPA vs. PPARA, and THRSPA vs. PPARGC1A) and metabolic enzymes (i.e., DIO3 vs. DIO2). Furthermore, we identified 9 of the 15 known members of the Kruppel-like transcription factor (KLF) family  as DEGs (KLF2, KLF3, KLF6, KLF9, KLF10, KLF11, KLF12, KLF13 and KLF15) from our 14 pairwise contrasts (Additional file 3). In particular, KLF15 seems to act as a homeorhetic transcription factor capable of “switching between lipogenesis and gluconeogenesis during fasting” . This dynamic adaptation in energy metabolism is achieved by KLF15 forming complexes with LXR-RXR heterodimers on the promoter of mouse SREBF1c, thereby inhibiting transcription of this major lipogenic transcription factor. The ability of KLFs to either activate or repress gene transcription depending upon metabolic demands was recently reviewed in an article that precisely coined their dynamic action as “crippling and uncrippling metabolic pathways” . The present examination of liver transcriptomes clearly shows interactions, interdependence and interconnectivity of multiple upstream regulators and their direct target genes in regulating lipolysis in late embryos and lipogenesis in fully-fed D3 and D9 hatchling chicks.
The NAD+-dependent deacetylase SIRT1 is involved in alternative energy utilization, like β-oxidation of lipid and gluconeogenesis [56,57,58], which are critical metabolic processes in liver of late embryos. Furthermore, SIRT1 enhances thyroid (T3) control over expression of several lipogenic genes, including CPT1A, PPARA, PPARGC1A, PDK4, PCK1 and SREBP-1c . Another transcription factor up-regulated in liver of embryos was prospero-related homeobox 1 (PROX1), which is important for hepatic embryogenesis and a negative regulator of triglyceride synthesis via mTOR signaling . Thus, the metabolism of the late chicken embryo is under the intricate control of multiple transcription factors and directed at utilization of yolk lipids that amass in liver. The retained yolk sac provides nutrients and lipid-derived energy until the hatchling eats its first meal and becomes a self-sufficient endotherm. In the present study, we included D1 hatchlings, which had not been fed and still dependent on yolk-lipid and swallowed albumen to reach and maintain endothermy. Full and untethered expression of multiple transcription factors and their target genes that contribute to lipogenesis occurs only after the hatchling chick is fully-fed (D3 and D9) and the residual yolk sac depleted.
Another interesting TF identified in the four contrasts of embryos and hatchlings was SERTA domain containing 2 (SERTAD2; also known as TRIP-Br2), which was expressed at highest levels (log2 ratio = + 1.75) in the E16 vs. D1 contrast. The differential expression of SERTAD2 was found in three additional contrasts of embryos against hatchlings (E vs. H, E16 vs. D3 and E18 vs. D1). This newly-discovered transcriptional co-regulator is highly expressed in visceral fat of obese humans; furthermore, SERTAD2 appears to control lipolysis, oxidative metabolism and thermogenesis in mouse models [61,62,63]. Knockout of the SERTAD2 gene in mice prevents diet-induced obesity, insulin resistance and inflammatory responses in visceral fat. Our discovery of up-regulated SERTAD2 in liver of E16 and E18 chick embryos suggest that this transcriptional co-regulator could contribute to transcriptional control of lipolysis and oxidative metabolism in the lipid-laden liver of late chick embryos via interactions with known lipolytic TFs (i.e., PPARA, PPARGC1A, SIRT1 and NR1H4).
Our present transcriptional study has also revealed acquisition of a competent blood clotting system in E16 and E18 embryos, which exhibit abundant differential expression of multiple coagulation factors and collagen genes (see Figs. 6a and 10b; Additional files 4 and 6). The deacetylase SIRT1 appears to control essential biological pathways like blood coagulation through selective acetylation of clotting factors and the closely related acute-phase response proteins . SIRT1 also enhances PPARA-PPARGC1A driven lipolysis and lipid depletion of adipocytes, which in turn activates gluconeogenesis a vital metabolic response of the embryo. The largest number of coagulation genes was found in the E18 vs. D3 contrast, where 10 clotting factors and 7 collagen genes were expressed higher in the E18 embryo. Similarly, 12 blood clotting factors were found in the E20 vs. D3 contrast, where 9 coagulation and 3 collagen genes were up-regulated in the liver of E20 embryos. In addition, fibrinogen alpha (FGA) was over-expressed in E18 and E20 embryos (see Additional file 3). The transcriptional analysis of the chicken yolk sac in late chicken embryos (E13-E21) revealed peak abundance of fibrinogen (β and γ subunits; FGB and FGG) genes at E19 . Presumably, a competent blood coagulation system is crucial to prevent excessive blood loss in the embryo during absorption of its residual yolk sac into the visceral cavity, destruction of the chorioallantois, and emergence of the E21 embryo as a hatchling chick .
Hepatic expression of the feather keratin gene FKER
Interestingly, FKER was the most DEG found in liver of E16 embryos (see Additional file 3). In fact, seven distinct transcripts represented on the Affymetrix Chicken Genome Array correspond to feather keratin or feather keratin-like genes (Table 8), which reside at multiple loci on the chicken genome (GGA1, GGA25 and GGA27) [65,66,67]. Five of these over-expressed FKER transcripts map to GGA25, while FK1 (BI064513) resides on GGA1 and XM_424568 is located on GGA27. At first glance, abundant expression of FKER in liver of the E16 embryo seems difficult to explain. However, we originally discovered the exceptional, and very high, abundance of FKER in liver of E15 embryos from Leghorn (egg-type) chickens carrying the riboflavin binding-protein deficiency mutation rd/rd, with or without riboflavin rescue . We proposed that the abundant hepatic expression of FKER exclusively in E15 embryos could be related to massive engorgement of the embryo’s liver with yolk lipids, which normally occurs at this stage of embryonic development. In fact, an early study described the presence of lipids in the embryo’s developing cutaneous feather, which was likely associated with “pre-keratin”, where lipids function as the scaffold for feather keratin formation . In the mouse, cytoskeletal keratin 8 (Krt8) is required for structure and integrity of hepatocytes . Furthermore, ablation of the Krt8 in mice results in oxidative stress and excessive accumulation of lipid and protein metabolites . Examination of the pigeon’s lactating crop transcriptome has revealed up-regulation of several β-keratin genes, including FKER, which are ultimately involved in production of lipid-laden crop milk . The unique synthesis of crop milk in Columbiformes involves prolactin-mediated cornification and lipid synthesis of crop epithelial cells, which are desquamated and sloughed off to form “crop milk”, which is fed to the altricial squabs. This transcriptional study of pigeon crop milk production clearly shows over-expression of multiple β-keratin genes, which are associated with de novo lipid synthesis and cornification of crop epithelium. Of particular interest, all of the differentially-expressed lipogenic genes found in the lactating pigeon crop sac were identified as DEGs in the present study, including the highly-expressed lipogenic transcription factor THRSPA . The present descriptive study of the liver transcriptome during the embryo-to-hatchling transition validates our original discovery of abrupt over-expression of FKER in E15 Leghorn (egg-type) embryos , which coincides with massive accumulation of yolk lipid in liver that fuels the final embryonic growth phase and emergence of the hatchling chick. Additional support for the ectopic expression of hepatic FKER transcripts in E16 Leghorn embryos is provided by our initial observation of differential expression of a cDNA target (pgf1n.pk001.j5; GenBank BI064513) on a low-density 3.2 K chicken microarray (NCBI GEO Platform GPL1742) . This cDNA clone (pgf1n.pk001.j5) was sequenced from a normalized adipose cDNA library derived from embryo and hatchling chickens. Furthermore, our recent transcriptional study of liver from rd/rd Leghorn embryos, riboflavin-deficient and riboflavin rescued, utilized the Arizona 20.7 K chicken oligo array (NCBI GEO Platform GPL6049), where three 70-mer oligo targets [Roslin Institute Gallus gallus (RIGG) oligo; RIGG10897, RIGG14163 and RIGG14953] had 19-times greater abundance in liver of riboflavin-rescued embryos than riboflavin-deficient embryos between E13 and E15. Thus, three independent chicken microarray platforms and multiple expressed-sequence tag (EST) sequences in GenBank clearly support abundant “ectopic” expression of FKER in liver of E15 and E16 embryos from both egg-type (Leghorn)  and, presently, meat-type (broiler) chicken breeds. However, the biological function of FKER in the lipid-laden liver of late embryos and the relationship between hepatic expression and cutaneous expression of FKER in embryos will require further definitive study.
The present analysis of the embryo-to-hatchling transition in meat-type chickens provides the first detailed view of the choreography of innate and dynamic transcriptional responses made by embryos and newly hatched chicks during the critical metabolic jump from lipid-laden ectotherm to free-living self-sufficient endotherm, respectively. The metabolism of late (E16-E20) embryos is dominated by β-oxidation of yolk lipids, glycolysis and gluconeogenesis, which utilizes the residual albumen and other proteins in allantoic fluid ingested by the embryo just prior to hatching. These metabolic processes are precisely controlled by multiple and interacting transcription factors (PPARA, PPARGC1A, NR1H4 and SIRT1, SERTAD2, KLF11, KLF13 and KLF15) and scores of metabolic enzymes, transporters, phosphatases and kinases. The metabolites derived from these processes interact with nuclear receptors, G-coupled protein receptors (GCPRs), transmembrane receptors and extracellular factors. Highly-expressed hepatic DIO3 inactivates thyroid hormones in the embryo, while its expression sharply declines in hatchling chicks. Although hepatic expression of DIO1 was rather constant during the peri-hatch period, another thyroid hormone-activating enzyme DIO2 reached peak abundance at two ages, E20 and D9. Upon pipping through the egg shell, the hatching chick immediately converts to the pro-hormone T4 into metabolically-active T3, which supports achievement of thermodynamic freedom in hatchlings and synergetic activation of THRSPA, the major lipogenic transcription factor in chickens. Furthermore, T3 bound to its receptor (THRB) forms heterodimers with other nuclear receptors, which then control programmed expression of major lipogenic enzymes (FASN, ME1, SCD, ELOV2, ELOV25, ELOV26, LSS, AGPAT2, ACACA, CYP7A, DHCR7 and DHCR724). In addition to THRSPA, other transcriptional regulators (PPARG, PPARD, LPIN1, KLF11 and SREBF2) also interact and contribute to the sharp increase in lipogenesis observed in liver of D3 and D9 hatchlings.
Our present observation of the extraordinary expression of multiple FKER transcripts in liver of E16 embryos from meat-type chickens validates our original discovery of highly-expressed FKER gene in liver of e15 Leghorn (egg-type) chickens. We have proposed that hepatic expression of FKER in e15-e16 embryos could be an adaptive response to the coincident and massive accumulation of yolk lipids. This idea is supported by the exceptional expression of multiple β-keratin transcripts, including FKER, in the lactating crop sac of the pigeon , which accumulates and synthesizes lipids in production of “crop milk” that nourishes the altricial squabs. Many of the DEGs involved in de novo synthesis of triglycerides in the pigeon’s lactating crop sac are also found in liver of hatchling chicks, including the major lipogenic transcription factor THRSPA. Nevertheless, the biological function of FKER expressed in liver of e15-e16 embryos and the relationship between hepatic FKER and cutaneous FKER remain to be elucidated by definitive study. The present study provides new insight into dynamic interaction of multiple transcription factors and their direct target genes that provide homeorhetric regulation of metabolism during the abrupt embryo-to-hatchling transition of the domestic chicken, Gallus gallus.
Average cycle time
Embryonic day 16
False discovery rate
Gene Expression Omnibus
- GGA :
Gallus gallus chromosome
General linear mode
Ingenuity Pathway Analysis
Linear Models for Microarray Data
Least significant difference
Minimum information about a microarray experiment
National Center for Biotechnology Information
Quantitative reverse transcriptase polymerase chain reaction
RNA integrity number
Statistical Analysis System
Thyroid hormone responsive Spot 14 protein, alpha
Decuypere E, Nouwen EJ, Kuhn ER, Geers R, Michels H. Iodohormones in the serum of chick embryos and post-hatching chickens as influenced by incubation temperature. Relationship with the hatching process and thermogenesis. Ann Biol Anim Bioch Biophys. 1979;19(6):1713–23.
Speake BK, Murray AMB, Noble RC. Transport and transformations of yolk lipids during development of the avian embryo. Prog Lipid Res. 1998;37:1–32.
Feast M, Noble RC, Speake BK, Ferguson MWJ. The effect of temporary reductions in incubation temperature on growth characteristics and lipid utilisation in the chick embryo. J Anat. 1998;193:383–90.
Noy Y, Sklan D. Nutrient use in chicks during the first week posthatch. Poult Sci. 2002;81:391–9.
Noy Y, Sklan D. Energy utilization in newly hatched chicks. Poult Sci. 1999;78:1750–6.
Sklan D. Fat and carbohydrate use in posthatch chicks. Poult Sci. 2003;82:117–22.
Uni Z, Tako E, Gal-Garber O, Sklan D. Morphological, molecular, and functional changes in the chicken small intestine of the late-term embryo. Poult Sci. 2003;82:1747–54.
Yadgary L, Wong EA, Uni Z. Temporal transcriptome analysis of the chicken embryo yolk sac. BMC Genomics. 2014;15:690.
Goodridge AG, Garay A, Silpananta P. Regulation of lipogenesis and the total activities of lipogenic enzymes in a primary culture of hepatocytes from prenatal and early postnatal chicks. J Biol Chem. 1974;249:1469–75.
Morris SM, Winberry LK, Fisch JE, Back DW, Goodridge AG. Developmental and nutritional regulation of the messenger RNAs for fatty acid synthase, malic enzyme and albumin in the livers of embryonic and newly hatched chicks. Mol Cell Biochem. 1984;64:63–8. https://doi.org/10.1007/BF00420929.
Goodridge AG, Fantozzi DA, Klautky SA, Ma XJ, Roncero C, Salati LM. Nutritional and hormonal regulation of genes for lipogenic enzymes. Proc Nutr Soc. 1991;50:115–22.
Bigot K, Mignon-Grasteau S, Picard M, Tesseraud S. Effects of delayed feed intake on body, intestine and muscle development in neonate broilers. Poult Sci. 2003;82:781–8.
Halevy O, Geyra A, Barak M, Uni Z, Sklan D. Early posthatch starvation decreases satellite cell proliferation and skeletal muscle growth in chicks. J Nutr. 2000;13:858–64.
Decuypere E, Kúhn ER, Clijmans B, Nouwen EJ, Michels H. Prenatal peripheral monodeiodination in the chick embryo. Gen Comp Endocrinol. 1982;47:15–7.
Darras VM, Visser TJ, Berghman LR, Kuhn ER. Ontogeny of type I and type III deiodinase activities in embryonic and posthatch chicks: relationship with changes in plasma triiodothyronine and growth hormone levels. Comp Biochem Physiol. 1992;103A:131–6.
De Groef B, Grommen SVH, Darras VM. Hatching the cleidoic egg: the role of thyroid hormones. Front Endocrinol. 2013;4:63.
Van der Geyten S, Van der Eynde I, Segers IB, Kuhn ER, Darras VM. Differential expression of iodothyronine deiodinases in chicken tissues during the last week of embryonic development. Gen Comp Endocrinol. 2002;128:65–73.
Van der Geyten S, Sanders JP, Kaptein E, Darras VM, Kuhn ER, Leonard JL, Visser TJ. Expression of chicken hepatic type I and type III iodothyronine deiodinases during embryonic development. Endocrinology. 1997;138:5144–52.
Cogburn LA, Wang X, Carré W, Rejto L, Porter TE, Aggrey SE, Simon J. Systems-wide chicken DNA microarrays, gene expression profiling and discovery of functional genes. Poult Sci. 2003;82:939–51.
R Core Team: R: A language and environment for statistical computing. https://www.R-project.org/; Accessed 11 Dec 2017.
Bioconductor: Open source software for bioinformatics. https://bioconductor.org/. Accessed 12 Mar 2018.
Wettenhall JM, Simpson KM, Satterley K, Smyth GK. affylmGUI: a graphical user interface for linear modeling of single channel microarray data. Bioinformatics. 2006;22:897–9.
Irizarry RA, Hobbs B, Collin F, Beazer-Barclay YD, Antonellis KJ, Scherf U, Speed TP. Exploration, normalization, and summaries of high density oligonucleotide array probe level data. Biostat. 2003;4:249–64.
Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK. Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43:e47.
Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J Royal Stat Soc (Series B). 1995;57:289–300.
Carlson M: chicken.db: Affymetrix chicken annotation data (chip chicken). R package version 3.2.3.; 2016.
Huibers RP, de Jong M, Dekter RW, Van den Ackerveken G. Disease-specific expression of host genes during downy mildew infection of Arabidopsis. Mole Plant-Microbe Interact. 2009;22:1104–15.
Ingenuity Pathway Analysis. https://www.qiagenbioinformatics.com/products/ingenuity-pathway-analysis/. Accessed 15 Nov 2017.
Biogazelle qbase+ software. https://www.qbaseplus.com/. Accessed 20 July 2016.
Cogburn LA, Wang X, Carré W, Rejto L, Aggrey SE, Duclos MJ, Simon J, Porter TE. Functional genomics in chickens: development of integrated-systems microarrays for transcriptional profiling and discovery of regulatory pathways. Comp Funct Genom. 2004;5:253–61.
Cogburn LA, Porter TE, Duclos MJ, Simon J, Burgess SC, Zhu JJ, Cheng HH, Dodgson JB, Burnside J. Functional genomics of the chicken--a model organism. Poult Sci. 2007;86:2059–94.
Wang X, Carré W, Saxton AM, Cogburn LA. Manipulation of thyroid status and/or GH injection alters hepatic gene expression in the juvenile chicken. Cytogenet Genome Res. 2007;117:174–88.
Resnyk CW, Carré W, Wang X, Porter TE, Simon J, Le Bihan-Duval E, Duclos MJ, Aggrey SE, Cogburn LA. Transcriptional analysis of abdominal fat in genetically fat and lean chickens reveals adipokines, lipogenic genes and a link between hemostasis and leanness. BMC Genomics. 2013;14:557.
Resnyk CW, Carré W, Wang X, Porter TE, Simon J, Le Bihan-Duval E, Duclos MJ, Aggrey SE, Cogburn LA. 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. BMC Genomics. 2017;18:626.
Cogburn LA, Smarsh DN, Wang X, Trakooljul N, Carré W, White HB 3rd. Transcriptional profiling of liver in riboflavin-deficient chicken embryos explains impaired lipid utilization, energy depletion, massive hemorrhaging, and delayed feathering. BMC Genomics. 2018;19:177.
Liu HC, Towle HC. Functional synergism between multiple thyroid hormone response elements regulates hepatic expression of the rat S14 gene. Mol Endocrinol. 1994;8:1021–37.
Jump DB, Thele AP, Mater MK. Functional interaction between sterol regulatory element-binding protein-1c, nuclear factor Y, and 3,5,3′-triiodothyronine nuclear receptors. J Bio Chem. 2001;276:34419–27.
MacDougald OA, Jump DB. Localization of an adipocyte-specific retinoic acid response domain controlling S14 gene transcription. Biochem Biophys Res Comm. 1992;188:470–6.
Breuker C, Moreau A, Lakhal L, Tamasi V, Parmentier Y, Meyer U, Maurel P, Lumbroso S, Vilarem MJ, Pascussi JM. Hepatic expression of thyroid hormone-responsive spot 14 protein is regulated by constitutive androstane receptor (NR1I3). Endocrinology. 2010;151:1653–61.
Ren J, Xu N, Zheng H, Tian W, Li H, Li Z, Wang Y, Tian Y, Kang X, Liu X. Expression of thyroid hormone responsive SPOT 14 gene is regulated by estrogen in chicken (Gallus gallus). Sci Rep. 2017;7:10243.
Towle HC, Kaytor EN, Shih H-M. Regulation of the expression of lipogenic enzyme genes by carbohydrate. Annu Rev Nutr. 1997;17:405–33.
Jump DB, Clarke SD, MacDougald O, Thelen A. Polyunsaturated fatty acids inhibit S14 gene transcription in rat liver and cultured hepatocytes. Proc Natl Acad Sci U S A. 1993;90:8454–8.
Mater MK, Thelen AP, Jump DB. Arachidonic acid and PGE2 regulation of hepatic lipogenic gene expression. J Lipid Res. 1999;40:1045–52.
Aipoalani DL, O'Callaghan BL, Mashek DG, Mariash CN, Towle HC. Overlapping roles of the glucose-responsive genes, S14 and S14R, in hepatic lipogenesis. Endocrinology. 2010;151:2071–7.
Wang X, Carré W, Zhou H, Lamont SJ, Cogburn LA. Duplicated spot 14 genes in the chicken: characterization and identification of polymorphisms associated with abdominal fat traits. Gene. 2004;332:79–88.
Resnyk CW, Chen C, Huang H, Wu CH, Simon J, Le Bihan-Duval E, Duclos MJ, Cogburn LA. RNA-Seq analysis of abdominal fat in genetically fat and lean chickens highlights a divergence in expression of genes controlling adiposity, hemostasis, and lipid metabolism. PLoS One. 2015;10:e0139549.
Cao ZP, Wang SZ, Wang QG, Wang YX, Li H. Association of Spot14α gene polymorphisms with body weight in the chicken. Poult Sci. 2007;86:1873–80.
d'Andre Hirwa C, Yan W, Wallace P, Nie Q, Luo C, Li H, Shen X, Sun L, Tang J, Li W, Zhu X, Yang G, Zhang X. Effects of the thyroid hormone responsive spot 14α gene on chicken growth and fat traits. Poult Sci. 2010;89:1981–91.
Zhan K, Hou ZC, Li HF, Xu GY, Zhao R, Yang N. Molecular cloning and expression of the duplicated thyroid hormone responsive spot 14 (THRSP) genes in ducks. Poult Sci. 2006;85:1746–54.
Su S, Zhu H, Li Q, Xie Z. Molecular cloning and sequence analysis of spot 14α in geese. Brit Poult Sci. 2009;50:459–66.
Wu J, Wang C, Li S, Li S, Wang W, Li J, Chi Y, Yang H, Xi K, Zhou Y, Dong C, Wang F, Xu G, Yang J, Gustafsson J, Guan Y. Thyroid hormone-responsive SPOT 14 homolog promotes hepatic lipogenesis, and its expression is regulated by liver x receptor α through a sterol regulatory element-binding protein -1c-dependent mechanism in mice. Hepatology. 2013;58:617–28.
Kalaany NY, Mangelsdorf DJ. LXRS and FXR: the yin and yang of cholesterol and fat metabolism. Annu Rev Physiol. 2006;68:159–91.
Garruti G, Wang HH, Bonfrate L, de Bari O, Wang DQH, Portincasa P. A pleiotropic role for the orphan nuclear receptor small heterodimer partner in lipid homeostasis and metabolic pathways. J Lipids. 2012;2012:304292.
Pollak NM, Hoffman M, Goldberg IJ, Drosatos K. Krüppel-like factors: Crippling and uncrippling metabolic pathways. JACC: Basic Translat Sci. 2018;3:132–56.
Takeuchi Y, Yahagi N, Aita Y, Murayama Y, Sawada Y, Piao X, Toya N, Oya Y, Shikama A, Takarada A, Masuda Y, Nishi M, Kubota M, Izumida Y, Yamamoto T, Sekiya M, Matsuzaka T, Nakagawa Y, Urayama O, Kawakami Y, Iizuka Y, Gotoda T, Itaka K, Kataoka K, Nagai R, Kadowaki T, Yamada N, Lu Y, Jain MK, Shimano H. KLF15 enables rapid switching between lipogenesis and gluconeogenesis during fasting. Cell Rep. 2016;16:2373–86.
Picard F, Kurtev M, Chung N, Topark-Ngarm A, Senawong T, de Oliveira RM, Leid M, McBurney MW, Guarente L. Sirt1 promotes fat mobilization in white adipocytes by repressing PPAR-γ. Nature. 2004;429:771–6.
Chang HC, Guarente L. SIRT1 and other sirtuins in metabolism. Trends Endocrinol Metab. 2014;25:138–45.
Simmons GE, Pruitt WM, Pruitt K. Diverse roles of SIRT1 in cancer biology and lipid metabolism. Intl JMol Sci. 2015;16:950–65.
Thakran S, Sharma P, Attia RR, Hori RT, Deng X, Elam MB, Park EA. Role of sirtuin 1 in the regulation of hepatic gene expression by thyroid hormone. J Biol Chem. 2013;288:807–18.
Kwon S, Jeon JS, Kim SB, Hong YK, Ahn C, Sung JS, Choi I. Rapamycin up-regulates triglycerides in hepatocytes by down-regulating Prox1. Lipids Health Dis. 2016;15:41.
Liew CW, Boucher J, Cheong JK, Vernochet C, Koh HJ, Mallol C, Townsend K, Langin D, Kawamori D, Hu J, Tseng YH, Hellerstein MK, Farmer SR, Goodyear L, Doria A, Blüher M, Hsu SIH, Kulkarni RN. Ablation of TRIP-Br2, a regulator of fat lipolysis, thermogenesis and oxidative metabolism, prevents diet-induced obesity and insulin resistance. Nat Med. 2013;19:217.
Qiang G, Kong HW, Fang D, McCann M, Yang X, Du G, Blüher M, Zhu J, Liew CW. The obesity-induced transcriptional regulator TRIP-Br2 mediates visceral fat endoplasmic reticulum stress-induced inflammation. Nat Commun. 2016;7:11378.
Qiang G, Kong HW, Gil V, Liew CW. Transcription regulator TRIP-Br2 mediates ER stress-induced brown adipocytes dysfunction. Sci Reps. 2017;7:40215.
Moran JET. Nutrition of the developing embryo and hatchling. Poult Sci. 2007;86:1043–9.
Presland RB, Gregg K, Molloy PL, Morris CP, Crocker LA, Rogers GE. Avian keratin genes I. a molecular analysis of the structure and expression of a group of feather keratin genes. J Mol Biol. 1989;209:549–59.
Wu P, Ng CS, Yan J, Lai YC, Chen CK, Lai YT, Wu SM, Chen JJ, Luo W, Widelitz RB, Li WH, Chuong C-M. Topographical mapping of α- and β-keratins on developing chicken skin integuments: functional interaction and evolutionary perspectives. Proc Natl Acad Sci U S A. 2015;112:E6770–9.
Bhattacharjee MJ, Yu CP, Lin JJ, Ng CS, Wang TY, Lin HH, Li WH. Regulatory divergence among beta-keratin genes during bird evolution. Mol Biol Evol. 2016;33:2769–80.
Bell E, Thathachari YT. Development of feather keratin during embryogenesis of the chick. J Cell Biol. 1963;16:215.
Loranger A, Duclos S, Grenier A, Price J, Wilson-Heiner M, Baribault H, Marceau N. Simple epithelium keratins are required for maintenance of hepatocyte integrity. Am J Pathol. 1997;151:1673–83.
Zhou Q, Ji X, Chen L, Greenberg HB, Lu SC, Omary MB. Keratin mutation primes mouse liver to oxidative injury. Hepatology. 2005;41:517–25.
Gillespie MJ, Crowley TM, Haring VR, Wilson SL, Harper JA, Payne JS, Green D, Monaghan P, Donald JA, Nicholas KR, Moore RJ. Transcriptome analysis of pigeon milk production--role of cornification and triglyceride synthesis genes. BMC Genomics. 2013;14:169.
The contribution of undergraduate research students (Bret Glass, Jessica Hall and Danielle Smarsh) toward completion of this project is gratefully acknowledged. We appreciate the support of Professor Karl Steiner, former Director of the Delaware Biotechnology Institute, for intramural funding and use of the Microarray Core Facility.
This research was supported by: the United States Department of Agriculture, Cooperative State Research, Education, and Extension Service, Initiative for Future Agriculture and Food Systems (USDA-CREES-IFAFS) Grant 00–52,100–9614 to LAC; the USDA-National Research Initiative (NRI) #2005–35206-15288 to LAC; and by the Delaware INBRE program (Grant P20GM103446) and institutional resources of the Center for Bioinformatics and Computational Biology at University of Delaware to CC, HH and CHW. The cost of this article’s publication was provided by the Department of Animal and Food Sciences, University of Delaware.
Availability of data and materials
The minimum information about a microarray experiment (MIAME)-compliant microarray data set obtained from this study are available from the NCBI GEO database under accession number GSE109451 [https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE109451]. Additional data generated or analyzed during this study are included within this published article and its supplementary files. Multiple worksheets in Additional file 3 provide lists of 14 pair-wise contrasts of DE gene sets embryos and hatchling chicks from statistical analysis of the microarray data.
Ethics approval and consent to participate
The care, maintenance and experimental use of embryos and hatchling chicks were performed in accordance with the United States Department of Agriculture (USDA) guidelines on the use of agricultural animals in research and approved by the University of Delaware Agricultural Animal Care and Use Committee (ACUC).
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Experimental design of microarray hybridizations. A Microsoft Excel file containing a single work sheet “Array Hybridization Design” which describes the hybridization scheme used for the 24 Affymetrix Chicken Genome Arrays. (XLSX 12 kb)
Primers used for qRT-PCR analysis. A Microsoft Excel file containing a single work sheet “qRT-PCR Primer Information”. This table provides the chicken gene symbol, GenBank ID number, 5′-3′ sequence for forward and reverse primers, and amplicon size (bp) for each gene used for qRT-PCR analysis. (XLSX 12 kb)
Annotated gene list for 14 pairwise contrasts of liver transcriptomes in embryos and hatchlings across six ages during the peri-hatch period (E16-D9). The last worksheet provides the list of 1005 commonly-shared AR-DEGs identified by the Venn diagram shown in Fig. 3b. Each worksheet provides a list of the “Analysis Ready” (AR)-DEGs (FDR adjusted P ≤ 0.05; ±0.75 cutoff) for each contrast. (XLSX 905 kb)
Canonical and biological pathways identified by the contrast of all embryos vs. all hatchlings. The AR-DEGs from the E vs. H contrast were assigned by Ingenuity Pathway Analysis (IPA) to canonical and functional and pathways. Fisher’s Exact Test is used in IPA to determine probability (P ≤ 0.05) that AR-DEGs belong to a particular canonical pathway or biological function accrued in the Ingenuity Knowledge Base. (XLSX 53 kb)
Canonical and biological pathways identified by IPA in the E18 vs. E20 contrast. The AR-DEGs from the E18 vs. E20 contrast were assigned by IPA to overrepresented canonical and functional pathways. IPA uses Fisher’s Exact Test to determine probability (P ≤ 0.05) of each AR-DEG belonging to a particular canonical pathway or biological function. (XLSX 36 kb)
Canonical-biological pathways identified by IPA in the contrast of E18 embryos vs. D3 hatchlings. The AR-DEGs identified in the E18 vs. D3 contrast were assigned by IPA to overrepresented canonical and functional and pathways. IPA uses Fisher’s Exact Test to determine probability (P ≤ 0.05) of each AR-DEG belonging to a particular canonical pathway or biological function. (XLSX 44 kb)
Canonical and biological pathways found by IPA in the D1 vs. D9 hatchling contrast. The AR-DEGs identified in the D1 vs. D9 contrast were assigned by IPA to overrepresented canonical and functional and pathways. IPA uses Fisher’s Exact Test to determine probability (P ≤ 0.05) of each AR-DEG belonging to a particular canonical pathway or biological function. Each worksheet provides a list of AR-DEGs assigned to particular canonical pathways or biological functions, annotated with gene symbol, Entrez name, log2 ratio, cellular compartment, type and Entrez protein ID. (XLSX 31 kb)
Pearson’s correlation analysis of gene expression levels of 15 “candidate” DEGs determined by genome microarray and verified by qRT-PCR analyses. A Microsoft Excel file containing the average normalized expression levels (log2 ratios) of 15 candidate DEGs as determined by genome microarray analysis and verified by qRT-PCR analysis. The Pearson’s correlation coefficient (r = 0.9287; P ≤ 0.01, 11 degrees of freedom) indicates a highly significant correlation between gene expression levels obtained from both microarray and qRT-PCR analyses. (XLSX 16 kb)
About this article
Cite this article
Cogburn, L.A., Trakooljul, N., Chen, C. et al. Transcriptional profiling of liver during the critical embryo-to-hatchling transition period in the chicken (Gallus gallus). BMC Genomics 19, 695 (2018). https://doi.org/10.1186/s12864-018-5080-4