Analysis of Molecular Cytogenetic Alteration in Rhabdomyosarcoma by Array Comparative Genomic Hybridization

Rhabdomyosarcoma (RMS) is the most common pediatric soft tissue sarcoma with poor prognosis. The genetic etiology of RMS remains largely unclear underlying its development and progression. To reveal novel genes more precisely and new therapeutic targets associated with RMS, we used high-resolution array comparative genomic hybridization (aCGH) to explore tumor-associated copy number variations (CNVs) and genes in RMS. We confirmed several important genes by quantitative real-time polymerase chain reaction (QRT-PCR). We then performed bioinformatics-based functional enrichment analysis for genes located in the genomic regions with CNVs. In addition, we identified miRNAs located in the corresponding amplification and deletion regions and performed miRNA functional enrichment analysis. aCGH analyses revealed that all RMS showed specific gains and losses. The amplification regions were 12q13.12, 12q13.3, and 12q13.3–q14.1. The deletion regions were 1p21.1, 2q14.1, 5q13.2, 9p12, and 9q12. The recurrent regions with gains were 12q13.3, 12q13.3–q14.1, 12q14.1, and 17q25.1. The recurrent regions with losses were 9p12–p11.2, 10q11.21–q11.22, 14q32.33, 16p11.2, and 22q11.1. The mean mRNA level of GLI1 in RMS was 6.61-fold higher than that in controls (p = 0.0477) by QRT-PCR. Meanwhile, the mean mRNA level of GEFT in RMS samples was 3.92-fold higher than that in controls (p = 0.0354). Bioinformatic analysis showed that genes were enriched in functions such as immunoglobulin domain, induction of apoptosis, and defensin. Proto-oncogene functions were involved in alveolar RMS. miRNAs that located in the amplified regions in RMS tend to be enriched in oncogenic activity (miR-24 and miR-27a). In conclusion, this study identified a number of CNVs in RMS and functional analyses showed enrichment for genes and miRNAs located in these CNVs regions. These findings may potentially help the identification of novel biomarkers and/or drug targets implicated in diagnosis of and targeted therapy for RMS.


Introduction
Rhabdomyosarcoma (RMS) is the most common soft tissue sarcoma in children, which has several subtypes including the more aggressive alveolar RMS (ARMS), the more prevalent embryonal RMS (ERMS), and the rare adult variant pleomorphic RMS (PRMS) [1]. Tumorigenesis for some RMSs is recognized, for example, the majority of ARMS tumors (about 85%) are characterized by recurrent translocation between genes encoding for transcription factors FKHR with either PAX3 or PAX7 [2]. The complete genetic etiology underlying RMS development and progression remains unclear.
Array comparative genomic hybridization (aCGH) is a technique that was developed for high-resolution, genome-wide screening of segmental genomic copy number variations [3,4]. aCGH allows for comprehensive interrogation of hundreds of genomic loci for DNA copy number gains and losses. For the large amount of data generated by high-resolution aCGH, in order to avoid random events of no biologic significance, researchers could deal with the data using various methods, for example GISTIC and waviCGH [5,6]. DNA copy number changes are common in cancer, and lead to altered expression and function of genes residing within the affected region of the genome. Identification of regions with copy number aberrations, as well as the genes involved, offers a basis for a better understanding of cancer development to provide improved tools for clinical management of cancer, such as new diagnostics and therapeutic targets [7]. Thus, detection of genomic imbalances and identification of these genes can elucidate RMS biology and help identify novel potential biomarkers and targets for clinical therapy.
Traditionally, microarray-based, high-throughput experiments (such as aCGH) produce massive gene lists without consideration of potential relationships among these genes. The gene-by-gene approach often lacks a coherent picture of disease-related pathologic interactions. Bioinformatics has attracted increasing interest in potential gene discovery. For an uploaded gene list, the DAVID bioinformatics resources [8] provide typical gene term enrichment analysis and tools that allow users to condense large gene lists into gene functional groups, visualize many-genes-tomany-terms relationships, categorize redundant and heterogeneous terms into groups, search for interesting and related genes or terms, dynamically view genes from their lists on biopathways, and other functions.
In addition to protein-coding genetic factors, microRNAs (miRNAs) are emerging as key non-protein-coding factors that affect the regulation of gene expression. Increasing evidence suggests that miRNAs participate in nearly all important biological processes, and miRNA dysfunctions are associated with various diseases [9]. Analyses of several human cancers have identified miRNA signatures associated with initiation, progression, diagnosis, or prognosis of tumors [10]. In the present study, high-resolution aCGH was used to identify the potential alterations that were involved in RMS pathogenesis. Genes and miRNAs that located in the altered genomic regions were identified. Finally, tools of DAVID [8] and TAM [11] were used to perform functional enrichment analysis for the identified genes and miRNAs, respectively.

Ethics Statement
Written informed consent was obtained from all participating patients before enrollment in the study. This study was approved by the institutional ethics committee at the First Affiliated Hospital of Shihezi University School of Medicine and conducted in accordance with the ethical guidelines of the Declaration of Helsinki. aCGH Isolation of genomic DNA (gDNA) from tumor tissues was completed using QIAamp DNA FFPE tissue kit following manufacturer protocols (Qiagen, Germany). The gDNA from the cell lines was isolated using the DNeasy blood and tissue kit (Qiagen, Germany). aCGH experiments were performed using standard NimbleGen protocols (NimbleGen Arrays User's Guide: CGH Analysis v5.1). We used pooled male and female reference gDNA provided by NimbleGen for comparison of male and female patient DNA samples. Tumor DNA fragments and digested references were labeled with Cy3 and Cy5, respectively. DNA was combined with 40 mL of diluted Cy3-random nonamers and water to a total volume of 80 mL, denatured at 98uC for 10 min, and immediately cooled on ice. This solution was then combined with 20 mL of dNTP/Klenow Master Mix, which was mixed well by pipetting up and down, and incubated for 2 h at 37uC. After clean-up, reference and tumor DNA probes were mixed and hybridized onto a Roche NimbleGen 36720 K probe platform for 48 h at 42uC. This platform has genome-wide probe spacing at approximately every 2509 bp. After hybridization, array slides were washed and dried in a NimbleGen Microarray Dryer.

Data processing and analysis
The arrays were scanned using MS200 scanner (NimbleGen) with 2 mm resolution, and fluorescent intensity data was extracted with NimbleScan 2.6 software. The hybridization controls (STC, Sample Tracking Controls) which were used to confirm that the correct samples were hybridized to each array.
For each spot on the array, the log 2 ratios of the Cy3-labeled test sample versus Cy-5 reference sample were computed. Before normalization and segmentation analysis, spatial correction was applied, which corrected position-dependent non-uniformity of signals across the array. Specifically, locally weighted polynomial regression (LOESS) is used to adjust signal intensities based on X, Y feature position [12]. Normalization was then performed using the q-spline method [13], compensated for inherent differences in signal between the two dyes, followed by segmentation using the CNVs calling algorithm segMNT [14]. The segMNT algorithm identified CNVs using a dynamic programming process that minimizes the squared error relative to the segment means, which showed better performance than the DNA copy algorithm [15]. The segments with |mean log 2 ratio| $0.25 and at least 5 consecutive probes were retained, |log 2 ratio|,0.25 represented ''unchanged''. Mean log 2 ratios of all probes in a chromosome region $0.25 were classified as genomic gains, and mean log 2 ratios $1.0 were classified as amplification. Meanwhile, mean log 2 ratios #-0.25 were regarded as losses, and mean log 2 ratios #-1.0 were regarded as deletions [16].

Real-time quantitative polymerase chain reaction (QRT-PCR)
QRT-PCR was used to validate the results of aCGH. The total RNA was extracted from RMS tissues and 14 cases normal muscle tissues as controls using RNeasy FFPE Kit (QIAGEN), and the total RNAs from all the samples were treated with DNAse I, and transcribed to single-stranded cDNA by reverse transcription using QuantiTect Reverse Transcription Kit (QIAGEN). A house keeping mRNA, ACTB (Hs_ACTB_2_SGQuantiTect Primer Assay, QT01680476), was assessed in all samples. The normal muscle tissues were used as control. The Primers GLI1 (QT00060501) and GEFT (QT00202916) genes were from QuantiTect Primer Assays (QIAGEN). The reaction was carried out on the ABI 7500 Real-Time PCR thermocycler (Applied Biosystems) using Quantifast SYBR Green PCR Kit (QIAGEN). The following thermal cycling program was applied: 5 min at 95uC, 40 cycles of 10 s at 95uC and 30 s at 60uC. Data were normalized for ACTB expression using comparative threshold cycle method. All PCRs were done in triplicates. Cycle threshold (Ct), the fractional cycle number at which the amount of amplified target reached a fixed threshold, was determined. DCt values were calculated by subtracting the ACTB Ct values from the target gene Ct values (DCt = Ct (GLI1 or GEFT gene in RMS/normal muscle sample) -Ct (ACTB gene in RMS/normal muscle sample)). Expression level was determined as 2-DCt.

Genomic map of the aberrant regions in chromosomes
A genomic map of the aberrant regions was created using Circos [17], which is a software package for displaying genomic data. Circos is a command-line application written in Perl, which can be easily deployed on any system for which Perl is available (http://www.cpan.org/ports/). Inputs were GFF-style data files and Apache-like configuration files, both of which can be easily generated by automated tools.

DAVID analysis
DAVID (which can be freely accessed at http://david.abcc. ncifcrf.gov/) is a web-based online bioinformatics resource that provides tools for the functional interpretation of large lists of genes/proteins. Amplification or deletion genes were subjected to separate cluster analyses. Each gene set was entered into the DAVID functional annotation clustering tool, which generated clusters of genes based on the similarity of the functional terms assigned to each gene. The clusters were then ranked according to scores of each term, with the higher ranked clusters selected for analysis. Within each cluster, the lowest P value (P value ,0.05) was selected as a representative functional term.

TAM analysis
The genomic coordinate data of miRNAs were downloaded from the miRBase (Release 19). The miRNAs were mapped to their corresponding amplification and deletion regions by an inhouse Java program. TAM (http://cmbi.bjmu.edu.cn/tam) was used to identify the enriched functions for the miRNAs within the above regions. Within each cluster, the lowest P value (P value ,0.05) was selected as a representative functional term.

Statistical analysis
SPSS software package (Version 17, Chicago, IL) was used for statistical analyses. Independent sample t test was used to evaluate differences in mRNA expression of GLI1 or GEFT between RMS groups and normal muscle tissues. Differences with a p value of ,0.05 were considered statistically significant.

Genomic profiles of RMS
aCGH analysis was carried out to identify genomic alterations in 20 RMS cases, and every RMS tumor displayed copy number changes with gains and losses. Figure 1 shows genomic map of the aberrant regions in human RMS chromosomes. The recurrent region was defined as a region that had a frequency of over 50% of cases. Based on this definition, recurrent regions of gain were 12q13.3, 12q13.3-q14.1, 12q14.1, and 17q25.1. Recurrent regions of loss were 2q14.1, 9p12-p11.2, 9q12, 14q32.33, 16p11.2, and 22q11.1 (Table 2). Then, we listed the frequency of over 20% of cases in regional amplification or deletion on chromosomes ( Table 3). The regions of amplification were

Recurrent copy number alterations in ERMS and ARMS detected by aCGH
We analyzed copy number alterations by examining chromosomal changes in ERMS and ARMS tumors. Table 4

Recurrent copy number alterations in RMS cell lines detected by aCGH
Chromosome imbalance was detected in the RMS cell line by aCGH. Figure 3A and B show the genomic maps of the aberrant regions in PLA-802 (ARMS) and RD (ERMS), respectively. As shown in Figure 3, frequent chromosomal changes were observed

GLI1 and GEFT mRNA Expression in RMS tissues and cell lines by QRT-PCR
We confirmed the overexpression of GLI1 mRNA in RMS by using QRT-PCR. We compared mRNA expression levels of GLI1 in 26 tumor specimens and GEFT in 33 tumor specimens to normal muscle tissues using real time PCR. To accurately quantify mRNA expression of GLI1 and GEFT, ACTB was similarly amplified as an internal control to normalize the results. As a whole, the mean mRNA level of GLI1 in RMS samples was 6.61fold higher compared with those in normal muscle tissues, as shown in Figure 4A  The mean mRNA level of GLI1 in PLA-802 and RD was lower compared with those in normal muscle tissues, being 0.0006-fold and 0.0076-fold, respectively. The mean mRNA level of GEFT in PLA-802 and RD was lower compared with those in normal muscle tissues, being 0.0015-fold and 0.03-fold, respectively.

Functional annotation clustering in RMS
Given that many genes are biologically related, grouping these highly connected genes by network analysis may reveal underlying functional processes in a manner complementary to standard differential expression analyses. We used DAVID functional annotation clustering to allow biological interpretation in a group level and analysis of the internal relationships among the clustered terms. Figure 5 listed the enrichment values associated with certain categories in RMS. It showed that many gene-enriched functional regions were involved in RMS. The representative amplification genes were related to the immunoglobulin domain, Rho-GAP domain, and induction of apoptosis. The representative deletion genes were related to defensin, amylase activity, wound healing. The functions of regions and genes are listed in Table 5.

Functional annotation clustering in ARMS and ERMS
We analyzed amplification and deletion regions genes using DAVID in ARMS and ERMS, respectively. Figures 6 and 7 listed the enrichment values associated with certain categories in ARMS and ERMS, respectively. Numerously enriched functions of genes were found within the amplification regions, but they differed between ARMS and ERMS. In ARMS, for example, enriched functions of genes within the amplification regions included cell cycle process and proto-oncogene. Functional annotation clustering amplification of the cell cycle process included CYP27B1, MDM2, CDK4, and high mobility group AT-hook 2 (HMGA2). Functional annotation clustering amplification of proto-oncogene included GLI1, MDM2, CDK4, HMGA2, MET, and DDIT3 (Table 6). In ERMS, enriched functions of genes within the amplification regions included immunoglobulin-like, IgG binding, and induction of apoptosis. Enriched functions of genes were observed within the deletion regions, include defensin, and wound healing, in the two types of RMS. The correlations of these genes and RMS tumorigenesis were previously unknown, and some could have a function in tumorigenesis processes.

miRNA functional enrichment analysis
Enriched miRNA functions were analyzed for the upregulated and downregulated miRNA in RMS by TAM. The upregulation of onco-miRNA, cell cycle-related miRNA, and muscle development miRNA were associated with RMS, as shown in Figure 8. The regulation of muscle development miRNAs included miR-24, miR-27a, and miR-331. A subset of onco-miRNAs (miR-24, miR-27a, and miR-146b) was associated with RMS (Table 7). No significant results were obtained for the downregulated miRNA in RMS.

Discussion
RMS, the most common soft tissue sarcoma in children, likely results from an imbalance in the proliferation and differentiation of precursor cells during the skeletal myogenesis program. Despite improved understanding of the molecular pathogenesis of RMS in recent decades, patient outcomes remain poor. To increase the accuracy of RMS outcome prediction, efficient molecular markers are needed. An increasing number of evidence shows that gene amplification or deletion is often involved in tumorigenesis and/or tumor progression. Correlations between genomic copy number and gene expression levels have been indicated [18,19,20].
In the current study, high-resolution aCGH was used to provide accurate molecular information on the pathogenesis of RMS. Table 5. Functional annotation cluster of the genes that located in the chromosome regions of amplification and deletion in 20 RMS cases. Amplification  immunoglobulin domain  FGFR2, IGHG1, CD244, IGHG3, IGHG4, VPREB1, LRIG3, LY9, IGHM, PVRL4, UNC5A, IGHA1, FCGR3A,  FCGR3B, IFNGR1, NCR2, F11R, BTNL8, MPZ, BTNL9, FLT4, BTNL3, SLAMF7, PALLD, CD84, FCGR2C, LRRN1,  TREML4, TREML3,  Combined with DAVID, we determined the potential relationships of these massive genes, and improved these genes from biological angles and biological interpretation in a network context. Only a few studies have reported chromosomal changes in frozen RMS or cell lines by aCGH. However, the resolution and number of genes covered by these aCGH chips vary substantially. Using frozen tissues and cell lines as materials, we summarize the results in Table 8. In our study, we used FFPE archival tissues as materials to efficiently detect chromosomal changes in RMS by high-resolution aCGH technique.

Change Annotation Cluster Genes
From Table 8, we found that many probes only included a few genes and regions. The focused regions and genes in previous studies include 12q13.3-q14.1, 8p11.2-11.2, and CDK4, MYCN, GLI, MDM2, FGFR1, and FGFR4, respectively [21,22,23,24,25,26]. Most of them exists frequent gains and amplifications. Using an aCGH platform to examine a specific subset of 26 frozen ERMS samples, Vera et al. have found that these tumors share a common genomic program with a high frequency of gains at 12q13.3 (about 50%) in ERMS [26]. In this study, we have observed high frequencies of gains at 12q13.3-q14.1 in RMS (about 60%), in ERMS (60%), and in ARMS (66.67%), respectively. The above regions include genes such as GLI1, GEFT, OS9, CDK4, PIP5K2C, and CYP27B1. Edoardo et al. indicated that overexpression of the CDK4 and MYCN genes is involved in RMS tumorigenesis, and CENPF, DTL, MYC, EYA2, and FGFR1 are functionally relevant [23]. Daniel et al. showed that the frequency of many specific amplifications and gains (CDK4 and MYCN) significantly varies between fusion genepositive ARMS and fusion gene-negative ARMS and ERMS, and CDK4 exhibits a high frequency of amplifications and gains in fusion gene-positive ARMS [25]. However, we found that the frequencies of CDK4 amplifications in ERMS (3/10, 30%) were similar to those in ARMS (3/9, 33.3%). One potential reason for the difference could be resulted from differences in probe resolution, sample quantity, and ethnicity. GLI1, as well as two other members of the GLI family, is a nuclear mediator of the Hh signaling pathway. Treatment with small-molecule Hh signaling inhibitors inhibits cell proliferation in the ERMS cell lines, which suggested that GLI1 could be an effective therapeutic target for ERMS [27]. Betulinic acid induces apoptosis and inhibits Hh signaling in RMS [28]. GLI1 is expressed at significantly higher levels in ERMS and fusion genenegative ARMS compared with those in fusion gene-positive ARMS. Targeted inhibition of Hh signaling can be an effective strategy for the development of future RMS treatments and prevention [29]. In this study, we found a high frequency of gain and amplification of GLI1, being 60% and 25%, respectively. By QRT-PCR, we found that the mean mRNA level of GLI1 in RMS samples is higher than that in normal muscle tissues. These findings suggest that GLI1 may play an important role in the pathogenesis of RMS, promising GLI1 as a potential candidate biomarker in treatment of RMS.
GEFT was identified as a novel Rho-family-specific guanine nucleotide exchange factors (GEFs), which is highly expressed in brain, heart, and skeletal muscles [30]. Vera et al. demonstrated that GEFT is amplified in ERMS [26], in keeping with our results as shown in RMS. Although the function of GEFT in tumors is unknown, many members of Rho GEFs exhibit increased abundance or activity in human tumors, and potentially affect cancer progression [31,32,33]. In this study, we found that the mean mRNA level of GEFT in RMS samples is higher than that in normal muscle tissues, suggesting that the involvement of GEFT in RMS pathogenesis could be a potential candidate biomarker in the treatment of RMS. To the best of our knowledge, the present investigation is the first study to report on the mRNA expression of GEFT in RMS samples. The specific role of GEFT in RMS needs further research.
MDM2 is an E3 ubiquitin ligase that regulates the protein level of p53 through ubiquitin-dependent degradation [34]. MDM2 is overexpressed in RMS cases [35]. Mitsuru et al. showed that inhibiting the MDM2-p53 pathway with a small-molecule antagonist of MDM2 suppresses tumor growth and induces death of RMS cells, and can be a potential therapy for patients with RMS [36]. Ziad A. et al reported that co-amplification of the CDK4 gene with MDM2 and GLI occurs in human sarcomas, including RMS, Ewing's sarcoma, osteosarcoma, and undifferentiated sarcoma [37]. In the present study, we found that these genes were highly expressed in RMS. Thus, these results are support by the published data.
Furthermore, functional annotation clustering of the expression data readily distinguished genomic alterations in RMS. The high resolution of aCGH combined with human genome database helped in determining the possible target genes present in amplification or deletion regions. By functional annotation clustering, we found many genes involved in tumorigenisis. Many of the abnormalities observed in this study encompass the hallmark chromosomal changes reported in the literature.
By functional annotation clustering, we found that the genes MET, HMGA2, CDK4, MDM2, GLI1, and DDIT3 played a role of proto-oncogene in ARMS. MET oncogene is a unique receptor tyrosine kinase (RTK) located on chromosome 7p, and is activated via its natural ligand hepatocyte growth factor/scatter factor (HGF/SF) ligand. Sandra et al. previously observed MET gain in one ARMS cell line by aCGH [21]. In this study, we have also found MET amplification in RMS. Some results showed that MET possibly has a function not only in ARMS, which carries a dominant genetic lesion in an upstream transcription factor, but   also in ERMS, in which the molecular mechanisms responsible for PAX3/PAX7 upregulation are more elusive [38,39,40]. The data of Riccardo indicated that MET may be necessary for RMS maintenance, and MET-directed therapies may be effective in the treatment of RMS [41]. Francesca et al. showed that MET is widely expressed in ARMS and ERMS at high levels in isolated marrow-infiltrating tumor cells. High levels of expression are associated with unfavorable clinical features, such as tumor marrow involvement [42]. The results of Hou possibly support the function of MET in the development and progression of RMS, and the inhibitor of MET can be an effective targeting therapy reagent for RMS, especially ARMS [43]. So, MET could play an important role in RMS and be expected to become a molecular therapeutic target. HMGA2 is an important regulator of cell growth and differentiation. HMGA2 is expressed during embryogenesis, but is absent or presented at low levels in terminally differentiated tissues. Overexpression of HMGA2 is associated with aggressive tumor growth, early metastasis, and poor prognosis [44,45,46,47]. Yang et al. showed that the expression of HMGA1 and HMGA2 in RMS with relapse or metastasis is higher than that in RMS without relapse or metastasis. Thus, the overexpression of HMGA1 and HMGA2 may be involved in the carcinogenesis and progression of RMS, and these two genes may also be prognostic indicators of the tumor and provide a new basis for targeted therapy [48]. Another report showed that HMGA2 is required for the proliferation and survival of ERMS cells both in vitro and in vivo [49]. In our study, the amplification of HMGA2 could have a co-effect function in RMS carcinogenesis.
In this study, we found that the gene deletion frequency of AMY2A was 20% in RMS, and we also found that AMY2A gene deletion could be involved in ERMS by functional annotation clustering. The AMY2A gene, which codes for salivary and pancreatic amylases, is located at 1p21.1 region. A homozygous deletion of the AMY2A gene is possibly involved in the pathogenesis of gastric carcinoma. AMY2A gene may have potential for tumor suppression in gastric carcinoma [50]. Mohammad showed that gene losses of AMY2A, TGFA, and REG1B in uterine leiomyosarcoma may be responsible for secondary changes that affect the progression and proliferation of the tumor [51]. At present, there are no studies that report AMY2A gene has role in RMS. It remains to be further studied whether AMY2A plays a role of tumor suppressor gene in RMS.
The participation of miRNAs in the pathogenesis of human cancer development has been suggested because dysfunction of specific miRNAs is associated with some cancers [52]. Multiple miRNAs are linked to oncogenes and tumor suppressor genes, including the Ras proto-oncogene, anti-apoptotic gene BCL2, potent p53 tumor suppressor gene, and MET oncogene. miRNAs act as tumor suppressors when they repress oncogenic genes, but act as oncogenes when they downregulate tumor suppressors [53]. In RMS, miR-1/206 suppresses MET expression and functions as a potent tumor suppressor in MET-overexpressing tumors [54]. TAM analysis showed upregulation of onco-miRNA (miR-24, miR-27a, and miR-146b). ALK4, MAPK14, and CDKN2A have been shown to be the target genes of multiple miRNAs, including miRNA-24 [55,56,57]. The CDKN2A gene is inactivated in numerous human tumors. Deletion of CDKN2A can lead to uncontrolled cell proliferation. Petra et al. reported that CDKN2A is an early event in urinary bladder transitional cell carcinoma [58]. These findings suggested the involvement of miRNAs in RMS.
In summary, this study has identified a number of differential changes in RMS-associated genetic alterations using aCGH and revealed several genes that may be candidate molecular targets for RMS. Taken together, the altered pathways may interact with one another in the induction of apoptosis, cell cycle, proto-oncogene, and amylase activity, all of which may ultimately contribute to the development and progression of RMS. Many of these implicated genes may be responsible for changes that affect tumor progression and proliferation. The findings presented here warrant further studies to investigate the pathophysiological functions of these candidate genes in RMS.

Author Contributions
Conceived and designed the experiments: FL. Performed the experiments: CL DL. Analyzed the data: CL. Contributed reagents/materials/analysis tools: JJ JH WZ YC XC YQ HZ WZ. Wrote the paper: CL.