Expression of Runx2 and Runx3 in articular cartilage during development of surgically induced OA
To examine Runx2, Runx3, and Sox9 expression in articular cartilage during OA development, we prepared an 8-week-old mouse OA model by surgical resection of the medial collateral ligament and medial meniscus (medial model).39 Immunohistochemistry of normal cartilage showed that Runx2 protein was predominantly located in the DZ of preoperative articular cartilage, while Runx3 and Sox9 was found in all layers (Fig. 1a). After surgery, Runx3 and Sox9 started to decrease in surgically treated knee joints, accompanied by Col2a1. Runx2 was continuously detected in the DZ. Immunohistochemistry revealed enhanced expression of Mmp13, especially 2 weeks after surgery and predominantly in the DZ of articular cartilage at 4, 6, and 8 weeks after surgery (Fig. 1a). We additionally examined changes in Runx2, Runx3, and Sox9 expression in articular cartilage with aging. In the knee joints of WT mice, their expression decreased over time and was markedly weakened at 12 months of age (Fig. 1b).
We performed in vitro analyses to examine roles of Runx2 and Runx3 under inflammation using IL-1β, which is widely employed to model chondrocyte degradation.14,40-42 Primary chondrocytes43 were exposed to IL-1β at different concentrations for 24 h (Fig 1c); or 1 ng/mL IL-1β for 0−48 h (Fig 1d). Runx3 was downregulated by IL-1β exposure in dose- (Fig 1c) and time-dependent manners (Fig 1d), but levels of Runx2 mRNA were unchanged following IL-1β exposure (Fig 1c, d). Conversely, exposure to IL-1β elicited dose- and time-dependent upregulation of Mmp13 and suppression of Sox9, followed by a decrease of Col2a1 (Fig 1c, d). These results indicate that exposure of in vitro chondrocyte to 1 ng/mL IL-1β for 24 h reproduces similar expression patterns of marker genes and proteins observed in vivo surgical OA models, including stable expression of Runx2, reduced Runx3, Sox9 and Col2a1, and increased Mmp13.
Runx3 knockout accelerated OA development through Prg4 and Acan suppression
To reveal roles of Runx3 in whole or SFZ articular cartilage after skeletal growth, we compared Runx3fl/fl littermates44 (R3-Cntl) with Col2a1-CreERT2;Runx3fl/fl (R3-cKOER)45 and Prg4-CreERT2;Runx3fl/fl (R3-pKOER)46 mice. We prepared three distinct OA models (Supplementary Fig. 1). As a result of the medial model,39 OA development was significantly accelerated in R3-cKOER (Fig. 2a, b) and R3-pKOER (Fig. 2c, d) knee joints compared with R3-Cntl joints. In the articular cartilage of sham-operated R3-cKOER joints, Runx3-positive cells decreased by 70% (Fig 2b). Positive areas of Acan and Prg4 decreased with Runx3 knockout (Fig. 2b). In the articular cartilage of sham-operated R3-pKOER joints, Runx3-positive cells decreased by 40% and the positive area of Prg4 was decreased (Fig. 2d). Additionally, eight weeks after destabilizing the medial meniscus (DMM)47 as a minor injury-induced OA model, OA development was also significantly accelerated in 24-week-old R3-pKOER mice compared with R3-Cntl mice (Supplementary Fig. 2). In the aging model (Supplementary Fig. 1c), OA development was unchanged in R3-cKOER mice compared with R3-Cntl mice (Fig. 2e, f). In the articular cartilage of 18-month-old joints, the rates of Runx3-positive cells were not different for either genotype (Fig. 2f), probably because Runx3 expression decreased with aging as shown in Fig.1b. Expression of the marker proteins in both mice was similar as well (Fig. 2f). By Runx3 knockout in the whole articular cartilage, chondrocyte apoptosis was not significantly changed in the sham side, the surgical OA model, or with aging, as well as by Runx3 knockout in the SFZ articular cartilage of the sham side (Supplementary Fig. 3a-d).
We then examined expression of OA associated genes in SFZ and DZ cells from Runx3 knockout and control mice. To enhance knockout efficiency, we prepared Runx3fl/fl (R3-Cntl) and Col2a1-Cre;Runx3fl/fl (R3-cKO) mice, which displayed normal skeletal development (Supplementary Fig. 4a-g). Prg4, was decreased in R3-cKO SFZ cells and Acan was also decreased in R3-cKO DZ cells (Fig. 2g). Other genes were not changed in SFZ cells or DZ cells (Fig. 2g). We then prepared SFZ and DZ cells from R3-pKOER and R3-Cntl littermates, and treated them with water-soluble tamoxifen. Only Prg4 was decreased in Runx3-deficient SFZ cells (Supplementary Fig. 5). We further performed Runx3 overexpression in Runx3 knockout cells using adenoviral transduction. Downregulation of Prg4 in R3-cKO SFZ cells and Acan in R3-cKO DZ cells was recovered to control levels by Runx3 overexpression (Fig. 2h, i). Additionally, we performed dimethylmethylene blue (DMMB) assays48 to investigate the catabolic effect of Runx3. Runx3 knockout did not significantly change the released glycosaminoglycan (GAG) content (Fig. 2j).
Development of surgically induced OA was inhibited by heterozygous Runx2 knockout, but accelerated by homozygous knockout
To investigate Runx2 loss-of-function in articular cartilage after skeletal growth in vivo, we generated a mutant mouse strain expressing a floxed Runx2 allele (Supplementary Fig. 6a). Alizarin red and Alcian blue staining of skeletal preparations was performed in E18.5 embryos. CAG-Cre;Runx2fl/lf mice exhibited a complete absence of mineralized matrix in developing skeletal bones, similar to global knockout mice (Supplementary Fig. 6b).25,26 Hypoplastic clavicles and open fontanelles were observed in CAG-Cre;Runx2fl/+ mice, which indicates a cleidocranial dysplasia phenotype, as reported in general heterozygous knockout mice (Supplementary Fig. 6b).19
We next mated Runx2fl/fl mice with Col2a1-CreERT2 mice45 and generated Runx2fl/+(R2-Hetero Cntl), Col2a1-CreERT2;Runx2fl/+(R2-Hetero cKO), Runx2fl/fl (R2-Homo Cntl), and Col2a1-CreERT2;Runx2fl/fl (R2-Homo cKO) male littermates. In the medial model,39 OA development was significantly inhibited in knee joints of R2-Hetero cKO mice and accelerated in R2-Homo cKO mice compared with their respective controls (Fig. 3a). In the DMM model,47 OA development was significantly inhibited in R2-Hetero cKO knee joints, similar to previous reports,31,32 but accelerated in R2-Homo cKO mice compared with R2-Cntl joints (Fig. 3a).
Immunohistochemistry showed that Runx2-positive cells were reduced by approximately 40% in R2-Hetero cKO cartilage and 80% in R2-Homo cKO cartilage of sham and OA joints (Fig. 3b, c), respectively. Protein levels of Sox9 and Col2a1 were unchanged in the sham joints of both genotypes (Fig. 3b), but Col2a1 was significantly suppressed in OA joints of R2-Homo cKO mice (Fig. 3c). Mmp13 was significantly reduced both in sham and OA joints of R2-Hetero and R2-Homo cKO mice (Fig. 3b, c). In the aging model, we further compared OA development between R2-Homo Cntl and R2-Homo cKO mice. At 18 months of age, OA development with aging was unchanged in R2-Homo cKO mice compared with R2-Cntl mice, as well as immunohistochemistry of Runx2, Sox9, Col2a1, or Mmp13 (Supplementary Fig. 7a, b).
To examine chondrocyte apoptosis, we performed TdT-mediated dUTP nick end labeling (TUNEL) staining of medial and DMM model knee joints. Percentages of TUNEL-positive cells were unchanged by Runx2 deletion in both models (Supplementary Fig. 8a). Conversely, TUNEL staining revealed upregulated chondrocyte apoptosis in the sham joints of R2-Homo KO mice (Supplementary Fig. 8b). Cell numbers in articular cartilage were decreased in R2-Homo cKO mice (Supplementary Fig. 8b). To examine short-term effects of Runx2 knockout, we injected tamoxifen into 7-week-old male littermate mice of the four genotypes daily for five days, and sacrificed 1 week after injection without any surgery. Chondrocyte apoptosis was enhanced in the articular cartilage of R2-Homo cKO mice with decreased cell numbers (Supplementary Fig. 8c).
Col2a1 expression was suppressed in homozygous Runx2-knockout chondrocytes under inflammatory conditions
Surgical induction differently altered OA development in R2-Hetero and Homo cKO mice, but had no effect on aging-model mice, indicating that Runx2 may be involved in cartilage degeneration by severe inflammation rather than aging with mild inflammation. To reveal roles of Runx2 in chondrocytes under inflammation, we performed ex vivo and in vitro analyses. We employed IL-1β because it is one of the essential cytokines responsible for human OA, increased in surgical OA models, and has been widely used in ex vivo and in vitro experiments for OA.49-51 Safranin-O staining of femoral heads cultured for 2 weeks with or without IL-1β displayed consistent decreases in proteoglycans for each genotype (Fig. 4a). Immunohistochemistry confirmed efficient knockdown of Runx2 in R2-Hetero and R2-Homo cKO femoral heads (Fig. 4a). Col2a1 protein levels were unchanged in normal cultures for both genotypes, but decreased following IL-1β exposure, especially in the R2-Homo cKO group (Fig. 4a). qRT-PCR using mRNA from the homogenized femoral heads confirmed that Runx2 was decreased according to genotype, but unchanged by IL-1β exposure (Fig. 4b). Mmp13 increased in all genotypes following IL-1β exposure, and the expression pattern of Mmp13 between genotypes was similar to that of Runx2 (Fig. 4b). Sox9 and Col2a1was decreased by IL-1β and unchanged between genotypes with a marked decrease in Col2a1 in R2-Homo cKO femoral heads following IL-1β exposure (Fig. 4b).
We next examined temporal changes in mRNA levels of marker genes in chondrocytes from the four genotypes exposed to 1 ng/mL IL-1β from 0−48 h. mRNA levels of Runx2 were efficiently decreased in R2-Hetero and R2-Homo cKO chondrocytes, and stable in each genotype following IL-1β exposure (Fig. 4c). Mmp13 increased with time and was suppressed in R2-Hetero and R2-Homo cKO chondrocytes compared with their respective controls at each time point. Sox9 gradually decreased, and there was no difference between genotypes at any point examined. Col2a1 was significantly suppressed in R2-Homo cKO more than 24 h after IL-1β exposure (Fig. 4c).
Genome-wide analysis of Runx2 and Runx3 association profiles in chondrocytes determined by ChIP-seq and RNA-seq
To investigate mechanisms underlying chondrocyte regulation by Runx3 and Runx2, we planned ChIP-seq with an anti-FLAG antibody. For Runx3, we prepared SFZ cells from WT mice transfected with the FLAG-tagged Runx3 expression vector. For the ChIP-seq, 21,730 raw peaks met the peak calling criterion and approximately 50% of all peaks mapped to an interval between ±50 and 500 kb from the transcriptional start sites (TSS) (Fig. 5a). Genomic Regions Enrichment of Annotations Tool (GREAT) Gene Ontology (GO) analysis52 indicated that the extracellular structure organization and collagen fibril organization terms were the most significantly enriched in the gene set (Fig. 5a). De novo motif analysis of the top 1,000 specific peaks using MEME-ChIP53 identified a previously predicted Runx motif, TG(T/C)GG(T/C) (Fig. 5b). ChIP-seq data around genes encoding Prg4 and Acan showed several peaks for Runx3 binding (Supplementary Fig. 9). Relative luciferase activities of the regions around Prg4 and Acan were increased by Runx3 overexpression (Supplementary Fig. 9). To further investigate the alterations of gene expression profiles by Runx3 knockout, we performed RNA-seq using SFZ and DZ cells from R3-cKO and R3-Cntl mice (Supplementary Fig. 10, Supplementary Table 1−4). Among 18 genes up- or downregulated by more than 2-fold in the R3-cKO SFZ cells, nine were the top one-third peak nearest genes as shown in the ChIP-seq, including Prg4 and Mmp9 (Fig. 5c).
For Runx2, we performed ChIP-seq in accordance with previous reports,12,54 using chondrocytes treated with vehicle control (primary chondrocytes) or exposed to 1 ng/mL IL-1β (inflamed chondrocytes) derived from Runx2-FLAG mice (Hojo H., et al. in prep.). For the ChIP-seq, 37,163 raw peaks in primary chondrocytes and 14,571 raw peaks in inflamed chondrocytes met the peak calling criteria (Fig. 5d). Peak distributions between primary chondrocytes and inflamed chondrocytes were similar. In both groups, a striking enrichment was observed around the TSS; approximately 24% of all peaks from two ChIP-seq data are within ± 500 bp of the TSS, even though this region represents only 0.001% of the genome52 (Fig. 5d). GREAT GO analysis52 identified “collagen fibril organization” as the fifth and second most significantly enriched term in the gene sets for primary chondrocytes and inflamed chondrocytes, respectively (Fig. 5d). We hypothesized that there were functional differences between the TSS-associated dataset and dataset excluding data ± 500 bp from the TSS, similar to Sox9.12 GREAT GO analyses showed cell cycle associated terms in TSS-associated Runx2 datasets of both primary and inflamed chondrocytes and terms for a Runx2-regulated skeletal program in the dataset excluding data ± 500 bp from TSS (Supplementary Fig. 11). A Venn diagram of GREAT GO analyses indicated that Runx2 was associated with collagen fibril organization in chondrocytes with or without inflammation (Fig. 5d).
We next performed de novo motif analysis of the top 1,000 specific peaks using MEME-ChIP.53 As expected, primary motifs were similar to the previously predicted Runx motif, TG(T/C)GGT, even in the dataset excluding data ± 500 bp from the TSS (Fig. 5e). Interestingly, we observed that poly-A sequence, reported as a SOX9 consensus sequences,55 was enriched in the TSS-associated Runx2 dataset. To identify the functional relationship between Sox9 and Runx2, we further analyzed our present dataset combined with previous Sox9 and related ChIP-seq datasets.12 The TG(T/C)GGT motif was highly enriched at the predicted center of Runx2 ChIP-seq peaks in primary and inflamed chondrocytes. Interestingly, the TG(T/C)GGT motif was slightly enriched at the predicted center of Sox9 ChIP-seq peaks (Fig. 5f). A poly-A sequence was also recovered, but there was no centering within the Runx2 and Sox9 peaks (Fig. 5f). Thus, the poly-A motif is unlikely to be the preferred primary site for Runx2 or Sox9 engagement in chondrocytes.
To test the integration of multiple regulatory inputs through Runx2- and Sox9-directed enhancer modules, we analyzed peak intensity associations. There were strong associations between Sox9-12 and Runx2-peak regions, which were notably enhanced in inflamed chondrocytes (Fig. 5g). Clear associations of Runx2-peak regions with enhancer signatures were evident, specifically: (1) bi-modal patterns of H3K4 dimethylation (H3K4me2) peaks at Runx2-peak center regions, which indicates both promoters and putative enhancers;56 (2) peaks of H3K27 acetylation (H3K27Ac) flanking Runx2 peaks, which indicate open chromatin,57,58 were more evident in inflamed chondrocytes compared with primary chondrocytes; and (3) RNA polymerase II peaks associated with Runx2-peak regions, consistent with active enhancers59 (Fig. 5d). Together, these results indicate that Runx2 enhances transcription of collagen fibrils by binding TG(T/C)GGT target sequences near Sox9 binding site, both in primary chondrocytes and inflamed chondrocytes, rather than interacting with a poly-A box within target chondrocyte enhancers.
To further investigate the general functions of Runx2 in chondrocytes, we performed RNA-seq analyses of R2-Homo Cntl and cKO chondrocytes with or without IL-1β exposure (Supplementary Fig. 12 and Supplementary Table 5 and 6). Principal component analysis highlighted the effect of Runx2 knockout under inflammation induced by IL-1β. As chondrocyte apoptosis was enhanced in R2-Homo cKO mice (Supplementary Fig. 8b, c), we examined apoptosis-related genes. Among the top 200 genes downregulated by Runx2 knockout in primary chondrocytes, eleven genes (including Ihh, Igf1, and Wnt11) were associated with “negative regulation of apoptosis” (Fig. 5h). Among the top 200 genes downregulated by Runx2 knockout in inflamed chondrocytes, ten (including Col2a1) matched with the top 500 genes nearest to Runx2 peaks detected by ChIP-seq (Fig. 5i).
Transcriptional regulation of Col2a1 and Mmp13 by Runx2 and Sox9
We mapped Runx2-FLAG, Sox9-FLAG,60 and active histone marks (H3K4me2 or H3K27Ac in chondrocytes12 or osteoblasts54 from 1-day-old mice) in ChIP-seq data using CisGenome browser. Sox9 and Runx2 sequence data showed similar patterns in introns 1 and 6 of Col2a1; specifically, minor peaks at intron 1 and major peaks at intron 6 (Fig. 6a). These peaks are identical to the functional enhancers of Col2a1 containing Sox9 motifs.8 Magnified views indicated that the peak centers of Runx2 and Sox9, containing their consensus motifs, were 300−400 bp apart (Supplementary Fig. 13). Based on the CisGenome browser, we prepared luciferase reporter vectors containing the Runx2 motifs around Col2a1 and Mmp13 genes shown in Fig. 6a. Relative luciferase activity of Col2a1 intron 1 and 6 fragments decreased in chondrocytes exposed to ≥ 1 ng/mL IL-1β (Fig. 6b), which exhibited decreased Sox9 mRNA while that of Runx2 remained stable (Fig. 1c, d). However, luciferase activities of Col2a1 and of Mmp13 TSS upstream regions were unchanged compared with controls (Fig. 6b). We next examined the effect of Sox9 on activity of each reporter. Sox9 overexpression increased activities of reporters containing the enhancers in Col2a1 introns 1 and 6 in a dose-dependent manner (Fig. 6c), but did not affect activity of the region upstream of the Col2a1 TSS and decreased activity of the region upstream of the Mmp13 TSS (Fig. 6c). Additionally, we co-transfected Runx2 and Sox9 (Fig. 6d). Notably, when the amount of Sox9 was decreased, luciferase activities of the Col2a1 intron 6 enhancer were increased by Runx2 co-transfection (Fig. 6d). In contrast, when Sox9 was overexpressed, luciferase activities of the Col2a1 intron 1 enhancer were not affected by Runx2 transfection (Fig. 6d). Activity of the region upstream of the Col2a1 TSS was increased by Runx2 transfection without Sox9 transfection (Fig. 6d).
We further performed luciferase assays using chondrocytes from the four genotypes with or without IL-1β exposure. Relative luciferase activities of Col2a1 intron 1 and 6 enhancers were decreased in inflamed chondrocytes compared with primary chondrocytes for each genotype (Fig. 6e), similar to the results shown in Fig. 6b. Additionally, luciferase activities of Col2a1 intron 1 and 6 enhancers, and the region upstream of the Col2a1 TSS, were suppressed in inflamed chondrocytes of R2-Homo cKO mice compared with other genotypes (Fig. 6e). Compared with controls, Runx2 heterozygous or homozygous knockout decreased the activities of reporters containing the region upstream of the Mmp13 TSS in both primary and inflamed chondrocytes (Fig. 6d). Taken together, the present results indicate that Runx2 could activate Col2a1 transcription through the enhancer in intron 6, in addition to the enhancer in intron 1 and region upstream of the Col2a1 TSS.
Enhanced expression of Runx3 inhibited surgically induced OA development
Finally, we investigated effects of adenoviral Runx2 or Runx3 overexpression in vitro and in vivo. Transduced of a Runx2 adenovirus (Ad-Runx2) at different doses to WT chondrocytes increased Runx2 expression in a dose-dependent manner (Fig. 7a). Although mRNA levels of Sox9 and Col2a1 were unchanged, Mmp13 was upregulated in Runx2 overexpressing chondrocytes (Fig. 7a).
When we transduced Runx3 adenovirus at different doses to primary chondrocytes from WT mice, Prg4 was increased in a dose-dependent manner (Fig. 7b). Acan and Sox9 were also increased by Runx3 overexpression (Fig. 7b). Catabolic factors were not changed, except for slightly increased Hif2a (Fig. 7b). We next introduced GFP or Runx3 adenoviral vectors into knee joints of WT mice that received surgical induction for the medial model at eight weeks of age. Runx3 expression was efficiently enhanced by the intra-articular injection of adenovirus, resulting in significant suppression of OA progression (Fig. 7c−d).