Salamander limb regeneration is dependent upon tissue interactions that are local to the amputation site. Communication among limb epidermis, peripheral nerves, and mesenchyme coordinate cell migration, cell proliferation, and tissue patterning to generate a blastema, which will form missing limb structures. An outstanding question is how cross-talk between these tissues gives rise to the regeneration blastema. To identify genes associated with epidermis-nerve-mesenchymal interactions during limb regeneration, we examined histological and transcriptional changes during the first week following injury in the wound epidermis and subjacent cells between three injury types; 1) a flank wound on the side of the animal that will not regenerate a limb, 2) a denervated limb that will not regenerate a limb, and 3) an innervated limb that will regenerate a limb. Early, histological and transcriptional changes were similar between the injury types, presumably because a common wound-healing program is employed across anatomical locations. However, some transcripts were enriched in limbs compared to the flank and are associated with vertebrate limb development. Many of these genes were activated before blastema outgrowth and expressed in specific tissue types including the epidermis, peripheral nerve, and mesenchyme. We also identified a relatively small group of transcripts that were more highly expressed in innervated limbs versus denervated limbs. These transcripts encode for proteins involved in myelination of peripheral nerves, epidermal cell function, and proliferation of mesenchymal cells. Overall, our study identifies limb-specific and nerve-dependent genes that are upstream of regenerative growth, and thus promising candidates for the regulation of blastema formation.
All animals regenerate some of their tissues by physiological turnover, yet only a select few can regenerate appendages. Vertebrates accomplish this feat by generating a blastema, a mass of lineage-restricted progenitor cells at the end of an amputation stump (Kragl et al., 2009). The cellular and molecular processes that coordinate blastema formation are poorly understood, likely because it is a complex process, requiring tissues of an anatomically complex amputation stump to coordinate wound healing, progenitor cell recruitment, cell proliferation, and tissue patterning. A major hurdle towards understanding appendage regeneration is to identify the necessary processes for regeneration and the molecular mechanisms by which these processes regulate blastema formation. For example, some cellular processes like inflammation, epidermal migration, and cell proliferation are common to all injury types, so it is necessary to devise experimental strategies that can distinguish pathways specific to general injury processes from those required for appendage regeneration.
The blastema of a regenerating salamander limb is a classic paradigm for studying appendage regeneration because it is an accessible experimental system that regenerates a morphologically complex structure. A critical tissue interaction that is necessary for blastema formation occurs between the wound epithelium (WE), which forms from rapid migration of adjacent epidermis, and the underlying mesenchymal stump cells. Blastema formation is inhibited if the WE is disrupted, either by suturing full thickness skin over the amputation stump (Mescher, 1976; Tassava and Garling, 1979), irradiation (Thornton, 1958), surgical removal (Thornton, 1957), or implantation of the limb stump into the body cavity (Goss, 1956) or dorsal fin to disrupt epidermal migration (Stocum and Dearlove, 1972). The WE gradually thickens after amputation to generate the apical epithelial cap (AEC). The AEC is a signaling center which supports mesenchymal cell proliferation (Boilly and Albert, 1990; Globus and Vethamany-Globus, 1985), promotes tissue histolysis (Singer and Salpeter, 1961), and regulates cell migration (Thornton, 1960b; Thornton, 1960a; Thornton and Steen, 1962; Thornton and Thornton, 1965). Molecules expressed in the AEC include; the transcription factors msx2, dlx3, id2, id3, hes1, sp9 (Satoh et al., 2008), the secreted signaling molecules wnt5a, wnt5b (Ghosh et al., 2008), fgf1, fgf2, fgf8, fgf10, the extracellular matrix molecules collagen type XII, collagen type IV, lamb1, and the enzymes mmp3/10b, mmp9 (Campbell and Crews, 2008). However, these studies were not designed to uncover molecules that are vital to the early function of the WE or that regulate its relationship with peripheral nerves and blastemal cells.
Innervation of the limb stump is necessary for regeneration. Transecting the spinal nerves that innervate the forelimb just prior to or shortly after limb amputation will block blastema formation (Singer, 1952). It is unknown why denervation leads to loss of regeneration in axolotls, but the favored hypothesis is that nerves provide trophic factors that support cell proliferation of the blastema, which is lost upon denervation. Several factors have been proposed as the trophic factor (Dungan et al., 2002; Globus et al., 1991; Mescher et al., 1997; Mullen et al., 1996; Satoh et al., 2008; Wang et al., 2000), but none have been clearly demonstrated to be the factor (Stocum, 2011). Furthermore, it is unclear if the nerve exerts its effects on the mesenchyme, epidermis, or Schwann cells. Part of the difficulty in identifying the exact mechanism of the nerve's influence is because nerve fibers quickly invade throughout the distal mesenchyme and wound epidermis after amputation (Singer, 1949; Taban, 1949; Thornton, 1954). Limbs containing only motor nerves that do not innervate the epidermis can regenerate, demonstrating that direct innervation of the epidermis is not necessary for limb regeneration (Sidman and Singer, 1960; Thornton, 1960b). Taken together, nerves need to invade the amputation stump to support cell proliferation, but the exact relationship between the nerves with the mesenchyme and epidermis is unclear. Identification of the downstream targets of the nerve during limb regeneration may reveal the relationship of the nerve with the amputation stump.
Although we do not yet understand the relationship between the epidermis, nerve, and mesenchymal cells at the molecular level, it is clear that the function of each tissue is dependent upon the presence of the other tissues and these interactions take place locally at the amputation plane. For this reason, the goal of this study was to identify a specific set of genes that are expressed in the WE and cells directly beneath the WE of the amputated limb stump early after injury. Our strategy was to use a custom microarray platform developed for the axolotl (Huggins et al., 2012) to compare gene expression differences over the first week after injury between an injured tissue that will not form a limb (flank wound), an example of aborted limb regeneration (denervated limb), and a regenerating limb (Fig. 1A–C). Using statistical approaches to identify significant transcript abundance differences between regenerating limbs and non-regenerating tissues (Fig. 1D), we were able to identify a regeneration-specific gene expression profile. A flank injury located outside the limb field was chosen because site-specific differences in the skin are known to mediate limb and flank morphology (Rinn et al., 2006; Satoh et al., 2007; Tank, 1984; Tank, 1987). Denervated limbs were chosen because they are an example of aborted limb regeneration – they present limb-specific gene expression patterns but do not generate a blastema. Lastly, gene expression during the first seven days after injury was chosen to identify genes expressed at the onset of AEC formation and blastema cell proliferation. In the following paragraphs we describe histological changes that take place over the first seven days in each injury type and then overlay transcriptional patterns of limb-enriched and nerve-dependent genes in the amputation stump.
Experimental design of the microarray analysis.
Histology and BrdU analysis of injured innervated limbs, denervated limbs, and flank wounds
We histologically characterized normal innervated limbs (NL), denervated limbs (DL), and flank wounds (FW) over the first seven days post injury (dpi) in order to examine differences between each injury response at the cellular level (Fig. 2). Masson's Trichrome staining revealed that the structure of the uninjured skin in NL, DL, and FW were similar with one another (data not shown). Uninjured epithelium consisted of an outer apical layer of epithelial cells, an interstitial layer of mucous secreting Leydig cells interspersed with keratinocytes, and a basal layer of germinative basal keratinocytes (Fox, 1986; Kelly, 1966). The underlying uninjured dermis consisted of mucous and granular glands interspersed with a loose network of fibroblasts that overlies muscle (Seifert et al., 2012).
Histology of NL, DL, and FW.
For all three injury types, the wound re-epithelialized within 24 hours after injury by migration of surrounding epidermis, generating a WE comprised of Leydig cells and keratinocytes (Fig. 2A–F). Underneath the WE was an accumulation of plasma and blood cells, with more blood in NL and DL versus FW, possibly because amputation severed major vasculature in the limb (Fig. 2B,D,F). Additionally, the FW was almost exclusively composed of muscle, while the limbs included bone, peripheral nerves, vasculature and muscle. In all three cases, the WE appeared to behave similarly during the first 24 hours after injury, although the extent of the hemostatic response and complexity of the underlying tissue is greater in the amputated limb compared to the flank.
By 7 dpi, the WE had thickened in all three injury types, but a distinct mound of epidermal cells was apparent in the middle of NL and DL WE, which was not present in FW (Fig. 2G–L). The epidermal mound may represent the maturation of the WE into the AEC, which was likely due to continuous cell migration from the wound margins rather than cell proliferation within the WE because BrdU-positive cells were evident at the margins of NL and DL (Fig. 3A,D), but not within the center of the WE (Fig. 3B,E) (Chalkley, 1954; Hay and Fischman, 1961). The WE of the FW was in direct contact with the subjacent muscle with little muscle dedifferentiation (Fig. 2K,L). In contrast, muscle fibers and peripheral nerves were becoming disorganized due to degeneration in both NL and DL (Fig. 2G,I). Limbs also had more plasma, red blood cells, and inflammatory cells. All together, these results suggest similar processes were taking place in each injury type, but tissue histolysis was more complete by seven days after injury in limbs versus FW.
BrdU staining of sections of injured limbs.
Taken together our data show that both limbs exhibited a hemostatic response and were histologically similar, containing a thickened WE, osteoclasts surrounding the bone, and histolysed tissues (Fig. 2A,C,G,I). The WE was lying directly over the bone in DL, while cells were present between the WE and bone in NL, suggesting that blastema growth was beginning within NL. BrdU analysis showed that DNA synthesis was taking place within the WE margin and mesenchyme (Fig. 3A,C,D,F) of both NL and DL, which is in accordance with previous studies showing that mesenchymal cells and epidermal cells enter S phase and divide in both denervated and innervated limbs (Maden, 1978). Loss of cell cycling in the mesenchyme of denervated limbs likely takes place after 7 dpi in the large-sized animals used in this study. These findings demonstrate that the time frame chosen for our study encompassed blastema formation rather than blastema outgrowth.
Commonly changed genes following injury
In order to characterize transcription during regeneration, transcript abundances were estimated from total RNA collected from the WE and a few subjacent cells from all three treatments (Fig. 1A–C). A total of 6684 probe sets yielded expression estimates that differed significantly as a function of RNA source and sample time (supplementary material Table S1). These genes were parsed to identify similarities and differences between injury types. First, probe sets that changed significantly from baseline to 1 dpi, 1 dpi to 3 dpi, or from 3 dpi to 7 dpi in each injury were identified to examine the commonalities between the injuries. We found that transcription was more similar between NL, DL, and FW than it was different with 1840 genes up-regulated and 1667 genes down-regulated in all three injuries (Fig. 4A,B). This high degree of similarity suggested the presence of a general wound-healing response regardless of whether or not a limb will regenerate. The list of genes that presented higher transcript abundances above baseline was significantly enriched for genes that annotate to gene ontology terms associated with processes known to take place during mammalian skin wound healing including immune system response (n = 187), macrophage activation (n = 36), and response to stimulus (n = 124) (supplementary material Table S2). Down-regulated genes belonged to ontology categories including lipid metabolic process (n = 110), chromosome segregation (n = 30), metabolic process (n = 604), and response to stress (n = 48). These results suggest that many of the same processes that take place during mammalian wound healing also occur in the axolotl following injury. Indeed, our histological analysis supports this result as well as an in-depth study on flank wound healing in the axolotl, which demonstrated that inflammation and a hemostatic response occurs in axolotls, but is dampened compared to mammals (Seifert et al., 2012). Although the set of genes common to all three injuries provides insight about wound closure, inflammation, and immunity, we focus in the following paragraphs on gene expression patterns that associate specifically with limb regeneration.
Summary of differentially regulated genes during limb regeneration.
Limb-enriched gene expression patterns
Overall, 2941 probe sets measured higher transcript abundances above baseline during the first seven days in NL (supplementary material Table S1). This list was filtered to identify genes with higher transcript levels in NL versus FW at 0, 1, 3, or 7 dpi; this yielded a list of 570 injury-induced, limb-enriched genes, of which 377 annotate to unique presumptive human orthologs (Fig. 4C; supplementary material Table S1). Annotation of gene ontologies for these probe sets identified a substantial number of genes involved in biological processes linked to limb development including developmental process (n = 125), mesoderm development (n = 70), and ectoderm development (n = 70) (supplementary material Table S3). In fact, mutations in 31 of these unique limb-enriched genes manifest human or mouse developmental limb defects when genetically disrupted, strongly suggesting that these genes play pivotal roles in vertebrate limb formation (Table 1). Furthermore, comparing our list to other genomic screens of limb regeneration (Campbell et al., 2011; Monaghan et al., 2009) we identified 73 genes that were commonly identified as highly expressed in amputated limbs (supplementary material Table S1).
Key regulatory genes involved in signaling pathways known to be necessary for limb development and limb regeneration were found in the limb-enriched list including genes integral to β-catenin-independent Wnt/planar cell polarity signaling (prickle1, prickle2, wnt5a, fzd2, fzd8, and ror2) (Stoick-Cooper et al., 2007), retinoic acid signaling (aldh1a3, crabp1, crabp2, and rdh10) (Blum and Begemann, 2012), insulin growth factor signaling (ctgf [igfbp8], cyr61 [igfbp10], igfbp2, igfbp3, htra1, and kazald1 [igfbp-rP10]) (Chablais and Jazwinska, 2010), FGF signaling (dusp6, fgfr1, and pdlim7) (Lee et al., 2009), and BMP signaling (bmp2, id3, bmp2r) (Guimond et al., 2010). Overall, this list supports the hypothesis that some gene expression programs used in development are re-deployed during limb regeneration (Muneoka and Sassoon, 1992). Surprisingly, genes associated with limb patterning and growth were up-regulated before considerable increases in cell proliferation and blastemal outgrowth, suggesting that patterning of the limb blastema may occur in parallel or prior to blastema growth.
Limb-enriched and nerve-dependent gene expression patterns
To identify nerve-dependent genes, limb-enriched genes were filtered to identify probe sets that measured higher transcript levels in NL versus DL at 1, 3, or 7 dpi. This strategy identified a short list of 56 genes (41 unique transcripts with presumptive human orthologs) that were up-regulated after injury, had higher transcript abundance in NL versus FW, and had higher transcript abundance in NL versus DL (Fig. 4C; supplementary material Table S1). This list was significantly enriched for genes that annotate to developmental process (n = 19), ectoderm development (n = 6), cell cycle (n = 8), and neurological system process (n = 7) ontology terms (supplementary material Table S4). Further annotation through PubMed searches showed that 14 of the 40 genes are important in epithelial function (Table 2), including genes important in maintaining the structure of the epithelia (krt8, kera, krt15, cldn19, col29a1, eppk1, and tgm1), epithelial cell growth factors (ereg), and transcription factors involved in keratinocyte growth and differentiation (zfp36l2, ifit5), suggesting that these genes are necessary for the maturation of the WE into the AEC. We also identified 9 out of 40 genes that are highly expressed in the peripheral nervous system of vertebrates including four genes that are highly abundant in peripheral nervous system myelin Schwann cells (mbp, pmp22, gldn, and mpz) (Table 2). Overall, this list of genes suggests that denervation affects maturation of the WE and behavior of Schwann cells within the first week of regeneration.
Localization of up-regulated, limb-enriched, and nerve-dependent transcripts
In situ hybridization was used to localize mRNA expression of limb-enriched and nerve-dependent genes at 7 dpi. The results show considerable variation in the location of transcripts among WE keratinocytes, blastema cells, and peripheral nerves (Fig. 5). For example, a putative S-adenosylmethionine-dependent methyltransferase (axo23458-r) was expressed in keratinocytes of the WE and not the underlying mesenchyme of 7 dpi limbs (Fig. 5A,B). This methyltransferase-like gene was highly expressed in each injury type at 1 dpi, but expression was sustained at higher levels in NL and DL at 7 dpi (Fig. 5C). This gene was identified in other genomic screens of limb regeneration (supplementary material Table S1), making it a promising candidate for its involvement in WE function after injury.
In situ hybridization (ISH) of limb-enriched genes in NL 7 dpi limbs.
A transcript highly similar to human krt5 was expressed in WE keratinocytes and the underlying mesenchyme and cartilage cells (Fig. 5D,E). In mice, krt5 and its binding partner krt15, are markers for salivary gland epithelial progenitor cells, which show decreased cell proliferation and krt5 expression upon removal of parasympathetic innervation (Knox et al., 2010). In our study, krt15 was limb-enriched and nerve-dependent, and krt5 was eight times higher in NL versus DL at 7 dpi (Fig. 5F), although highly variable estimates among replicates yielded a p-value below our statistical cutoff (P = 0.019). This strongly suggests that krt5 and krt15 are limb-specific and nerve-dependent gene candidates. In support of these results, a newt type II cytokeratin that is highly similar to our presumptive krt5 (blastn; 83% identical), is transcribed in the mesenchyme and WE of regenerating newt limbs (Ferretti et al., 1991; Ferretti and Ghosh, 1997). Together, these data suggest important roles for krt5 and krt15 in the blastema and more generally, nerve-mesenchyme-epidermis interactions that typify a normal regenerative response.
In situ hybridization also showed that genes associated with retinoic acid signaling were expressed in regenerating limbs (Fig. 5G–L). Aldh1a3, a retinaldehyde dehydrogenase that synthesizes retinoic acid during development and adulthood, was exclusively expressed in a subset of cells within peripheral nerve bundles 7 dpi (Fig. 5G,H) and did not rise above baseline levels in FW (Fig. 5I), likely because nerve bundles are only present in the limb samples. We also found that crabp1, an intracellular retinoic acid binding protein that regulates RA nuclear signaling, was up-regulated from baseline exclusively in the mesenchyme of the limb blastema and was both limb-specific and nerve-dependent (Fig. 5J–L). These expression patterns may explain why RA is necessary for appendage regeneration (Blum and Begemann, 2012; Maden, 1998) and can re-specify pattern in the regenerating axolotl limb (Maden, 1982). Further investigation is needed to identify whether aldh1a3 expressing cells are producing RA, signal to crabp1-expressing blastemal cells, and if this process is necessary for regeneration. Overall, our ISH analyses show that the genes identified in our study are expressed in three tissues that mediate blastema formation; the epidermis, mesenchyme, and peripheral nerve. It also suggests that each of these tissues is affected by denervation prior to blastema formation.
Differential expression between innervated and denervated limbs
Previous studies have shown that gene expression between innervated and denervated limbs is often quantitatively changed rather than absolutely (Monaghan et al., 2009). In order to address this possibility, differentially regulated genes were identified between innervated and denervated limbs regardless of expression changes from baseline at 1, 3, and 7 dpi. Comparing NL and DL at 1 dpi identified a small list of 25 unique genes (supplementary material Table S1) that presented higher transcript abundances in NL and were significantly enriched for genes that annotate to developmental, neurological, and systems process ontology terms (Table 3). This list includes genes associated with microvascular morphogenesis (krit1), blood coagulation (f5), retinoic acid binding (crabp1), extracellular matrix structure (col29a1), myelin synthesis and structure (mbp, pmp22, mpz), axon guidance (reln, homer1), and axon development (gldn). The genes that presented higher transcript abundances in denervated limbs included three genes, xdh, alox12b, and alox15b, that enriched 2 ontology terms, respiratory electron chain transport and generation of precursor metabolites and energy. Other genes in this list are predicted to function in the regulation of apoptosis (btg3), neurogenesis and pluropotency (rbbp9), and neurodegeneration (yars). These results show that denervation significantly alters transcription as early as 24 hours post amputation and we note that several of these differentially expressed genes were also found in our lists of limb-enriched and nerve-dependent genes.
More genes were differentially expressed at Day 3 than Day 1 (n = 52 unique genes) (supplementary material Table S1) and these enriched different biological process terms between innervated and denervated limbs (Table 3). Genes that were expressed more highly in innervated limbs enriched carbohydrate metabolism, transport, hematopoeisis, and B-cell immunity biological process terms. In addition, this list included genes associated with: (1) Schwann cells and neurons (gfpt2, tuba1a, glul, lnx1, marveld2), (2) extracellular matrix structure and synthesis (ugdh, slc23d2), (3) regulation of epithelial-mesenchymal transition (EMT) (fam3c, hmga2), (4) regulation of skeletal development (wsb1, tpp3), (5) FGF-signaling of angiogenesis (cav1) and limb development (pdlim7), and (6) regulation of epidermal cell differentiation and proliferation (tgm1, ovol2, lmo7, ehf, ereg, sorbs3, eppk1). The results show that a diverse group of developmentally important genes are differentially regulated between innervated and denervated limbs by 3 dpi.
The largest number of differentially expressed genes was discovered for 7 dpi (n = 103 unique genes) (supplementary material Table S1), and again, these enriched different biological process terms between innervated and denervated limbs (Table 3). Approximately 50% of the 103 genes that were expressed more highly in innervated limbs annotated to cell cycle and mitosis-related gene ontologies. These genes encode regulators of the cell cycle (ccna2, ccnb1, ccnb2, ccnb3, cdc2, tk1, ube2c, uhrf1), chromosome condensation and DNA repair (smc2, smc4, ncpag, c16orf75, pcna), genome replication (mcm2, mcm3, mcm4, mcm6, mcm7), and chromosome segregation (fbxo5, mad2l1, aurka, aurkb, aspm, kif11). In addition, 36 genes significantly enriched the nucleic acids metabolic process term, and the developmental process term was also enriched with genes associated with cell proliferation and differentiation (bcn1, zfhx4, krt15, krt8, crabp1, lingo1, lnx1, hmgb3, tk1). This strong signature of cell proliferation and DNA synthesis was missing in the list of genes that were significantly up-regulated >1.5 fold from baseline levels in NL (supplementary material Table S2). This suggests that cell cycle components are affected in a quantitative manner by denervation and they only begin to increase above baseline levels at 7 dpi in NL. This is in support of our histological analysis (Fig. 3), which showed that cell proliferation was present in both DL and NL at 7dpi, suggesting that proliferation dynamics are just beginning to diverge at this time (Fig. 3). Overall, the results show that denervation has a major effect on the transcription of proliferation-associated genes that are likely required for blastema growth.
Technical and biological replication of microarray results
To validate and extend the Affy microarray results, expression values were estimated using the nCounter platform. This analysis used 24 RNA samples from the Affymetrix analysis plus a newly generated set of 24 biological replicates. Custom Nanostring capture probes were designed for 50 genes (supplementary material Table S5) and fold change estimates were obtained between Day 0, 1, and 7 time points for innervated and denervated limbs. The correlation of fold change for the technical replicate samples was uniformly high across times and treatments (r = 0.95–0.97) (supplementary material Table S6). Thus, the Affy and nCounter platforms yielded precise estimates of fold change when the same RNA samples were processed. Precise estimates of fold change were also obtained between the Affy and nCounter platforms for biological replicates, and also between the two sets of 24 samples processed on the nCounter instrument. However, relative to the high correlation between technical replicates, the correlations between biological replicates were relatively lower for the Day1 and Day7 comparisons (r = 0.80–0.87) and lower still for all comparisons to Day 0 (r = 0.59–0.72). These results suggest that more variation is present between biological replicates than technical replicates and that the abundance of transcripts is most variable among samples that were collected at the time of limb amputation. This may be because animals were at different stages of the molting cycle at the time of collection. Overall, replication of the microarray results was high using the nCounter platform, demonstrating the reliability of each platform.
A recent comparison of transcription between innervated and denervated limbs of the Mexican axolotl provided the first global, transcriptional description of the limb regeneration program (Monaghan et al., 2009). That study used a small format microarray (∼4500 probes with 3271 presumptive human orthologs) to detail gene expression of whole blastemas at 5 days and 14 days after limb amputation. However, a more comprehensive analysis of gene expression was needed in concert with an earlier and more precise tissue-sampling scheme to thoroughly investigate the transcriptomics of blastema formation. To this end, we investigated transcription within axolotl epithelium and subjacent cells during the first week of limb regeneration with the primary goal of identifying a core set of genes that are likely to be necessary for limb regeneration. To meet this goal, we devised a strategy that allowed us to subtract out genes common to all injury responses as well as to identify genes that are uniquely expressed in limbs. We further selected genes specific to limbs that regenerate (NL) rather than regress (DL) to identify genes associated with blastema formation. Our study is the most detailed molecular analysis of limb regeneration to date and is the first to identify genes specific to the limb regeneration process by comparing the general wound healing response outside a limb field. Overall, the genes identified here will be useful as tissue specific markers for regenerating limbs and candidates for regulating blastema formation.
At both histological and transcriptional levels, we show that the initial injury response is similar between NL, DL, and FW. The time to re-epithelialization was within 1 dpi and many of the same genes were differentially regulated in NL, DL, and FW. Interestingly, many of these injury-response genes are similarly regulated in mammals, suggesting conservation of some aspects of wound healing among tetrapods.
Previous studies have shown that there are fundamental differences between limb skin and flank skin. When forelimb skin is replaced with grafts of flank skin in newts and axolotls, limb regeneration is defective (Tank, 1984; Tank, 1987). We compared gene expression differences between NL and FW and identified transcripts with limb-specific expression patterns. For example, the homeobox-containing transcription factor, emx2, was highly expressed in uninjured limb skin and was up-regulated after injury only in limb samples (supplementary material Table S1). Mouse emx2 null mutants fail to form a scapula during development (Table 1) (Capellini et al., 2010; Pellegrini et al., 2001) and newt emx2 is expressed in a graded proximodistal manner mainly in the epidermis of regenerating newt limbs (Beauchemin et al., 1998). Emx2 and other limb-specific genes identified in our study (supplementary material Table S1) may regulate limb-specific patterning events during regeneration.
Our analysis also identified a connection between salamander limb-enriched genes and orthologs that are associated with limb deformities in mammals (Table 1). For example, numerous genes involved in the Wnt/Planar cell polarity (PCP) signaling pathway were up-regulated in DL and NL, but remained at baseline levels in FW. Activation of PCP signaling by Wnt5 ligand through Vangl and Ror2 activation regulates limb bud elongation during mammalian development (Gao et al., 2011) and Wnt5a activity is necessary for axolotl limb regeneration (Ghosh et al., 2008). WNT/PCP signaling is thought to stabilize cellular polarity in epithelium of developing limbs, organize directional cell migration, and regulate directional cell proliferation (Wang et al., 2011). Overall, it is clear that activation of Wnt signaling through Wnt5a is necessary for limb outgrowth, but the key problem is to identify the property of salamander limbs that allows this pathway to re-activate after injury while not being induced after a flank injury. It is possible that sustained expression of genes like emx2 into post-embryonic and larval stages allows accessibility of this important signaling pathway in adult axolotls.
Genes involved in other important signaling pathways were also up-regulated specifically in limb samples. For example, genes associated with retinoic acid (RA) signaling were dynamically expressed in limbs after injury. Retinoic acid is an important signaling molecule involved in the development and regeneration of limbs; disruption of this pathway disrupts limb formation (Blum and Begemann, 2012; Kikuchi et al., 2011; Maden, 1998; Maden, 2007). We found that crabp1 was only up-regulated in NL and was expressed exclusively in the limb mesenchyme. In contrast, we found that crabp2 was up-regulated in NL, DL, and FW at 7 dpi. Our findings are in accordance with previous studies showing that CRABP protein is up-regulated during regeneration, although it is unclear whether these studies were detecting CRABP1 or CRABP2 (Maden et al., 1989; McCormick et al., 1988). CRABPs are intracellular RA binding proteins that are thought to shuttle RA to the nucleus to regulate RA-mediated transcription, which may explain why we observe their expression during limb regeneration (Noy, 2000). We also found that aldh1a3 and rdh10, enzymes involved in the synthesis of RA during development, were up-regulated in NL and DL at 1 dpi and aldh1a3 was expressed specifically in cells resembling perineural fibroblasts in peripheral nerve bundles. Altogether, our data suggest that RA signaling is a dynamic process during limb regeneration and identifies the genes that may mediate the necessity of RA during epimorphic regeneration (Blum and Begemann, 2012; Kikuchi et al., 2011).
Beyond signaling pathways, structural proteins showed very specific transcriptional profiles in regenerating tissues. Numerous keratins (krt5, krt8, krt15, and krt13) and keratin-associated molecules (eppk1, tgm1, kera) were up-regulated after injury and were enriched in limbs. Furthermore, some genes like krt8 and krt15 were highly nerve-dependent. Keratins are components of intermediate filaments that protect the structural integrity of cells, but have recently been implicated in other cellular processes including cell motility, cell signaling, cell growth, and cancer metastasis (Karantza, 2011; Windoffer et al., 2011). Although previous studies in newts have identified keratins NvKII, krt8, and krt18 in mesenchymal and WE cells during limb regeneration (Ferretti et al., 1991; Ferretti and Ghosh, 1997) and knockdown of krt8 and krt18 in newt blastemal cells in vitro decreased DNA synthesis (Corcoran and Ferretti, 1997), our understanding of these proteins during regeneration remains poor. Functional testing is necessary to determine if keratin proteins play solely a supportive, structural role during regeneration or whether they are mediating cell signaling to promote growth or patterning. Together, the highly limb-specific and nerve-dependent expression patterns of the keratin genes strongly suggest that they are integral to the formation of the blastema.
Other limb-enriched genes were more quantitatively different than FW rather than being expressed exclusively in the limb. For example, two possible salamander-specific genes, sodefrin-like (axo22108-r) and methyltransferase-like (axo23458-r), were up-regulated in NL, DL, and FW, exclusively in the epidermis (Fig. 5A,B; data not shown), but expression was only maintained in NL and DL. This suggests that these molecules are not limb-specific, although sustained expression in the limb WE may impose some necessary function to the WE during limb regeneration. Regardless, the fact that these genes seem to be unique to salamanders (Campbell et al., 2011) and show strong and specific expression in the WE warrants further functional studies.
A surprising result was the observation that myelin-associated genes were up-regulated and both limb-specific and nerve-dependent. Myelinated peripheral nerves permeate throughout the uninjured limb, but only naked sensory nerve fibers are found in uninjured epidermis of animals (Boulais and Misery, 2008). Hence, our tissue collection scheme did not sample myelinated nerve fibers in uninjured samples, yet injured NL and DL samples contained transected nerve bundles located just proximal to the WE. This likely explains why mRNA levels of myelin-associated genes increased above baseline in NL at 1 dpi. The fact that myelin-associated gene mRNA did not increase in DL suggests that expression of these genes was lost following denervation. A similar phenomenon takes place in mammals, where peripheral nerve fiber transection down-regulates expression of myelin-associated genes in distal Schwann cells (Hall, 2005). This result is interesting because it suggests that Schwann cells are affected early after denervation, which may have detrimental effects on downstream blastema formation. In newts, the protein Anterior Gradient 2 is expressed in Schwann cells after limb amputation and supplemental Anterior Gradient 2 can partially rescue regeneration in the denervated state (Kumar et al., 2007). Others have shown that denervation in axolotls induces peripheral nerves to become inhibitory to limb regeneration, suggesting that they may secrete inhibitory factors (Irvin and Tassava, 1998; Tassava and Olsen-Winner, 2003). It will be critical in future experiments to determine if the response of Schwann cells to denervation is the cause of a loss of blastema formation.
The proliferation of blastema cells is known to be a target of the nerve during limb regeneration (Stocum, 2011). In order to increase our sensitivity for identifying proliferation-associated genes during regeneration, we directly compared NL to DL without comparing samples to baseline or FW. This analysis showed that by 7 dpi, approximately 50% of the genes that were higher in NL versus DL were associated with the cell cycle, supporting the notion that the cell cycle is the primary target of denervation. Most of these genes were only different at 7 dpi, suggesting that our study identified the genes likely upstream of the cell proliferation effect of denervation. This result highlights that the limb-enriched and nerve-enriched genes we identified in our study are excellent candidates for regulating the increase in cell proliferation that is characteristic to limb regeneration. Overall, our study used a focused approach to identify the genes that are likely necessary for limb regeneration and showed that many of these genes are expressed in specific tissues and before considerable outgrowth takes place in the limb. The identification of these genes is an important advance in our ability to tease apart the cellular and molecular mechanisms that drive regeneration and will be a useful resource for regeneration researchers that may be looking for specific genes to analyze during early blastema formation.
Materials and Methods
Animals and surgical procedures
Axolotls were obtained from the Ambystoma Genetic Stock Center, Lexington, KY and raised to 7–10 cm snout to vent length. Animal care and use procedures were approved by the University of Florida IACUC (Application Number 201101534). Denervations were performed by anesthetizing animals in 0.01% benzocaine, making a small incision at the shoulder to expose the brachial nerves entering each forelimb, and severing the nerve bundles with surgical scissors. Limb amputations were performed at the mid-stylopod and the humerus was trimmed to make the amputation plane flush. Full thickness excisional wounds were performed along the flank of anesthetized animals using a 4 mm biopsy punch tool.
Histology and BrdU analysis
Tissues were processed for paraffin embedding, sectioned at 5 µm, and stained according to previous methods (Seifert et al., 2012) except that limb samples were decalcified in 10% EDTA for 3 days with daily changes at 4°C before histological processing. For DNA synthesis analysis, animals were injected with bromodeoxyuridine (BrdU) (conc. = 100 mg/g) 24 hours before tissue collection, harvested 24 hours later, and processed for paraffin embedding. After sectioning, sections were de-paraffinized, blocked for endogenous peroxidase activity in 3% H2O2 in methanol for 10 mins, rehydrated, treated for antigen retrieval in pH 6.0 sodium citrate buffer in a microwave for 25 mins, rinsed in water, incubated in 37°C 2N HCl for 15 mins, rinsed thoroughly in water, rinsed with TBS, blocked with rabbit serum, blocked for endogenous avidin and biotin, incubated with primary antibody rat anti-BrdU (1:500, Accurate Scientific), washed, incubated with biotinylated secondary anti-rat (1:400, Vector Scientific), washed and visualized using Vector ABC horseradish peroxidase and DAB reagents according to manufacturer's instructions. Tissue sections were counterstained with Hematoxylin (Vector).
Tissue collection for microarray analysis
One day prior to limb amputation, 12 axolotls were anesthetized and their forelimbs denervated. Approximately 24 hours later, these same axolotls with denervated limbs and 16 additional axolotls with innervated limbs were anesthetized and administered amputations at the mid-stylopod of both forelimbs (Fig. 1A,B). The epidermis adjacent to the amputation plane was taken from the arms of each of four individuals that were not denervated the day before; these served as Day 0 samples for the innervated and denervated limbs. To obtain a sufficient amount of RNA for microarray analysis, both forelimb samples from each individual were pooled to yield independent, replicate samples. An additional 16 axolotls were then anesthetized and a full thickness excisional wound was performed along the flank of each animal (Fig. 1C). After 1, 3, and 7 days post injury, the wound epithelium was removed and any cells that were adhered to the epithelium were included in the sample. Two tissues were pooled from each individual to obtain 4 replicate samples for each time point and tissue type.
RNA isolation and microarray analysis
Total RNA was isolated from all 44 samples and each was processed for hybridization to 44 independent and custom A. mexicanum (Amby_002) Affymetrix GeneChips (Huggins et al., 2012). Microarray results from FW samples are summarized in a different manuscript (Seifert et al., 2012), but the entire dataset can be found at the Gene Expression Omnibus (Accession number GSE37198). This Amby_002 GeneChip contains approximately ∼20,000 perfect match probesets. The probesets were designed using A. mexicanum expressed sequence tag contigs from Sal-Site (Smith et al., 2005). The GeneChips were processed by the University of Kentucky Microarrray Core Facility and expression values were extracted using RMA (Irizarry et al., 2003) and Affymetrix Expression Console software. The resulting data were subjected to one way Analysis of Variance using JMP Genomics version 4.1 and statistical estimates were defined to make 16 comparisons between groups (Fig. 1D). A gene was identified as differentially expressed if it passed a false discovery rate of ≤ 0.05 and had a fold change of ≥ 2, or passed a false discovery rate ≤ 0.001 and a fold change of ≥ 1.5. Genes that were defined as differentially expressed were analyzed further using pair-wise comparisons at a significance threshold of P<0.003 and a fold change cutoff of 1.5 for the comparison of interest. This threshold was determined using a Bonferroni correction to adjust for 16 pairwise comparisons at an alpha level of 0.05. The R package, VennDiagram was utilized to generate Venn diagrams (Chen and Boutros, 2011). Significant genes were annotated with gene ontology information from Panther (http://www.pantherdb.org) and gene lists were compiled and compared to identify biological processes that were statistically over-represented. For all analyses, the 11,131 probesets on the Ambystoma GeneChip that could be mapped to human orthologs in the Panther database were used to generate expected values (i.e., as the background). The count threshold was set to three and the significance threshold was set to P<0.05. The lists of significant biological process terms were manually inspected to remove redundant terms.
Cloning and RNA probe production
Axolotl genes were cloned using gene-specific primers designed using sequences collected from the Ambystoma Gene Collection (Smith et al., 2005). Total RNA was isolated from 7 dpi limb tissue and used to make cDNA template (iScript; BioRad). Genes were amplified as follows: krt5 primers were 5′ GAG GGA GCA GGT TCT GTG AG 3′ and 5′ ATC ACC CAG CCA GAA GAA TG 3′; aldh1a3 primers were 5′ CCT GCA TTG TGT TTG CTG AC 3′ and 5′ TGT CAG AGC CGG ATA ATT CA 3′; crabp1 primers were 5′ AGG AGT CCC CTG ACT TGG AG 3′ and 5′ TGC CAC CAC AAA TGA TGA GT 3′. PCR products were gel isolated, cloned into pGEM-T Easy Vectors (Promega), and sequence verified. methyltransferase-like primers 5′ TAA TAC GAC TCA CTA TAG GGA GAC AGC TCT GTG GAT CTG GTC A 3′ and 5′ ATT TAG GTG ACA CTA TAG AAG AGT CTC TAA GGT GCG GCT TGT T 3′ were used to make a PCR template that was used to generate a digoxygenin-labelled RNA probe using a Roche RNA labeling kit.
In situ hybridization
Limbs were collected 7 dpi and fixed overnight in 4% PFA at 4°C, mounted in optimal cutting medium, sectioned at 20 µm, dried for two hours, and processed for in situ hybridizations on the same day according to previously published methods (David Parichy, personal communication). Proteinase k treatment consisted of 10 minutes at 10 µg/ml concentration. Probe concentration was 0.5 µg/ml in hybridization solution at 55°C overnight. Anti-DIG antibody was incubated at 1:5000 dilution at 4°C overnight. Stained sections were mounted in 80% glycerol and images captured on a Nikon Eclipse 6600 upright compound microscope using a Cool-Snap Pro true color camera.
Nanostring nCounter development and analysis
Transcript abundance estimates obtained from the Ambystoma Affy GeneChip were compared to estimates obtained from the Nanostring nCounter System. The nCounter is a moderate throughput gene expression analysis instrument that estimates the number of RNA transcripts from samples of total RNA or lysed tissues. Nanostring staff designed capture probes for 50 genes from the Affy GeneChip (supplementary material Table S5) and processed 48 RNA samples. Twenty-four of the RNA samples corresponded to the same replicate Day 0, Day 1, and Day 7 RNA samples that were used in the Affymetrix experiment. The second group of 24 samples corresponded to a new set of replicate D0, D1, and D7 samples. The count data for all genes were normalized to the counts of two capture probes that were consistently expressed across the innervated and denervated treatments. Technical correlation of gene expression was examined between the Affymetrix and Nanostring platforms by calculating Pearson's correlation coefficient (r) across all 50 genes for fold change estimates obtained using the same RNA samples. Biological correlation of gene expression was examined between the Affymetrix and Nanostring platforms, and between the two sets of replicates processed on the Nanostring platform. Again, Pearson's correlation coefficient (r) was calculated across all 50 genes for fold change estimates obtained between different RNA samples.
The axolotls were obtained from the Ambystoma Genetic Stock Center at the University of Kentucky, which is funded by the National Science Foundation [DBI-0951484 to S.R.V.]. The Ambystoma GeneChip was generated under National Institute of Health [R24-RR016344 to S.R.V.] and Army Research Office [W911NF-09-1-0305 to S.R.V.]. The work was supported by National Institute of Health [RC2-NS069480 to S.R.V. and M.M.] and funding from The Regeneration Project at University of Florida to S.R.V., M.M., D.M.G., and J.R.M.
The authors have no competing interests to declare.