Skip to main content
Download PDF

The role of DNA methylation in chondrogenesis of human iPSCs as a stable marker of cartilage quality

Clinical Epigenetics volume 16, Article number: 141 (2024) Cite this article

Abstract

Background

Lack of insight into factors that determine purity and quality of human iPSC (hiPSC)-derived neo-cartilage precludes applications of this powerful technology toward regenerative solutions in the clinical setting. Here, we set out to generate methylome-wide landscapes of hiPSC-derived neo-cartilages from different tissues-of-origin and integrated transcriptome-wide data to identify dissimilarities in set points of methylation with associated transcription and the respective pathways in which these genes act.

Methods

We applied in vitro chondrogenesis using hiPSCs generated from two different tissue sources: skin fibroblasts and articular cartilage. Upon differentiation toward chondrocytes, these are referred to as hFiCs and hCiC, respectively. Genome-wide DNA methylation and RNA sequencing datasets were generated of the hiPSC-derived neo-cartilages, and the epigenetically regulated transcriptome was compared to that of neo-cartilage deposited by human primary articular cartilage (hPAC).

Results

Methylome-wide landscapes of neo-cartilages of hiPSCs reprogrammed from two different somatic tissues were 85% similar to that of hPACs. By integration of transcriptome-wide data, differences in transcriptionally active CpGs between hCiC relative to hPAC were prioritized. Among the CpG-gene pairs lower expressed in hCiCs relative to hPACs, we identified genes such as MGP, GDF5, and CHAD enriched in closely related pathways and involved in cartilage development that likely mark phenotypic differences in chondrocyte states. Vice versa, among the CpG-gene pairs higher expressed, we identified genes such as KIF1A or NKX2-2 enriched in neurogenic pathways and likely reflecting off target differentiation.

Conclusions

We did not find significant variation between the neo-cartilages derived from hiPSCs of different tissue sources, suggesting that application of a robust differentiation protocol such as we applied here is more important as compared to the epigenetic memory of the cells of origin. Results of our study could be further exploited to improve quality, purity, and maturity of hiPSC-derived neo-cartilage matrix, ultimately to realize introduction of sustainable, hiPSC-derived neo-cartilage implantation into clinical practice.

Background

Osteoarthritis (OA) is an age-related, heterogeneous joint disease, with degeneration of articular cartilage as an important hallmark of the pathophysiological process [1]. Clinical manifestations of OA are synovitis, chronic pain, stiffness, and joint deformities with associated dysfunction of the joint and loss of mobility that severely hampers the daily activities of patients [2]. Despite the critical impact on patients, there are currently no effective treatments that either stop or reverse OA except for costly joint replacement surgery at end-stage disease. Drawback of such replacements is the relatively limited lifespan of the implants, frequently leading to postponed surgery of patients particularly between 50 and 70 years of age with concurrent prolonged chronic pain and absence from work.

A major hurdle for the development of therapy in OA is the deficient inherent repair capacity of chondrocytes, the single-cell type present in articular cartilage [3]. On the other hand, articular cartilage is considered immunotolerant and therefore eligible for regenerative treatment strategies such as implantation of neo-cartilage [4]. To date, neo-cartilage implantation has mainly relied on autologous cell sources such as chondrocytes or mesenchymal stem cells (MSCs). Indeed, we showed previously in vitro that differentiated human primary articular chondrocytes (hPAC) readily deposit neo-cartilage that exhibits a DNA methylation landscape with an almost identical (99% similarity) profile to autologous cartilage [5]. Upon expansion, chondrocytes are, however, prone to dedifferentiate and lose their characteristic articular phenotype, resulting in deposition of a non-specific, mechanically inferior extracellular matrix which is less suitable for application in larger defects [6]. Similarly, autologous MSCs as chondrocyte precursor source are subjected to limited proliferation capacity resulting in decreased chondrogenic potential. Even more important, it has been generally acknowledged that due to its developmental fate, MSC-derived neo-cartilage is prone to deposit hypertrophic, mineralized cartilage [7].

To accomplish strong, durable, neo-cartilage implants at large scale, human-induced pluripotent stem cells (hiPSCs) emerge as an exciting novel cell source. hiPSCs have infinite expansion capacity enabling the generation of large quantities of differentiated chondrocytes for cartilage regeneration [8], while their use is minimally invasive [9]. To accommodate hiPSCs as sustainable source for articular cartilage regeneration, differentiation protocols resulting in human-induced chondrocytes that readily deposit cartilaginous tissue have been developed [10, 11]. Nonetheless, lack of insight, hence solutions, into aberrant cell fate decisions during hiPSC chondrogenesis hence quality and purity of deposited neo-cartilages, precludes application of this powerful technology toward regenerative solutions in the clinical setting. Such cell fate decisions are particularly facilitated by deposited methylation at transcriptionally active CpG sites that evoke stable set points of transcription [12]. Moreover, it was shown that the majority of such lineage-specific methylation patterns remain stable throughout life [13]. Owing to the notion that such an epigenetic memory of the somatic cell types exists [14,15,16,17], hiPSCs derived from chondrocytes may differentiate more readily and with higher similarity to articular cartilage than hiPSCs derived from cells of other tissues.

In the current study, we set out to test whether the epigenetic memory of the somatic cell type used for generation of hiPSCs is retained as a unique DNA methylation signature. Hereto, we applied established robust stepwise chondrogenesis protocol via mesodermal lineage [11] using hiPSCs from two different tissue sources, namely skin (fibroblasts) and articular cartilage (chondrocytes). In vitro hPAC-derived neo-cartilage was used as the golden standard [5]. To assess similarities, the methylome-wide landscape of neo-cartilages derived from these different cell sources was generated and compared. To subsequently explore and biologically interpret differences in cell identity between hPACs and hiPSC-derived chondrocytes in the neo-cartilages, we set out to specifically study discordant set points of DNA methylation that likely act on expression of positional genes. Since this layer of molecular information is key in critical lineage decisions, hence cell fates [12] data generated in this study could pave the way for large-scale hiPSCs-derived neo-cartilage formation of superior quality.

Methods

Tissue culture and chondrogenesis

Cell culture of hiPSCs and primary articular chondrocytes

To analyze potential differences in chondrogenic differentiation potential among hiPSCs generated from different cell sources (skin fibroblasts and articular cartilage; Fig. 1A), five independent hiPSCs lines/clones were used in the current study. All hiPSC cell lines were generated by the LUMC hiPSC Hotel while applying transformation with polycistronic lentivirus to induce expression of the Yamanaka factors as outlined in detail by Dambrot et al. [18]. The two skin fibroblast-derived hiPSCs have been officially registered at the human pluripotent stem cell registry (hPSCREG) and are from healthy donors without known genetic diseases. Lines are from a female donor (LUMC0030iCTRL12 or hPSCreg line LUMC004-B) and from a male donor (LUMC0004iCTRL10 or hPSCreg line LUMC029-B). The three chondrocyte-derived hiPSCs were generated from macroscopically preserved chondrocytes of a female OA donor who underwent a joint replacement surgery from the RAAK study [19] and are referred to as LUMC0131iCTRL02, LUMC0131iCTRL04, and LUMC0131iCTRL05. The generation of the hiPSCs line was approved by the Leiden University Medical Centre ethical committee under P13.080. Cells were characterized according to pluripotent potential and spontaneous differentiation capacity by the hiPSC core Hotel [18] and were karyotyped after 15 passages in culture (Fig. S1). hiPSCs were maintained under standard conditions (37 °C, 5% CO2) in TeSR-plus medium (STEMCELL Technologies) on VitronectinXF-coated plates (STEMCELL Technologies). The medium was refreshed three times a week, and cells were passaged in aggregates using Gentle Cell Dissociation Reagent (STEMCELL Technologies) upon reaching approximately 80% confluency. hiPSCs were used for our study between passages 21 and 27. Human primary articular cartilage was collected from the hip of N = 10 different OA patients that underwent joint replacement surgery as part of the RAAK study [19]. Of note, classification of OA as macroscopically preserved or lesioned for collection, expansion, and differentiation of the primary articular chondrocytes has been previously described [5]. Ethical permission for the described studies was obtained from the appropriate medical ethical committee under protocol numbers P08.239 and P19.013. Written informed consent was obtained from all participants.

Fig. 1
figure 1

Generation of 3D neo-cartilage derived from hiPSCs sources relative to hPAC. A Schematic representation of the study. Samples consisting of neo-cartilage from human chondrocyte-derived iPSCs (hCiC, N = 14), neo-cartilage from human fibroblast-derived iPSCs (hFiC, N = 3), and neo-cartilage from human primary articular chondrocytes (hPAC, N = 10). B Timelines and times of collection for the analyses of stepwise differentiation from hiPSCs toward paraxial mesoderm lineage and chondrocyte-like cells

Full size image

hiPSC differentiation to human chondroprogenitor cells

Generation of chondroprogenitor cells was based on a protocol previously described [11]. When hiPSCs reached 60% confluence, the culture medium was switched to mesodermal differentiation (MD) medium, composed of IMDM GlutaMAX (IMDM; Thermo Fisher Scientific) and Ham’s F12 Nutrient Mix (F12; Sigma-Aldrich) with 1% chemically defined lipid concentrate (Gibco), 1% insulin/human transferrin/selenous (ITS + ; Corning), 0.5% penicillin–streptomycin (P/S; Gibco), and 450 μM 1-thioglycerol (Sigma-Aldrich). Before induction of anterior primitive streak (day 0), hiPSCs were washed with wash medium (IMDM/F12 and 0.5% P/S) and then fed with MD medium supplemented with activin A (30 ng/ml; Stemgent), 4 μM CHIR99021 (CHIR; Stemgent), and human fibroblast growth factor (20 ng/ml; FGF-2; R&D Systems) for 24 h. Subsequently, the cells were washed again with wash medium, and paraxial mesoderm was induced on day 1, by MD medium supplemented with 2 μM SB-505124 (Tocris), 3 μM CHIR, FGF-2 (20 ng/ml), and 4 μM dorsomorphin (Tocris) for 24 h. Before induction of early somite (day 2), cells were washed with wash medium, and then, cells were fed with MD medium supplemented with 2 μM SB-505124, 4 μM dorsomorphin, 1 μM C59 (Cellagen Technology), and 500 nM PD173074 (Tocris) for 24 h. Subsequently, cells were washed with wash medium, and for induction of sclerotome, cells (days 3 to 5) were fed daily with MD medium supplemented with 2 μM purmorphamine (Stemgent) and 1 μM C59. To induce chondroprogenitor cells (days 6 to 14), cells were washed briefly with wash medium and fed daily with MD medium supplemented with human bone morphogenetic protein 4 (BMP-4; 20 ng/ml; Miltenyi Biotec). Five independent differentiations were performed per clone, and two neo-cartilages were pooled for the RNA and DNA isolations. Figure 1B shows a schematic overview of the timelines and times of collection for the analyses.

Chondrogenic differentiation of chondroprogenitor cells

Monolayer cultured chondroprogenitor cell aggregates present at day 14 of the differentiation were washed with MD medium, dissociated using Gentle Cell Dissociation Reagent (STEMCELL Technologies), and centrifuged for 4 min at 1200 rpm. Cell aggregates (250,000 cells per pellet) were subsequently maintained in Dulbecco’s modified Eagle’s medium/F12 (Gibco), supplemented with 1% ITS + , 55 μM 2-mercaptoethanol (Gibco), 1% non-essential amino acids (Gibco), 0.5% P/S, l-ascorbate-2-phosphate (50 μg/ml; Sigma-Aldrich), l-proline (40 μg/ml; Sigma-Aldrich), and transforming growth factor–β1 (10 ng/ml; PeproTech) for 7/21/35 days while refreshing medium every 3 to 4 days.

Histology and immunohistochemistry

3D neo-cartilage was fixed in 4% formaldehyde overnight and stored in 70% ethanol at 4 °C. The pellets were embedded in paraffin and sectioned at 5 μm. After sectioning, slides were deparaffinized and rehydrated prior to histology or immunohistochemistry. Overall cellular and tissue structure was visualized with hematoxylin–eosin (HE) staining. Glycosaminoglycan depositions were visualized by staining with 1% Alcian Blue 8-GX (Sigma-Aldrich) and nuclear fast red staining (Sigma-Aldrich). Alcian blue staining was quantified by loading the images in Fiji and splitting the color channels. Subsequently, gray values were measured of three to five separate squares per pellet and corrected for the gray value of the background. To detect collagen type ll (1:100; ab34712, Abcam) immunohistochemistry was performed with 3-diaminobenzidine (DAB) solution (Sigma-Aldrich) and hematoxylin (Klinipath) as described before (5). Pixel intensity quantification was performed by Fiji, and surface area of the pellets was measured with the CellSens Dimension software (Olympus, Leiderdorp, The Netherlands).

RNA analyses

RNA isolation and RT-qPCR

For RNA isolations neo-cartilages at day 35 of differentiation were used. Hereto, two pellets were pooled, and isolation was performed as described previously [5]. Total mRNA (150 ng) was processed with a first strand cDNA kit according to the manufacturer’s protocol (Roche Applied Science). cDNA was further diluted five times, and preamplification with TaqMan preamp master mix (Thermo Fisher Scientific Inc.) was performed for genes of interest. Gene expression was measured with RT-qPCR, and average of the two biological replicates was determined as relative levels (− ΔCt values) using expression levels of glyceraldehyde 3-phosphate dehydrogenase (GAPDH) and Acidic ribosomal phosphoprotein P0 (ARP) as housekeeping genes. P < 0.05 were considered statistically significant. Quality control of the results was performed as described before [20].

RNA sequencing

RNA was extracted from chondrocytes isolated from neo-cartilages at day 14 (hPACs) and day 35 (hCiCs and hFiCs) of differentiation. Hereto, RNA was extracted with chloroform and purified using the RNeasy Mini Kit (QIAGEN). RNA sequencing (polyA enriched) of 22 total RNA samples (hCiC; n = 9, hFiC; n = 3, and hPAC; n = 10 derived neo-cartilage) was performed using Illumina NovaSeq6000 sequencing, according to the standard operating procedures based on the Illumina protocol for Paired-End Sequencing (Per sample ~ 6 Gb, 20 million Paired-End reads RNA-sequencing mapping and quality control (QC) was performed using the in house pipeline BioWDL (https://biowdl.github.io/RNA-seq/v4.0.0/index.html) developed by the SASC team at Leiden University Medical Center. This pipeline was used to process FastQ files using Picard.v2.23.2 and Samtools.v1.10, adapter clipping using cutadapt.v2.10. Moreover, QC was performed using FastQC.v0.11.9 and MultiQCv. 1.9. Furthermore, mapping was performed with STAR.v2.7.5a software [21] and expression quantification and transcript assembly using HTSeq-Count.v.0.12.4 [22]. The reads were aligned to the reference transcriptome Ensemble GRChr38 release 102. RNA-seq data were normalized using the DESeq2_v.1.30.0 R package [23]. The data were further transformed using the variance-stabilizing transforming (VST) method.

Differential expression analysis

Differential expression analysis was performed using the DESeq2 package v. 1.30.0 using R version 4.0.2. A general linear model (GLM) assuming a negative binomial distribution was applied followed by a Wald-test between hCiC and hPAC samples. Benjamini–Hochberg multiple testing-corrected P values with significance cut-off of 0.05 are reported as false discovery rate (FDR).

DNA analyses

DNA isolation

Snap frozen neo-cartilage was powdered using a Mixer Mill 200 (Retsch, Germany) with continuous liquid nitrogen cooling. DNA was isolated using the Wizard Genomic DNA Purification kit (Promega) according to the manufacturer's protocol.

Methylation analysis

DNA methylation was assessed in chondrocytes isolated from neo-cartilages at day 14 (hPACs) and day 35 (hCiCs and hFiCs) of differentiation. Hereto, the Illumina Infinium Methylation EPIC (850 K) BeadChip was used according to standard operating procedures based on the Illumina Infinium II Protocol. To analyze methylation array data (MethylationEPIC 850 k array), the MethylAid R script [24] with the default settings was used for data quality assessment. All samples showed a detected CpG above 95%. The minfi.v_1.36.0 R package [25] was used to pre-process the data. We removed any probe that failed in one or more samples (p < 0.01). Probe level intensities were quantile normalized across samples prior to calculation of the ß-values.

MethylToSNP v0.99.0 R package was used to filter SNPs. This method looked for patterns in methylation array data and identified methylation probes with SNP-like patterns. This method allows to remove outliers, which adds robustness to the analysis and is enabled by default. A confidence score was calculated to show how close the observed pattern of methylation beta values was to a canonical case of a SNP in a homozygously methylated CpG locus. Additionally, MethylToSNP can overlap the SNPs identified in methylation data with known SNPs from dbSNP. The probes that have shown to be cross-reactive (demonstrated to map to multiple places in the genome) were filtered out [26]. The probes that were overlapping with rare SNPs (probes in TFBS that showed extreme methylation pattern) were filtered out [27]. To minimize the unwanted variation within and between samples, we used the functional normalization method from the minfi.1.36.0 R package [28].

Differential methylation analysis

We run differential mean analysis using t-moderated statistics. Using the MEAL.1.20.3 R package pipeline, which, relies on the lmFit from limma R package (design model =  ~ phenotype). CpGs after Bonferroni correction P < 6.243109e-08 (0.05/800883) was considered significant.

Similarities

To identify similarities between the cells, we used the Jaccard method with the jaccard_0.1.0 R package. Visualization was done using the corrplot_0.84 R package.

Integration methylation and RNA-seq

To identify transcriptive active CpGs, we integrated significant CpGs with expression of genes using the MEAL_v1.20.3 R package. The ß-values of the methylation data were used. We adapted the CorrelationMethExprs function, which estimates the correlation between methylation and expression. For each CpG, a range was defined by the position of the CpG plus the flank parameter (250 kb upstream and 250 kb downstream). Only those expressed genes that were entirely in this range were selected. After multiple testing correction by false discovery rate (FDR) and considering r2 > 0.5, CpGs correlated with genes were shown.

Gene ontology enrichment analysis

Gene ontology (GO) enrichment analysis was performed using clusterProfiler v4.0.5 R package to identify GO terms. Bonferroni multiple testing-corrected P values with a significance cut-off of 0.05 are reported as FDR.

Statistical analysis

Statistical analyses were performed in R version 4.0.2 (2020–06-22) and using SPSS version 25 (IBM). Results were considered significant for FDR ≤ 0.05.

Results

Generation and characterization of 3D neo-cartilage derived from different hiPSC sources relative to hPACs

After confirming adequacy of hiPSCs morphology, pluripotency, spontaneous differentiation capacity, and a normal karyotype [18] (Fig. S1A–D), human skin fibroblast-derived iPSCs (two cell lines; from two healthy donors) and human cartilage-derived iPSCs (three independent clones hence hiPSC lines; from macroscopically preserved cartilage chondrocytes of a patient who underwent joint replacement surgery due to OA) were differentiated toward chondrocytes (Fig. 1). Cellular and tissue structure during formation of the hiPSCs-derived neo-cartilages in comparison with that of hPACs was visualized by hematoxylin and eosin (H&E) staining at day 0, 7, 21, and 35 of differentiation from chondroprogenitor cells. At day 35, HE staining of both human cartilage-derived iPSCs differentiated to chondrocyte (hCi-Chondrocyte: hCiC) and human fibroblast-derived iPSCs differentiated to Chondrocyte (hFi-Chondrocyte: hFiC) showed similar matrix structure as compared to hPAC controls (Fig. S2A). Nevertheless, until day 21 of chondrogenic differentiation of both hCiC and hFiC also formation of off-target tissues was detected in developing neo-cartilage pellets as reflected by the darker purple edge of the neo-cartilages. Staining with cytokeratin (CK) showed that some of the off-target cells were keratinocytes. These were largely decreased at day 21 and no longer apparent at day 35 of neo-cartilage generation (Fig. S2B–D).

To validate glycosaminoglycan deposition in hCiC and hFiC relative to hPACs, Alcian blue staining was performed and quantified (Fig. 2A, B). Neo-cartilage from both sources of hiPSCs showed increasing amounts of glycosaminoglycan deposition, with 23-fold increase in hCiC and 13-fold increase in hFiC from chondroprogenitor stage to day 35 of differentiation (Fig. 2B). Nonetheless, the final level of glycosaminoglycan deposition in hPACs-derived neo-cartilage was 30% higher compared to that in hCiC (P < 0.05) and 20% higher compared to that hFiC (P < 0.05; Fig. 2B). Similarly, deposition of cartilaginous extracellular matrix determined by immunohistochemistry of collagen 2 protein (COL2; Fig. 2C, D) showed a significant increase from chondroprogenitor to day 35 of neo-cartilage across hiPSCs sources. However, final expression levels of COL2 in neo-cartilage in hPAC was 10% higher relative to both hCiC (P = 0.05) and hFiC (P = 0.05).

Fig. 2
figure 2

Characterization of 3D neo-cartilage derived from hiPSCs sources relative to hPAC. A Representative histological staining images of the neo-cartilages by Alcian Blue staining showing glycosaminoglycan deposition. B Quantification of Alcian Blue histology of hCiC and hFiC, respectively, in different time-points of differentiation compared to hPAC. C Representative immunohistochemistry images of neo-cartilage stained with Collagen-2 antibody. D Quantification of Collagen-2 immunohistochemistry of hCiC and hFiC, respectively, in different time-points of differentiation process compared to hPAC. Scale bars at day 7: 50 μm, other scale bars: 100 μm. All data are presented as mean ± STD. Two-tailed Student’s t tests were used to compare 1) each group relative to chondroprogenitor cells (stars on top of the bars) and 2) each group relative to hPAC (indicated by lines). *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001 (N = 6–9 pellets each)

Full size image

Next, neo-cartilage derived from hFiC and hCiC in comparison with hPAC was characterized by targeted RT-qPCR gene expression of well-known chondrogenic (COL2A1, SOX9 and ACAN), fibrotic (COL1A1), and hypertrophic (COL10A1) cartilage markers at day 35 of differentiation. In line with the histology, this showed that hPACs as compared to hiPSC-derived chondrocytes express significant higher levels of chondrogenic markers COL2A1, SOX9 and ACAN, but also of hypertrophic marker COL10A1, and fibrotic marker COL1A1 (Fig. S3). However, the COL2A1:COL1A1 ratios appeared highly comparable among hPACs, hCiCs and hFiCs (Fig. S3). Together, it can be concluded that neo-cartilage from both hiPSC sources at day 35 of differentiation has high similarity to the hPAC control pellets and was therefore used for further study.

Methylome-wide landscape of 3D neo-cartilage derived from hiPSCs sources relative to hPAC

Founded by the hypothesis that the epigenetic memory of the hiPSCs cartilage tissue-of-origin could contribute to chondrogenic potential of hiPSCs during differentiation, genome-wide DNA methylation datasets were generated from hiPSC derived neo-cartilage at day 35 of differentiation and from hPACs at day 14 of differentiation. Following quality control (QC), we obtained robust methylation data of 800,883 CpGs in total. As shown in Fig. 3A–C, we observe high overall Jaccard similarity scores within and between the DNA methylome-wide landscape from neo-cartilages of different cell sources. More specifically, as outlined in Fig. 3D, Jaccard similarities of 96% and 94% were observed in the DNA methylome-wide landscape from neo-cartilages of independent differentiations of, respectively, hFiC (N = 3) or hCiCs (N = 13), reflecting high consistency in the hiPSC chondrogenic differentiation protocol. Also, a Jaccard similarity of 92% was observed in the DNA methylome-wide landscape of neo-cartilages generated from hPACs of N = 10 independent donors, indicating high consistency and similarity in the in vitro neo-cartilages deposited by hPACs isolated from heterogeneous preserved cartilages of OA patients. As summarized in Fig. 3E, the Jaccard similarity between the methylome of the neo-cartilage of the hCiC and hFiC cell sources was also high with 92%, indicating consistency in the hiPSC chondrogenic differentiation protocol between hiPSC clones even when generated from different cell origins. Finally, in Fig. 3E, we outlined the average Jaccard similarities of 86% and 84% between the methylome-wide landscape of the neo-cartilages deposited by hPACs and, respectively, hFiC and hCiC. Together these data indicate that the applied hiPSC chondrogenic differentiation protocol is highly consistent within (~ 95%) and between (~ 92%) hiPSC lines, and that it resulted in very similar neo-cartilage as deposited by primary hPACs (~ 85%). Our data suggest that the epigenetic memory of hCiC (hiPSCs from cartilage as tissue-of-origin) did not have an additional positive effect on the chondrogenic differentiations. Despite the overall high similarities, principal component (PC) analysis on the methylome-wide landscape showed three specific clusters according to cell sources. Herein, hCiCs and hFiCs separated from hPACs by particularly PC1 explaining 96% of the variance and hCiCs separated from hFiCs by PC2 explaining 2% of the variation (Fig. S4A). Furthermore, Euclidean clustering analyses of the methylation of 1,000 most variable CpGs confirmed separation of hPACs, hCiCs, and hFiCs (Fig. S4B). These data indicated that the methylome landscape of neo-cartilages has also specific characteristics that reflect differences in the cell identity of the hCiC, hFiC or hPAC-derived chondrocytes.

Fig. 3
figure 3

DNA methylation landscape of 3D neo-cartilage. A Plot showing the Jaccard similarity within hCiC (N = 13), within hPAC (N = 10), and between hCiC vs hPAC. B Plot showing the Jaccard similarity within hCiC (N = 13), within hFiC (N = 3), and between two sources of hiPSCs-derived neo-cartilage hCiC vs hFiC. C Plot showing the Jaccard similarity within hFiC (N = 3), within hPAC (N = 10), and between hFiC vs hPAC. D Summary of the average Jaccard similarity index within the hCiC, hFiC and hPAC cell sources upon independent consecutive rounds of differentiation. E Summary of the average Jaccard similarity index between the hCiC, hFiC, and hPAC cell sources

Full size image

Transcriptome-wide landscape of 3D neo-cartilage-derived hiPSC sources relative to hPAC

To characterize and biologically interpret differences in the methylome-wide landscapes of neo-cartilage deposited by hiPSC-derived chondroprogenitors relative to hPACs, we next set out to prioritize on differences in the methylome that are most likely transcriptionally active. Hereto, we performed RNA sequencing of neo-cartilages derived from hPACs (N = 7, day 14 of differentiation), hCiCs (N = 9, day 35 of differentiation), and hFiCs (N = 3, day 35 of differentiation). Prior to assessing the differences in methylation that are most likely transcriptionally active, we first performed an exploratory analyses of the normalized VST values of the transcriptome of hCiCs relative to hPACs. Expression data of hFiCs are provided in Table S1, however, because of the low number of samples (N = 3) and the high similarity between the two hiPSCs cell sources, the hFiC transcriptomic data were not included in further downstream analyses. As shown in Fig. S5, the average expression levels of all genes in hPAC (Ave Exp = 4.9) were slightly lower as compared to hCiC (Ave Exp = 5.4), while the range in expression levels in hPACs (between 1.8 and 16) was larger as compared to hCiC (between 2.0 and 14). Moreover, as shown by the histogram the hPACs showed a large number of genes that are very lowly expressed (Fig. S5).

Upon subsequently assessing relative gene expression levels defined as the standard deviation of z-scores with SD1 being lowest expressed genes (< SD2) and SD4 being highest expressed genes (> SD2) we could explore genes that mark authentic human primary chondrocytes expression patterns in hPACs relative to that in hCiCs (Table S1). The STRING protein–protein interaction network of top 50 of 919 genes that were highly expressed in hPAC (> SD2) and hCiC (> SD2) showed a highly interconnected dense network with significant gene enrichment in pathways involved in collagen fibril organization (GO:0030199, P = 9.3 × 10−8), extracellular matrix organization (GO:0030198, P = 2.7 × 10−8), and cell adhesion molecular binding (GO: 0050839, P = 6.3 × 10−8), represented by well-known cartilaginous genes such as FN1, COL6A1/2, TNC, COL2A1, and MALAT1 [29], confirming quality of neo-cartilage deposited by hiPSCs-derived chondroprogenitors (Fig. S6). On the other hand, the STRING protein–protein interaction network of the n = 69 genes that were highly expressed in hPAC (> SD2) but lowly expressed in hCiC (< SD2) showed a network with notable significant gene enrichment in pathways involved in skeletal system development (GO:0001501, P = 8.0 × 10−4), and Heparin binding (GO: 0008201, P = 5.9 × 10−5), represented by relevant cartilaginous genes such as GDF5, TGFBR3, and CILP (Fig. S7A). Finally, the STRING protein–protein interaction network of the n = 124 genes that were highly expressed in hCiCs (> SD2) but lowly expressed in hPAC (> SD2) showed significant gene enrichment in pathways involved in nervous system development (GO:0007399, P = 2.5 × 10−21), represented by highly interconnected genes such as SOX2, SOX11, KIF1A, and STMN4 (Fig. S7B).

Discordant aspects of integrated methylome and transcriptome patterns between 3D neo-cartilage from hCiC relative to hPAC

To identify differences in the methylome that are most likely transcriptionally active, we integrated differential transcriptome-wide RNA sequencing data of hPACs (N = 7) and hCiCs (N = 9) neo-cartilage to the differential methylome-wide data. In total N = 94,771 significant differentially methylated CpGs (DM, P < 0.05 after Bonferroni correction) and N = 7251 differentially expressed genes (DEG, FDR < 0.05) were identified. To identify herein CpG-gene pairs that, most likely, reflect aberrant methylation set points of gene expression in hCiCs, we subsequently selected CpG sites that had a singular mapping to a UCSC reference gene and had a significant high correlation (r > 0.5 and FDR < 0.05) between methylation and gene expression (Fig. S8). As such we prioritized 2378 discordant CpG-gene pairs that were highly interconnected and marked potential discordant set points of gene expression between hCiC and hPACs (Table S2). Among these discordant CpG-gene pairs, the direction of differential methylation at CpG sites and gene expression was positive in N = 722 (30%) and inverse in N = 1656 (70%; Table S2). In the Circos plot in Fig. 4A, the distribution of all discordant CpG-gene pairs across the genome is plotted while highlighting the top 15 most significant ones such as CpGs near COMP encoding cartilage oligomeric protein, MMP7 encoding metalloproteinase 7 and the long noncoding RNA HOTAIR. Moreover, to visualize extend of FDR significant levels of differential methylation (Fig. S9) and expression (Fig. S10) of discordant CpG-gene pairs across the genome, we plotted FDR significant levels as Manhattan plots. As shown in Figs. S9 and S10, multiple highly significant CpG-gene pairs were recognized. Among them, CpGs near OA risk genes [30] such as MGP encoding Matrix-Gla protein or TNC encoding Tenacin, and notable genes marking OA pathophysiology [31, 32] such as TNFRSF11B encoding TNF Receptor Superfamily Member 11b, P4HB encoding prolyl 4-hydroxylase subunit beta, SMOC2 encoding SPARC related modular calcium binding 2, and MMP3 encoding metalloprotease 3 (Table S2). Additionally, we noted among the highly significant differential CpG-gene pairs CTNNB1 encoding beta-catenin. Together, by prioritizing on highly significant differential methylation at likely transcriptionally active CpG, we showed that differences in cell identity between hPACs and hCiCs are characterized by notable genes involved in OA (etio)pathophysiology (Table S2).

Fig. 4
figure 4

Prioritized 2378 discordant CpG-gene pairs marking discordant set points of gene expression between hCiC and hPACs. A Circos plot showing the distribution of N = 2378 differential CpG-gene pairs across the genome. Labeled are N = 15 CpG-gene pairs with highest correlation between methylation and expression levels. The outer circle displays the gene expression. Blue shows downregulated genes, while red shows upregulated genes. The inside circle represents DNA methylation from CpGs that mapped in direct vicinity of DEGs, where blue is hypomethylated and red is hypermethylated. B Visualization of notable significant pathways enriched among 2378 CpG-gene pairs lower and higher expressed in hCiCs relative to hPACs showing the links between genes and biological processes by using GO and KEGG networks. C Notable significant pathway enriched for CpG-gene pairs that were higher expressed in hCiCs relative to hPACs showing the links between genes and biological processes by using GO and KEGG networks. Size of the dots represents number of genes linked to each term. The color of the dots represents the fold change of the differentially expressed gene. P adjusted method = FDR depicted in the histograms

Full size image

Pathway enrichment analyses of 2378 differential set points of gene expression between 3D neo-cartilage from hCiC relative to hPAC

To obtain insight into the pathways in which the CpG-gene pairs act that mark differences in cell identify between hPACs and hCiC, we next performed gene enrichment pathway analyses for genes lower (N = 1237, Table S3) and higher (N = 1141, Table S4) expressed in hCiCs relative to hPACS. Among the significant pathways of the lower expressed genes in hCiC relative to hPACs (Fig. 4B), we recognized closely related GO-terms such as ‘extracellular matrix organization’ (GO:0030198, FDR = 7.7 × 10−18) with diverse cartilage component genes such as ACAN, COMP, DCN, COL2A1, and COL11A1, ‘chondrocyte differentiation’ (GO:0002062, FDR = 3.9 × 10−7) with NPR2, GDF5, TWSG1, and SOX9. The latter gene also represented in the enriched pathway ‘cartilage development’ (GO:0051216, FDR = 3.1 × 10−4) with MGP, CHRDL2, and CHI3L1 as notable other genes. Also, in Table S3 highly significant enriched KEGG pathway ‘Focal adhesion’ (hsa04510, FDR = 1.1 × 10−9), with CHAD and multiple integrin genes such as ITGA5, ITGAV, ITGA9 and ITGA10, many of which were overlapping with the ‘extracellular matrix’ pathway. Together, these lower expressed genes in hCiC relative to hPACs likely mark immaturity and/or lower quality of hCiC neo-cartilage (Fig. 4B, Table S3). On the other hand, some of the lower expressed genes in hCiC relative to hPACs such as TNFRSF11B and COL1A2 concurrent with genes such as CLEC3A, and OMD acting in the ossification pathway (GO:0001503 ‘FDR = 2.1 × 10−4, Table S4) actually reflect an OA phenotypic state of the chondrocytes. As such, these could merely reflect pathologic changes in set points of gene expression of the hPAC derived neo-cartilage relative to that of hCiCs. This, since hPACs were isolated from preserved cartilage of elderly patients.

Among the enriched pathways of the higher expressed genes in hCiC relative to hPACs (Fig. 4C, Table S4), we particularly recognized closely related GO-terms involved in neurogenesis being ‘neuron differentiation’ (GO: GO:0030182, FDR = 6.2 × 10−7) with diverse neurogenic genes such as NEUROG2, SOX11, and WNT genes, the latter genes overlapping with the enriched pathway ‘cell fate commitment’ (GO:0045165, FDR = 8.3 × 10−3). Additionally, we observed a highly significant enrichment within the ‘Axon guidance’ pathway (GO:0007411, FDR = 6.0 × 10−9) with multiple semaphorin genes such as SEMA3D, SEMA4C, SEMA5B, and SEMA6B. Semaphorins represent a family of transmembrane and secreted neuron-guidance molecules. Furthermore, among the genes enriched in the Axon guidance pathway we recognized NGFR, encoding the nerve growth factor receptor, and LHX2 encoding LIM Homeobox 2 that acts as a transcriptional activator of neural cell types. This may reflect off-target during differentiation of hiPSCs toward neurogenic lineage.

Prioritization of differential methylated CpGs that mark stable set points of transcription and that could facilitate strategies to improve hCiC neo-cartilage formation

Next, we wanted to prioritize for most eligible differential set points of the transcriptionally active methylation that could facilitate targeted modifications and improve hiPSC derived neo-cartilage formation. Hereto we prioritized among the 2378 CpG-gene pairs those that had the largest fold change (FC) in gene expression, the largest FC in methylation, and had the largest correlation. Additionally, to select for sensitively modifiable target genes we removed 21 CpG-gene pairs that were discordant, yet were highly expressed in both hPACs and hCiCs (> 1 Standard deviation, Table S5). This resulted in 195 CpG-gene pairs with high differences in set points of the transcriptionally active methylation that could reflect differences in phenotypic chondrocyte states (Fig. 5, Table S6).

Fig. 5
figure 5

STRING protein–protein interaction network of 195 eligible discordant set points of gene expression that could facilitate strategies to improve hCiC neo-cartilage formation. A Significant STRING protein–protein interaction network (P < 1 × 10−16) among downregulated genes in hCiC derived neo-cartilage relative to hPAC (N = 140). B Significant STRING protein–protein interaction network (P = 8.1 × 10−3) among of genes that were significantly upregulated in hCiC derived neo-cartilage relative to hPAC (N = 55). The edges indicate both functional and physical protein associations, with selected genes of Fig. 6 indicated in bold. Line thickness indicates the strength of data support. Minimum required interaction score = medium confidence (0.400)

Full size image

We recognized 140 CpG-gene pairs (72%) with lower and 55 CpG-gene pairs (28%) with higher gene expression in hCiCs relative to hPACs. To visually explore these potential epigenetically modifiable target genes, we generated STRING protein–protein interaction networks. Figure 5A represents lower expressed genes in hCiCs relative to hPACs with highly significant enriched protein–protein interactions (P = 1.0 × 10−16). Notable articular cartilage genes such as PRG4, SERPINA1, ACAN, CHAD, NT5E, and MGP that could mark cartilage immaturity were highly interconnected. Vice versa, Fig. 5B represents higher expressed genes in hCiCs relative to hPACs, also with significant enriched protein–protein interaction (P = 8.1 × 10−3). Here, we recognize genes such as SEMA5B, KIF1A, NKX2-2, DDR1 and GPC2 involved in the generation and differentiation of neurons.

As proof of concept, we finally inspected the genomic region of four of these differential methylated CpGs that mark stable set points of transcription. These genes were selected for being highly FDR significant, central and highly connected in the network (Fig. 5), and being either low in hCiC-derived neo-cartilage relative to hPACs (skeletal development pathway, i.e., CHAD and MGP) or high in hCiC-derived neo-cartilage relative to hPACs (neurogenic pathway, i.e., KIF1A and NKX2-2). As shown in Fig. 6A the genomic regions of CHAD, MGP (lower expressed in hCiCs) and in Fig. 6B the regions of KIF1A, NKX2-2 (higher expressed in hCiCs) were plotted including the location of the CpG sites, methylation-transcription correlation data, previously identified transcription factor binding sites (TFBS), RNA-seq expression levels, and CpG methylation levels. Notable is the overlap between the prioritized likely transcriptional active CpG sites and the mapping of TFBSs which underscore the validity of the applied prioritization scheme toward likely transcriptionally active CpG sites.

Fig. 6
figure 6

Genomic plots of selected genes. Genomic plots of selected CpG-gene pairs eligible as epigenetically modifiable target genes with mapping of genes and CpG sites, Log2FC of methylation, correlation between methylation and expression, mapping of transcription factors binding sites (ChIP-seq ENCODE) as well as differences in expression and methylation between hCiCs and hPACs. A Genomic plots of cg07730609-CHAD and cg20441426-MGP with lower expressed in hCiCs relative to hPACs. B Genomic plots of cg25834415-KIF1A and cg23425348-NKX2-2 with higher expressed in hCiCs relative to hPACs

Full size image

Discussion

To optimize generation of sustainable neo-cartilage constructs, we characterized molecular landscape of hiPSC-derived neo-cartilages in comparison with those deposited by its primary counterpart, the human articular chondrocyte. Moreover, by using hiPSCs generated from two different tissue sources (skin fibroblasts and articular cartilage chondrocytes), we explored whether the epigenetic memory of articular chondrocytes could further facilitate in vitro chondrogenesis. However, the high similarity among neo-cartilages generated from different hiPSC lines indicated that the articular cartilage epigenetic memory retained in hCiC does not further improve the consistency and quality of the in vitro chondrogenesis [33,34,35]. By subsequently prioritizing on likely transcriptionally active methylation discordant between neo-cartilage from hCiC relative to hPAC, we identified relevant differences in phenotypic cell states between chondrocytes derived from hiPSCs and hPACs. Since this molecular level of information is known to be important in on/off target cell fate decisions [12], results of our study could be exploited to improve quality, purity, and maturity of hiPSC-derived neo-cartilage matrix further. Ultimately to realize the introduction of sustainable, hiPSC-derived neo-cartilage implantation into clinical practice.

Among the CpG-genes pairs lower expressed in hCiCs relative to hPACs, we identified genes that were enriched in pathways such as ‘extracellular matrix organization,’ ‘chondrocyte differentiation,’ and ‘cartilage development’ and likely mark immaturity and/or lower quality of hCiC neo-cartilage relative to hPACs. Among the identified discordant genes we recognized many robust OA risk genes that have critical roles in articular cartilage maintenance such as MGP, GDF5, or SERPINA1. Moreover, other genes such as the identified integrins and CHAD are known to interact with structural ECM molecules, as well as, with cells and thus play a role in cartilage homeostasis [36]. CTNNB1, among the highest significant differential CpG-gene pairs, is known to control, in interaction with SOX9, chondrocyte differentiation and could be an important marker of the reduced efficiency for hiPSC chondrogenesis [37]. On a different note, MGP is known to regulate extracellular calcium levels via high affinity to its γ-carboxyglutamic acid (Gla) residues and is a critical marker of the chondrogenic cell lineage [38]. More recently, MGP was recognized among the most robust and targetable OA risk genes by virtue of its vitamin K dependency [39, 40]. In this respect, it is tempting to speculate that enhancing MGP action by supplementation of vitamin K during hiPSC chondrogenic differentiation could be a potential alternative strategy to improve quality of deposited neo-cartilage.

We prioritized N = 195 discordant genes that were highly significant differentially methylated and expressed between hCiC and hPACs, with large effect sizes, and with high correlation between methylation and expression, as eligible targets to improve hiPSC chondrogenic cell fate. The genomic regions plotted for compelling CpG-pairs highlighted that the mapping of these CpG sites indeed coincided with high confidence TFBSs. This adds to the validity that our prioritization scheme, although primarily based on association, has indeed identified differential set points of transcriptionally active methylation hence potential epigenetically modifiable target genes. Nonetheless, a full exploration of the N = 195 prioritized genes is required to implement strategies that improve quality, purity and maturity of hiPSC-derived neo-cartilages. We envision that such an exploration requires high throughput experimental validation by methodologies such as CRISPR-dCAS9 activation and interference (CRISPRi/a) during hiPSC chondrogenic differentiation with single cell read-out and preferably followed by system biological approaches to model interactions. Eventually, methodologies such as dCas9-DNMT/TET methods that allow long-term beneficial changes in set points of methylation could be performed for ultimate validation.

Notable was also the lower expression of genes such as TNFRS11B, P4HA2, and COL1A2 [41] in hCiC relative to hPAC neo-cartilage. These genes are also known to be consistently upregulated with OA pathophysiology [42]. As such these differences in gene expression levels are likely due to the source of hPAC being harvested from aged preserved articular cartilage. On the other hand, it confirms that set points of transcriptionally active methylation are stable throughout the harvesting and in vitro chondrogenesis steps of hPACs [5]. We also identified discordant CpG-gene pairs that were upregulated in hCiC-derived neo-cartilage relative to hPAC. These genes were enriched in pathways closely related to neurogenesis (Table S4). These genes could be a reflection of off-target differentiation toward neurogenic lineage as outlined previously by Wu et al. [43] using hiPSC derived from fibroblasts. Alternatively, these genes in the hiPSC-derived neo-cartilage could be a reflection of the generation of a more fetal or immature cartilage like phenotype.

Although we provided differential gene expression data between hPACs and hCiCs in Table S2, we reported an exploratory transcriptome analyses in the result section. This was done to circumvent the reporting of statistically significant yet subtle differences in expression between hPACs and hCiCs at the mature neo-cartilage stage. Such transcriptome-wide differences per definition reflect the phenotypic cell state in mature hiPSC-derived neo-cartilages and are actually not driving cell fate decisions, the primary focus of our study.

To identify differential methylation at CpG sites that likely mark aberrant set points of gene expression, a key layer of molecular information that does facilitate critical lineage decisions hence cell fates (11), we integrated transcriptome data to the methylome data solely to extract biological relevant information to the differential methylated sites. Herein, we took the assumption that differences in methylation that affect expression of genes should follow specific rules of correlation between these two molecular levels of information.

Upon studying the discordant likely transcriptionally active methylation between neo-cartilage, we disregarded the hFiC datasets. This to assure that the comparative integrated methylome and transcriptome differential expression analyses were performed in powerful datasets generated from homogeneous cell populations. Additionally, the fact that our methylome-wide data of neo-cartilages showed comparable high similarities within and between cell sources, we anticipated that the relative small hFiC dataset has likely no effect on the extent, nor content, of the conclusions of our study. Hence, we are confident that the N = 195 targets are robust across different cell sources used to generate hiPSCs.

Limitation of the current study is that we have used 3 hiPSC clones of chondrocytes of only one donor and 1 hiPSC clone of skin fibroblasts of 2 donors, that could have influenced the robustness of our results. Nonetheless, as reported, we have generated multiple differentiations of available clones and showed high methylome-wide similarities within and among clones and cell sources confirming consistency of our differentiation protocol and validity of our conclusions. Moreover, our previous work supports the consistency of the chondrogenesis protocol across many different additional hiPSC lines [10, 20].

Conclusion

By applying an integrative multi-omics approach, we demonstrated high methylome-wide similarity of neo-cartilages between and within cell sources underscoring the consistency and quality of the applied step-wise differentiation protocol [10, 11]. Nonetheless, to realize hiPSC-derived regenerative treatments into clinical practice, in-depth insight into the purity and quality of neo-cartilage is required. Hereto, we also set out to identify discordant aspects of stable set points of gene expression between neo-cartilage from hCiC relative to hPAC. Altogether, our findings provide important insights into discordant aspects of stable set points of gene expression between neo-cartilage derived from hiPSCs relative to that of its primary counterpart, the human articular chondrocyte. These insights could be exploited to improve quality, purity, and maturity of hiPSC-derived neo-cartilage matrix further, ultimately to realize the introduction of sustainable, hiPSC-derived neo-cartilage implantation into clinical practice.

Data availability

The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

Abbreviations

hCiC:

Human cartilage-derived iPSCs differentiated to chondrocytes

hiPSC:

Human-induced pluripotent stem cells

hFiC:

Human fibroblast-derived iPSCs differentiated to chondrocytes

hPAC:

Human primary articular chondrocytes

OA:

Osteoarthritis

References

  1. Goldring MB, Marcu KB. Cartilage homeostasis in health and rheumatic diseases. Arthritis Res Ther. 2009;11(3):224.

    Article  PubMed  PubMed Central  Google Scholar 

  2. Woolf AD, Erwin J, March L. The need to address the burden of musculoskeletal conditions. Best Pract Res Clin Rheumatol. 2012;26(2):183–224.

    Article  PubMed  Google Scholar 

  3. Tuan RS, Chen AF, Klatt BA. Cartilage regeneration. J Am Acad Orthop Surg. 2013;21(5):303–11.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Kimura T, Yamashita A, Ozono K, Tsumaki N. Limited immunogenicity of human induced pluripotent stem cell-derived cartilages. Tissue Eng Part A. 2016;22(23–24):1367–75.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Bomer N, Den HW, Suchiman H, Houtman E, Slieker RC, Heijmans BT, et al. Neo-cartilage engineered from primary chondrocytes is epigenetically similar to autologous cartilage, in contrast to using mesenchymal stem cells. Osteoarth Cartil. 2016;24(8):1423–30.

    Article  CAS  Google Scholar 

  6. Kamaraj A, Kyriacou H, Seah KTM, Khan WS. Use of human induced pluripotent stem cells for cartilage regeneration in vitro and within chondral defect models of knee joint cartilage in vivo: a preferred reporting items for systematic reviews and meta-analyses systematic literature review. Cytotherapy. 2021;23(8):647–61.

    Article  PubMed  Google Scholar 

  7. Kretlow JD, Jin YQ, Liu W, Zhang WJ, Hong TH, Zhou G, et al. Donor age and cell passage affects differentiation potential of murine bone marrow-derived stem cells. BMC Cell Biol. 2008;9:60.

    Article  PubMed  PubMed Central  Google Scholar 

  8. Qu C, Puttonen KA, Lindeberg H, Ruponen M, Hovatta O, Koistinaho J, et al. Chondrogenic differentiation of human pluripotent stem cells in chondrocyte co-culture. Int J Biochem Cell Biol. 2013;45(8):1802–12.

    Article  CAS  PubMed  Google Scholar 

  9. Guha P, Morgan JW, Mostoslavsky G, Rodrigues NP, Boyd AS. Lack of immune response to differentiated cells derived from syngeneic induced pluripotent stem cells. Cell Stem Cell. 2013;12(4):407–12.

    Article  CAS  PubMed  Google Scholar 

  10. Rodríguez Ruiz A, Dicks A, Tuerlings M, Schepers K, van Pel M, Nelissen R, et al. Cartilage from human-induced pluripotent stem cells: comparison with neo-cartilage from chondrocytes and bone marrow mesenchymal stromal cells. Cell Tissue Res. 2021;386:309–20.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Adkar SS, Wu C-L, Willard VP, Dicks A, Ettyreddy A, Steward N, et al. Step-wise chondrogenesis of human induced pluripotent stem cells and purification via a reporter allele generated by CRISPR-Cas9 genome editing. Stem Cells. 2019;37(1):65–76.

    Article  CAS  PubMed  Google Scholar 

  12. Parry A, Rulands S, Reik W. Active turnover of DNA methylation during cell fate decisions. Nat Rev Genet. 2021;22(1):59–66.

    Article  CAS  PubMed  Google Scholar 

  13. Ziller MJ, Gu H, Müller F, Donaghey J, Tsai LT, Kohlbacher O, et al. Charting a dynamic DNA methylation landscape of the human genome. Nature. 2013;500(7463):477–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Doi A, Park IH, Wen B, Murakami P, Aryee MJ, Irizarry R, et al. Differential methylation of tissue- and cancer-specific CpG island shores distinguishes human induced pluripotent stem cells, embryonic stem cells and fibroblasts. Nat Genet. 2009;41(12):1350–3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Kim K, Doi A, Wen B, Ng K, Zhao R, Cahan P, et al. Epigenetic memory in induced pluripotent stem cells. Nature. 2010;467(7313):285–90.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Lister R, Pelizzola M, Kida YS, Hawkins RD, Nery JR, Hon G, et al. Hotspots of aberrant epigenomic reprogramming in human induced pluripotent stem cells. Nature. 2011;471(7336):68–73.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Ohi Y, Qin H, Hong C, Blouin L, Polo JM, Guo T, et al. Incomplete DNA methylation underlies a transcriptional memory of somatic cells in human iPS cells. Nat Cell Biol. 2011;13(5):541–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Dambrot C, van de Pas S, van Zijl L, Brändl B, Wang JW, Schalij MJ, et al. Polycistronic lentivirus induced pluripotent stem cells from skin biopsies after long term storage, blood outgrowth endothelial cells and cells from milk teeth. Differentiation. 2013;85(3):101–9.

    Article  CAS  PubMed  Google Scholar 

  19. Tuerlings M, van Hoolwerff M, Houtman E, Suchiman E, Lakenberg N, Mei H, et al. RNA sequencing reveals interacting key determinants of osteoarthritis acting in subchondral bone and articular cartilage: identification of IL11 and CHADL as attractive treatment targets. Arthritis Rheumatol. 2021;73(5):789–99.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. van Hoolwerff M, Rodriguez Ruiz A, Bouma M, Suchiman HED, Koning RI, Jost CR, et al. High-impact FN1 mutation decreases chondrogenic potential and affects cartilage deposition via decreased binding to collagen type II. Sci Adv. 2021;7(45):eabg8583.

    Article  PubMed  PubMed Central  Google Scholar 

  21. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29(1):15–21.

    Article  CAS  PubMed  Google Scholar 

  22. Anders S, Pyl PT, Huber W. HTSeq–a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31(2):166–9.

    Article  CAS  PubMed  Google Scholar 

  23. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.

    Article  PubMed  PubMed Central  Google Scholar 

  24. van Iterson M, Tobi EW, Slieker RC, den Hollander W, Luijk R, Slagboom PE, et al. MethylAid: visual and interactive quality control of large Illumina 450k datasets. Bioinformatics. 2014;30(23):3435–7.

    Article  PubMed  Google Scholar 

  25. Aryee MJ, Jaffe AE, Corrada-Bravo H, Ladd-Acosta C, Feinberg AP, Hansen KD, et al. Minfi: a flexible and comprehensive Bioconductor package for the analysis of Infinium DNA methylation microarrays. Bioinformatics. 2014;30(10):1363–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Chen YA, Lemire M, Choufani S, Butcher DT, Grafodatskaya D, Zanke BW, et al. Discovery of cross-reactive probes and polymorphic CpGs in the Illumina Infinium HumanMethylation450 microarray. Epigenetics. 2013;8(2):203–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Martin-Trujillo A, Patel N, Richter F, Jadhav B, Garg P, Morton SU, et al. Rare genetic variation at transcription factor binding sites modulates local DNA methylation profiles. PLoS Genet. 2020;16(11):e1009189.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Fortin JP, Triche TJ Jr, Hansen KD. Preprocessing, normalization and integration of the Illumina HumanMethylationEPIC array with minfi. Bioinformatics. 2017;33(4):558–60.

    Article  CAS  PubMed  Google Scholar 

  29. Zhang D, Xue J, Peng F. The regulatory activities of MALAT1 in the development of bone and cartilage diseases. Front Endocrinol (Lausanne). 2022;13:1054827.

    Article  PubMed  Google Scholar 

  30. Boer CG, Hatzikotoulas K, Southam L, Stefánsdóttir L, Zhang Y, Coutinho de Almeida R, et al. Deciphering osteoarthritis genetics across 826,690 individuals from 9 populations. Cell. 2021;184(18):4784–818.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Coutinho de Almeida R, Mahfouz A, Mei H, Houtman E, den Hollander W, Soul J, et al. Identification and characterization of two consistent osteoarthritis subtypes by transcriptome and clinical data integration. Rheumatology. 2020.

  32. Coutinho de Almeida R, Ramos YFM, Mahfouz A, Den Hollander W, Lakenberg N, Houtman E, et al. RNA sequencing data integration reveals an miRNA interactome of osteoarthritis cartilage. Ann Rheum Dis. 2019;78(2):270–7.

    Article  PubMed  Google Scholar 

  33. Kilpinen H, Goncalves A, Leha A, Afzal V, Alasoo K, Ashford S, et al. Common genetic variation drives molecular heterogeneity in human iPSCs. Nature. 2017;546(7658):370–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Kim K, Zhao R, Doi A, Ng K, Unternaehrer J, Cahan P, et al. Donor cell type can influence the epigenome and differentiation potential of human induced pluripotent stem cells. Nat Biotechnol. 2011;29(12):1117–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Nishizawa M, Chonabayashi K, Nomura M, Tanaka A, Nakamura M, Inagaki A, et al. Epigenetic variation between human induced pluripotent stem cell lines is an indicator of differentiation capacity. Cell Stem Cell. 2016;19(3):341–54.

    Article  CAS  PubMed  Google Scholar 

  36. Hessle L, Stordalen GA, Wenglén C, Petzold C, Tanner E, Brorson SH, et al. The skeletal phenotype of chondroadherin deficient mice. PLoS ONE. 2014;8(6):e63080.

    Article  PubMed  Google Scholar 

  37. Akiyama H, Lyons JP, Mori-Akiyama Y, Yang X, Zhang R, Zhang Z, et al. Interactions between Sox9 and beta-catenin control chondrocyte differentiation. Genes Dev. 2004;18(9):1072–87.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Luo G, D’Souza R, Hogue D, Karsenty G. The matrix Gla protein gene is a marker of the chondrogenesis cell lineage during mouse development. J Bone Miner Res. 1995;10(2):325–34.

    Article  CAS  PubMed  Google Scholar 

  39. Houtman E, de Almeida RC, Tuerlings M, Suchiman E, Broekhuis D, Nelissen RGHH, et al. Characterization of dynamic changes in Matrix Gla Protein (MGP) gene expression as function of genetic risk alleles, osteoarthritis relevant stimuli, and the vitamin K inhibitor warfarin. Osteoarthr Cartil. 2021;29:1193–202.

    Article  CAS  Google Scholar 

  40. Boer CG, Szilagyi I, Nguyen NL, Neogi T, Meulenbelt I, Ikram MA, et al. Vitamin K antagonist anticoagulant usage is associated with increased incidence and progression of osteoarthritis. Ann Rheum Dis. 2021;80(5):598–604.

    Article  CAS  PubMed  Google Scholar 

  41. Xu H, Ding C, Guo C, Xiang S, Wang Y, Luo B, et al. Suppression of CRLF1 promotes the chondrogenic differentiation of bone marrow-derived mesenchymal stem and protects cartilage tissue from damage in osteoarthritis via activation of miR-320. Mol Med. 2021;27(1):116.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Rodríguez Ruiz A, Tuerlings M, Das A, Coutinho de Almeida R, Suchiman HED, Nelissen R, et al. The role of TNFRSF11B in development of osteoarthritic cartilage. Rheumatology (Oxford). 2022;61(2):856–64.

    Article  PubMed  Google Scholar 

  43. Wu CL, Dicks A, Steward N, Tang R, Katz DB, Choi YR, et al. Single cell transcriptomic analysis of human pluripotent stem cell chondrogenesis. Nat Commun. 2021;12(1):362.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank Prof. Judith Bovee and all members of her team in the Dept. Pathology, Leiden University Medical Center, Leiden, The Netherlands. No personal data are included in this study. All authors gave consent for publication.

Funding

This study is partially supported by the partners of Regenerative Medicine Crossing Borders (RegMed XB), a public–private partnership that uses regenerative medicine strategies to cure common chronic diseases and by the Dutch Scientific Research council NWO/ZonMW VICI scheme (91816631/528). The funding bodies played no role in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.

Author information

Author notes

  1. Ghazaleh Hajmousa, Rodrigo Coutinho de Almeida and Niek Bloks have contributed equally to this work.

Authors and Affiliations

  1. Department of Biomedical Data Sciences, Section Molecular Epidemiology, Leiden University Medical Center, Post-zone S-05-P, P.O. Box 9600, 2300 RC, Leiden, The Netherlands

    Ghazaleh Hajmousa, Rodrigo Coutinho de Almeida, Niek Bloks, Alejandro Rodríguez Ruiz, Yolande F. M. Ramos & Ingrid Meulenbelt

  2. Department of Anatomy and Embryology and Human iPSC Hotel, 2333 ZA, Leiden, The Netherlands

    Marga Bouma & Christian Freund

  3. Department of Cell and Chemical Biology, Leiden University Medical Center, Leiden, The Netherlands

    Roderick Slieker & Thomas B. Kuipers

  4. Department of Orthopedics, Leiden University Medical Center, Leiden, The Netherlands

    Rob G. H. H. Nelissen

  5. Orthopaedic Biomechanics, Department of Biomedical Engineering, Eindhoven University of Technology, P.O. Box 513, 5600 MB, Eindhoven, The Netherlands

    Keita Ito

Authors
  1. Ghazaleh Hajmousa
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Rodrigo Coutinho de Almeida
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Niek Bloks
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Alejandro Rodríguez Ruiz
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Marga Bouma
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Roderick Slieker
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Thomas B. Kuipers
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Rob G. H. H. Nelissen
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. Keita Ito
    View author publications

    You can also search for this author in PubMed Google Scholar

  10. Christian Freund
    View author publications

    You can also search for this author in PubMed Google Scholar

  11. Yolande F. M. Ramos
    View author publications

    You can also search for this author in PubMed Google Scholar

  12. Ingrid Meulenbelt
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Conception and design were done by CF, YR, and IM. Collection and/or assembly of data was done by GH, ARR, MB, and RN. Data analysis and interpretation were done by GH, RCA, NB, RS, TK, KI, YR, and IM. Writing of the manuscript were done by GH, RCA, NB, YR, and IM. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Ingrid Meulenbelt.

Ethics declarations

Ethics approval and consent to participate

Ethical permission for the described studies was obtained from the medical ethical committee (METC) of the Leiden University Medical Center: 1) “Research Artrotisch Articulair Kraakbeen (RAAK)” (P08.239 and P19.013, approved, respectively, d.d. December 2, 2008 and May 20, 2019) 2) “Parapluprotocol: hiPSCs” (P13.080 approved d.d. July 2, 2014). Written informed consent was obtained from all participants.

Competing interests

The authors declare no competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Yolande F. M. Ramos and Ingrid Meulenbelt Shared last authorship.

Supplementary Information

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Hajmousa, G., de Almeida, R.C., Bloks, N. et al. The role of DNA methylation in chondrogenesis of human iPSCs as a stable marker of cartilage quality. Clin Epigenet 16, 141 (2024). https://doi.org/10.1186/s13148-024-01759-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13148-024-01759-y

Keywords

Clinical Epigenetics

ISSN: 1868-7083

Contact us