Transcriptomic analyses of regenerating adult feathers in chicken
BMC Genomics volume 16, Article number: 756 (2015)
Feathers have diverse forms with hierarchical branching patterns and are an excellent model for studying the development and evolution of morphological traits. The complex structure of feathers allows for various types of morphological changes to occur. The genetic basis of the structural differences between different parts of a feather and between different types of feather is a fundamental question in the study of feather diversity, yet there is only limited relevant information for gene expression during feather development.
We conducted transcriptomic analysis of five zones of feather morphologies from two feather types at different times during their regeneration after plucking. The expression profiles of genes associated with the development of feather structure were examined. We compared the gene expression patterns in different types of feathers and different portions of a feather and identified morphotype-specific gene expression patterns. Many candidate genes were identified for growth control, morphogenesis, or the differentiation of specific structures of different feather types.
This study laid the ground work for studying the evolutionary origin and diversification of feathers as abundant data were produced for the study of feather morphogenesis. It significantly increased our understanding of the complex molecular and cellular events in feather development processes and provided a foundation for future studies on the development of other skin appendages.
The genetic and developmental basis of morphological complexity is one of the most important issues in evolutionary biology [1, 2]. Avian feather provides an excellent system for studying the evolution and development of novel morphological traits because it has diverse forms [3–6], and the complex structure of feathers allows for various types of morphological changes to occur. Yet, feathers are homologous with the simpler scales of reptiles and could have evolved from a scale-like epidermal appendage of dinosaurian ancestors of birds [7–10].
Feathers have evolved to have different forms in color, morphology and mechanical properties not only among different bird species, but also among different body regions of a bird individual, giving us an excellent model to study the molecular basis of phenotypic variation of an important structure in a single species. The feather has been used as a model to study morphogenesis of skin appendages [11, 12]. Several candidate genes have been found to be involved in feather formation . For examples, barb and rachis are formed by a periodic invagination and regulated by BMP, NOG, SPRY, and FGF. Moreover, the basal branch pattern is formed by differential cell death and regulated by NCAM, SHH, and caspase. In addition, radial, bilateral symmetric, and asymmetric branching patterns are formed by modulating basal branching circuit that is involved the WNT3A gradient and SPRY. However, the genetic basis of feather variation is still largely unknown. A better understanding of the molecular dynamics associated with the process of feather growth will provide insight into the evolution of diverse feather structures.
The feather is also an excellent example of exaptation. Feathers initially might have evolved for heat regulation, but were then co-opted for display, and later co-opted for flight. These and other evolutionary novelties probably have occurred by changing the expression patterns of genes involved in feather development. The evolutionary co-option of plesiomorphic molecular signaling modules allows for the morphological innovations of feathers to originate and evolve [14, 15].
A hypothesis of morphological evolution postulates that form evolves largely through altering the expression of conserved genes . The molecular and developmental mechanisms that produce the diversification of feather are still poorly understood. The epithelium and the mesenchyme are two major components in feather follicles [16–18]. The epithelium includes both the epithelium enwrapping the mesenchyme and the feather wall epithelium that is connected with the interfollicular epidermis. The mesenchyme includes the dermal papilla and the pulp [19, 20]. The invagination of the multilayered epithelium in the ramogenic zone starts branching morphogenesis. The rachis is formed by fusion of barb ridges at the anterior end of the feather. The marginal plate in basal layer flanking each barb ridge and axial cells undergo apoptosis after the barbule plates are keratinized. The feather branches open in the more mature distal end after the apoptosis of feather sheath and pulp epithelium. Thanks to the feasibility of experimental manipulation and observation, feather regeneration can be analyzed in a comprehensive way and has been proposed to be a unique model for understanding organogenesis .
High-throughput sequencing technologies have been applied to characterize transcriptome architectures [21–26]. Systems biology study provides a new technology platform that can reveal molecular expression profiles associated with different morphological developments. Bioinformatic analyses are used to identify genes associated with feather and scale differences . These technologies and skills were used in this study.
The main goal of this study was to identify differentially expressed genes between different portions of feather using RNA-Seq. We characterized and quantified mRNAs that are expressed in the feather base during feather development in the domestic chicken. Feathers develop from the distal end to the proximal end in a temporal-spatial manner, thus providing an opportunity to analyze gene expression profiles associated with different zones of a feather (Fig. 1a) . Two zones of body feather and three zones of flight feather were selected to represent morphological, structural, and mechanical property differences in feathers (Additional file 1: Figure S1). We made four comparisons: 1) between pennaceous and plumulaceous portions of body feather for understanding how the morphological differences between two parts of a body feather arise; 2) between the pennaceous portions of body feather and flight feather for understanding differences in physical and mechanical properties; 3) between the distal pennaceous portion and the proximal pennaceous portion of flight feather for understanding how the morphological differences between two parts of a flight feather arise; 4) between the proximal pennaceous portion and the calamus of a flight feather for understanding how the barb and rachis are differentiated. These analyses shed light on the genetic basis of feather diversity.
Results and discussions
We used Illumina Hiseq 2000 to obtain five feather transcriptomes from the pennaceous and plumulaceous portions of body feather, the distal and proximal parts of flight feather, and the flight feather calamus (Additional file 1: Figure S1). Each sample was a pool of mRNA from two-three follicle epidermis of one individual. Three cDNA libraries with insert lengths ranging from 300 to 400 bp were constructed for each sample. The RNA-seq data had been used previously to study the expression pattern of α-and β-keratin genes . In this study we conducted a detailed analysis of the expression patterns for all genes (Additional file 2: Table S1).
To validate the expression data obtained by RNA-Seq, ten genes were selected from the differentially expressed protein-coding genes to perform real time quantitative PCR (RT-qPCR) assays. Reproducibility of the data was confirmed by the strong correlation between the values of gene expression obtained by RNA-Seq and RT-qPCR (Additional file 1: Figure S2). When the pattern of gene expression levels was compared, strong correlations (R 2) ranging from 0.833 to 0.998 between RT-qPCR and RNA-Seq platforms were observed for 90 % of the expressed genes exception for one sample (with R 2 = 0.618), confirming the high reproducibility of the data.
For those genes with a FPKM > 0.1 (FPKM = Fragments Per Kilobase of transcript per Million mapped reads), there were 12,608 genes expressed in all three biological replicates of at least one type of the feather epithelium in the total transcriptomes. Principal components analysis (PCA) showed that samples from the same group clustered together (Fig. 1b). This observation suggests that each sample harbors transcriptomic features that are unique to the feather types or regions. The genes identified in the GO enrichment analysis of the most abundant transcripts in these samples are involved in protein translation, reflecting the rapid production of a protein-made structure (Fig. 2). These results are consistent with the efficient biosynthesis of proteins in the feather follicles. Genes involved in developmental morphogenesis as well as cytoskeletons are also highly abundant in developing feather epithelia.
An aggregate total of 13,973 expressed genes was expressed (FPKM > 0.1) in at least one of the 15 samples (Additional file 4: Table S3), among which 9,638 genes were expressed in all 15 samples (Fig. 3). In general, the flight feather has more specifically expressed genes than the body feather. The GO enrichment analysis showed that the specifically expressed genes of the body feather are not significantly enriched for any known functions, whereas those of the flight feather are significantly enriched for several functions (Fig. 3). We also found that 49 genes previously identified to have undergone rapid evolution and/or positive selection in avian lineages  are expressed in all feather samples (Additional file 5: Table S4). Most of these rapidly evolving and/or positively selected genes are enriched for cytoskeleton and cell adhesion. These proteins may have evolved new functions or properties in feathers and it is worth further investigation.
Transcriptomic comparison between distal and proximal body feathers
The distal end of a feather forms earlier than the proximal end and the structure and shape of a body feather change dynamically from the distal end to the proximal end. The distal end of a body feather is mainly pennaceous while the barbs of the proximal end become plumulaceous (Additional file 1: Figure S1). Among the 957 differentially expressed protein-coding genes (DEGs), 223 were up-regulated, while 734 were down-regulated in the plumulaceous portion in comparison to the pennaceous portion (Fig. 4a, Additional file 6: Table S5). IPA (Ingenuity Pathway Analysis) canonical pathway analysis showed that several genes involved in semaphorin signaling in neurons (PLXNA1, NRP1, DPYSL3, MAPK1, CDK5) were differentially expressed between pennaceous and plumulaceous body feathers (Fig. 5, Additional file 10: Table S9). Semaphorin signaling is known to play an important role in intersomitic vessels, lung, and kidney branching morphogenesis [30, 31], but has not been reported to play any role in feather morphogenesis. A bone morphogenetic protein, BMP2, was predicted to be the upstream regulator for gene expression differences (Table 1). Level of BMP activity has been shown to determine barb ridge branching morphogenesis . Interactions between activators and inhibitors involving Sonic hedgehog (SHH) and BMP2 have been suggested to be involved in the formation of barb ridges in feathers .
Among the up-regulated genes in the pennaceous portion compared to the plumulaceous portion of the body feather, the expression of genes involved in protein translation, oxidative phosphorylation, pyrimidine metabolism, ribosomal subunit, contractile fiber part, Peptidyl-prolyl cis-trans isomerase activity, inorganic cation transmembrane transporter activity was increased significantly (Table 2). The up-regulation of genes in these groups may be due to the need of large quantities of protein synthesis in the pennaceous portion of the body feather.
Genes involved in amino acid transmembrane transporter activity, lysosome, tube development, enzyme linked receptor protein signaling pathway, carboxylic acid biosynthetic process, and plasma membrane were increased significantly when the growth of a body feather turns to the plumulaceous portion. Six genes involved in tube development (SP3, NRP1, EPAS1, TP63, BMPR-II, BMPR1A) and six genes involved in enzymes linked receptor protein signaling pathway (SMAD5, MADH2, BMPR1A, NRP1, BMPR-II, EGFR) were up-regulated. Among these genes, TP63, a transcription factor of the p53 family, is known to be essential for the development of epidermis and its derivatives in vertebrates [33, 34]. In situ hybridization studies in chickens have shown that TP63 is highly expressed in the apical ectodermal ridge (AER) of the limb buds, interdigital tissues, epithelium of branchial arches, and feather buds . Two receptors of BMPs were differentially expressed. Signaling via BMPRIA and BMPRIB is required to regulate intramembranous bone formation, chondrogenesis, and feather formation in chicken embryos . The antagonistic balance between noggin and BMP4 has been shown to play a critical role in feather branching, with BMP4 promoting rachis formation and barb fusion, and noggin enhancing rachis and barb branching . Epidermal growth factor (EGF) signaling is known to be required to pattern the feather array by promoting the interbud development .
Transcriptomic comparison between pennaceous body and flight feathers
Among the 1,287 DEGs between pennaceous body and flight feathers, 988 were up-regulated and 299 genes were down-regulated in the pennaceous body feather (Fig. 4b, Additional file 7: Table S6). IPA canonical pathway analysis showed that these DEGs included several genes involved in the Sertoli cell-Sertoli cell junction signaling (PVRL1, TJP1, TUBB3, CLDN3, CLDN4, CGN, TUBB4B, RAB8B, PRKAR1A, MAP3K1, TUBA1B, MAPK1, TUBA4A, TJP3) and in the regulation of the epithelial-mesenchymal transition pathway (FGFR2, FGFR3, EGFR, PDGFD, PARD6A, JAG1, CLDN3, SMAD2, FZD2, JAK1, PYGO2, HIF1A, MAPK1, WNT6, HMGA2) (Fig. 5, Additional file 11: Table S10).
Among the up-regulated genes in the pennaceous portion of the body feather, the expression of genes involved in translation, oxidative phosphorylation, pyrimidine metabolism, ribosomal subunit, ubiquinol-cytochrome-c reductase activity, and RNA polymerase was increased significantly (Table 2). The up-regulation of genes in these groups may be due to the need of large quantities of protein synthesis in the pennaceous portion of the body feather.
Genes involved in skeletal system development, tube development, melanin biosynthetic process, regulation of RNA metabolic process, drug metabolism, respiratory tube development, and lytic vacuole were expressed significantly higher in early-grow flight feather (Table 2). The genes in melanin biosynthetic process (TYR, TYRP1, PMEL) were up-regulated simply because the color is usually darker in the flight feather of this breed than in the contour feather of both TCC_L2 and white leghorn chickens. Tyrosinase (TYR) and tyrosinase-related protein 1 (TYRP1) are known to be involved in the feather pigment pattern formation . TYR and TYRP1 are found to be associated with melanic plumage color differences in chickens, Korean quails (Coturnix coturnix), ducks, geese, and pigeons [39–43]. The mutation of premelanosome protein (PEML) can cause hypopigmentation in chickens .
Several genes involved in tissue morphogenesis (PRKAR1A, LMO4, TP63, TWSG1, JAG1, FGFR2, ACVR1, CA2) and in the regulation of developmental process (SMAD5, JAG1, FGFR2, HIF1A, ACVR1, WNT7B) were upregulated in distal flight feather (Table 2). Jagged-1 (JAG1), a Notch ligand, is involved in the orientation of feather bud elongation . Three fibroblast growth factor (FGF) receptor genes, FGFR1, FGFR2 and FGFR3, have been suggested to be involved in feather morphogenesis .
Transcriptomic comparison between distal and proximal parts of flight feather
Toward the proximal end, the rachis gradually increases in width and eventually turns into the calamus. Among the 1,167 DEGs, 534 genes were up-regulated and 633 genes were down-regulated in the proximal flight feather in comparison to the distal flight feather (Fig. 4c, Additional file 8: Table S7). IPA canonical pathway analysis revealed several DEGs that were involved in WNT/β-catenin Signaling (JUN, DKK3, WNT5A, CDH3, CD44, CDH5, PPARD, WNT2B, SFRP4, SOX14, SFRP2, PPP2R2B, SOX7, FRZB, RARB) and in TGF-β Signaling (JUN, INHBA, SMAD2, SMAD1, RUNX3, MAPK11, INHBB) (Fig. 5, Additional file 12: Table S11). Other pathways basically overlap with the WNT/β-catenin signaling or the TGF-β signaling pathway.
Nine genes involved in skeletal system development (MGP, GLI1, WWOX, IGFBP5, GJA5, SHH, SOX14, SMAD1, CBFB) were increased in expression in the early-grow (distal) flight feather (Table 2). BMP4 and matrix gla protein (MGP) are considered an activating and an inhibitory morphogen, respectively, and their interaction is important for vascular branching . MGP may promote rachis and barb branching in feather. SHH is a secreted protein expressed in the epidermis that is involved in the mitogenic and morphogenetic processes throughout feather development [32, 48–52]. The interactions between SHH and BMP2 signaling during feather barb ridge morphogenesis may be critical for the initial formation of a meristic pattern of barb ridges and the variation in barb morphogenesis in feathers . The activation of the SHH signaling pathway leads to the expression of the transcription factor glioma-associated oncogene 1 (GLI1), a SHH-targeted mediator . Insulin-like growth factor binding protein 5 (IGFBP5) is expressed in human hair follicle dermal papilla and plays a specific role in the local modulation of IGF action during the hair growth cycle .
Nine genes involved in tube development (WNT5A, LMO4, LIPA, SHH, EDNRA, RARB, GJA5, CRH, GLI1) were increased in expression in the early-growth flight feather (Table 2). The expression levels of WNT ligands such as WNT5A/WNT5B/WNT6 were found to be high in the feather epithelium and pulp compared to dermal papillae . WNT5A is involved in non-canonical pathways but its downstream signaling events are not known yet. LIM domain-only protein 4 (LMO4) is expressed in mouse hair follicles, especially in the sebaceous glands, undifferentiated bulb cells, and the outer epithelial root sheath . Retinoic acid receptor beta (RARB) is a receptor of retinoic acid which regulates cell proliferation, differentiation, and morphogenesis and is involved in the feather-bud formation . Gap junction alpha-5 protein (GJA5), also known as connexin 40 (CX40), is an integral membrane protein that oligomerizes to form intercellular channels that are clustered at gap junctions which are present in supportive cells located in the vicinity of barbule cells . Corticotropin-releasing hormone (CRH) peptides modulate human hair growth/cycling [60, 61].
Many genes involved in extracellular region and cell adhesion were up-regulated significantly in the middle-grow flight feather (Table 2). Cell adhesion molecules (CAMs) may regulate feather morphogenesis by constraining cell motion and forming borders. Several adhesion molecules, including L-CAM, N-CAM, integrin, tenascin, as well as proteoglycan, are involved in feather development [62–64]. Tenascin-C has been shown to evolve rapidly in avian lineages .
Many collagen genes were up-regulated in the middle-grow feather portion compared to the early-growth flight feather portion (Table 2). The orientation of collagen fibers in the feather buds may promote feather growth by creating a gradient of stiffness, thus triggering the pressure sensitive growth factors . Collagen types I and III, and fibronectin are known to be involved in feather morphogenesis in the chick embryo . Matrix metalloproteinases (MMPs) and their inhibitors are important in tissue development remodeling for the formation of feather follicles such as epithelium invagination and mesenchymal cell proliferation . Several collagens and a MMP expressed in feathers have been found to evolve rapidly in a previous study . Dickkopf-related protein 2 (DKK2), which presumably encodes a WNT signaling inhibitor, regulates feather regeneration in the dermal papillae . The expression of CD44, which is also known to evolve rapidly in avian lineages , correlates with the onset of keratinocyte stratification and mesenchymal maturation into fibrous dermis in fetal human skin . Tissue inhibitor of metalloproteinases-3 (TIMP3) is expressed in epithelial outer root sheath cells of growing hair follicles of human fetus . Cysteine-rich secretory protein 1 (CRISP1) is expressed in murine hair follicles and down-regulated in mice overexpressing a homeobox gene HOXC13 .
State-dependent signaling by Cav1.2 regulates hair follicle stem cell function by regulating the production of the bulge-derived BMP inhibitor follistatin-like1 (FSTL1), derepressing stem cell quiescence [71, 72]. Expression of muscle-related genes are known to be enriched in the feather dermal papilla, including ACTG2 (smooth muscle actin, gamma 2), ACTA2 (smooth muscle actin, alpha 2), Desmin, MYH11 (myosin heavy chain11), MYL4 (myosin light chain4), MYL9 (myosin light chain 9), MYLK (myosin light chain kinase), etc. . Our results showed that genes involved in smooth muscle contraction, such as ADRA1A, ITPR3, PLA2G4A, ACTG2, ARHGEF12, PLCB4, PLA2G12A, ITPR2, PLA2G10, RAMP2, and PLA2G4C are differentially expressed. CLR/RAMP2-overexpressing mice revealed a defined phenotype with thinning of the hair during postnatal development .
Transcriptomic comparison between proximal flight feather and calamus
Among the 702 DEGs, 263 genes were up-regulated and 404 genes were down-regulated in the proximal flight feather in comparison to the calamus (Fig. 4d, Additional file 9: Table S8). IPA canonical pathway analysis showed that several genes involved in the TGF-β signaling (INHBA, RUNX3, PMEPA1, RUNX2, INHBB), the Sertoli cell-Sertoli cell junction signaling (TUBB3, TUBA1B, CLDN4, TJP3, JAM3, ACTN1, MTMR2) and the germ cell-Sertoli cell junction (CDH2, TUBB3, GSN, TUBA1B, ACTN1, MTMR2) signaling were differentially expressed. Other pathways basically overlap with the TGF-β signaling pathway (Fig. 5, Additional file 13: Table S12).
Compared to the close proximal part of the flight feather, the calamus expressed significantly more genes involved in extracellular matrix and cell adhesion (Table 2). The calamus can basically be seen as the rachis of the flight feather. The genes involved in extracellular matrix and cell adhesion may be required for making a tougher feather structure.
The molecular mechanism of feather branching morphogenesis
Major signaling pathways are involved in feather branching morphogenesis, including the Wnt/β-catenin, SHH/BMP, and Notch pathways [11, 18, 19, 45, 48, 56, 74, 75]. Besides feathers, epithelial tissues such as the vascular system, kidney, lung, and mammary gland arise through branching morphogenesis of a pre-existing epithelial structure [13, 76–78]. Common morphological stages and a similar set of developmental regulations are shared by these tissues. The spatial and temporal controls of branching are controlled by developmental decisions requiring regulation of cell proliferation, apoptosis, invasiveness, and cell motility. Similar molecular mechanisms could exist for the epithelial branching program, even though the feather is an evolutionary novel tissue. Key branching morphogenetic molecules include central signaling molecules such as BMPs, TGF-β, FGF, and MMPs [13, 76–78]. Our study supports the previous findings that temporal and spatial variation of BMP signals are critical for generating branching differences between pennaceous and plumulaceous body feathers because genes involved in BMP signaling were significantly upregulated in the plumulaceous portions.
Genes involved in axon guidance (MYL4, CDK5, SEMA4B, PRKAR1A, NFATC3, PLXNA1, PLXNB2, ARPC4, NRP1, MAPK1, WNT6, PRKCI, RASA1, ECE2, SEMA5A) from IPA pathway analysis (Fig. 5), especially those in semaphorin signaling, were differentially expressed between the pennaceous and plumulaceous portions of body feather. This observation suggests that they are recruited in feather development and play a critical role in controlling the morphological differences in feathers, and perhaps are involved in changing the extracellular environment for signals that instruct the cell of the barb plate which direction to grow by affecting the cytoskeleton. The differential expression profile of these genes among different feather types suggests that they are involved in critical guidance cues during feather morphogenesis, although functional studies remain to be demonstrated.
Genes involved in Sertoli cell-Sertoli cell junction signaling and germ cell-Sertoli cell junction are recruited in feather development. Differentiating barb/barbule cells have been found to have many adhesion junctions, some gap junctions and fewer tight junctions during early stages of feather development . The cytological details on the type of cell junctions present in barb/barbules of feathers are poorly known. Our study provides the data to characterize the types of cell junctions, and their molecular nature that are critical in feather morphogenesis.
Although the hair and the feather are not homologous, they share many pathways. Hair follicle morphogenesis also depends on WNT, SHH, NOTCH, BMP and other signaling pathways that interact between epithelial and mesenchymal cells. However, as hairs have no branching structures, the genes that are involved in feather morphogenesis but not in hair development may participate in generating branching structures. A detailed comparison of transcriptomes between feathers and hairs may reveal the molecular mechanism shared and distinct between these two types of keratinized skin appendages. The genes involved in vessel and tube development are differentially expressed in feather epithelium, suggesting a role in regulating the morphology of feather branching. We found that genes involved in developing the vessel, tube and kidney were enriched. They may have been co-opted to develop an evolutionary novelty. The origin and diversification of a novel structure may not require the evolution of new gene or gene duplication because existing genes can be recruited to have new expression pattern and regulation. The genes with similar functions can be reused to construct a new network.
This study has significantly increased our understanding of the expression profiles of feather related genes. We examined the expression profiles of genes associated with the development of feather structure and compared the gene expression patterns in different types of feathers and different portions of a feather to advance our understanding of the molecular mechanisms of feather growth and the molecular basis of variation in feather structure. Our results are a valuable resource for understanding the molecular mechanisms of avian feather development. This study produced abundant data for the analysis of gene expression during feather morphogenesis. Morphotype-specifically expressed genes were identified from five zones of feather filament epithelia. Some identified genes may be associated with the growth control during feather regeneration, the formation of special branching structures, or barb differentiation themselves. The present study provides a basis for future study of the complex molecular and cellular events during feather development.
All the animals used in this study were processed following the approved protocol of Institutional Animal Care and Use Committees of the National Chung Hsing University (Taichung, Taiwan). For total RNA extraction, we used the Taiwan County Chicken (TCC_L2) breed chicken for wing flight feathers and white leghorn for body contour feathers. TCC_L2 and white leghorn chicken contour feathers are different in color but highly similar in morphology and structure.
Total RNA Isolation and RNA-seq
We collected regenerating pennaceous and plumulaceous portions of body contour feathers, the distal and proximal portions of primary flight feathers, and the calamus of primary flight feathers. Total RNA was isolated from early or late grow fresh feather epithelial tissues corresponding, respectively, to the distal and proximal part of a feather (Additional file 1: Figure S3), which was dissected from the follicle tissue and separated from the mesenchyme in Calcium-Magnesium Free Saline (CMFS 2X) on ice . White leghorn chickens were used for body contour feathers to avoid melanin contamination, which is difficult to remove and can inhibit essential enzymatic reactions for RNA-seq [80, 81]. RNA-seq and analysis of paired-end reads were performed as described in Ng et al. 2014 . Reads were mapped onto the chicken genome assembly ICGSC Gallus_gallus-4.0 (GCA_000002315.2).
Validation by real-time quantitative PCR
A total of 2 μg RNA of each sample was reverse transcribed with MultiScribe Reverse Transcriptase (Thermo Fisher Scientific, Waltham, MA) into cDNA for both Reverse Transcription PCR (RT-PCR) and Quantitative Reverse Transcription PCR (qRT-PCR) reactions. Total RNA was incubated with RT enzymes at 25 °C for 10 min prior to the RT reaction. RT reactions were performed at 37 °C for 2 h followed by the inactivation of RT enzyme at 85 °C for 10 s. For RT-PCR, 1 μl of 10× diluted cDNA was amplified by Fast Start Taq DNA polymerase (Roche Applied Science, Penzberg, Germany) in a total of 10 ul reaction. For RT-qPCR, 1 μl of 10× diluted cDNA products was quantified with 2 × SYBR Green Master Mix (Kapa Biosystems, Wilmington, MA) in a total of 10 ul reaction and performed on a Roche LightCycler 480 Instrument II. All the data were analyzed by the HTC1 software (Roche Applied Science). The 2-ΔΔCt method was used to calculate relative expression levels . The cycling parameters of RT-qPCR were as follows: 95 °C for 3 min, then 40 cycles of 95 °C for 10 s, and annealing for 20 s. Gene names and primer sequences are shown in Additional file 3: Table S2. Each sample was analyzed in duplicates, and gene expression levels were normalized against the corresponding TATA-binding protein (TBP) expression level.
Prior to statistical analyses performed with R v2.15.3 (R Development Core Team 2011), raw read counts were normalized by Fragments Per Kilobase of transcript per Million mapped reads (FPKM). Principal component analysis (PCA) was performed on the covariance matrix f using a custom R script based on the “prcomp” R package.
Identification of differentially expressed genes
We used the non-parametric method to identify differentially expressed genes (DEGs) between two samples . Here, we set the q value (differentially expression probability) in the method to be 0.75 (this value is equivalent to an odd of 3:1, i.e., the gene is three times more likely to be differentially expressed than it is not) and require at least a 2-fold change in FPKM between the two samples.
Gene functional annotation, canonical pathway and upstream regulator analyses
Functional annotations of gene loci were compared with the complete genome using annotations from the Database for Annotation, Visualization, and Discovery (DAVID), which uses fuzzy clustering to group genes into functionally related classes based on the similarity of their annotations [84, 85]. Pathway analyses of differentially-expressed genes were carried out using the Ingenuity Pathways Analysis software (IPA; Ingenuity Systems, www.ingenuity.com). Each gene identifier was mapped to its corresponding gene object in the Ingenuity Pathways Knowledge Base. A canonical pathways analysis was generated to identify the pathways from the IPA library that were most significant. Fischer’s exact test was employed to calculate the p-value which determines the probability that each biological function or/and canonical pathway is due to chance alone. The Upstream Analysis section of the core analysis was used to determine which upstream regulators were associated with the observed differently expressed genes.
Availability of supporting data
The full data sets have been submitted to NCBI Sequence Read Archive (SRA) under accession nos. SRX528281, SRX528834, SRX528843, SRX529337, SRX529339, SRX529353-SRX529362. Bioproject: PRJNA245063.
Fragments Per Kilobase of transcript per Million mapped reads
Real time quantitative PCR
Differentially expressed genes
Principal components analysis
Database for Annotation, Visualization, and Discovery
Ingenuity pathways analysis
Losos JB, Arnold SJ, Bejerano G, Brodie 3rd ED, Hibbett D, Hoekstra HE, et al. Evolutionary biology for the 21st century. PLoS Biol. 2013;11(1):e1001466.
Carroll SB. Evo-devo and an expanding evolutionary synthesis: a genetic theory of morphological evolution. Cell. 2008;134(1):25–36.
Widelitz RB, Jiang TX, Yu M, Shen T, Shen JY, Wu P, et al. Molecular biology of feather morphogenesis: a testable model for evo-devo research. J Exp Zool B Mol Dev Evol. 2003;298(1):109–22.
Prum RO, Brush AH. The evolutionary origin and diversification of feathers. Q Rev Biol. 2002;77(3):261–95.
Chen CF, Foley J, Tang PC, Li A, Jiang TX, Wu P, et al. Development, regeneration, and evolution of feathers. Annu Rev Anim Biosci. 2015;3:169–95.
Lucas AM, Stettenheim PR. Avian anatomy. Integument. In: Agriculture Handbook 362. Washington, D.C: U.S. Dept. Agric; 1972.
Sawyer RH, Knapp LW. Avian skin development and the evolutionary origin of feathers. J Exp Zool B Mol Dev Evol. 2003;298(1):57–72.
Prum RO. Development and evolutionary origin of feathers. J Exp Zool. 1999;285(4):291–306.
Chuong CM, Chodankar R, Widelitz RB, Jiang TX. Evo-devo of feathers and scales: building complex epithelial appendages. Curr Opin Genes Dev. 2000;10(4):449–56.
Xu X, Zhou Z, Dudley R, Mackem S, Chuong CM, Erickson GM, et al. An integrative approach to understanding bird origins. Science. 2014;346(6215):1253293.
Lin SJ, Wideliz RB, Yue Z, Li A, Wu X, Jiang TX, et al. Feather regeneration as a model for organogenesis. Develop Growth Differ. 2013;55(1):139–48.
Wu P, Hou L, Plikus M, Hughes M, Scehnet J, Suksaweang S, et al. Evo-Devo of amniote integuments and appendages. Int J Dev Biol. 2004;48(2–3):249–70.
Chuong CM, Bhat R, Widelitz RB, Bissell MJ. SnapShot: Branching Morphogenesis. Cell. 2014;158(5):1212–2. e1211.
Prum RO. Evolution of the morphological innovations of feathers. J Exp Zool B Mol Dev Evol. 2005;304(6):570–9.
True JR, Carroll SB. Gene co-option in physiological and morphological evolution. Annu Rev Cell Dev Biol. 2002;18:53–80.
Lillie FR, Juhn M. Physiology of development of the feather. II. General principles of development with special reference to the after-feather. Physiol Zool. 1938;11:434–48.
Yue Z, Jiang TX, Widelitz RB, Chuong CM. Mapping stem cell activities in the feather follicle. Nature. 2005;438(7070):1026–9.
Yu M, Yue Z, Wu P, Wu DY, Mayer JA, Medina M, et al. The biology of feather follicles. Int J Dev Biol. 2004;48(2–3):181–91.
Yu M, Wu P, Widelitz RB, Chuong CM. The morphogenesis of feathers. Nature. 2002;420(6913):308–12.
Lillie FR. Physiology of development of the feather. III. Growth of the mesodermal constituents and blood circulation in the pulp. Physiol Zool. 1940;13:143–75.
Mutz KO, Heilkenbrinker A, Lonne M, Walter JG, Stahl F. Transcriptome analysis using next-generation sequencing. Curr Opin Biotechnol. 2013;24(1):22–30.
Costa V, Angelini C, De Feis I, Ciccodicola A. Uncovering the complexity of transcriptomes with RNA-Seq. J Biomed Biotechnol. 2010;2010:853916.
Wang Z, Gerstein M, Snyder M. RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009;10(1):57–63.
Oshlack A, Robinson MD, Young MD. From RNA-seq reads to differential expression results. Genome Biol. 2010;11(12):220.
Ozsolak F, Milos PM. RNA sequencing: advances, challenges and opportunities. Nat Rev Genet. 2011;12(2):87–98.
Kassahn KS, Waddell N, Grimmond SM. Sequencing transcriptomes in toto. Integr Biol (Camb). 2011;3(5):522–8.
Chang KW, Huang NA, Liu IH, Wang YH, Wu P, Tseng YT, et al. Emergence of differentially regulated pathways associated with the development of regional specificity in chicken skin. BMC Genomics. 2015;16:22.
Ng CS, Wu P, Fan WL, Yan J, Chen CK, Lai YT, et al. Genomic organization, transcriptomic analysis, and functional characterization of avian alpha- and beta-keratins in diverse feather forms. Genome Biol Evol. 2014;6(9):2258–73.
Nam K, Mugal C, Nabholz B, Schielzeth H, Wolf JB, Backstrom N, et al. Molecular evolution of genes in avian genomes. Genome Biol. 2010;11(6):R68.
Reidy K, Tufro A. Semaphorins in kidney development and disease: modulators of ureteric bud branching, vascular morphogenesis, and podocyte-endothelial crosstalk. Pediatr Nephrol. 2011;26(9):1407–12.
Horowitz A, Simons M. Branching morphogenesis. Circ Res. 2008;103(8):784–95.
Harris MP, Fallon JF, Prum RO. Shh-Bmp2 signaling module and the evolutionary origin and diversification of feathers. J Exp Zool. 2002;294(2):160–76.
Mikkola ML. p63 in skin appendage development. Cell Cycle. 2007;6(3):285–90.
Koster MI, Roop DR. p63 and epithelial appendage development. Differentiation. 2004;72(8):364–70.
Yasue A, Tao H, Nohno T, Moriyama K, Noji S, Ohuchi H. Cloning and expression of the chick p63 gene. Mech Dev. 2001;100(1):105–8.
Ashique AM, Fu K, Richman JM. Signalling via type IA and type IB bone morphogenetic protein receptors (BMPR) regulates intramembranous bone formation, chondrogenesis and feather formation in the chicken embryo. Int J Dev Biol. 2002;46(2):243–53.
Atit R, Conlon RA, Niswander L. EGF signaling patterns the feather array by promoting the interbud fate. Dev Cell. 2003;4(2):231–40.
Lin SJ, Foley J, Jiang TX, Yeh CY, Wu P, Foley A, et al. Topology of feather melanocyte progenitor niche allows complex pigment patterns to emerge. Science. 2013;340(6139):1442–5.
Li S, Wang C, Yu W, Zhao S, Gong Y. Identification of genes related to white and black plumage formation by RNA-Seq from white and black feather bulbs in ducks. PLoS One. 2012;7(5):e36592.
Liu WB, Chen SR, Zheng JX, Qu LJ, Xu GY, Yang N. Developmental phenotypic-genotypic associations of tyrosinase and melanocortin 1 receptor genes with changing profiles in chicken plumage pigmentation. Poult Sci. 2010;89(6):1110–4.
Xu Y, Zhang XH, Pang YZ. Association of Tyrosinase (TYR) and Tyrosinase-related Protein 1 (TYRP1) with Melanic Plumage Color in Korean Quails (Coturnix coturnix). Asian-Australas J Anim Sci. 2013;26(11):1518–22.
Domyan ET, Guernsey MW, Kronenberg Z, Krishnan S, Boissy RE, Vickrey AI, et al. Epistatic and combinatorial effects of pigmentary gene mutations in the domestic pigeon. Curr Biol. 2014;24(4):459–64.
Wang Y, Li SM, Huang J, Chen SY, Liu YP. Mutations of TYR and MITF Genes are Associated with Plumage Colour Phenotypes in Geese. Asian-Australasian J Anim Sci. 2014;27(6):778–83.
Kerje S, Sharma P, Gunnarsson U, Kim H, Bagchi S, Fredriksson R, et al. The Dominant white, Dun and Smoky color variants in chicken are associated with insertion/deletion polymorphisms in the PMEL17 gene. Genetics. 2004;168(3):1507–18.
Li A, Chen M, Jiang TX, Wu P, Nie Q, Widelitz R, et al. Shaping organs by a wingless-int/Notch/nonmuscle myosin module which orients feather bud elongation. Proc Natl Acad Sci U S A. 2013;110(16):E1452–1461.
Noji S, Koyama E, Myokai F, Nohno T, Ohuchi H, Nishikawa K, et al. Differential expression of three chick FGF receptor genes, FGFR1, FGFR2 and FGFR3, in limb and feather development. Prog Clin Biol Res. 1993;383B:645–54.
Yao Y, Nowak S, Yochelis A, Garfinkel A, Bostrom KI. Matrix GLA protein, an inhibitory morphogen in pulmonary vascular development. J Biol Chem. 2007;282(41):30131–42.
Lin CM, Jiang TX, Widelitz RB, Chuong CM. Molecular signaling in feather morphogenesis. Curr Opin Cell Biol. 2006;18(6):730–41.
Ting-Berreth SA, Chuong CM. Sonic Hedgehog in feather morphogenesis: induction of mesenchymal condensation and association with cell death. Dev Dyn. 1996;207(2):157–70.
Morgan BA, Orkin RW, Noramly S, Perez A. Stage-specific effects of sonic hedgehog expression in the epidermis. Dev Biol. 1998;201(1):1–12.
McKinnell IW, Makarenkova H, de Curtis I, Turmaine M, Patel K. EphA4, RhoB and the molecular development of feather buds are maintained by the integrity of the actin cytoskeleton. Dev Biol. 2004;270(1):94–105.
Nohno T, Kawakami Y, Ohuchi H, Fujiwara A, Yoshioka H, Noji S. Involvement of the Sonic hedgehog gene in chick feather formation. Biochem Biophys Res Commun. 1995;206(1):33–9.
Harris MP, Williamson S, Fallon JF, Meinhardt H, Prum RO. Molecular evidence for an activator-inhibitor mechanism in development of embryonic feather branching. Proc Natl Acad Sci U S A. 2005;102(33):11734–9.
Fisher CE, Howie SEM. Shh and Gli Signalling in Development. New York, NY: Springer; 2007.
Batch JA, Mercuri FA, Werther GA. Identification and localization of insulin-like growth factor-binding protein (IGFBP) messenger RNAs in human hair follicle dermal papilla. J Invest Dermatol. 1996;106(3):471–5.
Chu Q, Cai L, Fu Y, Chen X, Yan Z, Lin X, et al. Dkk2/Frzb in the dermal papillae regulates feather regeneration. Dev Biol. 2014;387(2):167–78.
Sum EY, O'Reilly LA, Jonas N, Lindeman GJ, Visvader JE. The LIM domain protein Lmo4 is highly expressed in proliferating mouse epithelial tissues. J Histochem Cytochem. 2005;53(4):475–86.
Obinata A, Akimoto Y. Effects of retinoic acid and Gbx1 on feather-bud formation and epidermal transdifferentiation in chick embryonic cultured dorsal skin. Dev Dyn. 2012;241(9):1405–12.
Alibardi L. Gap and tight junctions in the formation of feather branches: A descriptive ultrastructural study. Ann Anat. 2010;192(4):251–8.
Ito N, Ito T, Kromminga A, Bettermann A, Takigawa M, Kees F, et al. Human hair follicles display a functional equivalent of the hypothalamic-pituitary-adrenal axis and synthesize cortisol. FASEB J. 2005;19(10):1332–4.
Kauser S, Slominski A, Wei ET, Tobin DJ. Modulation of the human hair follicle pigmentary unit by corticotropin-releasing hormone and urocortin peptides. FASEB J. 2006;20(7):882–95.
Chuong CM, Chen HM, Jiang TX, Chia J. Adhesion molecules in skin development: morphogenesis of feather and hair. Ann N Y Acad Sci. 1991;642:263–80.
Chuong CM, Edelman GM. Expression of cell-adhesion molecules in embryonic induction. II. Morphogenesis of adult feathers. J Cell Biol. 1985;101(3):1027–43.
Chuong CM, Edelman GM. Expression of cell-adhesion molecules in embryonic induction. I. Morphogenesis of nestling feathers. J Cell Biol. 1985;101(3):1009–26.
Paetzold H, Goepfert C, Huber G, Hoenig E, Portner R, Schilling AF, et al. The development of the collagen fibre network in tissue-engineered cartilage constructs in vivo. Engineered cartilage reorganises fibre network. Eur Cell Mater. 2012;23:209–21.
Mauger A, Demarchez M, Herbage D, Grimaud JA, Druguet M, Hartmann D, et al. Immunofluorescent localization of collagen types I and III, and of fibronectin during feather morphogenesis in the chick embryo. Dev Biol. 1982;94(1):93–105.
Jiang TX, Tuan TL, Wu P, Widelitz RB, Chuong CM. From buds to follicles: matrix metalloproteinases in developmental tissue remodeling during feather morphogenesis. Differentiation. 2011;81(5):307–14.
Tuhkanen AL, Agren UM, Tammi MI, Tammi RH. CD44 expression marks the onset of keratinocyte stratification and mesenchymal maturation into fibrous dermis in fetal human skin. J Histochem Cytochem. 1999;47(12):1617–24.
Airola K, Ahonen M, Johansson N, Heikkila P, Kere J, Kahari VM, et al. Human TIMP-3 is expressed during fetal development, hair growth cycle, and cancer progression. J Histochem Cytochem. 1998;46(4):437–47.
Peterson RL, Tkatchenko TV, Pruett ND, Potter CS, Jacobs DF, Awgulewitsch A. Epididymal cysteine-rich secretory protein 1 encoding gene is expressed in murine hair follicles and downregulated in mice overexpressing Hoxc13. J Investig Dermatol Symp Proc. 2005;10(3):238–42.
Yucel G, Altindag B, Gomez-Ospina N, Rana A, Panagiotakos G, Lara MF, et al. State-dependent signaling by Cav1.2 regulates hair follicle stem cell function. Genes Dev. 2013;27(11):1217–22.
Coulombe PA, Caterina MJ. The incidental pore: CaV1.2 and stem cell activation in quiescent hair follicles. Genes Dev. 2013;27(12):1315–7.
Ittner LM, Husmann K, Muff R, Gotz J, Bonneick S, Suter U, et al. Defective hair formation in calcitonin-like receptor transgenic mice. Exp Dermatol. 2004;13(9):582–2.
Yue Z, Jiang TX, Wu P, Widelitz RB, Chuong CM. Sprouty/FGF signaling regulates the proximal-distal feather morphology and the size of dermal papillae. Dev Biol. 2012;372(1):45–54.
Suksaweang S, Jiang TX, Roybal P, Chuong CM, Widelitz R. Roles of EphB3/ephrin-B1 in feather morphogenesis. Int J Dev Biol. 2012;56(9):719–28.
Pohl M, Stuart RO, Sakurai H, Nigam SK. Branching morphogenesis during kidney development. Annu Rev Physiol. 2000;62:595–620.
Iber D, Menshykau D. The control of branching morphogenesis. Open Biol. 2013;3(9):130088.
Metzger RJ, Krasnow MA. Genetic control of branching morphogenesis. Science. 1999;284(5420):1635–9.
Chuong CM. Skin morphogenesis. Embryonic chicken skin explant cultures. Methods Mol Biol. 2000;136:101–6.
Eckhart L, Bach J, Ban J, Tschachler E. Melanin binds reversibly to thermostable DNA polymerase and inhibits its activity. Biochem Biophys Res Commun. 2000;271(3):726–30.
Giambernardi TA, Rodeck U, Klebe RJ. Bovine serum albumin reverses inhibition of RT-PCR by melanin. BioTechniques. 1998;25(4):564–6.
Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(−Delta Delta C(T)) Method. Methods. 2001;25(4):402–8.
Tarazona S, Garcia-Alcalde F, Dopazo J, Ferrer A, Conesa A. Differential expression in RNA-seq: a matter of depth. Genome Res. 2011;21(12):2213–23.
Huang DW, Sherman BT, Lempicki RA. Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists. Nucleic Acids Res. 2009;37(1):1–13.
Huang DW, Sherman BT, Tan Q, Collins JR, Alvord WG, Roayaei J, et al. The DAVID Gene Functional Classification Tool: a novel biological module-centric algorithm to functionally analyze large gene lists. Genome Biol. 2007;8(9):R183.
The authors thank Drs. Meng-Shin Shiao, Meng-Ru Ho, Hsu-Chen Cheng, Pin-Chi Tang, Chua-Ti Ting, and Tzi-Yuan Wang for help. We thank two anonymous reviewers for their constructive comments. This study was supported by Ministry of Science and Technology (formerly known as National Science Council), Taiwan (MOST 104-2621-B-001-003-MY3 and 99-2321-B-001-041-MY2) and Academia Sinica, Taiwan. It was also supported by the National Institute of Arthritis and Musculoskeletal and Skin Diseases of the National Institutes of Health under Award Number AR 47364 and 60306. C.S.N. was supported by a postdoctoral fellowship from Academia Sinica, Taiwan. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.
The authors declare that they have no competing interests.
CNS designed experiments. CNS, CKC, WLF analyzed data. SMW, JJC, YTL, CTM performed laboratory experiments. MYJL, DRC, ZSL, KJY performed library construction and sequencing. CNS wrote the manuscript with help from other authors. PW assisted CNS on the developmental timing of feather follicle analyses. YAS, TCT, CFC provided chickens for feather samples. CMC and WHL conceived the project and help interpret the significance. WHL managed the project and revised the manuscript. All authors read and approved the final version of the manuscript.
The feather samples used in this study. Figure S2. qPCR validation of 10 genes with biological replicates. Figure S3. The feather samples used for RNA extraction. (PDF 270 kb)
Summary of the 15 feather epithelial transcriptomes. (XLSX 15 kb)
Gene name and primer sequences used in RNA-seq validation. (XLSX 12 kb)
Transcriptome Expression Data. Table of mapped reads to Galgal4 transcripts for all 15 data sets. FPKM (Fragments per kilobase of exon per million fragments mapped): normalized transcript abundance values for each gene in the indicated tissues. (CSV 1314 kb)
Positively selected and/or rapid evolving genes in avian lineages  that are expressed in all feather samples. (XLSX 11 kb)
cEB vs. cLB differentially expressed gene set (shown in FPKM and fold change). (XLSX 81 kb)
cEB vs. cEF differentially expressed gene set (shown in FPKM and fold change). (XLSX 94 kb)
cEF vs. cMF differentially expressed gene set (shown in FPKM and fold change). (XLSX 86 kb)
cMF vs. cLF differentially expressed gene set (shown in FPKM and fold change). (XLSX 62 kb)
Canonical pathways for cEB vs. cLB differentially expressed gene set. (XLS 50 kb)
Canonical pathways for cEB vs. cEF differentially expressed gene set. (XLS 35 kb)
Canonical pathways for cEF vs. cMF differentially expressed gene set. (XLS 46 kb)
Canonical pathways for cMF vs. cLF differentially expressed gene set. (XLS 30 kb)
About this article
Cite this article
Ng, C.S., Chen, CK., Fan, WL. et al. Transcriptomic analyses of regenerating adult feathers in chicken. BMC Genomics 16, 756 (2015). https://doi.org/10.1186/s12864-015-1966-6