- Research article
- Open Access
Digital gene expression tag profiling of bat digits provides robust candidates contributing to wing formation
BMC Genomicsvolume 11, Article number: 619 (2010)
As the only truly flying mammals, bats use their unique wing - consisting of four elongated digits (digits II-V) connected by membranes - to power their flight. In addition to the elongated digits II-V, the forelimb contains one shorter digit (digit I) that is morphologically similar to the hindlimb digits. Here, we capitalized on the morphological variation among the bat forelimb digits to investigate the molecular mechanisms underlying digit elongation and wing formation. Using next generation sequencing technology, we performed digital gene expression tag profiling (DGE-tag profiling) of developing digits in a pooled sample of two Myotis ricketti and validated our sequencing results using real-time quantitative PCR (RT-qPCR) of gene expression in the developing digits of two Hipposideros armiger.
Among hundreds of genes exhibiting significant differences in expression between the short and long digits, we highlight 14 genes most related to digit elongation. These genes include two Tbx genes (Tbx3 and Tbx15), five BMP pathway genes (Bmp3, RGMB, Smad1, Smad4 and Nog), four Homeobox genes (Hoxd8, Hoxd9, Hoxa1 and Satb1), and three other genes (Twist1, Tmeff2 and Enpp2) related to digit malformations or cell proliferation. In addition, our results suggest that Tbx4 and Pitx2 contribute to the morphological similarity and five genes (Acta1, Tnnc2, Atp2a1, Hrc and Myoz1) contribute to the functional similarity between the thumb and hindlimb digits.
Results of this study not only implicate many developmental genes as robust candidates underlying digit elongation and wing formation in bats, but also provide a better understanding of the genes involved in autopodial development in general.
Bats are a speciose group which constitute approximately 20% of extant mammals, and are the only mammalian group that have powered flight . Unlike the bird wing, made up of an elongated radius and ulna and modified wrist and hand bones that support the flight feathers, bat wings are primarily shaped by elongated digits (digits II-V) that support the wing membranes or dactylopatagia (Figure 1A) . This digit elongation is extreme. For their body length, bats have the longest relative digit lengths with the longest digit (digit III) exceeding the length of the body. However, not all of the forelimb digits are elongated, digit I - the thumb - maintains a short morphology similar to the hindlimb digits and is used to grasp surface or move on the ground (Figure 1). This unique morphology of the bat forelimb digits has sparked interest and motivated investigation by many developmental and evolutionary biologists. Despite these efforts, the mechanisms underlying development of these two distinct digit morphologies (i.e., short and long) of the bat hand remains unresolved. In this study, we capitalize on variation in digit development morphology within the forelimb and between the fore- and hindlimb in bats to identify candidate genes associated with development of distinct digit morphologies. Using DGE-tag profiling (described below) of the thumb, the elongated forelimb digits and the hindlimb digits we identified 14 candidate genes likely associated with digit elongation and seven candidate genes likely associated with the morphological and functional similarity between the thumb and hindlimb digits.
Elongation of the posterior digits (digits II-V) in the bat forelimb is initiated at an early stage when digit anlagen initially appear, embryonic Stage 15 . By embryonic Stage 17, the elongation becomes obvious with lengthening of forelimb digits II-V greatly exceeding that of the thumb and hindlimb digits [3, 4]. This rapid elongation of the posterior forelimb digits is thought to results from modification in expression of key developmental genes [5, 6]. In fact, several experimental and empirical studies have found differences in expression of limb and digit associated genes between mouse and bat autopodia. First, a mouse transgenic experiment that replaced an enhancer of the mouse Prx1 - a transcription factor that regulates long bone growth in mouse - with the orthologous putative enhancer sequence from a bat resulted first in an increase in Prx1 expression and second in a small, but significant elongation of the mouse forelimb digits [7, 8]. Second, Shh, a key factor in limb development, exhibits a unique second phase of signaling late in development of the bat wing that is not present in the development mouse autopodia . Next, other developmental genes including Fgf8 and Bmp2/4 are up regulated in bat forelimb digits compared to bat hindlimb or mouse digits [10, 11]. Finally, the distribution Hoxd13 expression is more posteriorly restricted in the bat forelimb relative to mouse fore- and hindlimbs and bat hindlimbs during late embryonic stages [12, 13]. While these results indicate that differences in expression of Hoxd13, Fgf8, Bmp2/4 and Prx1 during forelimb digit development may contribute to bat digit elongation, the complex developmental changes required to generate the bat wing likely require numerous molecular changes . Thus, further examinations of gene expression during development of the bat digits are critical to identify the genes associated with digit elongation in the posterior forelimb digits and those responsible for the morphological and functional similarities between the short thumb and hindlimb digits.
Here, using a DGE-tag profiling approach, we compared patterns of gene expression among three cDNA libraries formed from distinct regions of the fetal autopodia (the thumb, the four posterior forelimb digits and the all five hindlimb digits) in Rickett's big-footed bat, Myotis ricketti. The method of DGE-tag profiling combines the SAGE principle with the massively parallel sequencing technology to study genome-wide gene expression profiles . DGE-tag profiling has three main advantages compared to the conventional microarray technology. First, it is high-throughput sequencing. As a result, data from a single lane of a flow cell (consisting of eight lanes allowing for eight samples) is enough to analyze of genome wide gene expression including weakly expressed genes which cannot be assessed using microarray [15–17]. Next, this type of sequencing data is highly replicable and accurate such that technical replicates have little variation . Finally, sequencing can be performed without a priori sequence information, which is required in microarrays . These features of DGE-tag profiling aided our identification of genes associated with digit morphology in bats. Overall, our examination of differentially expressed genes revealed candidate genes associated with both digit elongation as well as morphological and functional similarities between the thumb and hindlimb digits. Finally, our results revealed genes associated with digit position in both limbs.
All fetuses obtained were in the "Fetal Stage" which is the last prenatal stage according to bat embryonic staging systems . At this stage, the overall appearance of the embryo is similar to the former stages, but the body size as well as the digit length increases rapidly [3, 18]. The genital tubercle had become a vagina or a penis allowing determination of specimen sex. All the specimens used in this study had relatively weak pigmentation indicating an early Fetal Stage as overall body pigmentation increases as fetal development progresses [3, 18]. Claws were sharp and keratinized in the thumb and hindlimb digits. Crown rump length (CRL) and forearm length (FA) for both species were recorded in Table 1.
Statistics of tag sequencing
In total, 4377261, 4650734 and 7514943 reads (not including the adapter tags) were sequenced from Library Hand DI, Library Hand DII-V and Library Foot, respectively. Library Hand DI contains 375701 distinct tags, Library Hand DII-V contains 483043, and Library Foot contains 532850. These distinct tags and their genomic frequency as well as the raw data were deposited in NCBI Gene Expression Omnibus (GEO) database with the accession number [GEO: GSE20038]. Most of the distinct tags (over 62%) appeared only once (frequency = 1); however, a small number of distinct tags (less than 1.7%) with frequency higher than 100 made up over 63% of all the tags in all three libraries (Additional file 1, Table S1). When sequencing depths reach 1,000,000 total tags, the number of distinct tags discovered (especially those with frequency >1) drops dramatically in all the three libraries (Figure 2). From that point, increasing sequence depth results in a slow and stable accumulation of new distinct tags indicating that the sequencing has reached saturation (Figure 2).
Differences in expression of genes among the three libraries are indicated in Figure 3. In total, 7080, 8205 and 8967 human genes and 4742, 5666 and 6379 mouse genes were mapped using the tags from the three libraries Hand DI, Hand DII-V and Foot. In total, 38.5% human genes (10490 genes) and 24.6% mouse genes (7843 genes) were mapped using the tags from the three libraries (Additional file 2, Table S2 and additional file 3, Table S3). Some tags mapped to human genes but not mouse genes or vice versa.
DGE-tag profiling of marker gene expression
The expression patterns of five marker genes (Tbx5, Pitx1, Meis2, Hoxd10 and Krt17) in the three DGE-tag profiling libraries are consistent with previously reported patterns of expression at embryonic earlier stages (reviewed in the Discussion; Table 2). We found that expression of Tbx5 was significantly higher in forelimb digits than the hindlimb digits (p < 0.001) and Pitx1 expression was significantly lower in forelimb than in hindlimb digits (p < 0.001). Similarly, as expected, we found that expression of Meis2 and Hoxd10 were significantly higher in Library Hand DII-V relative to the other two libraries (p < 0.001) and Krt17 was significantly lower in Library Hand DII-V relative to the other two libraries (p < 0.001).
Differences in gene expression across the libraries
Additional files 2 and 3 list the results of DEGseq analysis for all the homologs of human (Additional file 2, Table S2) and mouse genes (Additional file 3, Table S3). For the comparison between elongated digits and the short digits, we found that 286 homologs of human genes and 89 homologs of mouse genes were up regulated in Library Hand DII-V compared with other two libraries (p < 0.001, Figure 4) and 97 homologs of human genes and 75 homologs of mouse genes were down regulated in Library Hand DII-V compared with other two libraries (p < 0.001, Figure 4). For the comparison between forelimb digits and the hindlimb digits, we found that 187 homologs of human genes and 163 homologs of mouse genes were down regulated in Library Foot compared with other two libraries (p < 0.001, Figure 4) and 198 homologs of human genes and 84 homologs of mouse genes were up regulated in Library Foot compared with other two libraries (p < 0.001, Figure 4).
Among the differentially expressed genes, we found 14 genes likely related to digit elongation based on results of previous studies (see Discussion). The expression of Tbx3, Tbx15, Bmp3, Rgmb, Smad1, Smad4, Hoxd8, Hoxd9, Satb1, Hoxa1, Twist1, Tmeff2 and Enpp2 were significantly higher in Library Hand DII-V compared with other two libraries, and the expression of Nog was significantly lower in Library Hand DII-V (Table 2). Moreover, our results hightlighted 7 genes were likely related to the morphological and functional similarity between the thumb and the hindlimb digits. The levels of Tbx4, Pitx2, Acta1, Tnnc2, Atp2a1, Hrc and Myoz1 were significantly higher in Libraries Hand DI and Foot than in Library Hand DII-V (Table 2).
Real-time quantitative PCR (RT-qPCR)
To validate the results of DGE-tag profiling, we performed RT-qPCR of 14 genes which were differentially expressed between libraries (Table 2). All 14 genes exhibited the same overall expression pattern using the distinct methods (Figure 5). Meis2, Hoxd10, Tbx3, Tbx15, Bmp3, Rgmb, Hoxa1, Twist1, Tmeff2 and Enpp2 showed significantly higher expression in the posterior forelimb digits (HDII-V) compared with the thumb (HDI) and all five hindlimb digits (Foot). Krt17, Tbx4, Acta1 and Tnnc2 showed significantly higher expression in both the thumb (HDI) and all five hindlimb digits (Foot) than in the posterior forelimb digits (HDII-V). The fold change in expression among the samples was often lower in RT-qPCR as compared to DGE-tag profiling; however, this is not surprising given the differences between the two methods in sensitivity of gene detection as well as the fact that the methods were performed in different species.
Massively parallel sequencing is a powerful technology enabling extensive investigation of gene expression profiles. In this study, we obtained over 16.5 million sequences of expressed tags from three bat autopodial libraries, including the thumb (Hand DI), the posterior forelimb digits (Hand DII-V), and all five hindlimb digits (Foot). Using these data, we compared gene expression profiles among the libraries to identify candidate genes associated with elongation of the posterior forelimb digits (Hand DII-V). In addition, we examined gene expression patterns associated with the short digit morphology of the thumb and hindlimb digits as well as gene expression patterns associated with the fore- versus the hindlimb digits. The results of DGE-tag profiling were validated by examining gene expression in a subset of genes using RT-qPCR.
To assess the reliability of our gene expression profiles, we first examined expression of several genes with previously established expression patterns. We compared expression of Tbx5 and Pitx1, fore- and hindlimb marker genes, in the forelimb and hindlimb digits (Libraries Hand DI and Hand DII-V and Library Foot). Consistent with previous studies, we found that expression of Tbx5 was significantly higher in the forelimb digits and Pitx1 expression was significantly higher in the hindlimb digits although the embryonic stages used here are much later than those used in previous studies [19, 20]. In addition, we compared expression of Meis2 and Hoxd10 in the fore- and hindlimb digits. Previously published microarray data revealed high levels of expression of these two genes in bat autopodia at Stages 16 and 17 . We found that expression of both genes, Meis2 and Hoxd10, remained significantly higher in Library Hand DII-V in the Fetal Stage (Table 2). Finally, Krt17, a gene expressed in nails and important for the development of nails, exhibited significantly higher expression in the thumb and hindlimb digits . Because only the thumb (Hand DI) and the hindlimb digits have claws, we would expect higher expression of Krt17 in the thumb and hindlimb digits relative to the posterior forelimb digits (Hand DII-V) that make up the wing.
In addition to comparing our gene expression profiles to genes with known expression patterns, we further validated our results by performing RT-qPCR using 14 genes which were identified as being differentially expressed between elongated digits and short digits in DGE-tag profiling. With the goal of identifying genes related to wing formation across bat taxa, we collected samples from a different species which belongs to the other suborder of Chiroptera (H. armiger). We selected 14 genes likely important for wing formation or for the similarity between the thumb and foot digits (see below). All the 14 genes exhibited the same overall expression pattern in the RT-qPCR as in the DGE-tag profiling (Figure 5). This finding not only provides strong support for the DGE-tag profiling method used in this study, but also maintenance of gene expression patterns across these two phylogenetically distant species provides strong support for the role of these genes in digit elongation in the posterior forelimb digits or morphological and functional similarities between the thumb and hindlimb digits.
Overall, comparisons of gene expression profiles between digit morphologies and limbs identified hundreds of differentially expressed genes. Several interesting patterns have emerged from this data. Specifically, we highlight 21 genes likely related to wing formation or to morphological and functional similarities between thumb and hindlimb digits (Table 2). First, we found 14 genes that are likely associated with digit elongation in bats - two Tbx genes (Tbx3 and Tbx15), five genes from the BMP pathway (Bmp3, Rgmb, Smad1, Smad4 and Nog), four Homeobox genes (Hoxd8, Hoxd9, Satb1 and Hoxa1) and three other gene (Twist1, Tmeff2 and Enpp2) related to either digit malformation or cell proliferation. Next, we identified seven genes (Tbx4, Pitx2, Acta1, Tnnc2, Atp2a1, Hrc and Myoz1) that are likely associated with the morphological and functional similarities between the thumb and hindlimb digits. Expression patterns of these genes in the three distinct digit libraries as well as known function of these genes are described in detail below.
Genes associated with digit elongation
In general, Tbx3 and Tbx15 belong to T-box gene family and are involved in skeletogenesis [23, 24]. In digit formation specifically, Tbx3 is thought to modulate posterior digit identity through Shh and BMP signaling, and Tbx15 insufficiency causes Cousin Syndrome which includes congenital dwarfism, moderate brachydactyly (or shortening of the digits) and leg shortening [25–27]. Mutants of Tbx15 have reduced bone size and changes in bone shape in the forelimb skeleton indicating an important role in the skeletal development of forelimbs .
Bmp3, Rgmb, Smad1, Smad4 and Nog are all components of the BMP pathway, a pathway critical for bone formation . Bmp3 belongs to the transforming growth factor-beta (TGFβ) superfamily and stimulates cell proliferation and differentiation in a concentration dependent fashion . Overexpression of Bmp3 in the chick wing bud has be shown to result in an increase in cell proliferation resulting in lengthening of skeletal elements . Rgmb, also named as Dragon, enhances BMP signaling by directly binding to BMP proteins including Bmp2 . Overexpression of Bmp2 and Bmp4 in the developing chick limb leads to a dramatic increase in the volume of cartilage elements . In addition, Bmp2/4 has been shown to be up regulated in the bat forelimb embryonic digits relative to both mouse forelimb digits and bat hindlimb digits . This up-regulation is correlated with increased rate of cartilage proliferation . Beyond up-regulation of BMPs themselves, Smad1 mediates BMP signaling. Specifically, in its phosphorylated state Smad1 forms a complex with Smad4 and then accumulates in the nucleus. This SMAD complex regulates transcription of several target genes . Embryonic manipulation of Smad1 has confirmed that its expression at the phalanx-forming region influences digit identity in the chick limb . In addition, cartilage-specific knockout of Smad1 and Smad4 reduce chondrocyte proliferation and increase apoptosis resulting in defects of limb elements in the mouse [34, 35]. While all the components of the BMP pathway mentioned above stimulate cartilage formation, Nog is an antagonist of BMPs and as such inhibits cartilage formation . Mutations of Nog disturb the balance of BMP signaling leading to brachydactyly . Thus, our results - exhibiting up regulation of Bmp3, Rgmb, Smad1, Smad4 and down regulation of Nog in the elongated bat forelimb digits - suggest a role for BMP signaling in the lengthening of the posterior bat forelimb digits.
Hoxd8, Hoxd9, Hoxa1 and Satb1 are homeobox genes which specify the anterior-posterior axis and have been implicated as digit identity determining genes . Although Hoxd8 and Hoxd9 mainly contribute to proximal limb development and are only mildly expressed in the developing digits, we found that Hoxd8 and Hoxd9 were highly expressed in bat wing digits II-V relative to the thumb and hindlimb digits . Our results also indicate that the homeobox gene Satb1 and Hoxa1 were highly expressed in bat wing digits II-V relative to the thumb and hindlimb digits. These two genes are highly expressed in cancer cells promoting tumor growth [40, 41]. To our knowledge, there are no studies examining expression of Satb1 and Hoxa1 in limb development making them interesting genes for future investigation.
Similarly to Satb1 and Hoxa1, Tmeff2 and Enpp2 are up regulated in carcinomas [42, 43]. Previous study showed Tmeff2 contributes to cell proliferation. Enpp2 which is expressed in precartilaginous condensations can be strongly activated by Hoxa13 and Hoxd13 which are crucial factors for digit identity and development . Finally, knowledge of the role of Twist1 and the digit development comes from pathological studies . Mutations Twist1 have been found in patients with Saethre-Chotzen syndrome which is associated with digit malformation . Overexpression of Twist1 inhibits osteoblast differentiation and prevents premature fusion of skeletons .
Genes associated with similarities between the thumb and hindlimb digits
In addition to genes related to bat digit elongation, we found several genes that may contribute to the morphological and functional similarities between the bat thumb and hindlimb digits. First, Tbx4 and Pitx2 are primarily associated with hindlimb development. Tbx4 has also been shown to regulate bone size. Specifically, in the developing mouse forelimb, overexpression of these genes results in shortening of the forelimb skeletons [47, 48]. Next, we found several genes that are likely associated with the distinct musculature and function of the short digits. Because the clawed thumb of bats is used to cling to the cave while roosting and the clawed hindlimb digits are used for hanging upside down these digits should be richer in skeletal muscles than the elongated digits. In fact, expression of several skeletal muscle associated genes - Acta1, Tnnc2, Atp2a1, Hrc and Myoz1 - was significantly higher in the thumb and hindlimb digits then in the elongated posterior forelimb digits .
Using genome-wide sequencing method we have identified many genes differentially expressed in the developing digits of the bat autopodia. Several of these genes are highly related to digit formation and thus likely play important roles in the development of bat wing. Our results provide robust candidate genes for future functional tests and molecular evolution comparison with other mammals to fully understand the mechanisms of wing formation and the evolution of bat flight.
Animal collection and breeding
All procedures were carried out in accordance with the Policy on the Care and Use of Animals, approved by the Ethical Committee, East China Normal University. Myotis ricketti, Rickett's big-footed bat, belongs to suborder Yangochiroptera and is widely distributed in China . Their feet - which are used to forage on fishes - are larger than those of other bat species (Figure 1A) . As in many species of bats, copulation in Myotis ricketti occurs in autumn and spermatozoa are stored in the female's uterus and oviducts through the winter hibernation. Ovulation and subsequent fertilization begin in the spring with parturition occurring in summer . Because gravid M. ricketti could not be found during the summer breeding season, we capitalized on sperm storage by females and captured 10 hibernating female M. ricketti from a cave (39°42'N, 115°43'E) in Beijing on February 13, 2008. In addition to Myotis ricketti, we captured three gravid Hipposideros armiger from a cave (30°20'N, 117°50'E) in Anhui province of China on April 24, 2009 to use for quantitative real-time PCR confirmation sequencing results. H. armiger, the great leaf-nosed bat, belongs to the other suborder Yinpterochiroptera and is widely distributed in subtropical and tropical zones of Asia .
After capture, bats were kept in a large cage (1.2 × 1 × 1 m) covered with wire netting allowing them to hang. The cage was kept in a dark, temperature controlled room. The temperature was maintained between 18 and 24°C. Water was freely available and meal worms mixed with a powdered multivitamin and calcium tablets were provided daily.
M. ricketti in the captive colony were euthanized by decapitation on May 28 and June 6, 2008 (Table 1). H. armiger were euthanized on May 20, 2009 (Table 1). Fetuses from the gravid females were put into ice cold PBS. The thumb, the remaining forelimb digits (metacarpi, phalanges and dactylopatagium) and all five hindlimb digits (metatarsi and phalages) were dissected from the fetuses (Figure 1B and 1C) and stored in liquid nitrogen until use. All other components of the fetus were fixed overnight in Bouin's fluid, washed with several changes of 70% ethyl alcohol and stored at room temperature. Crown rump length (CRL) and forearm length (FA) were measured for further estimation of the embryonic stages.
Digital gene expression tag profiling (DGE-tag profiling)
To reduce sex-specific gene expression patterns, we pooled the samples from a male and a female fetus of M. ricketti (Bat M1 and Bat M2 respectively, see Table 1). Tag libraries were prepared using the DGE-Tag Profiling for NlaIII Sample Prep kit from Illumina according to the manufacturer's instructions. In brief, total RNA was extracted from forelimb digit I, forelimb digits II-V and all five hindlimb digits using TRIzol® Reagent (Invitrogen, Cat. No. 15596026, US), treated with DNase I (Roche, Cat. No. 04716728001, Germany) and converted to double stranded cDNA using Oligo dT beads. Subsequently, the cDNA samples were digested using the restriction enzyme NlaIII, which recognizes and cuts the most 3' "CATG". Once digested into fragments, cDNA samples were ligated to Illumina specific adapter A, which contains a recognition site for enzyme MmeI. The Illumina specific adapter B was ligated after MmeI digestion (Illumina, Cat. No. FC-102-1005, US). Using this technique, three DGE libraries (Libraries Hand DI, Hand DII-V and Foot) from the three samples (Figure 1B and 1C) were constructed such that each molecule in the library is a 21 bp tag derived from a single transcript with Illumina specific adapters attached to both ends. Next, the 21 bp tags were enriched using PCR primers that anneal to the adaptors and clusters were created on a flow cell using a Solexa Cluster Station following the manufacturer's instructions (Illumina, Cat. No. FC-103-1004, US). Finally, massively parallel sequencing-by-synthesis was performed on the Illumina Genome Analyzer. Image analysis, base calling, extraction of 21 bp tags, and tag counting were performed using the Illumina pipeline. The number of times that a unique tag sequence is detected represents the quantitative expression level of the corresponding transcript in the tissue. The raw data, including tag sequences and counts, were deposited in NCBI Gene Expression Omnibus (GEO) database with the accession number [GEO: GSE20038]
In silico analysis
We analyzed saturation of each of the three libraries by examining the addition of new distinct tags associated with the increase of sequencing depth (frequency = 1, > 1 and all distinct tags). Because little gene annotation exists for the bat genome, all distinct tags were aligned against the well annotated NCBI human and mouse Refseq cDNA sequences using SOAP2 [54, 55]. Taking into account the differences between species, one base pair mismatches were allowed during the mapping of tags to cDNA sequences. Tags that matched to more than one gene were not used in analyzing differential expression among the libraries. The identified cDNA sequences were mapped to NCBI Entrez GeneID. Differential gene expression analysis was performed using DEGseq software . We compared libraries Hand DII-V with Hand DI and Foot to identify differentially expressed genes between elongated digits and short digits. We also compared libraries Hand DII-V and Hand DI with Foot to identify differentially expressed genes between the forelimb and hindlimb digits.
Real-time quantitative PCR (RT-qPCR)
Total RNA were extracted from forelimb digit I, forelimb digits II-V and all five hindlimb digits of two H. armiger (Bat H1 and H2, Table 1), using TRIzol® Reagent (Invitrogen, Cat. No. 15596026, US) and treated with DNase I (Roche, Cat. No. 04716728001, Germany). cDNA synthesis and qPCR was performed using the SYBR® PrimeScript® RT-PCR Kit (TaKaRa, Cat. No. DRR063A, Japan) on a Applied Biosystems 7300 Real-Time PCR System (ABI, US). RT-qPCR for each H. armiger (Bat H1 and H2, Table 1) was performed in duplicate. β-actin was taken as a reference gene and it was amplified by a pair of previously published primers . We designed 14 pairs of degenerate primers to amplify 14 target genes (Meis2, Hoxd10, Tbx3, Tbx15, Bmp3, Rgmb, Twist1, Hoxa1, Tmeff2, Enpp2, Krt17, Tbx4, Acta1 and Tnnc2). PCR products of each gene were ligated into a pMD®19-T vector (TaKaRa, Cat. No. D102B, Japan), cloned, and sequenced using the Big Dye Terminator kit on an ABI 3730 DNA sequencer (Applied Biosystems, US). Using the obtained sequences, we designed gene specific primers for each target gene for qPCR. 2-ΔΔCT method was used to quantify gene expression, data are expressed as mean ± SD . Student's t-test was performed to examine the difference of gene expression between two different samples. Sequences of forward and reverse primers for normal PCR and RT-qPCR of each gene are shown in additional file 4, Table S4. The partial CDS of the genes were deposited at GenBank under accession number HM777024 - HM777038.
Crown rump length
Hunter P: The nature of flight - The molecules and mechanics of flight in animals. Embo Reports. 2007, 8 (9): 811-813. 10.1038/sj.embor.7401050.
Hedenstrom A, Johansson LC, Spedding GR: Bird or bat: comparing airframe design and flight performance. Bioinspir Biomim. 2009, 4 (1): 015001-10.1088/1748-3182/4/1/015001.
Cretekos CJ, Weatherbee SD, Chen CH, Badwaik NK, Niswander L, Behringer RR, Rasweiler JJ: Embryonic staging system for the short-tailed fruit bat, Carollia perspicillata, a model organism for the mammalian order Chiroptera, based upon timed pregnancies in captive-bred animals. Dev Dyn. 2005, 233 (3): 721-738. 10.1002/dvdy.20400.
Hockman D, Mason MK, Jacobs DS, Illing N: The Role of Early Development in Mammalian Limb Diversification: A Descriptive Comparison of Early Limb Development Between the Natal Long-Fingered Bat (Miniopterus natalensis) and the Mouse (Mus musculus). Developmental Dynamics. 2009, 238 (4): 965-979. 10.1002/dvdy.21896.
Sears KE: Molecular determinants of bat wing development. Cells Tissues Organs. 2008, 187 (1): 6-12. 10.1159/000109959.
Cooper KL, Tabin CJ: Understanding of bat wing evolution takes flight. Gene Dev. 2008, 22 (2): 121-124. 10.1101/gad.1639108.
Behringer RR, Rasweiler JJt, Chen CH, Cretekos CJ: Genetic regulation of Mammalian diversity. Cold Spring Harb Symp Quant Biol. 2009, 74: 297-302.
Cretekos CJ, Wang Y, Green ED, Martin JF, Rasweiler JJ, Behringer RR, Progra NCS: Regulatory divergence modifies limb length between mammals. Genes & Development. 2008, 22 (2): 141-151.
Hockman D, Cretekos CJ, Mason MK, Behringer RR, Jacobs DS, Illing N: A second wave of Sonic hedgehog expression during the development of the bat limb. Proc Natl Acad Sci USA. 2008, 105 (44): 16982-16987. 10.1073/pnas.0805308105.
Cretekos CJ, Deng JM, Green ED, Rasweiler JJ, Behringer RR: Isolation, genomic structure and developmental expression of Fgf8 in the short-tailed fruit bat, Carollia perspicillata. Int J Dev Biol. 2007, 51 (4): 333-338. 10.1387/ijdb.062257cc.
Sears KE, Behringer RR, Rasweiler JJ, Niswander LA: Development of bat flight: morphologic and molecular evolution of bat wing digits. Proc Natl Acad Sci USA. 2006, 103 (17): 6581-6586. 10.1073/pnas.0509716103.
Chen CH, Cretekos CJ, Rasweiler JJ, Behringer RR: Hoxd13 expression in the developing limbs of the short-tailed fruit bat, Carollia perspicillata. Evol Dev. 2005, 7 (2): 130-141. 10.1111/j.1525-142X.2005.05015.x.
Ray R, Capecchi M: An examination of the Chiropteran HoxD locus from an evolutionary perspective. Evol Dev. 2008, 10 (6): 657-670. 10.1111/j.1525-142X.2008.00279.x.
Velculescu VE, Zhang L, Vogelstein B, Kinzler KW: Serial analysis of gene expression. Science. 1995, 270 (5235): 484-487. 10.1126/science.270.5235.484.
't Hoen PA, Ariyurek Y, Thygesen HH, Vreugdenhil E, Vossen RH, de Menezes RX, Boer JM, van Ommen GJ, den Dunnen JT: Deep sequencing-based expression analysis shows major advances in robustness, resolution and inter-lab portability over five microarray platforms. Nucleic Acids Res. 2008, 36 (21): e141-10.1093/nar/gkn705.
Marioni JC, Mason CE, Mane SM, Stephens M, Gilad Y: RNA-seq: an assessment of technical reproducibility and comparison with gene expression arrays. Genome Res. 2008, 18 (9): 1509-1517. 10.1101/gr.079558.108.
Feng L, Liu H, Liu Y, Lu Z, Guo G, Guo S, Zheng H, Gao Y, Cheng S, Wang J: Power of deep sequencing and agilent microarray for gene expression profiling study. Mol Biotechnol. 45 (2): 101-110. 10.1007/s12033-010-9249-6.
Wang Z, Han N, Racey PA, Ru B, He G: A comparative study of prenatal development in Miniopterus schreibersii fuliginosus, Hipposideros armiger and H. pratti. BMC Dev Biol. 2010, 10 (1): 10-10.1186/1471-213X-10-10.
Logan M, Tabin CJ: Role of Pitx1 upstream of Tbx4 in specification of hindlimb identity. Science. 1999, 283 (5408): 1736-1739. 10.1126/science.283.5408.1736.
Rodriguez-Esteban C, Tsukui T, Yonei S, Magallon J, Tamura K, Izpisua Belmonte JC: The T-box genes Tbx4 and Tbx5 regulate limb outgrowth and identity. Nature. 1999, 398 (6730): 814-818. 10.1038/19769.
Mason M, Hockman D, Jacobs D, Illing N: Differences in the wing and hindlimb transcriptomes of the natal long-fingered bat, Miniopterus natalensis, during embryonic development. Mech Dev. 2009, 126: S17-03. 10.1016/j.mod.2009.06.1015.
Mclean WHI, Rugg EL, Lunny DP, Morley SM, Lane EB, Swensson O, Doppinghepenstal PJC, Griffiths WAD, Eady RAJ, Higgins C, et al: Keratin-16 and Keratin-17 Mutations Cause Pachyonychia-Congenita. Nature Genetics. 1995, 9 (3): 273-278. 10.1038/ng0395-273.
King M, Arnold JS, Shanske A, Morrow BE: T-genes and limb bud development. Am J Med Genet A. 2006, 140 (13): 1407-1413.
Singh MK, Petry M, Haenig B, Lescher B, Leitges M, Kispert A: The T-box transcription factor Tbx15 is required for skeletal development. Mech Dev. 2005, 122 (2): 131-144. 10.1016/j.mod.2004.10.011.
Nissim S, Allard P, Bandyopadhyay A, Harfe BD, Tabin CJ: Characterization of a novel ectodermal signaling center regulating Tbx2 and Shh in the vertebrate limb. Dev Biol. 2007, 304 (1): 9-21. 10.1016/j.ydbio.2006.12.010.
Suzuki T, Takeuchi J, Koshiba-Takeuchi K, Ogura T: Tbx genes specify posterior digit identity through Shh and BMP signaling. Dev Cell. 2004, 6 (1): 43-53. 10.1016/S1534-5807(03)00401-5.
Lausch E, Hermanns P, Farin HF, Alanay Y, Unger S, Nikkel S, Steinwender C, Scherer G, Spranger J, Zabel B, et al: TBX15 mutations cause craniofacial dysmorphism, hypoplasia of scapula and pelvis, and short stature in Cousin syndrome. Am J Hum Genet. 2008, 83 (5): 649-655. 10.1016/j.ajhg.2008.10.011.
Anderson GJ, Darshan D: Small-molecule dissection of BMP signaling. Nat Chem Biol. 2008, 4 (1): 15-16. 10.1038/nchembio0108-15.
Carrington JL, Chen P, Yanagishita M, Reddi AH: Osteogenin (bone morphogenetic protein-3) stimulates cartilage formation by chick limb bud cells in vitro. Dev Biol. 1991, 146 (2): 406-415. 10.1016/0012-1606(91)90242-U.
Gamer LW, Ho V, Cox K, Rosen V: Expression and function of BMP3 during chick limb development. Dev Dyn. 2008, 237 (6): 1691-1698. 10.1002/dvdy.21561.
Samad TA, Rebbapragada A, Bell E, Zhang Y, Sidis Y, Jeong SJ, Campagna JA, Perusini S, Fabrizio DA, Schneyer AL, et al: DRAGON, a bone morphogenetic protein co-receptor. J Biol Chem. 2005, 280 (14): 14122-14129. 10.1074/jbc.M410034200.
Duprez D, Bell EJ, Richardson MK, Archer CW, Wolpert L, Brickell PM, Francis-West PH: Overexpression of BMP-2 and BMP-4 alters the size and shape of developing skeletal elements in the chick limb. Mech Dev. 1996, 57 (2): 145-157. 10.1016/0925-4773(96)00540-0.
Suzuki T, Hasso SM, Fallon JF: Unique SMAD1/5/8 activity at the phalanx-forming region determines digit identity. Proc Natl Acad Sci USA. 2008, 105 (11): 4185-4190. 10.1073/pnas.0707899105.
Zhang J, Tan X, Li W, Wang Y, Wang J, Cheng X, Yang X: Smad4 is required for the normal organization of the cartilage growth plate. Dev Biol. 2005, 284 (2): 311-322. 10.1016/j.ydbio.2005.05.036.
Retting KN, Song B, Yoon BS, Lyons KM: BMP canonical Smad signaling through Smad1 and Smad5 is required for endochondral bone formation. Development. 2009, 136 (7): 1093-1104. 10.1242/dev.029926.
Groppe J, Greenwald J, Wiater E, Rodriguez-Leon J, Economides AN, Kwiatkowski W, Affolter M, Vale WW, Belmonte JC, Choe S: Structural basis of BMP signalling inhibition by the cystine knot protein Noggin. Nature. 2002, 420 (6916): 636-642. 10.1038/nature01245.
Lehmann K, Seemann P, Silan F, Goecke TO, Irgang S, Kjaer KW, Kjaergaard S, Mahoney MJ, Morlot S, Reissner C, et al: A new subtype of brachydactyly type B caused by point mutations in the bone morphogenetic protein antagonist NOGGIN. Am J Hum Genet. 2007, 81 (2): 388-396. 10.1086/519697.
Montavon T, Le Garrec JF, Kerszberg M, Duboule D: Modeling Hox gene regulation in digits: reverse collinearity and the molecular origin of thumbness. Genes Dev. 2008, 22 (3): 346-359. 10.1101/gad.1631708.
Zakany J, Duboule D: The role of Hox genes during vertebrate limb development. Curr Opin Genet Dev. 2007, 17 (4): 359-366. 10.1016/j.gde.2007.05.011.
Zhang X, Emerald BS, Mukhina S, Mohankumar KM, Kraemer A, Yap AS, Gluckman PD, Lee KO, Lobie PE: HOXA1 is required for E-cadherin-dependent anchorage-independent survival of human mammary carcinoma cells. J Biol Chem. 2006, 281 (10): 6471-6481. 10.1074/jbc.M512666200.
Han HJ, Russo J, Kohwi Y, Kohwi-Shigematsu T: SATB1 reprogrammes gene expression to promote breast tumour growth and metastasis. Nature. 2008, 452 (7184): 187-193. 10.1038/nature06781.
Glynne-Jones E, Harper ME, Seery LT, James R, Anglin I, Morgan HE, Taylor KM, Gee JM, Nicholson RI: TENB2, a proteoglycan identified in prostate cancer that is associated with disease progression and androgen independence. Int J Cancer. 2001, 94 (2): 178-184. 10.1002/ijc.1450.
Nouh MA, Wu XX, Okazoe H, Tsunemori H, Haba R, Abou-Zeid AM, Saleem MD, Inui M, Sugimoto M, Aoki J, et al: Expression of autotaxin and acylglycerol kinase in prostate cancer: association with cancer development and progression. Cancer Sci. 2009, 100 (9): 1631-1638. 10.1111/j.1349-7006.2009.01234.x.
Williams TM, Williams ME, Kuick R, Misek D, McDonagh K, Hanash S, Innis JW: Candidate downstream regulated genes of HOX group 13 transcription factors with and without monomeric DNA binding capability. Dev Biol. 2005, 279 (2): 462-480. 10.1016/j.ydbio.2004.12.015.
Firulli BA, Krawchuk D, Centonze VE, Vargesson N, Virshup DM, Conway SJ, Cserjesi P, Laufer E, Firulli AB: Altered Twist1 and Hand2 dimerization is associated with Saethre-Chotzen syndrome and limb abnormalities. Nat Genet. 2005, 37 (4): 373-381. 10.1038/ng1525.
Bialek P, Kern B, Yang X, Schrock M, Sosic D, Hong N, Wu H, Yu K, Ornitz DM, Olson EN, et al: A twist code determines the onset of osteoblast differentiation. Dev Cell. 2004, 6 (3): 423-435. 10.1016/S1534-5807(04)00058-9.
Marcil A, Dumontier E, Chamberland M, Camper SA, Drouin J: Pitx1 and Pitx2 are required for development of hindlimb buds. Development. 2003, 130 (1): 45-55. 10.1242/dev.00192.
Holmberg J, Ingner G, Johansson C, Leander P, Hjalt TA: PITX2 gain-of-function induced defects in mouse forelimb development. BMC Dev Biol. 2008, 8: 25-10.1186/1471-213X-8-25.
Welle S, Brooks AI, Delehanty JM, Needler N, Bhatt K, Shah B, Thornton CA: Skeletal muscle gene expression profiles in 20-29 year old and 65-71 year old women. Exp Gerontol. 2004, 39 (3): 369-377. 10.1016/j.exger.2003.11.011.
Ma J, Dai Q, Zhang SY, Shen JX, Liang B: Distribution of Ricketti's big-footed bat (Myotis ricketti). Sichuan J Zool. 2003, 22 (3): 155-156.
Ma J, Zhang J, Liang B, Zhang L, Zhang S, Metzner W: Dietary characteristics of Myotis ricketti in Beijing, North China. J Mammal. 2006, 87: 339-344. 10.1644/05-MAMM-A-183R1.1.
Wang Z, Liang B, Racey PA, Wang Y, Zhang S: Sperm storage, delayed ovulation and menstruation of the female Rickett's big-footed bat (Myotis ricketti). Zool Stud. 2008, 47: 215-221.
Simmons NB: Order Chiroptera. Mammal species of the world: a taxonomic and geographic reference. Edited by: Wilson DE, Reeder DM. 2005, Baltimore: Johns Hopkins University, 312-529.
Li RQ, Yu C, Li YR, Lam TW, Yiu SM, Kristiansen K, Wang J: SOAP2: an improved ultrafast tool for short read alignment. Bioinformatics. 2009, 25 (15): 1966-1967. 10.1093/bioinformatics/btp336.
Pruitt KD, Maglott DR: RefSeq and LocusLink: NCBI gene-centered resources. Nucleic Acids Res. 2001, 29 (1): 137-140. 10.1093/nar/29.1.137.
Wang L, Feng Z, Wang X, Zhang X: DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010, 26 (1): 136-138. 10.1093/bioinformatics/btp612.
Chen J, Sun M, Liang B, Xu A, Zhang S, Wu D: Cloning and expression of PDK4, FOXO1A and DYRK1A from the hibernating greater horseshoe bat (Rhinolophus ferrumequinum). Comp Biochem Physiol B Biochem Mol Biol. 2007, 146 (2): 166-171. 10.1016/j.cbpb.2006.10.095.
Schmittgen TD, Livak KJ: Analyzing real-time PCR data by the comparative C(T) method. Nat Protoc. 2008, 3 (6): 1101-1108. 10.1038/nprot.2008.73.
We thank Junpeng Zhang, Yinan Wang and Qian Yao for bat capture and care. This work was funded by a grant awarded to S. Zhang under the Key Construction Program of the National "985" Project and "211" Project.
ZW participated in design of the study. ZW and DD carried out the data analysis and wrote the manuscript. ZW, BR and TG performed the experiments. RLY contributed to the final manuscript preparation. NH fed and examined the bats. SZ conceived of the study, participated in its design, and supervised the work. All authors read and approved the final manuscript.