Objectives To elucidate the functional epigenomic landscape of articular cartilage in osteoarthritis (OA) affected knee and hip joints in relation to gene expression.
Methods Using Illumina Infinium HumanMethylation450 BeadChip arrays, genome-wide DNA methylation was measured in 31 preserved and lesioned cartilage sample pairs (14 knees and 17 hips) from patients who underwent a total joint replacement due to primary OA. Using previously published genome-wide expression data of 33 pairs of cartilage samples, of which 13 pairs were overlapping with the current methylation dataset, we assessed gene expression differences in differentially methylated regions (DMRs).
Results Principal component analysis of the methylation data revealed distinct clustering of knee and hip samples, irrespective of OA pathophysiology. A total of 6272 CpG dinucleotides were differentially methylated between the two joints, comprising a total of 357 DMRs containing 1817 CpGs and 245 unique genes. Enrichment analysis of genes proximal of the DMRs revealed significant enrichment for developmental pathways and homeobox (HOX) genes. Subsequent transcriptomic analysis of DMR genes exposed distinct knee and hip expression patterns.
Conclusions Our findings reveal consistent DMRs between knee and hip articular cartilage that marked transcriptomic differences among HOX genes, which were not reflecting the temporal sequential HOX expression pattern during development. This implies distinct mechanisms for maintaining cartilage integrity in adulthood, thereby contributing to our understanding of cartilage homeostasis and future tissue regeneration approaches.
- Gene Polymorphism
Statistics from Altmetric.com
Articular cartilage (AC) is a highly specialised and characteristic tissue in all synovial joints at the ends of longitudinal bones. Its main function is to facilitate protection of subchondral bone against heavy loads, while maintaining smooth locomotor function of the articular joint.1 AC of the load bearing knee and hip joints is morphologically similar, both in health and disease.2–4 Furthermore, histological assessment and expression profiling of preserved and osteoarthritic AC have revealed generic processes and pathways to be involved in osteoarthritis (OA) pathophysiology, independent of the affected joint.5–7 Nevertheless, epidemiological studies and genome-wide approaches have respectively shown distinct prevalence patterns and genetic risk factors for OA at different joints.8 ,9 Moreover, although pathway analysis of gene expression data has revealed transcriptomic commonalities between knee and hip AC in OA, individual gene expression differences have been reported.7
In general, tissue identity is marked by the epigenetic landscape of respective cells and is, among others, reflected in the DNA methylation profile.10 DNA methylation, in which the cytosine residue in cytosine-phosphate-guanine dinucleotides (CpGs) acquires a methyl group, is known to regulate gene expression upon environmental changes such as age and disease. Overall differences in the methylome on the tissue level are commonly reflected in differentially methylated regions (DMRs), while single CpGs do usually not harbour this property10 ,11 and possibly only mark environmental, stochastic or individual differences. These observations raise the question whether knee and hip AC are either epigenetically distinct or similar tissues. Although some differences on the epigenetic level between knee and hip joints have been reported,12 it is currently unknown whether the observed differences in DNA methylation have any functional properties in terms of regulating expression of putative joint specific genes.
In the current study, pairwise preserved and lesioned AC from knee and hip joints was sampled from patients undergoing joint replacement surgery due to primary OA. Genome-wide DNA methylation was measured to assess communalities and discrepancies of the AC methylome in knee and hip joints and with respect to preserved and lesioned AC. Furthermore, we have subsequently combined epigenomic and transcriptomic data to gain a functional understanding of the observed methylation differences. To our knowledge, this is the first study in which highly similar tissues, being AC from either knee or hip joints, are compared comprehensively on the epigenomic and transcriptomic level.
Materials and methods
The RAAK cohort and sampling
Ethical approval was obtained from the medical ethics committee of the LUMC (P08.239) and informed consent was obtained from all participants.6 Participant details are listed in online supplementary table S1 (for sampling details see online supplementary methods and13). Macroscopically preserved as well as macroscopically lesioned cartilage was sampled from patients who underwent a total joint replacement due to primary OA of either the knee (N=14) or hip (N=17). From an additional 3 knee and 3 hip joints, healthy cartilage was sampled.
DNA was isolated using the Promega Wizard Genomic DNA Purification kit according to the manufacturer's protocol.
DNA was bisulfite treated using the ZymoResearch EZ DNA Methylation kit. DNA methylation was assessed using Illumina Infinium HumanMethylation450 BeadChips. All methylation values are reported as fractions between 0 and 1, commonly known as the β value. For additional details see online supplementary methods.
Normalised expression data were downloaded from GEO (GSE57218).13
Principal component analysis (PCA) and statistical procedures were carried out in R-3.0.2. All analyses were corrected for technical covariates as well as sex, disease status, age and Body Mass Index (BMI). A random effect for patient ID was included to correct for putative correlations between preserved and OA affected AC from the same joint. Additional details are listed in the supplementary methods.
Knee and hip AC show distinct methylation profiles
Genome-wide DNA methylation profiling was performed in all samples, consisting of both macroscopically preserved and lesioned AC derived from 14 knee and 17 hip joints. By means of PCA we observed two distinct clusters of samples reflecting the joint type from which the cartilage was sampled (figure 1A). Although in knee samples preserved AC tended to cluster apart from lesioned AC, clustering by joint type was largely independent of the OA affection status (figure 1C–F). This was further emphasised when 6 truly healthy samples (3 knees, 3 hips) were included in the PCA, which showed that irrespective of age and disease, samples clustered according to their joint type (figure 1B). Next, to elucidate the specific CpGs driving this distinct clustering by joint type, we fitted a linear mixed model to identify the specific CpGs that were differentially methylated between hip and knee cartilage, while correcting for sex, age, BMI and OA affection status. After adjustment for multiple testing (Benjamini–Hochberg), 6272 CpGs were significantly differentially methylated between knee and hip AC (p<0.05) by at least 0.1 β, covering a total of 2726 unique genes (see online supplementary table S2).
DMRs among HOX containing genes
Next, we applied a sliding window algorithm10 to distinguish inherent tissue differences in the methylation data from possible environmental, stochastic or individual differences and thereby observed 357 DMRs, consisting of 1817 CpGs and 245 unique genes (see online supplementary table S3). Pathway analysis revealed significant enrichment among the constructed DMRs mostly for developmental pathways (such as limb development and skeletal system morphogenesis) and, more specifically, homeodomain containing genes (see online supplementary table S4). Rather strikingly in this respect is the presence of 42 DMRs in all four canonical homeobox (HOX) clusters, comprising over 10% of the observed DMRs. Here again we observed no major distinction between OA affected joints and healthy joints, confirming that the DMRs are highly joint specific (see online supplementary figure S1). Visual inspection of representative CpGs in DMRs ratified the distinct and consistent differences in methylation between the two joint types, while this was less evident between preserved and OA affected AC (figure 2).
Putative functionality of DMRs in adult AC
Finally, to investigate the putative functionality of the observed DMRs, in terms of respective mRNA expression of proximal genes, we assessed the expression patterns of the 245 DMR associates genes (see online supplementary table S3). Previously, gene expression was quantified for 33 pairs of preserved and lesioned knee (N=11) and hip (N=22) AC, of which 13 pairs were overlapping with the methylation data (GSE57218).13 Hereby, differential expression analysis of the entire GSE57218 dataset revealed that independent of the OA affection status, 28 out of 245 genes were differentially expressed between knee and hip AC (table 1). Respectively, 6 and 11 genes were only expressed in either knee or hip AC, while 11 genes were expressed significantly different between knee and hip AC. Among the 28 differentially expressed genes, genes from all four HOX clusters were present, as well as multiple HOX containing cofactors (such as PITX1, MEIS2, DLX5 and IRX3).
In the current study, we report on differences in the epigenetic landscapes between knee and hip AC. Based on the entire DNA methylation landscape, knee and hip AC show distinct epigenomic profiles independent of the tissue disease state (figure 1). Subsequent indepth analysis of the CpG dinucleotides conferring these distinct profiles (figure 2, see online supplementary table S2) revealed significant enrichment for developmental genes such as the canonical homeotic clusters and HOX cofactors (see online supplementary table S4). Furthermore, integration of epigenomic and transcriptomic data revealed significant differences in expression among these enriched loci between knee and hip AC, mediated by tissue specific DMRs (table 1).
In order to identify generic and subsequently functional joint related changes, as opposed to possible stochastic, environmental or individual related differences, we have constructed DMRs, since DMRs are known to consistently reflect the tissue of origin.10 ,11 We here report on the fact that despite the morphological and functional similarities between knee and hip AC,1–4 they contain inherently distinct cellular phenotypes based on their functional epigenomic landscape. Nonetheless, our data show that although methylation profiles at DMRs are highly tissue specific, they do not necessarily correlate to gene expression, as only a minority of DMRs appear to be associated with joint specific gene expression differences (28 out of 245 genes). This emphasises the need for comprehensive integration of multiple levels of genome-wide data, such as transcriptomics, for the interpretation of epigenomic studies in OA.
Although specific HOX gene functions remain partly elusive, in part due to complex interactions with HOX cofactors, increasingly more developmentally distinct functions are being ascribed to the various HOX genes.14 ,15 In adult tissues, however, regulation of HOX gene expression and their respective function remains largely unknown and is likely tissue specific.16 ,17 Multiple studies have reported on distinct expression patterns of homeotic genes and related HOX cofactors in adult tissues, reflecting the collinear embryonic HOX code.16–18 However, up to date not much is known about AC in this respect. If the observed differences in this study are due to retainment of the embryonic HOX code, i.e. the spatiotemporal expression pattern during development, then across the four canonical HOX clusters similar differences in methylation between knee and hip AC were to be expected. Here, however, functional differences in DNA methylation across the HOX clusters were observed, while the embryonic HOX code or colinearity were absent, as is reflected by unique knee and hip methylation patterns observed across the four HOX clusters (see online supplementary figure S1). This observation suggests specific functional roles for the basal HOX transcription factors and likely marks differences in cellular identity between chondrocytes residing in either knee or hip AC.
In the field of tissue engineering, major efforts are made to understand cartilage homeostasis, thereby contributing to the development of novel therapeutic approaches for treatment of degenerative joint diseases, including OA.19 The results presented in the current study, suggesting differences in cellular identity between chondrocytes residing in either knee or hip AC independent of OA pathophysiology, could putatively have implications for future regenerative approaches. As HOX genes are crucially involved in AC development,20 ,21 the observed epigenomic and transcriptomic differences in this study could indicate that directing articular chondrocytes into extracellular matrix production and/or active remodelling of damaged AC could hypothetically be achieved only via distinct mechanisms, depending on the joint type of a cartilage lesion. Moreover, it has been shown that expression of certain sets of HOX genes regulate the regenerative propensity of neural crest cells,22 presumptively indicating that chondrocytes originating from either knee or hip AC exhibit unequal regenerative capacities.
To our knowledge, we are the first to report on inherent differences between knee and hip AC by virtue of the joint specific epigenetically regulated transcriptomic landscape of HOX gene clusters and related cofactors. Recent studies on genome-wide DNA methylation in AC have focused on either the knee12 ,23 or hip12 ,24 joints and more primarily on comparing OA affected with control tissues. Nevertheless, Rushton et al12 did report epigenetic differences between knee and hip AC at specific CpGs. However, they have only analysed single CpGs as opposed to DMRs, did not report on the specific HOX loci nor for that matter did they study the functionality of the reported differences in terms of transcriptomic regulation. We have here consequently shown that especially the last is crucial to comprehend the results from epigenetic studies.
Although we here have jointly analysed preserved and lesioned cartilage originating from the same joint, as shown in figure 2 there exists only little epigenetic variation within sample pairs, while the joint specific differences are markedly present. This is even further pronounced when non-OA samples, derived from healthy joints, cluster tightly according to their joint specific epigenome (figure 1 and see online supplementary figures S1 and S2). Furthermore, putative confounding due to pooling is corrected for in all statistical analyses by including a random effect for patient ID. Of note, we do neither rule out nor disregard the presence of relevant gene specific epigenetic differences between preserved and lesioned AC. In the current study, however, we have restricted our focus on inherent epigenetic tissue differences between knee and hip joints.
In conclusion, we have observed consistent DMRs between knee and hip joints among HOX domain containing genes, both in the four canonical homeotic clusters as well as HOX cofactors. They were found to mark differential expression of genes residing in or near these DMRs in AC of knee and hip joints. The different methylation profiles of knee and hip AC likely mark distinct cellular identities, which could have relevant implications for the field of AC tissue engineering. Together, these findings contribute to our understanding of cartilage homeostasis and future repair strategies.
We thank all the participants of the RAAK study.
This web only file has been produced by the BMJ Publishing Group from an electronic file supplied by the author(s) and has not been edited for content.
Files in this Data Supplement:
Handling editor Tore K Kvien
Contributors All authors meet the ICMJE recommendations. More specifically: WdH, YFMR, SDB, NB, RvdB, NL and WJdD were involved in sampling and wet-lab work. WdH performed all analyses and dry-lab work. BJD and RGHHN were involved in patient communication and subject acquisition. WdH, YFMR, NB, PES and IM were involved in study design development and data interpretation. All authors have critically reviewed the manuscript.
Funding The RAAK studies were supported by the Leiden University Medical Centre, the Dutch Arthritis Association (DAA 101-402 and Reumafonds LRR) and the Centre of Medical System Biology and Netherlands Consortium for Healthy Aging, both in the framework of the Netherlands Genomics Initiative (NGI). Furthermore, we acknowledge support by TreatOA and IDEAL, which are funded by the European Union's Seventh Framework Program (FP7/2007-2011) under respective grant agreement nos. 200800 and 259679. The funders had no role in study design, data collection and analysis, decision to publish or preparation of the manuscript.
Competing interests None.
Patient consent Obtained.
Ethics approval LUMC Medical Ethics Committee (P08.239).
Provenance and peer review Not commissioned; externally peer reviewed.
If you wish to reuse any or all of this article please use the link below which will take you to the Copyright Clearance Center’s RightsLink service. You will be able to get a quick price and instant permission to reuse the content in many different ways.