Comparative Analysis of Gene Expression Profiles of Human Dental Fluorosis and Kashin-Beck Disease

To explore the pathologies of Kashin-Beck disease (KBD) and KBD accompanied with dental fluorosis (DF), we conducted a comparative analysis of gene expression profiles. 12 subjects were recruited, including 4 KBD patients, 4 patients with KBD and DF and 4 healthy subjects. Genome-wide expression profiles from their peripheral blood mononuclear cells were evaluated by customized oligonucleotide microarray. R programming software was used for the microarray data analysis followed by functional enrichment analysis through KOBAS. Several potential biomarkers were identified, and quantitative real-time reverse transcription–polymerase chain reaction (qRT-PCR) was used for their validation. In this study, 28 genes and 8 genes were found to be up- and down-regulated respectively in KBD patients compared with health subjects. In patients with KBD and DF, we obtained 10 up-regulated and 3 down-regulated genes compared with health controls. Strikingly, no differential expression gene (DEG) was identified between the two groups of patients. A total of 10 overlaps (DUSP2, KLRF1, SRP19, KLRC3, CD69, SIK1, ITGA4, ID3, HSPA1A, GPR18) were obtained between DEGs of patients with KBD and patients with KBD and DF. They play important roles in metabolism, differentiation, apoptosis and bone-development. The relative abundance of 8 DEGs, i.e. FCRL6, KLRC3, CXCR4, CD93, CLK1, GPR18, SRP19 and KLRF1, were further confirmed by qRT-PCR analysis.

Microarray Hybridization. Total RNA was reverse-transcribed into complementary DNA (cDNA), and then transcribed into cRNA and labled with Cy-Dye using Amino Allyl MessageAmp aRNA Kit (Ambion) following the manufacturer's instructions. Thereafter, 0.5 μg of each labeled cRNA was purified separately and then mixed with hybridization buffer before being applied on the microarray. The hybridization solution was prepared with the In Situ Hybridization Plus kit (Agilent Technologies), and hybridization was performed in the hybridization chamber (Gene-Machines, San Carlos, CA, USA). Conditions of hybridization and washing were in accordance with the manufacturer's recommendations (Agilent Technologies). Quantitative Real-Time Reverse Transcription PCR. Total RNA was prepared for qRT-PCR. These RNA samples were transformed into complementary DNA (cDNA) using Superscript II reverse transcriptase (Invitrogen, Carlsbad, CA) and random primers. qRT-PCR was operated using the ABI 7500 Real-Time PCR system (Applied Biosys-tems, Foster City, CA) according to the manufacturer's specification.

Results
Microarray Data Analysis. The GeneChip ® PrimeView ™ Human Gene Expression Array provides comprehensive coverage of the human genome in a cartridge array. Figure 1A illustrates the overall expression profiles in all samples. Comparable expression levels were obtained after normalization, which should be suitable for the following analysis. With the specified thresholds, 28 genes at higher level and 8 genes at lower level in KBD patients compared with controls were obtained. Tables 1 and 2 shows the up-and down-regulated genes and their enriched functions respectively. For patients with KBD and DF, we identified 11 up-regulated and 7 down-regulated genes compared with healthy controls. Tables 3 and 4 is the up-and down-regulated genes and their enriched functions respectively. Figure 1B,C illustrates the heatmap of DEGs in KBD patients and KBD and DF patients respectively in which green and red represents low and high expression level. A total of 10 overlaps, including 7 up-regulated (DUSP2, KLRF1, SRP19, KLRC3, CD69, SIK1 and ITGA4) and 3 down-regulated (ID3, HSPA1A and GPR18) genes, were identified between the two lists of DEGs. Strikingly, no gene was found to be significantly differential expression between patients with KBD and patients with DF.
qRT-PCR Validation of Microarray Data. 8 genes were further verified by qRT-PCR. The results of qRT-PCR experiment were consistent with microarray analysis. According to qRT-PCR results, the expression levels of KLRC3, KLRF1, SRP19 and CXCR4 were higher in KBD and KBD with DF than controls, while expression levels of CLK1 and GPR18 were lower in both KBD and KBD with DF samples compared with healthy controls (shown in Fig. 2). Besides, FCRL6 was at higher level and CD93 was at lower level only in patients with KBD and DF.

Discussion
KBD and DF are complex diseases that determined by genetic to a large extent. In this study, we conducted a comparative analysis of gene expression profiles for patients with KBD, KBD and DF, and healthy controls using Affymetrix PrimeView ™ Human Gene Expression Array. qRT-PCR was used to validate the oligonucleotide array data. Based on the gene ontology enrichment analysis, we divided these genes into different categories, including metabolism, apoptosis, cytoskeleton, signal transduction and bone development-related genes. Cartilage damage is the main characteristics of pathological changes in KBD, including necrosis in deeper articular cartilage, excessive apoptosis of chondrocytes, extracellular matrix degradation and so on [13][14][15] . Endemic fluorosis is a chronic systemic diseases, characterized with lumbocrural pain, anchylosis, bone-deformity. It is necessary to study bone  Table 3. List of genes differentially expressed in KBD with DF (up-regulated genes).
Scientific RepoRts | (2018) 8:170 | DOI:10.1038/s41598-017-18519-z development-related genes of KBD and DF. Significant differences in gene expression pattern were observed between patients with KBD and healthy controls, as well as between patients with KBD and DF and healthy controls. Ten genes (7 up-regulated and 3 down-regulated) were found to be consistently differentially expressed in KBD and KBD with DF samples, which contains several bone development-related genes, such as DUSP2, ITGA4, ID3, GPR18, and they should provide valuable information for further understanding of KBD and DF. In mammalian cells, the dual-specificity phosphatase (DUSP) family is responsible for the dephosphorylation of threonine and tyrosine residues. Hamamura et al. showed that Dusp2 could suppress inflammation in antibody-induced arthritis in a mouse model through down-regulating inflammatory signs 16 . Besides, DUSP2 is involved in response to oxidative stress and apoptotic signaling, which play important roles in the development of KBD. Yin et al. reported that DUSP2 transcription was induced in response to oxidative stress, causing p53-dependent apoptosis 17 . Moreover, DUSP2 involved in the process of salvianolic acid a (SAA) effects rat cardiomyocytes by mediating regulation of the ERK1/2/JNK pathway 18 .
It was reported that the expression of integrin was associated with the osteoarthritis severity, especially ITGA4 (integrinα4). Becerril M et al. declared that ITGA4 played an important role in the loss of proteoglycans and clusters formation at OA early stages 19 . Proteoglycans, the main component of extracellular matrix of cartilage, were associated with articular cartilage metabolism in patients with KBD 1,20 . Consistent with the previous studies, we also identified ITGA4 as an important biomarker in the pathology of KBD.
Inhibitor of DNA Binding 3 (ID3), a transcription factor involved in the development of T cell and growth inhibition of a B cell progenitors, plays an important role in controlling cell cycle progression [21][22][23] . Thornemo et al. reported that ID3 is important for chondrocyte differentiation and ID proteins are expressed in a lot of cell types and decrease in various cell lines during differentiation 24 . Here, ID3 was also proved to be down-regulated in KBD, as well as KBD and DF samples, which should indicate its roles in the progression of KBD and DF.GPR18, one of the orphan G protein-coupled receptors, has been found to be a receptor for endogenous lipid neurotransmitters. Ramos et al. reported that GPR18 was differentially expressed in osteoarthritis patients 25 . Takenouchi et al. declared that GPR18 involved in the regulation of apoptosis 26 , and apoptosis plays an important role on pathological process of KBD, so GPR18 might contribute the development of KBD.
In summary, we conducted a comparative analysis of gene expression profiles to explore the common pathogenesis and the underlying molecular functions between KBD and DF. Significant differences in gene expression pattern were found between KBD, KBD with DF samples and healthy controls. Our results should provide novel insights for further study of the molecular mechanism of KBD and DF. While, further studies should be conducted to confirm our findings.  Table 4. List of genes differentially expressed in KBD with DF (down-regulated genes).