RSPO2 defines a distinct undifferentiated progenitor in the tendon/ligament and suppresses ectopic ossification.

Ectopic endochondral ossification in the tendon/ligament is caused by repetitive mechanical overload or inflammation. Tendon stem/progenitor cells (TSPCs) contribute to tissue repair, and some express lubricin [proteoglycan 4 (PRG4)]. However, the mechanisms of ectopic ossification and association of TSPCs are not yet known. Here, we investigated the characteristics of Prg4-positive (+) cells and identified that R-spondin 2 (RSPO2), a WNT activator, is specifically expressed in a distinct Prg4+ TSPC cluster. The Rspo2+ cluster was characterized as mostly undifferentiated, and RSPO2 overexpression suppressed ectopic ossification in a mouse Achilles tendon puncture model via chondrogenic differentiation suppression. RSPO2 expression levels in patients with ossification of the posterior longitudinal ligament were lower than those in spondylosis patients, and RSPO2 protein suppressed chondrogenic differentiation of human ligament cells. RSPO2 was induced by inflammatory stimulation and mechanical loading via nuclear factor κB. Rspo2+ cells may contribute to tendon/ligament homeostasis under pathogenic conditions.

INTRODUCTION

Tendons and ligaments are connective tissues that are essential for joint stability and motion. They are characterized as dense fibers that consist of collagen fibrils (). Tendons and ligaments connect muscle to bone and bone to bone, respectively. Ectopic ossification of the tendon is common in populations with high risk of tendon injury such as athletes and workers with repetitive tendon overuse (). Although the underlying mechanisms have not been fully elucidated, it is thought that a failed injury repair process that begins with inflammation is involved in ectopic ossification of the tendon (). Histologically, the process is usually followed by endochondral ossification, which begins with chondrocyte differentiation, hypertrophy, and replacement by osteoblasts (). Similar pathological conditions occur in the ligaments. Ectopic ossification is often observed in the ligaments, e.g., ossification of the posterior longitudinal ligament (OPLL) of the spine and ossification of the yellow ligament (, ). A previous genome-wide association study (GWAS) identified several susceptibility loci that are located near the HAO1, RSPO2, EIF3E, EMC2, and CCDC91 genes (). Among them, R-spondin 2 (RSPO2) is associated with the regulation of endochondral ossification ().

Because there is no effective therapy for the abovementioned diseases, several studies have aimed to determine the mechanisms of tendon/ligament homeostasis, degeneration, and ossification. Recently, several markers for tenocytes and tendon stem/progenitor cells (TSPCs) have been reported. Regarding differentiated tendon cells, transcription factors including scleraxis (SCX), mohawk (MKX), and glycoprotein tenomodulin (TNMD) are well-established markers (). Regarding TSPCs, a recent study using single-cell transcriptomics identified the tubulin polymerization-promoting protein family member 3 (TPPP3) as a novel marker for TSPC (). The cell population that expresses both Tppp3 and platelet-derived growth factor receptor-alpha (Pdgfra), a representative marker for mesenchymal stem cells (MSCs), is defined as TSPC (). Tppp3/Pdgfra-positive (+) cells can generate new tenocytes and undergo self-renew upon injury (); however, the role of these cells in tendon/ligament homeostasis and their association with pathological conditions, such as enthesopathy and ligament ossification, remain unknown.

Proteoglycan 4 (PRG4), also known as lubricin, is an extracellular matrix protein that is expressed in superficial cells of the articular cartilage and synovium, ligaments, and tendons (, ). It plays essential roles in the boundary lubrication of synovial joints and tendons (). In humans, loss-of-function mutations in PRG4 leads to camptodactyly-arthropathy-coxa vara-pericarditis syndrome (). Cell tracking analyses using Prg4-Cre mice demonstrated that Prg4 cells in the superficial layer of articular cartilage serve as a progenitor population for chondrocytes in deeper layers (). We recently reported another role of lubricin, in which it contributes to joint homeostasis by suppressing differentiation of the superficial cells and the ectopic chondrogenic differentiation of synovial cells (). Prg4 is expressed in tendons and ligaments and has been reported to be expressed in Tppp3 cells (); however, Prg4-expressing cells in these tissues have not yet been characterized. Considering the progenitor properties of Prg4 cells and the inhibitory effects of PRG4 against cell differentiation in the cartilage and synovium, we hypothesized that PRG4 is also associated with the differentiation of progenitors in connective tissues such as tendons and ligaments and that Prg4 cells in those locations have distinct features.

Here, we describe the identification and characterization of an Rspo2 cluster within the Prg4 cell cluster of a mouse model of tendon ectopic ossification. We examined cell populations during ectopic ossification of a mouse tendon model by single-cell RNA sequencing (scRNA-seq) using Prg4-Cre and R26-tdTomato mice. We characterized the Rspo2 cluster and analyzed interactions between the identified clusters. We compared our findings with deposited data from another mouse heterotopic tendon ossification model and further examined the role of the RSPO2 protein in the Prg4 clusters in tendon homeostasis using tissue-specific transgenic mice. We then investigated the expression and role of RSPO2 in human ligament cells and surgical specimens by RNA-seq and immunostaining. Last, we show that RSPO2 is induced by mechanical loading and nuclear factor κB (NF-κB).

RESULTS

Ectopic endochondral ossification occurs in a mouse Achilles tendon puncture model

To examine cell populations during ectopic endochondral ossification in vivo, we used a mouse Achilles tendon puncture (ATP) model, in which ectopic endochondral ossification in the tendon was produced by minimal invasion (). As previously described, mouse Achilles tendons were punctured with a 27-gauge needle at five different sites. Micro–computed tomography (μCT) detected ossification between 8 and 12 weeks after ATP induction (Fig. 1A). Before the ossification, expression of sex-determining region Y-box transcription factor 9 (SOX9) and collagen type 2 (COL2) began to increase at 1 week, and COL2-positive cell rates remained higher than those in the sham group throughout (Fig. 1B). COL2 was detected around cuboid- or spherical-shaped cells (Fig. 1B). PRG4-positive cells tended to increase at 1 week, but this was not significant compared with the sham group (Fig. 1B). We also confirmed no obvious changes in the sham tendons during the experimental period (fig. S1, A and B). These data indicate that ectopic chondrogenic differentiation begins within 1 week after the ATP.

Fig. 1.

Ectopic endochondral ossification in a mouse ATP model.

(A) μCT images of the tendons in a time-course analysis after ATP. The right panel indicates the bone volume of the tendons. Sham mice received only skin puncture and were analyzed 12 weeks after surgery. (B) Immunohistochemistry of SOX9, COL2, PRG4, and RSPO2 in a time-course analysis after ATP. Right panels indicate positive cell rates. All data are expressed as dot plots and means ± SD (n = 8 mice per group). *P < 0.05, **P < 0.01, and ****P < 0.0001 versus sham, one-way analysis of variance (ANOVA) and Tukey’s multiple comparisons test.

Characteristics of the TSPC population change after invasion

To optimize the protocol for scRNA-seq, we confirmed the profile of Prg4-expressing cells using deposited scRNA-seq data of a heterotopic ossification model of the mouse Achilles tendon, which was induced by burn and tenotomy (GSE126060) (). The uniform manifold approximation and projection (UMAP) plotting of the scRNA-seq data at days 0 and 7 identified nine distinct cell clusters (fig. S2A). On the basis of the expression distribution of genes related to each differentiation pathway, two clusters, 1 and 2, were characterized as mesenchymal cells (fig. S2, A to C). In the secondary analyses of the mesenchymal cell clusters, seven clusters were identified (fig. S3A). In the normal mouse Achilles tendon at day 0, most cells were divided into cluster 2 or 6. Cluster 2, which highly expresses Tppp3 and Pdgfa, was characterized as TSPCs (fig. S3, B and C). Cluster 6 expressed Prg4 and Tnmd, markers of tenocytes (fig. S3C). Notably, at day 7, the number of cluster 2 cells was markedly decreased and several Pdgfra clusters emerged (fig. S3B). In addition to cluster 6, Prg4 cells were detected in cluster 7 at day 7 (fig. S3, B and C). Cluster 7 expressed some progenitor markers including Pdgfra, Cd151, and Bglap, instead of tenogenic or tenocyte markers (fig. S3C). Together, in addition to the Prg4 expression level, the characteristics of Prg4 cells were altered in response to invasion.

Rspo2 is specifically expressed in a subpopulation of the Prg4+ progenitor cluster

Next, we prepared Prg4-Cre mice to analyze the Prg4 cell population by scRNA-seq. To examine the postinvasive changes of Prg4 expression and Prg4 cells, we injected tamoxifen twice into Prg4-Cre mice 3 and 4 days after ATP induction and performed scRNA-seq using cells obtained from the Achilles tendons 1 week after the ATP induction. Histological analysis was used to confirm that tdTomato expression was not present in the tendon without tamoxifen injection (fig. S4A) and that Prg4 and tdTomato were coexpressed after tamoxifen injection (fig. S4B). The UMAP using the Seurat package () showed 10 distinct clusters (Fig. 2A and fig. S5). Clusters 0, 2, and 5 were estimated as undifferentiated by CytoTRACE () (Fig. 2B) and highly expressed undifferentiated mesenchymal cell markers (fig. S5). We then performed secondary analyses of the mesenchymal cell clusters and identified seven different cell clusters (Fig. 2C). Prg4 was highly expressed in clusters 6 and 7, and tdTomato was only detected in these two clusters (Fig. 2D). Notably, Rspo2 was specifically expressed in cluster 6 (Fig. 2D). RSPO family members are secreted proteins that can activate canonical WNT signaling (). The WNT pathways play essential roles in mesenchymal cells during skeletal development (). RSPO2 has been identified as an OPLL susceptible gene (), and it suppresses the early stage of the endochondral ossification process (). In contrast to cluster 6, another Prg4+ cluster, cluster 7, did not express Rspo2 (Fig. 2D). Cluster 7 expressed marker genes for chondrogenic and endochondral ossification, such as Col2a1, Acan, and Sox9, indicating that the Rspo2 cluster can be characterized as a chondrogenic progenitor population (Fig. 2D). In the deposited scRNA-seq data of GSE126060, Rspo2 was specifically expressed in one of the Prg4+ clusters (cluster7; fig. S3C). Notably, the Rspo2 cluster was rarely detected at day 0 (fig. S3, B and C).

Fig. 2.

scRNA-seq of Achilles tendons from Prg4-Cre ATP mice.

(A) Uniform manifold approximation and projection (UMAP) of 10 clusters from the Achilles tendon. Overall, 7157 cells were analyzed from five mice. Clusters 0, 2, and 5, enclosed in a dotted line, are characterized as mesenchymal cells. (B) Predicted differentiation stages for each cluster using CytoTRACE. (C) UMAP from the mesenchymal cell clusters [clusters 0, 2, and 5 in (A)]. (D) Violin plots of marker genes for various cell types in mesenchymal cell clusters.

To confirm Rspo2 expression in Prg4+ cells, we collected tdTomato+ cells from the Achilles tendons of Prg4-Cre mice by fluorescence-activated cell sorting (FACS). We used cells derived from the Achilles tendons of CAG-Cre;R26-tdTomato and wild-type (WT) mice as FACS-positive and FACS-negative controls, respectively. We collected tdTomato+ and tdTomato− cells from the Achilles tendons of Prg4-Cre mice 1 week after ATP induction followed by subsequent tamoxifen administration, according to the threshold values using positive and negative controls (fig. S6A). Prg4 and Rspo2 were clearly and specifically expressed in the tdTomato+ cells (fig. S6B). We then compared tdTomato+ cells from the Achilles tendons of Prg4-Cre mice 1 week after ATP or sham surgery. Notably, tdTomato+ cells from sham mice expressed Tppp3, but mRNA levels of Rspo2, Pdgfra, and Procr were much lower (fig. S6C). In contrast, tdTomato+ cells from ATP mice highly expressed Rspo2, as well as Pdgfra, Tppp3, and Procr (fig. S6C). These data were consistent with the deposited scRNA-seq data of GSE126060 at days 0 and 7 (fig. S3, B and C).

We performed immunofluorescence to examine RSPO2 expression in the Achilles tendons of Prg4-Cre mice, which received tamoxifen administration after ATP induction. tdTomato+ cells, but not RSPO2+ or COL2+ cells, were present in sham tendons (Fig. 3 and fig. S7). RSPO2+tdTomato+ and COL2+tdTomato+ cells were observed in the center of tendons at 1 and 2 weeks after ATP (Fig. 3 and fig. S7). RSPO2+tdTomato+ cells, but not COL2+ cells, were present in the tendon sheath (Fig. 3). In addition, we injected tamoxifen into Prg4-Cre mice 1 week before ATP induction. Unlike the expansion of tdTomato+ cells in the postinjury tamoxifen-injected mice (Fig. 3), few tdTomato+ cells were present, and RSPO2+tdTomato+ cells were not observed in the preinjury tamoxifen-injected mice (fig. S8).

Fig. 3.

Immunofluorescence of tdTomato+ cells in the tendons of the Prg4-Cre ATP model.

Immunofluorescence of RSPO2, tdTomato, and COL2 in the center or sheath of Prg4-Cre Achilles tendons 1 and 2 weeks after ATP induction. Tamoxifen was injected twice, 3 and 4 days after ATP induction. Eight mice were analyzed, and representative images are shown. Positive cell rates are shown in fig. S7.

The Rspo2 cluster is characterized as a cluster of undifferentiated progenitors

In addition to Rspo2 and Prg4, cluster 6 cells expressed stem cell markers such as Pdgfra, Tppp3, Itga6, Tspan15, and Procr (Fig. 2D) (, –). In addition, Bglap-positive cells have been reported to be involved in tendon repair (). Tppp3 cells, characterized as TSPCs, were detected in clusters 1, 3, and 6 (Fig. 2D). Compared with cluster 6, clusters 1 and 3 predominantly expressed Col15a1, Fst, and Runx2 (Fig. 2D). Col15a1 is an indicator for the differentiation of MSCs toward osteogenic differentiation (). Follistatin, coded by Fst, is involved in osteogenesis (). Runx2 is a master transcription factor for osteoblast differentiation (–). Scx, a master transcription factor for tenogenesis, was highly expressed in clusters 0, 2, 4, and 5 (Fig. 2D). Of these, cluster 5 abundantly expressed Mkx and Tnmd, indicating that the Scx cluster can be characterized as a tenogenic progenitor population (Fig. 2D).

The Gene Ontology (GO) pathways enriched in the seven clusters were analyzed by clusterProfiler () and are shown in Fig. 4A. “Extracellular structure organization,” “extracellular matrix organization,” “external encapsulating structure organization,” and “collagen-containing extracellular matrix” were up-regulated terms in clusters 0, 1, 2, 5, and 7 (Fig. 4A). Conversely, “ossification” was detected only in cluster 7 (Fig. 4A). GO pathways related to cell cycle regulation and proliferation were detected in cluster 4 (Fig. 4A). The GO pathways related with specific differentiation tendencies were not up-regulated in cluster 6 (Fig. 4A).

Fig. 4.

GO and pseudotime analyses of the scRNA-seq data of the tendons from Prg4-Cre ATP mice.

(A) Representative GO pathway enrichment of the 11 clusters in Fig. 2C using clusterProfiler. (B) Diffusion map pseudotime analysis of the mesenchymal cell clusters in Fig. 2C. dpt, diffusion pseudotime. (C) Marker gene plotting in pseudotime by Monocle2. (D) Trajectory analysis of the mesenchymal cell clusters by Monocle2. Characterization of each cluster is as follows: 0, 2, and 4, tenogenic; 1, 3, and 6, undifferentiated progenitor; 5, tenocyte; 7, chondrogenic, as shown in Fig. 2C.

Next, we performed pseudotime analysis. In diffusion map pseudotime (), cluster 6 was clearly located in the most undifferentiated position (Fig. 4B). Rspo2 and other stem cell markers were highly expressed in the undifferentiated state during the pseudotime course estimated by Monocle2 () (Fig. 4C). On the basis of gene expression distribution patterns, a trajectory starting from the Rspo2 cluster (cluster 6) was created (Fig. 4D). In the trajectory, the Tppp3 clusters (clusters 1 and 3) were downstream of the Rspo2 cluster, and the tenogenic and chondro/osteogenic clusters were located in a different downstream trajectory (Fig. 4D). These data indicate that Rspo2 cells are undifferentiated progenitors.

RSPO2 inhibits ectopic endochondral ossification in tendons

To estimate the role of the Rspo2 cluster and the RSPO2 protein in tendons, we evaluated cell-cell interrelationships between the Rspo2 cluster and the chondrogenic, tenogenic, Tppp3, and Rspo2 clusters using scTensor (tables S1 to S4) (). Leucine-rich repeat containing G protein–coupled receptors 4 to 6 (LGR4 to LGR6) are known as receptors for Rspo2 (). Notably, Rspo2-Lgr6 was ranked fourth in the ligand-receptor pairs from the Rspo2 cluster to the chondrogenic progenitor (table S1 and fig. S9). Considering previous findings that showed that RSPO2 suppresses the early stage of chondrogenic differentiation (), it was hypothesized that the Rspo2 cluster suppresses differentiation of the chondrogenic cluster via RSPO2 and other secreted factors.

To examine roles of RSPO2 in tendons in vivo, we generated Cre-dependent transgenic mice (CAG-EGFP-Rspo2; Fig. 5A). We confirmed the transgenic expression of Rspo2 in the Achilles tendon of Prg4-Cre mice by measuring the RSPO2 protein levels (Fig. 5B). Bone mass of ectopic ossification areas was significantly decreased in Prg4-Cre mice compared with Prg4-Cre mice (Fig. 5C). mRNA levels of Sox9 and Col2a1 were significantly reduced in Prg4-Cre tendons (Fig. 5D). Immunofluorescent analysis of Prg4-Cre tendons demonstrated the efficient expression of RSPO2 protein in the tdTomato+ cells and the suppression of COL2 expression (Fig. 5E).

Fig. 5.

Rspo2 overexpression suppresses ectopic ossification in a mouse ATP model.

(A) Schematic of a transgene for Cre-dependent Rspo2 overexpression (CAG-EGFP-Rspo2 mouse). (B) Protein levels of RSPO2 in the Achilles tendon of 9-week-old Prg4-Cre and Prg4-Cre mice that received tamoxifen injection at 8 weeks of age. n = 3 mice per group. (C) Ectopic ossification in the Achilles tendon of 20-week-old Prg4-Cre and Prg4-Cre mice, which received tamoxifen injection and ATP at 8 weeks of age. Left panels are representative images of μCT, and right panels show quantified bone volume in the Achilles tendons. n = 14 mice per group. (D) mRNA levels of Col2a1, Sox9, Prg4, and Rspo2 in the Achilles tendons of Prg4-Cre and Prg4-Cre mice 1 week after ATP induction. n = 6 mice per group. (E) Immunofluorescence of RSPO2, tdTomato, and COL2 in the Achilles tendons of Prg4-Cre and Prg4-Cre mice 1 week after ATP induction. All data are expressed as dot plots and means ± SD. *P < 0.05, ***P < 0.001, and ****P < 0.0001, unpaired t test.

We examined the effects of an antibody against RSPO2 injected intraperitoneally into WT mice 3 days after ATP induction. Bone mass was significantly decreased in RSPO2 antibody–injected mice compared with vehicle-injected mice (fig. S10A). mRNA levels of Sox9 and Col2a1 were significantly increased in the tendons of RSPO2 antibody–injected mice (fig. S10B). All these data indicate that RSPO2 protein inhibits ectopic endochondral ossification via suppression of chondrogenic differentiation.

RSPO2 suppresses chondrogenic differentiation of human ligament cells

We next investigated the expression and role of RSPO2 in ligament cells. We performed RNA-seq of human primary ligament cells that had undergone chondrogenic differentiation (fig. S11, A to C). The heatmap shows the increase of many chondrogenic transcription factors and cartilage matrix genes, such as SOX9, SOX6, SOX5, RUNX3, COL2A1, COL11A1, COL11A2, and ACAN (Fig. 6A). GO analysis indicated an up-regulation of processes such as “skeletal system development,” “ossification,” “extracellular matrix organization,” “extracellular structure organization,” and “tissue development” (Fig. 6B). Quantitative reverse transcription polymerase chain reaction (qRT-PCR) confirmed the down-regulation of RSPO2, accompanied by a decrease of AXIN2, a representative downstream molecule of canonical WNT pathway, and the MSC markers, ENG (CD105) and NT5E (CD73; Fig. 6C).

Fig. 6.

RSPO2 expression decreases during the chondrogenic differentiation of human ligament cells.

RNA-seq of human ligament cells following chondrogenic differentiation. Control (Cntl) and differentiation (Diff) groups (n = 3 each) were analyzed. (A) Heatmap of differentially expressed genes extracted by the terms “chondrocyte,” “development,” and “differentiation” in the GO biological process. The blue arrow indicates RSPO2. (B) GO biological processes altered by chondrogenic differentiation. (C) Validation of mRNA levels of RSPO2 and marker genes by qRT-PCR. Data are expressed as dot plots and means ± SD (n = 9 per group, three independent experiments using cells from three different donors). ***P < 0.001 and ****P < 0.0001, unpaired t test.

To compare these gene expression changes with those in standard chondrogenic differentiation, we reanalyzed deposited RNA-seq data of human MSCs that had undergone chondrogenic differentiation (). Similar to the differentiation of human ligament cells, the heatmap showed the increase of many chondrogenic transcription factors and cartilage matrix genes during human MSC differentiation (fig. S12). RSPO2 increased at day 1 and decreased thereafter (fig. S12). In the comparison of data from days 0 and 7, GO processes, such as “skeletal system development,” “ossification,” “extracellular matrix organization,” “extracellular structure organization,” and “tissue development,” were also up-regulated (fig. S13, A and B).

To examine the role of RSPO2 protein in the differentiation of ligament cells, we next performed in vitro experiments using human primary ligament cells. We induced chondrogenic differentiation of human ligament cells in the presence of recombinant human (rh) RSPO2 and LiCl, an activator of canonical WNT signaling, or a vehicle. Both rhRSPO2 and LiCl treatments reduced Alcian blue staining (fig. S14A) and the expression of COL2A1 and ACAN (fig. S14B), whereas AXIN2 was increased (fig. S14B). PRG4 expression was mildly decreased by rhRSPO2 and LiCl treatments (fig. S14B). Luciferase activity of TOPflash, a reporter of the canonical WNT pathway, was clearly up-regulated by both rhRSPO2 and LiCl (fig. S14C).

RSPO2 expression is down-regulated in patients with OPLL

Because RSPO2 is a known OPLL susceptibility gene, we next examined its expression in the spine. In a WT mouse model, RSPO2 protein was detected in the posterior longitudinal ligament and posterior anulus fibrosus (Fig. 7A). These areas were stained with X-gal in a TOPGAL reporter mouse to investigate canonical WNT signaling (Fig. 7B). We next examined the expression level of RSPO2 in human paraspinal ligaments. Because it was difficult to obtain posterior longitudinal ligaments from non-OPLL patients, we compared yellow ligaments from patients with cervical OPLL and cervical spondylotic myelopathy (CSM). Both the mRNA level of RSPO2 and the RSPO2+ cell rate were significantly lower in the OPLL ligaments than those in CSM ligaments (Fig. 7, C and D). RNA-seq of OPLL and CSM samples showed a decrease in RSPO2 expression, and ingenuity pathway analysis indicated a down-regulation of the canonical WNT pathway in the OPLL samples in comparison to CSM samples, although the Z score (−1.291) is less than the significance threshold (Fig. 7, E and F, and fig. S15, A to C). These findings are consistent with previous studies (, ) and the inhibitory effect of RSPO2 on endochondral ossification.

Fig. 7.

RSPO2 expression is down-regulated in patients with OPLL.

(A) Immunohistochemistry of Rspo2 in the cervical vertebrae of C57BL/6J mice. The inset box in the left panel indicates the location of the enlarged image shown in the right panel. (B) X-gal staining in the cervical vertebrae of TOPGAL mice. The inset box in the left panel indicates the location of the enlarged image shown in the right panel. (C) mRNA levels of RSPO2 in spinal ligaments from patients with OPLL (n = 21) and CSM (n = 17). (D) Immunohistochemistry of RSPO2 in spinal ligaments from patients with OPLL (n = 14) and CSM (n = 12). The two left panels show representative images; the right panel indicates the positive cell rates. (E) Heatmap of differentially expressed genes (false discovery rate < 0.1, |logFC| > 4) ordered by gene expression level from RNA-seq analysis of human spinal ligaments obtained from patients with OPLL and CSM (n = 2 per group). The blue arrow indicates RSPO2. (F) Ingenuity pathway analysis of the RNA-seq data. |Z score| > 2 was the significance threshold. All data are expressed as dot plots and means ± SD. **P < 0.01, unpaired t test.

We examined the location of RSPO2 in the posterior longitudinal ligaments resected from patients with OPLL and found a representative ligament tissue in which endochondral ossification was progressing (e.g., containing regenerated ligament, cartilage, and bone; fig. S16). RSPO2+ cells were detected in the boundary region between COL2+ and degenerated ligament areas but were not detected in the ossification lesions (fig. S16).

RSPO2 is induced by NF-κB signaling

Last, we analyzed the upstream pathway of RSPO2. We analyzed the levels of transcription factors that played regulatory roles in each of the identified mesenchymal cell clusters from Fig. 2C, using SCENIC () (fig. S17), and extracted those in the Rspo2 undifferentiated progenitor (cluster 6; Fig. 8A). FOXO1 (forkhead box O1) and CREB5 (cAMP responsive element binding protein 5) were reported to regulate PRG4 expression in articular cartilage (, ). Among the 14 transcription factors estimated to be involved in regulation of cluster 6, we focused on Rela, a representative transcription factor of NF-κB, because NF-κB is widely associated with various biological phenomena including responses to injury, tissue repair, and inflammation. Although the mRNA level of Rela does not correlate with NF-κB activity, the Rela expression level was relatively higher in the Rspo2 cluster compared with other clusters (Fig. 8, B and C). Human ligament cells cultured with interleukin-1β (IL-1β) showed a marked induction of RSPO2 (Fig. 8D); this induction was significantly reduced by Bay11-7085, an IκB kinase (IKK) inhibitor (Fig. 8E).

Fig. 8.

The NF-κB pathway induces RSPO2.

(A) Heatmap of upstream transcription factors that regulate the Rspo2+ cluster in Fig. 2C, extracted by SCENIC. Heatmap for all clusters is shown in fig. S17. Feature plot (B) and violin plot (C) showing Rela expression in the cluster depicted in Fig. 2C. (D) mRNA levels of RSPO2 in human ligament cells cultured in the presence of vehicle and various amounts of IL-1β. (E) mRNA levels of RSPO2 in human ligament cells cultured in the presence of vehicle and IL-1β (1 ng/ml), with or without IKK inhibitor, Bay11-7085 (Bay). (F) mRNA level of RSPO2 in human ligament cells before (BL) and after tensile stress loading (+) or without loading (−). (G) mRNA level of RSPO2 in human ligament cells after tensile stress loading with or without Bay. (H) Protein level of inhibitor of nuclear factor κBα (IκBα), phosphorylated IκBα (p-IκBα), RSPO2, and glyceraldehyde 3-phosphate dehydrogenase (GAPDH) in human ligament cells after tensile stress loading with or without Bay. (I) mRNA levels of RSPO2 in human ligament cells after tensile stress loading with or without Rac1 inhibitor, EHT1864 (EHT; left) or NSC23766 (NSC; right). All data are expressed as dot plots and means ± SD. n = 9 per group, three independent experiments using cells from three different donors (D to G and I); n = 3 per group, using cells from three different donors (H). *P < 0.05, **P < 0.01, and ****P < 0.0001, one-way ANOVA and Tukey’s multiple comparisons test.

In addition to inflammatory stimulation, we previously reported that NF-κB also has essential roles in mechanical overloading-induced cartilage degeneration; mechanical stress induces gremlin-1 via the Rac1–NF-κB pathway in chondrocytes, and induced gremlin-1 induces the excessive activation of the NF-κB pathway (). Because mechanical stress is associated with ectopic ossification of the tendon and ligament, we examined whether mechanical loading induced RSPO2 expression. Notably, RSPO2 was significantly induced by cyclic tensile stress loading (Fig. 8F), and this induction was also abolished by Bay11-7085 (Fig. 8G). RSPO2 was increased by cyclic tensile stress loading at the protein level, which was accompanied by increased inhibitor of nuclear factor κBα (IκBα) phosphorylation (Fig. 8H). The induction of RSPO2 by cyclic tensile stress loading was also diminished by two different Rac1 inhibitors, EHT1864 and NSC23766 (Fig. 8I). These data suggest that inflammation or mechanical loading on the tendon/ligament induces RSPO2 expression via the Rac1–NF-κB pathway, similarly to gremlin-1 induction in chondrocytes ().

DISCUSSION

The present study demonstrated that Rspo2 is coexpressed with Prg4 in a subset of Tppp3 progenitors in the mouse Achilles tendon. This population is characterized by undifferentiated progenitors rather than Tppp3 progenitors. Rspo2 overexpression in Prg4 cells suppressed ectopic endochondral ossification of the mouse Achilles tendon. Similarly, RSPO2 inhibited chondrogenic differentiation of human ligament cells, and its expression is decreased in patients with OPLL. Last, we identified NF-κB as an upstream regulator of RSPO2 and displayed that inflammation or mechanical loading induced RSPO2 via the Rac1–NF-κB pathway.

Nakajima et al. () performed a GWAS of OPLL in ~8000 individuals followed by a replication study using an additional ~7000 individuals. The GWAS revealed six susceptibility loci for OPLL (). The rs374810 locus, one of the most significant loci, is located 121 bp upstream of the transcription start site of the RSPO2 gene (). RSPO2 is a potent activator of canonical WNT signaling (). The activation of the canonical WNT pathway in immature chondrocytes suppresses maturation (, ). Nakajima et al. () also reported roles of RSPO2 protein in chondrocyte differentiation. In the chondrogenic cell line, ATDC5, overexpression of RSPO2 inhibited chondrocyte differentiation, while RSPO2-neutralizing antibody stimulated it (). In the present study, the Rspo2-Lgr6 pair was highly estimated to function between the Rspo2 undifferentiated progenitor and Rspo2 chondrogenic progenitor by the CCI (cell-cell interaction) analysis (table S1 and fig. S9). RSPO2 overexpression reduced ectopic ossification in the mouse Achilles tendon (Fig. 5), and rhRSPO2 treatment suppressed chondrogenic differentiation of human ligament cells (fig. S14). The anti-chondrogenic differentiation effects of RSPO2 protein can be explained by the activation of the canonical WNT signaling pathway. We did not perform scRNA-seq using human ligament cells in the current study. Further comprehensive studies will be necessary to determine the roles of the RSPO2+ subset in human ligaments and the pathophysiology of ligament ossification.

Staverosky et al. () first reported that Tppp3 is expressed in cells at the circumference of the developing tendon, where there are abundant progenitors. Recently, Harvey et al. () reported that the Tppp3 cell population can generate new tenocytes and can undergo self-renewal upon injury, indicating that expression of Tppp3 and Pdgfra defines TSPCs. In our analysis, the Tppp3 population consisted of heterogeneous cells, as shown by previous and present scRNA-seq data (Fig. 2D and fig. S3C). We found that Rspo2 is specifically expressed in a distinct cluster of Tppp3 cells. RSPO2 and LGRs are involved in the maintenance of stemness and self-renewal (, ). In addition to Rspo2, Tppp3, and Pdgfra, the Rspo2+ cluster expresses Itga6, Tspan15, and Procr (Fig. 2D), which are known stem cell and progenitor markers (–). The pseudotime analyses also support the undifferentiated status of the Rspo2cluster (Fig. 4, B to D). All these data consistently indicate that Rspo2 is expressed in a distinct undifferentiated progenitor cluster. However, the RSPO2+ cluster is rarely detected in normal tendons and emerges after invasion (Fig. 3 and fig. S3B). Because RSPO2 expression is increased by inflammation or mechanical loading through NF-κB, these stimulations may be involved in the induction of the RSPO2+ cluster. Considering that RSPO2 suppresses ectopic endochondral ossification after tendon injury, one of the major roles of the RSPO2+ cluster might be to support the appropriate regeneration of tendons. Although we did not determine the origin of the RSPO2+ population in this study, we showed that it does not derive from PRG4+ lineage cells (fig. S8). Some TSPCs in normal tendons may acquire the de novo expression of RSPO2 and PRG4 by invasion. In addition, we did not determine the fate of the RSPO2+ cluster. Although pseudotime analyses indicated the undifferentiated state of the Rspo2 cluster, we did not obtain data to demonstrate its subsequent differentiation.

In the articular cartilage, PRG4 is expressed in the superficial zone (). Kozhemyakina et al. () reported that Prg4 cells serve as a progenitor population for mature chondrocytes in deeper layers. They showed that Prg4 cells in the superficial zone expand into deeper regions in a time-course analysis with cell tracking using Prg4-Cre mice (). In addition, the PRG4 protein itself suppresses chondrocyte differentiation of superficial immature cells and synovial cells (). However, Prg4 cells in the articular cartilage have not been characterized by scRNA-seq, and their heterogeneity has not yet been revealed. The present study indicated the presence of two populations of Prg4 cells in the tendon/ligament: Rspo2 undifferentiated progenitors and Rspo2 chondrogenic progenitors. Although RSPO2 in the growth plate chondrocytes is involved in endochondral ossification during skeletal development () and RSPO2 derived from the synovium is involved in osteoarthritis development (, ), the expression and role of RSPO2 in the articular cartilage are unknown. We previously reported that the canonical WNT pathway is activated by mechanical loading in the superficial zone of articular cartilage and contributes to the transcriptional induction of Prg4 (). However, Prg4 was not induced by RSPO2 overexpression in the mouse Achilles tendon (Fig. 5D). Characteristics of Prg4 or Rspo2 cells and the role of RSPO2 may be quite different between the articular cartilage and the tendon/ligament.

Recently, we reported that mechanical overload induced gremlin-1 via the Rac1–NF-κB pathway in chondrocytes (). Consistently, we found that Rspo2 is induced by the same pathway, which was determined by estimation of the transcription factors associated with the Rspo2cluster using SCENIC (Fig. 8A). In addition to gremlin-1, Rspo2 is induced by IL-1β and mechanical loading, and its induction is disrupted by NF-κB inhibitors (Fig. 8). Because NF-κB is widely involved in inflammatory responses, Rspo2 may be induced by proinflammatory cytokines other than IL-1β during tissue injury and repair processes. Repetition of inflammation or mechanical overloading can cause micro-injury, degeneration, and consequently ectopic ossification of the tendon/ligament (). Together, the Rspo2 cluster likely contributes to appropriate maintenance or repair of the tendon/ligament under such pathological conditions by suppression of chondrogenic differentiation in progenitor cells. We did not show the biological significance of mechanical loading and NF-κB in the regulation of RSPO2 expression in this study. Furthermore, we did not examine other upstream pathways that might induce RSPO2. In addition to inflammation and mechanical loading, aging is another factor associated with this pathophysiology (). In the present study, we did not examine alterations in Rspo2 expression levels or in Rspo2 cells with aging. These issues should be investigated to fully determine how RSPO2 and Rspo2 cells are regulated against various conditions such as inflammation, overuse, and aging.

We identified an Rspo2 undifferentiated cluster in our scRNA-seq data from a mouse ATP model and confirmed this cluster using deposited scRNA-seq data of a different mouse Achilles tendon model. However, we could not perform scRNA-seq of ligaments of human samples because resection surgery of ossified tendons/ligaments is rarely performed, and therefore, appropriate human surgical specimens were not available. Moreover, we could not obtain a sufficient number of cells for scRNA-seq analysis of mouse ligaments. Considering that tendons and ligaments are not identical, the presence and role of Rspo2 cells in ligaments and in human samples should be examined in the future.

In conclusion, we found that Rspo2 is specifically expressed in one of two Prg4 clusters in the tendon/ligament. Rspo2 expression defines a distinct undifferentiated progenitor subpopulation of Tppp3 cells. RSPO2 protein suppresses ectopic endochondral ossification in the tendon/ligament. Moreover, Rspo2 is induced by inflammatory stimulation and mechanical loading via the NF-κB pathway. The present findings contribute to our understanding of tendon/ligament homeostasis and of the mechanisms underlying ectopic ossification that is observed in diseases, such as OPLL.

MATERIALS AND METHODS

Animals

All animal experiments were conducted in accordance with the guidelines described in the Guide for the Use and Care of Laboratory Animals of the Institute for Laboratory Animal Research (ILAR, 2011). Furthermore, animal experiments were approved by the Institutional Animal Care and Use Committee (IACUC) of The University of Tokyo (M-P17-091). All animals were housed in a facility under the supervision of the IACUC at 18° to 22°C with a 12-hour light/dark cycle and were allowed free access to food and water. Genotypes were analyzed in each experiment. For all in vitro and in vivo studies requiring WT mice, adult male (8 weeks old) C57BL/6J mice were purchased from Sankyo Lab Services (Tokyo, Japan). The comparative groups were maintained on a C57BL/6J background. TOPGAL mice, Prg4-Cre, and Rosa26-tdTomato mice were obtained from the Jackson Laboratory (Sacramento, CA). We generated Prg4-Cre mice by mating Prg4-Cre male mice with R26-tdTomato female mice. The primer sequences used for genotyping are shown in table S5. For the generation of CAG-EGFP-Rspo2 mice, the transgene shown in Fig. 5A was constructed. DNA purification and microinjection were performed in accordance with standard protocols.

Mouse ATP model

All surgeries were performed on 8-week-old male mice under general anesthesia with the aid of an operating microscope. The mice were anesthetized by intraperitoneal injection with a mixture of medetomidine (0.3 mg/kg; #0002949, Nippon Zenyaku Kogyo, Koriyama, Japan), midazolam (4.0 mg/kg; #169738, Astellas, Tokyo, Japan), and butorphanol (5.0 mg/kg; #016440, Meiji Seika Pharma, Tokyo, Japan). A 27-gauge needle was used for the puncture, which was performed percutaneously from the outside of the Achilles tendon body, and this process was repeated five times at different sites on the Achilles tendon body of each mouse, as previously described (). In the sham operation, the needle was inserted through the skin without touching the Achilles tendon (). For tamoxifen induction, mice were injected intraperitoneally with tamoxifen (10 mg/ml; #T006000, Toronto Research Chemicals, Toronto, Canada) in corn oil at a dose of 50 μg/g body weight. The injection was performed twice over consecutive days before or after ATP induction. For the administration of an antibody against RSPO2, ATP was induced in 8-week-old WT mice, which then received 10 μl of RSPO2 antibody solution or vehicle (10 μg/ml) [phosphate-buffered saline (PBS)] injected into the Achilles tendon every day for 20 days starting from the day after ATP induction.

Micro–computed tomography

CT images were acquired using an inspeXio SMX-100CT (Shimadzu Corporation, Kyoto, Japan). We defined a 2 mm by 2 mm by 4 mm measurement space 0.1 mm proximal to the Achilles tendon attachment site and validated the volume of ossification within the measurement space. The imaging conditions were as follows: voltage, 75 kV; current, 120 μA; voxel size, 0.006 mm per voxel; image size, 1024 by 1024 pixels; and slice thickness, 0.04 mm. TRI/3D-BON-FCS64 (RATOC System Engineering, Tokyo, Japan) was used for analysis.

Histological analyses and immunohistochemistry

Mice were euthanized by carbon dioxide (CO2) inhalation and perfusion-fixed in 10% buffered formalin from the left ventricle for 5 min. The ankles and Achilles tendons were then dissected, fixed in 10% buffered formalin for 24 hours, demineralized in 10% EDTA (pH 7.4) for 14 days, and paraffin-embedded, and then sagittal sections of 4-μm thickness were cut out. The specimens were used for hematoxylin and eosin staining, nuclear fast red staining, X-gal staining, and immunohistochemistry. Information on antibodies is provided in table S6. After reacting with horseradish peroxidase–conjugated antibodies, positive signal colors were generated using Histofine Simple stain DAB (3,3’-diaminobenzine) (#415172, Nichirei Biosciences, Tokyo, Japan). The positively stained areas were measured using ImageScope (Leica Microsystems, Wetzlar, Germany). Fluorescence staining was evaluated by confocal microscopy (LSM880, Carl Zeiss, Oberkochen, Germany) using Alexa Fluor (#A32731, #A32728, and #A11058, Invitrogen, Waltham, CA) as a secondary antibody followed by a 4′,6-diamidino-2-phenylindole–containing inclusion agent (#12745-74, Nacalai Tesque Inc., Kyoto, Japan).

Quantitative reverse transcription polymerase chain reaction

Total RNA was extracted using the Direct-zol RNA Kit (#R2062, Zymo Research, Irvine, CA) in accordance with the manufacturer’s protocol. We confirmed that A260/A280 = 1.8 to 2.0 in all RNA samples. Total RNA was reverse-transcribed into cDNA using Rever Tra Ace qPCR RT Master Mix (#FSQ-201, TOYOBO, Osaka, Japan). Then, the qRT-PCR was performed using THUNDERBIRD SYBR qPCR Mix (#FSQ-201, TOYOBO) and Thermal Cycler Dice Real-Time System III (Takara Bio, Kusatsu, Japan). Information on primers is given in table S7. The absolute standard curve method was used to measure the expression level, which was normalized to glyceraldehyde 3-phosphate dehydrogenase expression.

scRNA-seq using 10x genomics

ATP was induced in five Prg4-Cre mice at 8 weeks of age, followed by tamoxifen injection 3 and 4 days after ATP, and were then euthanized 1 week after the ATP. The Achilles tendons were harvested and minced in ice-cold PBS. The harvested tissues were digested in RPMI 1640 medium (#30264-85, Nacalai Tesque Inc.) containing 0.3% type 1 collagenase (#CLS1, Worthington, Columbus, OH) and 0.4% dispase II (#04942078001, Roche, Basel, Switzerland) at 37°C for 45 min with constant agitation. Cells were filtered through 100-, 70-, and 40-μm strainers, centrifuged, counted, and resuspended at a concentration of 1000 cells/μl. Cell viability was assessed using the trypan blue exclusion method with a LUNA-FL (Logos Biosystems, Anyang, South Korea) automated counter, and samples with viability greater than 80% were processed for further sequencing. Single-cell libraries were prepared using the Chromium Controller (10x Genomics, Pleasanton, CA) according to the protocol of Chromium Next GEM Single Cell 3′ Reagent Kits (v3.1 10x Genomics). Libraries were prepared using Chromium Controller (10x Genomics). Alignment, quantification, and aggregation of sample count matrices were performed using the 10x Genomics Cell Ranger pipeline (v.5.0.0), retaining an average of 61,308 reads per cell mapped read depth normalization.

Bioinformatics analysis of single-cell sequencing data

The Cell Ranger summary showed a cell count of 7157, with a median of 3242 genes per cell. We performed the downstream analysis steps using the Seurat R package (). Cells with less than 1000 genes per cell and a mitochondrial read content greater than 5% were filtered out. Downstream analysis steps for each sample type included normalization, identification of genes with high variability among single cells, scaling based on UMI (unique molecular identifier) count and batch effects, dimensionality reduction (UMAP), unsupervised clustering, and discovery of cell type–specific markers with different expression levels. The clustering data were further subclustered using CytoTRACE (). The subcluster data were subjected to trajectory analysis using diffusion map pseudotime () and Monocle2 (), and cell-cell interactions were analyzed using scTensor (). SCENIC () was used to estimate upstream transcription factors.

Flow cytometry

Cells were obtained from the Achilles tendons of Prg4-Cre mice, as described above. The cells were filtered through 100-, 70-, and 40-μm strainers; centrifuged; and suspended in PBS with 0.1% bovine serum albumin for analysis. tdTomato+ cells were sorted using Cell Sorter SH800 (Sony, Tokyo, Japan).

Human samples

For human samples, written informed consent was obtained from all individuals, and approval was obtained from The University of Tokyo Ethics Committee (Institutional Review Board 0622-12). All relevant codes of ethics were complied with. Samples were collected from patients who underwent cervical spine surgery for OPLL or CSM.

Primary culture of human ligament cells

The harvested human paraspinal ligament tissues were digested in RPMI 1640 medium (#30264-85, Nacalai Tesque Inc.) containing 0.3% type 1 collagenase (#CLS1, Worthington) and 0.4% dispase II (#04942078001, Roche) at 37°C for 45 min with constant agitation. Primary cells were then seeded into 24-well collagen-coated microplates (#4820-010, AGC Techno Glass, Haibara-gun, Japan) at a density of 1 × 105 cells per well. For chondrogenic differentiation, the cells were cultured with chondrogenic medium (Dulbecco’s modified Eagle’s medium and 1% ITS (Insulin-Transferrin-Selenium) + Premix; #354352, Corning, New York City, NY), l-ascorbic acid diphosphate (50 μg/ml; #49752, Sigma-Aldrich), and l-proline (40 μg/ml; #P0380, Sigma-Aldrich, St. Louis, MO) supplemented with bone morphogenetic protein-2 (100 ng/ml; #100-21, PeproTech) and transforming growth factor–β1 (10 ng/ml; #100-21, PeproTech) at 37°C in a humidified atmosphere containing 5% CO2. The medium was changed twice a week.

RNA sequencing

Total RNA was extracted from human cells. We confirmed that A260/A280 = 1.8 to 2.0, and the RNA integrity number was >7 in all RNA samples. In accordance with the manufacturer’s protocol, the sequence libraries for expression analysis of human surgical samples were prepared using the NEBNext Ultra II RNA Library Prep Kit for Illumina (New England Biolabs, Ipswich, MA). Sequencing was performed on an Illumina HiSeq 4000 System with two 150-bp paired-end reads (Veritas Genetics, Danvers, MA). The raw sequence data were filtered using Trimmomatic (version 0.39) to remove adapter sequences, low-quality reads, N-rich sequences, and reads less than 50 bp in length. The filtered data were aligned against the human reference genome (GRCh38.p13) using STAR (version 2.7.3a). In addition, gene expression counts and transcripts per million value were calculated using RSEM (version 1.3.3). The data for RNA-seq of the chondrogenic differentiation samples were analyzed using DNBseq (BGI, Shenzhen, China) with 50-bp single-end or 100-bp paired-end reads. The filtered data were aligned against the human reference genome (GRCh38.p12). Data analysis was performed by iDEP.91 (http://bioinformatics.sdstate.edu/idep/) and ingenuity pathway analysis (QIAGEN, Düsseldorf, Germany).

Cyclic tensile strain loading

Human ligament cells were seeded into silicon stretch chambers coated with fibronectin (#356008, Corning) at a density of 1 × 105 cells per chamber; each chamber had a culture surface of 2 cm by 2 cm. After 48 hours, cyclic tensile strain (0.5 Hz, 10% elongation) was applied for 24 hours using the STB-140 mechanical stretch system (Strex, Osaka, Japan) in a CO2 incubator. Control cells were seeded into the same chambers and cultured without cyclic tensile strain. Cells were treated with 50 nM BAY11-7085 (#S7352, Selleck Chemicals, Houston, TX), 10 μM EHT1864 (#S7482, Selleck), and NSC23766 (#SML0952, Merck, Darmstadt, Germany).

Western blotting

Cells were lysed with T-PER (Tissue Protein Extraction Reagent) (#78510, Thermo Fisher Scientific, Waltham, MA) supplemented with protease and phosphatase inhibitors (#25955-24 and #07575-51, Nacalai Tesque Inc.). After mixing the samples with 4× LDS (lithium dodecyl sulfate) sample buffer (#2178723, Thermo Fisher Scientific) and dithiothreitol (#14128-91, Nacalai Tesque Inc.), the samples were boiled at 70°C for 10 min. Next, 10 μg of protein was applied to each lane of a 4 to 12% Bolt Bis-Tris Plus precast polyacrylamide gel (#NW04122BOX, Thermo Fisher Scientific) and separated by electrophoresis. After blocking with Blocking One (#03953-95, Nacalai Tesque Inc.) for 30 min, the polyvinylidene difluoride membranes (#033-22453, Wako, Osaka, Japan) were probed with the corresponding antibodies overnight at 4°C or for 1 hour at room temperature. The immunoreactions were visualized using ChemiLumi One Super (#02230-30, Nacalai Tesque Inc.) and iBright 1000 (Thermo Fisher Scientific). Band signals were measured using iBright Analysis Software (Thermo Fisher Scientific). Information on antibodies is given in table S8.

Luciferase assay

The human chondrosarcoma cell line SW1353 and human embryonic kidney cells 293 were seeded in 24-well plates. Lipofectamine 3000 (#L3000015, Invitrogen) was used to cotransfect the pGL4.49 (luc2P/TCF-LEF-RE/Hygro) reporter vector (Promega, Madison, WI) and pRL-TK (Promega). After 48 hours, cells were treated with rhRSPO2 (#120-43, PeproTech, NJ, USA) or lithium chloride (#123-01162, Wako) for another 24 hours. Luciferase assays were performed using a PicaGene Dual SeaPansy Luminescence Kit (#PGD-S, Toyo Bnet, Tokyo, Japan) and a GloMax 96 Microplate Luminometer (Promega). Data were presented as the ratio of Fleffi activity to Renilla activity.

Statistical analyses

GraphPad Prism version 8.4.3 (GraphPad Software, San Diego, CA) was used for statistical analysis. Values are expressed as means ± SD. The Student’s unpaired two-tailed t test was applied to compare two groups. For multiple comparisons, statistical analysis was performed by one-way analysis of variance (ANOVA) and Tukey’s multiple comparisons test. The sample size was chosen on the basis of the fact that previous literature using similar methods had reported moderate effect sizes. P < 0.05 was considered significant.