Abstract
Clinical variables contribute to the severity of Kashin-Beck disease (KBD). However, it is unclear if there is a correlation between gene expression and clinical variables. Peripheral blood samples were collected from 100 patients with KBD and 100 healthy controls from KBD-endemic areas to identify differentially expressed genes in KBD. Correlation analysis and multiple logistic regression analysis were performed using gene expression and clinical parameters. Immunohistochemistry (IHC) was used to detect the expression of related proteins in articular cartilage tissues. Thirty-nine differentially expressed genes were identified in patients with KBD. Nine differentially expressed genes were correlated with the metacarpal length/metacarpal breadth index. FZD1 was identified as having statistical significance in establishing the regression model of clinical parameters and gene expression. FZD1 expression levels were remarkably reduced in patients with KBD. Our results indicate that FZD1 could be involved in the pathological process of phalanges tuberositas and brachydactylia and may provide new insight into the pathogenesis of articular cartilage destruction observed in patients with KBD.
1. Introduction
Kashin-Beck disease (KBD) is a result of lesions in articular cartilage and epiphyseal plate cartilage from necrosis in deep zone chondrocytes, as well as excessive apoptosis and dedifferentiation in chondrocytes, and an important cause of disability in northwestern China [1–3]. According to the clinical features of the disease, including tuberosity of the fingers and brachydactylia, deformed and enlarged joints, and limited motion of the joints, patients with KBD are classified as different grades (I, II, and III) [4, 5]. Although examinations of lesions in articular cartilage and epiphyseal plate cartilage have provided insights into the cellular and molecular mechanisms involved in the necrosis, excessive apoptosis, and dedifferentiation, further exploration is needed to improve the relationship between clinical variables and gene expression in the pathogenesis of KBD.
Some clinical parameters are connected with disease severity in patients with KBD. Patients with limitation of motion, arthralgia, and deformity of the limbs are considered to have advanced KBD [6]. It is possible that the clinical features that represent the KBD grade follow gene expression patterns. Modern molecular medicine has been performed in an exploration to identify genes expressed in the pathologic tissues of KBD. Although research results on pathogenesis in KBD are limited, there are plenty of experimental data in the studies on mechanisms involved in gene expression [7, 8].
However, there are few research studies focused on the relationship between gene expression pattern and clinical parameters in KBD. In this study, we have studied the association between clinical parameters and the gene expression signatures of patients with KBD to test the hypothesis that gene expression is related to patients’ clinical parameters. If clinical parameters are significant determinants of gene expression alterations in KBD, these results could provide a molecular insight into understanding the pathogenesis of cartilage injuries and the subsequent development of KBD.
2. Materials and Methods
The study was approved by the Ethical Committee of Xi’an Jiaotong University. All participants gave their informed consent by signing a document that had been carefully reviewed by the Ethical Committee of Xi’an Jiaotong University.
2.1. Patients and Tissue Samples
Random selection was used for enrolling subjects from Yong-shou and Lin-you counties in Shaanxi Province, two endemic areas of KBD with prevalence rates of 20.4% and 10.5%, respectively. The national diagnostic criteria for KBD in China (WS/T207-2010) was used to diagnose the patients. Healthy controls without clinical symptoms of KBD were collected and matched with KBD patients by gender and age. Subjects with other types of osteoarthropathy and other chronic diseases, such as cardiovascular disease, diabetes, and hypertension, were excluded. As mentioned above, 100 KBD patients (50 degree I and 50 degree II) and 100 healthy controls were selected for this study. The average age of patients was 57.9 years, ranging from 43 to 79 years, and the average age of the controls was 53.4 years, ranging from 40 to 77 years. Both patients and controls consisted of 44 males and 56 females. Peripheral blood was collected from each subject included in this study using heparinized Vacutainer® tubes (Becton Dickenson, San Jose, CA, USA). The Hemavet 950 (Drew Scientific, Oxford, CT, USA) was used for determining leukocyte cell numbers. The peripheral blood was centrifuged at 1500 × g for 20 min to separate peripheral blood mononuclear cells (PBMCs). The cell pellets were resuspended in Hanks’ balanced salt solution (Gibco BRL/Invitrogen, Carlsbad, CA, USA). The cell suspensions were layered over 5 mL of lympholyte-H (Cedarlane Labs, Hornby, BC, Canada) in a 15 mL Falcon tube and were centrifuged at 1500 × g for 40 min. The cells were stored in RNAlater® (Ambion Inc., Austin, TX, USA) until RNA isolation after rinsing twice with cold Hanks’ balanced salt solution.
According to the inclusion and exclusion criteria of the sample selection described above, the articular cartilage specimens were collected from five patients with KBD and five healthy controls (Table 1). The donors of the adult KBD samples were from the KBD-endemic area of Yong-shou County. Control sample donors were also from the KBD-endemic areas. The articular cartilage samples of the KBD patients and the healthy controls were collected from individuals who had undergone an arthroplasty in the knee or who had suffered an amputation due to an accident.
All of the articular cartilage samples were obtained within 6 h after surgery or patient death. After removing the skin and muscles from the fingers, we opened the joint capsules carefully and harvested articular cartilage. Chondrocytes were isolated as described previously [9, 10].
2.2. Clinical Parameter Questionnaire
A KBD questionnaire with a series of clinical variables and general information was completed by all KBD subjects in this study. The clinical variables used for this study included phalanges tuberositas, brachydactylia, the deformity of the elbow and knee, and dyskinesia in joints of the wrist, elbow, shoulder, knee, and ankle. General information included age, gender, BMI, metacarpal length, and metacarpal breadth. The Manouvrier’s Skelic Index, metacarpal length/metacarpal breadth index, and metacarpal length/height index were calculated based on the following equations:
2.3. Selection of the Candidate Genes
In this study, we reviewed two databases of previously published microarray analyses [7, 11] and selected 169 genes including 11 housekeeping genes as target genes to be measured in a custom-made microarray. From the previously published studies, genes with more than 2- or less than 0.5-fold change in gene expression level between KBD and healthy persons, with a value less than 0.05 and biologic function related to cartilage, were selected. Ultimately, the 158 genes identified by application of the selection criteria above were selected as target genes.
2.4. Microarray Hybridization
Total RNA was isolated from the PBMCs using TRIzol® reagent (Life Technologies Inc., Carlsbad, CA, USA) following the manufacturer’s recommended protocol. A high-resolution electrophoresis system (Agilent 2100 Bioanalyzer, Agilent Technologies, Palo Alto, CA, USA) was used for determining the quality and integrity of the extracted total RNA. The total RNA isolated from all patients with KBD and healthy controls were transcribed into complementary DNA (cDNA) and then were reverse-transcribed into cRNA and labelled with CyDye using the Amino Allyl Message Amp™ and RNA Amplification Kit (Ambion) according to the manufacturer’s instructions. The custom-made primer array contained 169 oligonucleotide probes representing the selected 169 human genes (manufactured by the National Engineering Research Center for Miniaturized Detection System in Xi’an, China). Microarray hybridization was performed following the recommended protocol [12].
2.5. Gene Expression Analysis
Expression ratios were calculated for each gene with more than 2-fold change or less than 0.5-fold change that was regarded as statistically significant and considered a differentially expressed gene. values were calculated using the value log ratio algorithm: in which is calculated as
represents twice the integral of the Gaussian distribution, with a mean value of 0 and variance of 0.5, and is the deviation of the log ratio from 0. This calculation gives the statistical significance of the log ratio for each feature (i.e., transcript levels) between the red and green channels. Only values less than 0.05 were regarded as significant.
2.6. Data Mining and Statistical Analysis
Statistics for comparison and correlation were analysed by SPSS18.0. Differences in categorical variables were evaluated with by Pearson’s Chi-square test. Differences of means were determined by one-way analysis of variance (ANOVA) for multiple comparisons. The t-test was applied to determine the difference between two groups. A normality test was first applied in continuous variables before any deeper comparison analyses. Nonparametric methods (e.g., Mann-Whitney U and Kruskal Wallis tests for pairwise comparisons) were applied when the data were not normally distributed. Correlation analysis was used between gene expression and indexes of clinical parameters, such as the Manouvrier’s Skelic Index, the metacarpal length/metacarpal breadth index, and the metacarpal length/height index. Univariate and multiple logistic regression analyses were used to identify the clinical sign-specific genes.
2.7. Immunohistochemical Localization Verification in Cartilage Tissues from Patients with KBD
Cartilage tissues were fixed with 4% paraformaldehyde for 24 hours after removal of the tissue and decalcified in 10% ethylenediaminetetraacetic acid (EDTA). The samples were dehydrated in a series of alcohol, cleared in xylene, and embedded in paraffin wax. Paraffin sections were cut into 5 μm sections, mounted on slides, and stored at room temperature until used for staining. The paraffin-embedded sections were deparaffinized with xylene and then rehydrated in decreasing graded ethanol. Endogenous peroxidase activity was blocked by 3% hydrogen peroxide for 10 min at room temperature, then the sections were washed with PBS and incubated in 10 M urea solution and trypsin at 37°C for 20 min to unmask the antigen. After blocking with 5% goat serum for 20 min at room temperature, anti-FZD1 (1 : 50 dilution), as well as negative control IgG, were applied on the sections, and the samples were incubated overnight at 4°C. After washing with PBS, sections were incubated using the SAP kit (Zhongshan Jinqiao, Guangzhou, China). The substrate 3, 3-diaminobenzidine (DAB) was added to stain the sections with haematoxylin counterstaining. Finally, the sections were dehydrated and mounted under alcohol-cleaned cover slips. The percentage of positive cells in different zones of cartilage was calculated using the Kuettner standard [13] for defining the zones. The localizations of FZD1 in each cartilage zone were assessed systematically by the rate of positively stained cells calculated by the equation below:
3. Results and Discussion
After checking whether the questionnaires were consistent with the inclusion and exclusion criterion by two well-trained research assistants, eighty questionnaires of patients with KBD were ultimately included in this analysis.
3.1. Univariate Analysis in Clinical Parameters and General Information
Univariate analysis was applied in clinical parameters and general information. There were significant differences between KBD patients with grade I and II disease in four items among nine clinical parameters (brachydactylia, ; deformity of the knee, ; and dyskinesia in joints of the elbow and knee, ; ) (Table 2), and only four items among fourteen general information items were significantly different between KBD patients with grade I and II disease (age, ; metacarpal length, ; blurred vision, ; and number of family members with KBD, ) (Tables 2 and 3).
3.2. Differentially Expressed Genes
39 differentially expressed genes, including 21 downregulated and 18 upregulated genes, were identified in the PBMCs of patients with KBD (Table 4). There were 12 differentially expressed genes in KBD grade I patients, and 92 differentially expressed genes were identified in KBD grade II patients. These genes belong to different functional pathways, including metabolism, transcription factors, cytokine factors, and signal transduction.
3.3. Correlation Analysis and Univariate Analysis in Gene Expression and Clinical Parameter Indexes
The correlation analysis was performed to identify the relationship between differentially expressed genes and clinical parameter indexes. We found that only nine genes among the 39 differentially expressed genes were correlated with the metacarpal length/metacarpal breadth index (Table 5). By using Pearson’s chi-squared test, the differential expression rate of all nine genes was identified as statistically significant between KBD grades I and II (Table 6), six genes were identified as statistically significant in phalanges tuberositas (Table 7), and eight genes were identified as statistically significant in brachydactylia (Table 8).
3.4. Multiple Logistic Regression Analysis of Clinical Symptoms (Phalanges Tuberositas and Brachydactylia) and Differentially Expressed Genes
We used the six and eight genes that were identified as having a statistically significant correlation with the clinical symptom by univariate regression analysis and multiple logistic regression analysis. Clinical signs (phalanges tuberositas and brachydactylia) were considered as dependent variables, and differentially expressed genes were considered as independent variables in multiple logistic regression analysis. Finally, we found that FZD1 and TNFSF11 genes were identified as statistically significant in establishing the regression model of phalanges tuberositas and gene expression (Table 9). The FZD1 gene was identified as a statistically significant factor in establishing the regression model of brachydactylia and gene expression (Table 10).
3.5. The Expression of FZD1 in Cartilage Tissue from Patients with KBD
To evaluate the expression level of FZD1 in patients with KBD, FZD1 in articular cartilage tissues was detected using immunohistochemistry (IHC). In cartilage, positive staining for FZD1 was found mainly in the cytomembrane and cytoplasm of the superficial and middle zones of articular cartilage from the adult control group. FZD1-positive staining was lower in the superficial and middle zones in the adult KBD samples than in the adult control samples (Figure 1).
This study first investigated the comprehensive gene expression of patients with KBD in relation to the patients’ clinical variables. The most important finding is that the expression levels of FZD1 and TNFSF11 were statistically significant in establishing the regression model of phalanges tuberositas, brachydactylia, and gene expression. Our study suggests that the molecular and biological effect of KBD could be linked with KBD clinical variables.
At present, the diagnosis and determination of KBD grades are primarily dependent on symptoms and changes observed by X-ray, but it is only effective for diagnosing the advanced cases of KBD with grades II and III. Thus, it is extremely difficult to correctly diagnose the early stages of KBD with grades I and II, leading to inevitable exacerbation of symptoms [5, 14]. In this study, we found that brachydactylia, deformity of the knee, and dyskinesia in joints of the elbow and knee were significantly different between KBD patients with grades I and II, which could provide the new indicators for diagnosis of KBD patients with grades I and II.
The Manouvrier’s Skelic Index, metacarpal length/metacarpal breadth index, and metacarpal length/height index were used for reflecting the severity of KBD. The patients with KBD included in this study mainly consisted of grade I and II patients, whose lesions mainly occurred on their hands and joints. Therefore, gene correlation analysis showed that TNFSF11, GDF5, CTSC, and FZD1 were strongly and significantly associated with metacarpal length/metacarpal breadth index regardless of the Manouvrier’s Skelic Index and metacarpal length/height index.
Differential expression rates of a set of genes (ATR, CSGALNACT, CTSC, COL1A1, FZD1, GDF5, SLC14A1, SSBP1, and TNFSF11) were associated with KBD grades I and II; these genes tended towards abnormal expression in grades II, and all of these genes were differentially expressed in patients with phalanges tuberositas and brachydactylia except CSGALNACT, COL1A1, and SSBP1. Some of the genes above have also been reported to be differentially expressed in chondrocytes from patients with KBD [1, 5, 10]. This observation is consistent with our analysis and suggests that the involvement of these genes in pathophysiological changes is associated with clinical variables in KBD.
Ageing and joint disease are characterized by disruption of the equilibrium within cartilage tissue, in which the synthesis of new matrix components is more than the loss of collagens and proteoglycans [15–17]. This imbalance between anabolic and catabolic processes may result in progressive cartilage degeneration. The Wnt signalling pathway, composed of several important components, such as frizzled, a cell membrane receptor [18–20], plays an important role in embryonic development, skeletal development, and the maintenance of cartilage homeostasis [21–24]. Cartilage homeostasis in healthy joints is maintained by the balance of synthesis and degradation of extracellular matrix (ECM). Cartilage undergoes destruction when this balance is lost. In this study, we found that FZD1 (a member of the frizzled family) was downregulated in the peripheral blood of KBD patients and the expression level of FZD1 was associated with phalanges tuberositas and brachydactylia, two clinical symptoms of KBD analysed by multiple logistic regression analysis in the current study. Additionally, IHC was performed to evaluate the expression of these genes in articular cartilage from adult and child patients with KBD. From the results of IHC, expression levels of FZD1 were significantly reduced in patients with KBD compared with normal controls. It suggested that the abnormal expression of FZD1 in the articular cartilage of patients with KBD results in disequilibrium of the anabolism and catabolism of ECM by abnormally activating the Wnt signalling pathway. This disequilibrium in ECM involved in the abnormal expression of FZD1 could be one of the reasons for the process of pathological change in phalanges tuberositas.
Apoptotic chondrocytes are increased in the middle zone of articular cartilage in adult KBD patients in comparison to healthy controls, with the exception of multiple focal chondronecrosis [6, 25]. It suggested that the signal for apoptosis has already occurred in the middle zone of articular cartilage in KBD children. In previous studies, Bcl-2, Bax, and Fas were expressed in enhanced amounts in the upper zone and middle zone, and iNos expression was more abundant in the entire articular cartilage [6]. The TNFSF gene encodes a member of the tumor necrosis factor (TNF) cytokine family which is a ligand for osteoprotegerin and functions as a key factor for osteoclast differentiation and activation. This protein was shown to activate antiapoptotic kinase AKT/PKB through a signalling complex involving SRC kinase and tumor necrosis factor receptor-associated factor (TRAF) 6, which indicated this protein may have a role in the regulation of cell apoptosis [26, 27]. In this study, we found that TNFSF11 was upregulated in the peripheral blood of KBD patients and the expression level of TNFSF11 was associated with phalanges tuberositas in multiple logistic regression analysis. These results suggest that the abnormal expression of relevant chondrocyte molecules triggered apoptosis in KBD cartilage. Abnormal chondrocyte apoptosis can not only affect the development of bone and cartilage but also lead to phalanges tuberositas.
4. Conclusions
In summary, our study first investigated the correlation between gene expression and clinical parameter indexes in KBD. The results from multiple logistic regression analysis and immunohistochemistry indicate that FZD1 could be involved in the process of pathological change in phalanges tuberositas/brachydactylia and may provide further new information contributing to explanations for the articular cartilage destruction and pathology observed in patients with KBD.
Data Availability
The gene expression data used to support the findings of this study can be searched at https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE59446.
Conflicts of Interest
The authors declare that they have no conflict of interest.
Acknowledgments
This work is supported by the General Financial Grants from the China Postdoctoral Science Foundation (2017M613153 and 2017M623197), the National Natural Science Foundation of China (81472924, 81803179, and 81803178), and the National Key R&D Program of China (2016YFE0119100); the Fundamental Research Funds for the Central Universities; and the Shaanxi Natural Science Basic Research Project-General Project (Youth) (2018JQ8025 and 2018JQ8012).