Integrative Analysis of Genome-Wide Association Studies and DNA Methylation Profile Identified Genetic Control Genes of DNA Methylation for Kashin-Beck Disease

Abstract

Objective

Epigenetic modifications of DNA are regarded as a crucial factor for understanding the molecular basis of complex phenotypes. This study aims to uncover insight into the epigenetic modifications for Kashin-Beck disease (KBD) by integrating genome-wide association studies (GWAS), methylation quantitative trait loci (meQTLs), and DNA methylation profiles data.

Design

The knee articular cartilages of 5 KBD patients and 5 healthy controls were collected for DNA methylation profiling, using Illumina Infinium HumanMethylation450 BeadChip. Mass spectrograph validation of identified differently methylated genes was conducted using independent samples of 4 KBD patients and 3 healthy controls, together with a previous sample of 2743 Han Chinese individuals of GWAS study for KBD and a study of 697 normal subjects for meQTLs annotation datasets. KBD GWAS single nucleotide polymorphisms (SNPs) and normal meQTLs SNPs were integrated with DNA methylation profiles of KBD articular cartilage to identify genetic control (GC) genes of DNA methylation for KBD. Quantitative polymerase chain reaction (qPCR) was performed to validate the mRNA expression of several identified candidate genes.

Results

A total of 162 CpG sites, 253 SNPs, and 123 GC genes for KBD were identified. Enrichment analysis detected 642 marked GO terms and 19 KEGG pathways (P < 0.05). Six potential key GC genes were conducted for qPCR experiment (ERG, MN1, MITF, WISP1, TRIO, and NOSTRIN).

Conclusions

The results suggest that GC genes of DNA methylation may lead to the erosion of cartilage in KBD, which may help us in understanding the epigenetic alteration of KBD.

Keywords

Kashin-Beck disease genetic control genes DNA methylation GWAS meQTLs

Introduction

Kashin-Beck disease (KBD) is an endemic chronic osteochondropathy, which is prevalent at the Eastern Siberia of Russia, North Korea, and the broad diagonal belt from northeast to southwest in China.¹ In 2015, 37.9 million people living at KBD prevalent counties are at the risk of KBD, and 0.57 million patients are suffered KBD in China, including 12,730 children younger than 13 years. KBD is mainly characterized as cartilage degeneration, cartilage matrix degradation, chondrocyte necrosis and apoptosis in growth plate cartilage and articular cartilage. The typical symptoms of KBD include serious joint pain, limited motion, and deformities of joints. According to the diagnosis criteria of KBD in China (WS/T 207-2010), KBD is clinically divided into 3 grades (grades I, II, and III). Recent studies demonstrated the importance of genetic factors in the development of KBD. Nearly half of KBD risk can be explained by genetic components.² Several susceptibility genes have been identified for KBD, such as HLA-DRB1, GPX4, and ITPR2.^3-5 Additionally, genome-wide RNA expression profiling of KBD articular cartilage observed altered RNA expression profiles and reported a group of differently expressed mRNA, microRNA, and lnRNA molecules for KBD.^6-8 However, the molecular mechanism of KBD remains unclear. A large part of grade III KBD patients finally need total knee replacement, resulting in heavy medical and economic burden to the patients.

Epigenetic modifications of DNA are regarded as a crucial factor for understanding the molecular basis of complex phenotypes.⁹ DNA methylation is one of the major epigenetic modifications that plays a pivotal role in gene expression regulation. Recent studies have demonstrated the implication of DNA methylation alteration in the development of arthritis.¹⁰ Jeffries et al.¹¹ performed a genome-wide DNA methylation profiling of cartilage in osteoarthritis (OA) patients and healthy controls, identifying multiple differently methylated genes and correlation between CpG methylation and histologic severity in OA. Another study identified 1,859 differentially methylated loci in the fibroblast-like synoviocytes (FLS) of rheumatoid arthritis (RA), indicating DNA methylation altering the gene expression levels of FLS and leading to the joint damage of RA.¹²

Moreover, recent studies observed that DNA methylation status is under genetic control. The regulatory genetic variants affecting DNA methylation levels are named as methylation quantitative trait loci (meQTLs).¹³ It has been identified that peripheral blood DNA methylation associated with the risk of some malignancies even causing mortality.¹⁴ Accordingly, a large amount of meQTLs have been detected and highlighted for complex diseases. In addition, it has been found that the disease-associated loci identified by GWAS were enriched in regulatory genetic variants.¹⁵ Therefore, integrative analysis of GWAS, meQTLs, and DNA methylation profiles data may uncover novel insight into the epigenetic modifications for complex diseases. In this study, we conducted a trans-omics integrative analysis of KBD GWAS and KBD DNA methylation profiles to detect novel candidate methylation genetic control (GC) genes and their corresponding pathways for KBD.

Methods

Ethics Approval and Consent of Participants

All subjects gave their informed consent for inclusion before they participated in the study. The study was conducted in accordance with the Declaration of Helsinki, and the protocol was approved by the Human Ethics Committee of Xi’an Jiaotong University (project identification code: 2015-237). A written informed consent was obtained from all the patients or relatives of donors.

Study Subjects

A total of 9 grade III KBD patients and 8 healthy control subjects, matched for age and sex were enrolled for genome-wide DNA methylation profiling (5 grade III KBD patients and 5 control subjects) and mass spectrograph validation (4 grade III KBD patients and 3 control subjects). All the KBD patients and healthy control subjects were Chinese Han. The clinical data of study subjects were obtained from nurse-administered questionnaires, including self-reported ethnicities, lifestyle characteristics, health statuses, family, and medical histories. KBD was diagnosed according to the diagnosis criteria of KBD in China (WS/T 207-2010). Healthy control cartilage specimens were collected from the subjects undergoing amputation caused by traffic accidents. The subjects with genetic bone and cartilage diseases, RA and other skeletal diseases, or arthritis were excluded. The total of 17 articular cartilage specimens were harvested from the same anatomic areas of femoral condyles of knee joints.

DNA Methylation Profiling of KBD and Healthy Controls

DNA was extracted from the cartilage specimens using QIAamp DNA Mini Kit (QIAGEN, Germany) following the standard protocol. Agarose gel electrophoresis was employed to evaluate the quality of extracted DNA. Illumina Infinium HumanMethylation450 BeadChip was conducted for genome-wide DNA methylation profiling of 5 KBD and 5 healthy controls. A total of 500 ng DNA was used for bisulfite conversion according to the standard protocol for the EZ DNA Methylation Kit (Zymo Research, USA). Following addition of whole-genome amplification, the samples were incubated overnight at 37°C. After hybridization to HumanMethylation450 array, the array was washed and stained following the instructions of Infinium HD Assay Methylation. The raw image intensities were obtained by IScan SQ scanner (iScan System, Illumina, USA).

Differentially methylated CpG sites were identified by empirical Bayes moderated t test implemented in limma package.¹⁶ Adjusted P value was calculated by Benjamini-Hochberg method. Significant CpG sites were defined as (1) adjusted P value (P_adj) ≤0.05; (2) |β-value difference (Δβ)| ≥ 0.2. For quality control, the CpG sites with missing values, or detection P values >0.05 in more than 90% of cartilage specimens were excluded.

Mass spectrograph validation was conducted using an independent cartilage sample set of 4 KBD patients and 3 healthy control subjects to assure the reliable outcomes of DNA methylation. EpiTect Bisulfate kit (Qiagen, German) was applied to process extracted DNA specimens. Sequenom MassARRAY platform was used for measuring the methylation level of selected CpG sites following standard protocol. Methprimer was used for polymerase chain reaction (PCR) primers design (http://epidesigner.com). The methylation data were generated by Epityper software version 1.0 (Sequenom, CA).

KBD GWAS Data

The GWAS data of KBD were drawn from our published study.⁴ In brief, a GWAS of KBD was conducted. Following the extreme phenotype sampling approach, 90 grade III KBD patients and 1,627 healthy control subjects were recruited for GWAS. Peripheral blood samples (5 mL) were obtained from each participant. E.Z.N.A Blood DNA Midi Kit (Omega Bio-Tek) or a Pure gene DNA isolation kit (Gentra Systems) was used to extract genomic DNA from peripheral blood leukocytes. Affymetrix Genome-Wide Human single nucleotide polymorphisms (SNPs) Array 6.0 conducted for genotyping. Data management and quality control curbed by Affymetrix GeneChip Command Console software. PLINK were used for genome-wide association study. The overage call rate of GWAS samples was 98.92%. In terms of quality control, 106,625 SNPs (call rates of <98.00%), 311,776 SNPs (minor allele frequencies [MAF] < 0.01) and 6,000 SNPs (P < 0.001 from Hardy-Weinberg equilibrium [HWE] testing) were excluded. Significant associations were identified at a traditional genome-wide threshold (P< 5.0 × 10⁻⁸). After quality control filtering, a total of 532,894 SNPs were retained for the GWAS. Detailed procedures could be found in our published study.⁴

meQTLs Annotation Datasets

meQTLs annotation datasets were drawn from a previous study.¹³ Briefly, the peripheral blood specimens of 697 normal subjects were collected. Next-generation sequencing was conducted to assay blood DNA methylation at approximately 4.5 million loci, and methylation measures at each locus were checked with 4.5 million SNPs to exhaustively screen for meQTLs.¹³ The meQTLs were divided into 3 groups, including 542,061 “local” meQTLs (SNP ≤ 1 Mb from methylation site, corresponding to 25,844 nearby genes), 2,761 “distant-same chromosome” meQTLs (SNP > 1 Mb from methylation site, corresponding to 1,327 nearby genes) and 286 “cross chromosome” meQTLs (corresponding to 156 nearby genes).¹³

mRNA Expression in KBD Chondrocytes

Quantitative real-time polymerase chain reaction (qPCR) was used to validate the functional relevance of identified candidate genes in KBD. Human articular chondrocytes were obtained from another independent samples of 3 grade III KBD patients and 3 normal controls with femoral neck fracture (FNF) for qPCR. The collected cartilages tissues were cut into 5 to 10 mm³ slices and incubated with 0.2% trypsin at 37°C for 10 to 20 minutes. After washing with phosphate buffered saline (PBS), the tissue slices were treated for 10 hours with 0.2% type II collagenase at 37°C in constant temperature shaking table. Then, the cells were isolated and cultured at 37°C in 5% CO₂ in Dulbecco’s modified Eagle medium (DMEM)/F12 supplemented with 10% fetal bovine serum (Sijiqing, Zhejiang Tianhang Biotechnology Company, Hangzhou, China) and 1% penicillin/streptomycin. First-passage chondrocytes of all the samples were used in this experiment.

Total RNA from human chondrocytes was isolated according to TRIzol protocol and reverse transcribed to cDNA using One step TaKaRa PrimeScript RT reagent Kit. GAPDH was reference gene. The mRNA expression of genes was detected by qPCR in triplicate analysis, which was performed by specific primers and TaKaRa SYBR Premix Ex Taq II (Tli RNaseH Plus) in the Bio-Rad CFX96 Touch real-time PCR detection system. The reactions were performed in a total volume of 25 µL, containing 2 µL cDNA, 1 µL forward primer (10 µM), 1 µL reverse primer (10 µM), SYBR Premix Ex Taq II (Tli RNaseH Plus) 12.5 µL and sterilized distilled water 8.5 µL. The primer pairs used for qPCR are listed in Table 1 . The cycle of threshold (CT) of each sample was averaged and normalized to that of GAPDH. The results were analyzed by comparative CT equation (2^−∆∆CT).

Table 1.

The Primary Pairs Used for Quantitative Polymerase Chain Reaction.

Gene	Forward	Reverse
ERG	5′-TTGCCAAACCTCTGTGCTGTG-3′	5′-GGATGGTTGGGCTCCGTCTA-3′
MN1	5′-AGAGTTTGCAGATCTGGGCTTTC-3′	5′-ACAACGTGGGTGTTGTATGGTGA-3′
MITF	5′-TCAGATGTGTGAATGGACGGTTTAG-3′	5′-TGATAGCATCACCAATGTTTCCAAG-3′
WISP1	5′-GGGAAGAAGTGTCTGGCTGTG-3′	5′-TGGGTTGATAGGAGCGTGTG-3′
TRIO	5′-TGCAGTCCACCCAGAGCAAC-3′	5′-GCTGGCCAGAATTTGAACGA-3′
NOSTRIN	5′-TCAAGCCAGGCAAGATGATGAG-3′	5′-TGTGTTGCCAGCATTTGAAGGTA-3′
GAPDH	5′-GCACCGTCAAGGCTGAGAAC-3′	5′-TGGTGAAGACGCCAGTGGA-3′

Statistical Analysis

According to SNP ID, the GWAS summary data of KBD was first aligned with the meQTLs annotation datasets to obtain the meQTLs and nearby genes showing association signals (GWAS P < 0.05) for KBD. The identified genes were then compared with the significantly differently methylated genes (detected by DNA methylation profiling) to detect overlapped genes, known as GC genes for KBD. Functional enrichment analysis of the GC genes was performed by ClusterProfiler.¹⁷ T test was performed to compare the mRNA expression levels of 6 candidate genes measured by qPCR.

Results

Genome-wide DNA methylation profiles identified 1,212 different methylated CpG sites and 611 corresponding genes in KBD patients (Padj<0.05). After aligning the GWAS summary data with the meQTLs annotation datasets, we obtained 5,224 “local” meQTLs SNPs (corresponding to 2,059 genes), 12 “distant-same chromosome” meQTL SNPs (corresponding to 4 genes) and 5 “cross chromosome” SNPs (corresponding to 3 genes). Further comparing the differently methylated genes with KBD associated meQTLs SNP lists, we detected 162 CpG sites, 253 SNPs as well as 123 differently methylated GC genes, including 59 hyper-methylated genes and 64 hypomethylated genes in KBD articular cartilage (Supplementary Table S1). All those identified genes are of “local” meQTLs SNPs. Neither “distant-same chromosome” nor “cross chromosome” meQTL SNPs have found corresponding genes related to the DNA methylation of KBD. Six potentially key genes from the identified 123 GC genes were detected by qPCR (ERG, MN1, MITF, WISP1, TRIO, and NOSTRIN). Statistically analysis observed that only the mRNA expressions of MITF (P = 0.0315, ratio = 0.659) and WISP1 (P = 0.0027, ratio = 2.465) are significant in KBD chondrocytes ( Fig. 1 ).

Figure 1.

The mRNA expression of ERG, MN1, MITF, NOSTRIN, TRIO, and WISP1, in patients with Kashin-Beck disease (KBD) compared with normal samples.

Gene ontology (GO) enrichment analysis detected 624 GO terms (P < 0.05), including 476 GO terms for Biological Process (BP), 37 for Cellular Component (CC), and 111 for Molecular Function (MF) (Supplementary Table S2). The top 10 GO terms are shown in Figure 2 . We observed that neuron migration (GO:00001764, P = 2.20 × 10⁻⁵) in BP, cell leading edge (GO:0031252, P = 0.00144) in CC, phosphoric ester hydrolase activity (GO:0042578, P = 0.000361) in MF were the most significant GO terms. Pathway enrichment analysis identified 19 KEGG pathways (P < 0.05), such as AMP-activated protein kinase (AMPK) signaling pathway (hsa04152, P = 0.006791), Adipocytokine signaling pathway (hsa04920, P = 0.0091) and the forkhead box O (FOXO) signaling pathway (hsa04068, P = 0.0496) ( Table 2 ).

Figure 2.

The top 10 gene ontology (GO) terms in Biological Process (BP), Cellular Component (CC), and Molecular Function (MF)

Table 2.

KEGG Enrichment Analysis Results.

ID	Description	P
hsa05202	Transcriptional misregulation in cancer	0.0010
hsa04930	Type II diabetes mellitus	0.0029
hsa04931	Insulin resistance	0.0044
hsa04152	AMPK signaling pathway	0.0068
hsa04920	Adipocytokine signaling pathway	0.0091
hsa00562	Inositol phosphate metabolism	0.0106
hsa00030	Pentose phosphate pathway	0.0150
hsa05410	Hypertrophic cardiomyopathy (HCM)	0.0150
hsa00051	Fructose and mannose metabolism	0.0180
hsa04211	Longevity regulating pathway	0.0181
hsa04022	cGMP-PKG signaling pathway	0.0187
hsa04925	Aldosterone synthesis and secretion	0.0221
hsa04960	Aldosterone-regulated sodium reabsorption	0.0224
hsa04070	Phosphatidylinositol signaling system	0.0239
hsa00290	Valine, leucine and isoleucine biosynthesis	0.0249
hsa04510	Focal adhesion	0.0356
hsa00280	Valine, leucine and isoleucine degradation	0.0363
hsa04270	Vascular smooth muscle contraction	0.0400
hsa04068	FOXO signaling pathway	0.0496

Discussion

In this study, a total of 123 GC genes were identified significantly for KBD, and they were only associated with “local” meQTLs SNPs. Showing that GC genes of methylation occur where SNP and CpG site within 1 Mb of each other in the same chromosome in KBD. A few previous studies have identified the correlation between meQTLs and expression QTLs (eQTLs) and found that only 4.8%¹⁸ and 11.8%¹⁹ were significant as both meQTL and eQTL. This could explain what we have found in the mRNA expression of KBD chondrocytes that only 2 out of 6 GC genes significantly expressed.

Among all the identified 123 genes, 48.0% were hypermethylated genes in KBD chondrocytes, in which, MITF was found as the top significant hypermethylated gene, with 3 methylated CpG sites, including cg18503031, cg13523819, and cg13151171. The significant site is at the locus on chromosome 3 and the associated SNP is rs6782870. It also showed significant lower mRNA expression in KBD chondrocytes than controls. MITF encodes melanogenesis associated transcription factor. Numerous studies have identified roles for MITF in oxidative stress response, cell growth, differentiation, miRNA regulation, and DNA methylation.^20-22 A model has illustrated that lower MITF levels favor cell invasion and suppress proliferation.²⁰ Strub et al.²³ recognized that knocking down MITF led to DNA damage, finally ending in cellular senescence. Additonally, Wang et al.²⁴ have identified oxidative stress mediates the initiation and progression of KBD and leads to necrosis of chondrocytes. Taken together, MITF may play a vital role in the progression of KBD, probably involving in the oxidative stress response, cell growth, and DNA damage in chondrocytes of KBD.

Meningioma 1 (MN1), known as transcriptional regulator, is another significant hypermethylated GC gene for KBD. The significant site is at the locus on chromosome 22 and the associated SNP is rs134133. Its methylated CpG site is cg01214458 in KBD chondrocytes. Prior study has shown that MN1 is a novel target in osteoblastic cells and overexpressed dramatically during differentiation of primary osteoblastic cells,²⁵ suggesting that MN1 may have pivotal roles in modulating osteoblast proliferation and inhibiting cell growth. KBD is a complicate disorder with skeletal growth and development failure, including osteopenia and osteoporosis.²⁶ Therefore, more biological research on MN1 might identify the mechanism of osteopenia and osteoporosis in KBD.

WISP1 is a GC gene with hypomethylation. Its methylated CpG site is cg10191240 in KBD chondrocytes. The significant site is at the locus on chromosome 8 with the associated SNP rs4736441. It showed pronounced overexpression in KBD chondrocytes. WISP1 is known as WNT1 inducible signaling pathway protein 1 and the fourth member of the CCN family (cysteine-rich angiogenic protein 61, connective tissue growth factor, nephroblastoma overexpressed). It is considered to have the functions of cell differentiation, cell proliferation, tissue regeneration, and the regulation of transduction.^27,28 Wnt signaling pathway plays an important role in embryonic development of cartilage and bone, and the expression of various Wnt ligands is strongly increased in the synovium in experimental OA.²⁹ It has also found that overexpression of WISP1 induces matrix metalloproteinase and aggrecanase expression and contributes to cartilage damage.²⁹ In addition, several studies have been identified the abnormal activation of Wnt signaling pathway in KBD.³⁰ In this study, DNA methylation profiling, GWAS, meQTLs, and mRNA expression have confirmed the hypomethylated and overexpression of WISP1 in KBD patients.

NOSTRIN, a GC gene with hypomethylation, named as nitric oxide synthase (NOS) trafficking, is well known for its relation to nitric oxide (NO). The significant site is at the locus on chromosome 2, with the associated SNP rs17808106, and the methylated CpG sites cg00174992 in KBD chondrocytes. NO regulates a variety of biologic processed such as neurotransmission, inflammatory response, and vascular homeostasis. It stimulates metalloproteinase production and activity, increases susceptibility to cell damage by reactive oxygen species (ROS), inhibits collagen and proteoglycan synthesis, and induces the apoptotic cascade.³¹ The enzyme system responsible for producing NO is NOS, which exists in three distinct, such as constitutive neuronal NOS (nNOS), inducible NOS (iNOS) and constitutive endothelial NOS (eNOS).³² NOSTRIN binds the enzyme responsible for NO production and inhibits production of NO by suppressing the activation of eNOS.³³ eNOS has been identified as acting a mechanistic role in the function of osteocytes, osteoblasts, and osteoclasts as well as chondrocytes.³¹ Li et al.³⁴ found that NO levels significantly increased in adult KBD patients. Together, we speculated NOSTRIN-eNOS-NO may participate in the cartilage damage of KBD.

Pathway enrichment analysis detected AMPK signaling pathway significant in KBD. AMPK is a serine threonine kinase that is highly conserved through evolution and acts as an intracellular sensor modulating the energy balance within the cell.³⁵ The signaling pathways accelerating AMPK activation have been link to several catabolic responses in diseases.³⁶ Zhou et al.³⁷ found AMPK deficiency in chondrocytes aggravated the instability-induced and aging-associated OA in adult mice. Additionally, the AMPK-mTOR signaling pathway is vital in autophagy, which contributed to the imbalance of articular cartilage homeostasis and chondroprotection, leading to cartilage erosion during OA.³⁸ Therefore, more studies are required to investigate the potential roles and mechanism of AMPK signaling pathway in the cartilage damage of KBD.

It is notable that the FOXO signaling pathway was enriched in KBD patients. The transcription factors of FOXO family regulate the expression of genes in cellular physiological events including apoptosis, cell-cycle control, glucose metabolism, oxidative stress resistance, and longevity.³⁹ FOXO family contains 4 members in mammals, FOXO1, FOXO3, FOXO4, and FOXO6. Akasaki et al.⁴⁰ identified FOXO factors is downregulated in human OA cartilage, contributing to increased apoptotic cell death, intracellular ROS generation as well as decreased antioxidant proteins and autophagy-related proteins in human chondrocytes. Increasing apoptotic chondrocytes, defective autophagy, and ROS have been demonstrated playing important roles in cartilage damage of KBD, which meet the functions of FOXO factors to some extent.

New candidate genes and pathways may provide new clues for understanding the mechanisms of KBD. However, this study only tested 6 of the identified GC genes in KBD chondrocytes. Therefore, it is worthy of more studies to illustrate all those GC genes and associated pathways in the evolution of KBD in future work.

Conclusions

Integrating GWAS data of KBD, meQTL of normal subjects and DNA methylation profiles of KBD, we have identified 123 GC genes significantly for KBD. Our results suggest that GC genes of DNA methylation may lead to the cartilage erosion in KBD. In order to interpret associations with genetic variants or disease phenotype, 123 GC genes may suggest mechanisms in the alteration of DNA methylation in KBD and may help us in understanding the epigenetic alteration of KBD.

Supplemental Material

Additional_Table_S1_5 – Supplemental material for Integrative Analysis of Genome-Wide Association Studies and DNA Methylation Profile Identified Genetic Control Genes of DNA Methylation for Kashin-Beck Disease

Supplemental material, Additional_Table_S1_5 for Integrative Analysis of Genome-Wide Association Studies and DNA Methylation Profile Identified Genetic Control Genes of DNA Methylation for Kashin-Beck Disease by Ping Li, Cuiyan Wu, Xiong Guo, Yan Wen, Li Liu, Xiao Liang, Yanan Du, Lu Zhang, Mei Ma, Shiqiang Cheng, Bolun Cheng, Sen Wang and Feng Zhang in CARTILAGE

Supplemental Material

Additional_Table_S2_GO_Enrichment_P_value_1 – Supplemental material for Integrative Analysis of Genome-Wide Association Studies and DNA Methylation Profile Identified Genetic Control Genes of DNA Methylation for Kashin-Beck Disease

Supplemental material, Additional_Table_S2_GO_Enrichment_P_value_1 for Integrative Analysis of Genome-Wide Association Studies and DNA Methylation Profile Identified Genetic Control Genes of DNA Methylation for Kashin-Beck Disease by Ping Li, Cuiyan Wu, Xiong Guo, Yan Wen, Li Liu, Xiao Liang, Yanan Du, Lu Zhang, Mei Ma, Shiqiang Cheng, Bolun Cheng, Sen Wang and Feng Zhang in CARTILAGE

Footnotes

Author Contributions

FZ takes responsibility for the integrity of the work as a whole, from inception to the finished manuscript. GX gave the conception, design and some critical suggestion to this study. PL, CYW and YW performed all the experiments and analysis. LL, XL, YND, LZ, MM, SQC, BLC, and SW partially performed some experiments. PL and CYW convinced the project and wrote the manuscript.

Acknowledgments and Funding

The work was supported by National Natural Science Foundation of China Grant Numbers 81472925, 81673112, and 81502766; Key Projects of International Cooperation among Governments in Scientific and Technological Innovation Grant Number 2016YFE0119100; and the Fundamental Research Funds for the Central Universities. We would like to thank the authors who generously shared their data as well as all study participants.

Declaration of Conflicting Interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Ethical Approval

The study was conducted in accordance with the Declaration of Helsinki, and the protocol was approved by the Human Ethics Committee of Xi’an Jiaotong University (project identification code: 2015-237).

Informed Consent

A written informed consent was obtained from all the patients or relatives of donors.

Trial Registration

Not applicable.

ORCID iD

Ping Li

Supplemental Material

Supplementary material for this article are available online.

References

Guo

Zhang

Ren

Lammi

. Recent advances in the research of an endemic osteochondropathy in China: Kashin-Beck disease. Osteoarthritis Cartilage. 2014;22(11):1774-83. doi:10.1016/j.joca.2014.07.023

Wen

Guo

Yang

Shen

Chen

, et al. Genetic association, mRNA and protein expression analysis identify ATG4C as a susceptibility gene for Kashin-Beck disease. Osteoarthritis Cartilage. 2017;25(2):281-6. doi:10.1016/j.joca.2016.09.019

Shi

Liu

Wang

Huang

Liu

, et al. Genetic variants in the HLA-DRB1 gene are associated with Kashin-Beck disease in the Tibetan population. Arthritis Rheum. 2011;63(11):3408-16. doi:10.1002/art.30526

Zhang

Wen

Guo

Zhang

Wang

Yang

, et al. Genome-wide association study identifies ITPR2 as a susceptibility gene for Kashin-Beck disease in Han Chinese. Arthritis Rheumatol. 2015;67(1_suppl):176-81. doi:10.1002/art.38898

Xiong

Zou

Song

Sun

, et al. Association study between polymorphisms in selenoprotein genes and susceptibility to Kashin-Beck disease. Osteoarthritis Cartilage. 2010;18(6):817-24. doi:10.1016/j.joca.2010.02.004

Wang

Guo

Duan

Zhang

, et al. Comparative analysis of gene expression profiles between the normal human cartilage and the one with endemic osteoarthritis. Osteoarthritis Cartilage. 2009;17(1_suppl):83-90. doi:10.1016/j.joca.2008.05.008

Wen

Xiao

Hao

Zhang

, et al. Comparison of microRNA expression profiles of Kashin-Beck disease, osteoarthritis and rheumatoid arthritis. Sci Rep. 2017;7(1_suppl):540. doi:10.1038/s41598-017-00522-z

Liu

Zhang

Shao

Yang

Ning

, et al. Long noncoding RNA expression profile reveals lncRNAs signature associated with extracellular matrix degradation in Kashin-Beck disease. Sci Rep. 2017;7(1_suppl):17553. doi:10.1038/s41598-017-17875-0

Petronis

. Epigenetics as a unifying principle in the aetiology of complex traits and diseases. Nature. 2010;465(7299):721-7. doi:10.1038/nature09230

10.

Fernández-Tajes

Soto-Hermida

Vázquez-Mosquera

Cortés-Pereira

Mosquera

Fernández-Moreno

, et al. Genome-wide DNA methylation analysis of articular chondrocytes reveals a cluster of osteoarthritic patients. Ann Rheum Dis. 2014;73(4):668-77. doi:10.1136/annrheumdis-2012-202783

11.

Jeffries

Donica

Baker

Stevenson

Annan

Humphrey

, et al. Genome-wide DNA methylation study identifies significant epigenomic changes in osteoarthritic cartilage. Arthritis Rheumatol. 2014;66(10):2804-15. doi:10.1002/art.38762

12.

Nakano

Whitaker

Boyle

Wang

Firestein

. DNA methylome signature in rheumatoid arthritis. Ann Rheum Dis. 2013;72(1_suppl):110-7. doi:10.1136/annrheumdis-2012-201526

13.

McClay

Shabalin

Dozmorov

Adkins

Kumar

Nerella

, et al. High density methylation QTL analysis in human blood via next-generation sequencing of the methylated genomic DNA fraction. Genome Biol. 2015;16:291. doi:10.1186/s13059-015-0842-7

14.

Zhang

Wilson

Heiss

Breitling

Saum

Schöttker

, et al. DNA methylation signatures in peripheral blood strongly predict all-cause mortality. Nat Commun. 2017;8:14617. doi:10.1038/ncomms14617

15.

Nicolae

Gamazon

Zhang

Duan

Dolan

Cox

. Trait-associated SNPs are more likely to be eQTLs: annotation to enhance discovery from GWAS. PLoS Genet. 2010;6(4):e1000888. doi:10.1371/journal.pgen.1000888

16.

Ritchie

Phipson

Law

Shi

, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47. doi:10.1093/nar/gkv007

17.

Wang

Han

. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012;16(5):284-7. doi:10.1089/omi.2011.0118

18.

Gibbs

van der Brug

Hernandez

Traynor

Nalls

Lai

, et al. Abundant quantitative trait loci exist for DNA methylation and gene expression in human brain. PLoS Genet. 2010;6(5):e1000952. doi:10.1371/journal.pgen.1000952

19.

Schulz

Ruppert

Herms

Wolf

Mirza-Schreiber

Stegle

, et al. Genome-wide mapping of genetic determinants influencing DNA methylation and gene expression in human hippocampus. Nat Commun. 2017;8(1_suppl):1511. doi:10.1038/s41467-017-01818-4

20.

Vachtenheim

Ondrušová

. Microphthalmia-associated transcription factor expression levels in melanoma cells contribute to cell invasion and proliferation. Exp Dermatol. 2015;24(7):481-4. doi:10.1111/exd.12724

21.

Lauss

Haq

Cirenajwis

Phung

Harbst

Staaf

, et al. Genome-wide DNA methylation analysis in melanoma reveals the importance of CpG methylation in MITF regulation. J Invest Dermatol. 2015;135(7):1820-8. doi:10.1038/jid.2015.61

22.

Haq

Shoag

Andreu-Perez

Yokoyama

Edelman

Rowe

, et al. Oncogenic BRAF regulates oxidative metabolism via PGC1α and MITF. Cancer Cell. 2013;23(3):302-15. doi:10.1016/j.ccr.2013.02.003

23.

Strub

Giuliano

Bonet

Keime

Kobi

, et al. Essential role of microphthalmia transcription factor for DNA replication, mitosis and genomic stability in melanoma. Oncogene. 2011;30(20):2319-32. doi:10.1038/onc.2010.612

24.

Wang

Wei

Luo

Cao

Yang

, et al. Oxidative stress and status of antioxidant enzymes in children with Kashin-Beck disease. Osteoarthritis Cartilage. 2013;21(11):1781-9. doi:10.1016/j.joca.2013.08.002

25.

Zhang

Dowd

Moore

Kranenburg

Meester-Smoor

Zwarthoff

, et al. Meningioma 1 is required for appropriate osteoblast proliferation, motility, differentiation, and function. J Biol Chem. 2009;284(27):18174-83. doi:10.1074/jbc.M109.001354

26.

Wen

Guo

Hao

Xiao

Wang

, et al. Integrative analysis of genome-wide association studies and gene expression profiles identified candidate genes for osteoporosis in Kashin-Beck disease patients. Osteoporos Int. 2016;27(3):1041-6. doi:10.1007/s00198-015-3364-y

27.

van den Bosch

Blom

Kram

Maeda

Sikka

Gabet

, et al. WISP1/CCN4 aggravates cartilage degeneration in experimental osteoarthritis. Osteoarthritis Cartilage. 2017;25(11):1900-11. doi:10.1016/j.joca.2017.07.012

28.

Maiese

. Picking a bone with WISP1 (CCN4): new strategies against degenerative joint disease. J Transl Science. 2016;1(3):83-5.

29.

Blom

Brockbank

van Lent

Van-Beuningen

Geurts

Takahashi

, et al. Involvement of the Wnt signaling pathway in experimental and human osteoarthritis: prominent role of Wnt-induced signaling protein 1. Arthritis Rheum. 2009;60(2):501-12.

30.

Zheng

Zhang

Hou

, et al. Abnormal expression of chondroitin sulphate N-acetylgalactosaminyltransferase 1 and Hapln-1 in cartilage with Kashin-Beck disease and primary osteoarthritis. Int Orthop. 2013;37(10):2051-9. doi:10.1007/s00264-013-1937-y

31.

Teixeira

Ischiropoulos

Leboy

Adams

Shapiro

. Nitric oxide-nitric oxide synthase regulates key maturational events during chondrocyte terminal differentiation. Bone. 2005;37(1_suppl):37-45. doi:10.1016/j.bone.2005.03.010

32.

Ricciardolo

Nijkamp

Folkerts

. Nitric oxide synthase (NOS) as therapeutic target for asthma and chronic obstructive pulmonary disease. Curr Drug Targets. 2006;7(6):721-35.

33.

Zimmermann

Opitz

Dedio

Renne

Muller-Esterl

Oess

. NOSTRIN: a protein modulating nitric oxide release and subcellular distribution of endothelial nitric oxide synthase. Proc Natl Acad Sci U S A. 2002;99(26):17167-72. doi:10.1073/pnas.252345399

34.

Guo

Wang

Geng

Kang

Wang

, et al. Serum hyaluronic acid, tumor necrosis factor -alpha, vascular endothelial growth factor, NO, and Se levels in adult patients with Kashin-Beck disease [in Chinese]. Nan Fang Yi Ke Da Xue Xue Bao. 2007;27(7):941-4.

35.

Yang

Lai

Bai

Jin

. Chondrocyte-specific ablation of AMPKα1 does not affect bone development or pathogenesis of osteoarthritis in mice. DNA Cell Biol. 2016;35(3):156-62.

36.

Sondergaard

Schultz

Madsen

Bay-Jensen

Kassem

Karsdal

. MAPKs are essential upstream signaling pathways in proteolytic cartilage degradation—divergence in pathways leading to aggrecanase and MMP-mediated articular cartilage degradation. Osteoarthritis Cartilage. 2010;18(3):279-88. doi:10.1016/j.joca.2009.11.005

37.

Zhou

Chen

, et al. AMPK deficiency in chondrocytes accelerated the progression of instability-induced and ageing-associated osteoarthritis in adult mice. Sci Rep. 2017;7:43245.

38.

Liu-Bryan

Terkeltaub

. Emerging regulators of the inflammatory process in osteoarthritis. Nat Rev Rheumatol. 2015;11(1_suppl):35-44. doi:10.1038/nrrheum.2014.162

39.

Daitoku

Fukamizu

. FOXO transcription factors in the regulatory networks of longevity. J Biochem. 2007;141(6):769-74. doi:10.1093/jb/mvm104

40.

Akasaki

Alvarez-Garcia

Saito

Caramés

Iwamoto

Lotz

. FoxO transcription factors support oxidative stress resistance in human chondrocytes. Arthritis Rheumatol. 2014;66(12):3349-58. doi:10.1002/art.38868

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.00 MB

0.21 MB

0.25 MB