Single-Cell RNA Sequencing of In Vitro Expanded Chondrocytes: MSC-Like Cells With No Evidence of Distinct Subsets

Abstract

Objective

To investigate the heterogeneity of in vitro expanded chondrocytes used for autologous chondrocyte implantation.

Methods

Human articular chondrocytes were expanded in vitro for 14 days, sorted into 86 single cells using fluorescence-activated cell sorting and subjected to single-cell RNA sequencing. Principal component, Cross R² hierarchical clustering, and differential gene expression analyses were used for data evaluation. Flow cytometry and single-cell RT-qPCR (reverse transcriptase quantitative polymerase chain reaction) was used to validate the results of the RNA sequencing data Polyclonal chondrocyte populations from the same donor were differentiated in vitro toward the osteogenic and adipogenic lineages.

Results

There was considerable variation in gene expression between individual cells, but we found no evidence for separate cell subpopulations based on principal component, hierarchical clustering, and differential gene expression analysis. Most of the cells expressed all the markers defining mesenchymal stem cells, and as polyclonal chondrocyte populations from the same donor were shown to differentiate into osteocytes and adipocytes in vitro, these cells formally qualify as mesenchymal stem cells.

Conclusions

In vitro expanded chondrocytes consist of one single population of cells with heterogeneity in gene expression between the cells. Dedifferentiated chondrocytes qualify as mesenchymal stem cells as they fulfill all the criteria suggested by the International Society for Cellular Therapy.

Keywords

single-cell RNA sequencing chondrocytes mesenchymal stem cells articular cartilage ACI

Introduction

Articular cartilage is complex, with 4 zones containing different amounts and types of extracellular matrix (ECM) and an uneven distribution of articular chondrocytes (ACs).¹ In the outermost layer, the superficial zone (SZ), the cells are flat and elongated and produce collagen fibers that are oriented in parallel to the articular surface. It has been suggested that the SZ contains progenitor or stem cells.² Below the SZ is the middle zone. Here the collagen fibrils are more randomly distributed and the ECM contains higher amounts of proteoglycans. In the deep zone (DZ), the collagen fibrils are aligned perpendicular to the joint surface and the cells are arranged into column-like structures. Below the DZ is the tide mark that marks the boundary between the DZ and the calicified zone (CZ) below. In the CZ, the cells are large and produce type X collagen, a non-fibrillar collagen. Here the ECM is calcified and integrated with the underlying subchondral bone.¹ Thus, the cells in the different cartilage zones are not the same based on cell morphology, gene expression, and secretion of ECM molecules.³

For more than 20 years autologous chondrocyte implantation (ACI) has been used to treat focal lesions of the articular cartilage.⁴ During an ACI procedure, a small piece of articular cartilage, containing cells from all zones, is harvested under arthroscopic guidance, and the ACs are isolated and expanded in culture before being transplanted into the focal lesion. This leads to clinical improvement for many patients, but the repair tissue does not only consist of the desired hyaline cartilage, but also of fibrocartilage or a mixture of fibrocartilage and hyaline cartilage.⁵ This is because ACs dedifferentiate in the course of in vitro culture, leading to synthesis of fibrocartilage-specific ECM molecules, such as COL1 and versican instead of hyaline cartilage-specific ECM molecules, such as COL2 and ACAN.^6,7 Fibrocartilage does not have the mechanical properties of articular cartilage, and in some patients the repair tissue will degrade over time resulting in recurrence of pain and the possibility of developing osteoarthritis (OA) later in life. Thus, there is a need for cells that make repair tissue consisting of hyaline cartilage rather than fibrocartilage.

Single-cell RNA sequencing has recently been used to identity subpopulation of cells, with distinct gene expression pattern in cell cultures and after isolation of cells from tissues.^8-10 During an ACI procedure millions of cells are transplanted into the lesion. This cell population is derived from cells that were initially isolated from all the cartilage zones. However, it is not currently known if the expanded cells consist of one homogenous group of cells or of several subpopulations of cells with varying potential for cartilage repair. Some studies suggest that human cartilage also contains stem cells or progenitor cells that may be particularly suitable in cell-based cartilage repair strategies.² In this study, we used single-cell RNA sequencing to analyze 86 individual ACs after 14 days of in vitro expansion. Principle component analysis (PCA) showed that the ACs clustered as one distinct group of cells. Clustering analyses did not reveal subsets among the ACs, although the variation in gene expression between individual ACs was considerable. Practically all the cells expressed the markers used to identify mesenchymal stem/progenitor cells (MSCs) and polyclonal populations of these cells showed differentiation potential typical of MSCs. The results presented here suggest that in vitro expanded ACs are heterogeneous with regard to gene expression, but that this heterogeneity must have explanations other than the presence of separate subsets within the ACs.

Materials and Methods

Isolation and Culture of Human ACs

ACs from femur were isolated from discarded OA-cartilage tissue after total knee replacement surgery and cultured for 14 days (passage 3) as previously described.¹¹ Cells from one donor were used in this study. The cartilage biopsy contained cells from all cartilage zones and only tissue with no macroscopic signs of OA was used. Briefly, the cartilage was minced into small pieces and digested with collagenase type XI (Sigma-Aldrich, St. Louis, MO) for 120 to 180 minutes and seeded out in culture medium consisting of DMEM/F12 (Thermo Fisher Scientific/Gibco, Waltham, MA) containing 10% human platelet lysate-rich plasma (hPLP), 100 units/mL penicillin, and 100 µg/mL streptomycin (Sigma-Aldrich). The donor provided written, informed consent. The study was approved by the Regional Committee for Medical Research Ethics, Southern Norway, Section A. All methods and experiments were performed in accordance with the relevant guidelines and regulations. hPLP was prepared as previously described with minor modifications.¹²

In Vitro Differentiation

Adipogenic Differentiation

Confluent cultures were incubated in DMEM/F12 containing 10% FBS (Thermo Fisher Scientific/Gibco), 0.5 µM 1-methyl-3 isobutylxanthine (Sigma-Aldrich), 1 µM dexamethasone (Galenpharma, GmbH, Kiel, Germany), 10 µg/mL insulin (Novo Nordisk, Copenhagen, Denmark), and 100 µM indomethacin (Dumex-Alpharma, Copenhagen, Denmark). Fresh induction medium was replaced every 3 days. After 4 weeks, differentiated cells were fixed with 4% formalin, washed with phosphate-buffered saline (PBS) without Ca²⁺ and Mg²⁺ (Gibco), and subsequently incubated for 10 minutes with Oil-Red O to visualize lipid droplets.

Osteogenic Differentiation

Cells were incubated at 3,000 cells per cm² in DMEM/F12 containing 10% FBS, 100 nM dexamethasone (Galenpharma), 10 mM β-glycerophosphate (Sigma-Aldrich), and 0.05 mM L-ascorbic acid-2-phosphate (Sigma-Aldrich). Fresh induction medium was replaced every 3 days. After 4 weeks, differentiated cells were fixed with 4% formalin and washed with PBS without Ca²⁺ and Mg²⁺ (Gibco). Mineralization of the extracellular matrix was visualized by staining with 40 mM Alizarin Red S, pH 4.2, for 5 minutes.

Fluorescence-Activated Sorting of One Cell per Well

ACs were trypsinzed using trypsin-EDTA (Sigma-Aldrich) and washed twice with PBS (Gibco). A total of 500,000 cells were resuspended in 0.5 mL PBS containing 0.5 µL of RNaseOUT Recombinant Ribonuclease Inhibitor (Life Technologies, Carlsbad, CA) and 2 µg/mL of Hoechst 34580 (Life Technologies). Single cells were then sorted into wells of a twin-tec PCR 96-well plate (Eppendorf, Hamburg, Germany) containing 5 µL of 1× TCL lysis buffer (Qiagen, Hilden, Germany) by using a FacsAria II instrument (BD Biosciences, Franklin Lakes, NJ) and the Diva 6 software (BD Biosciences, San Jose, CA). After sorting, the plate was quickly sealed, spun down at 1300 × g for 30 seconds, frozen on dry ice and stored at −80°C until RNA was isolated.

Smart-Seq2 Libraries

Smart-Seq2 version 1.01 (Illumina, San Diego, CA) libraries were prepared by the Broad Technology Labs and sequenced by the Broad Genomics Platform according to the smartseq2 protocol¹³ with minor modifications.¹⁴ Briefly, total RNA was purified using RNA-SPRI beads. Poly(A)+ mRNA was reverse transcribed to cDNA and then amplified. cDNA was subject to transposon-based fragmentation that used dual-indexing to barcode each fragment of each converted transcript with a combination of barcodes specific to each sample. Each cell was given its own combination of barcodes. Prior to sequencing the barcoded cDNA fragments were pooled. Sequencing (NextSeq500, Illumina) was carried out on one flowcell as paired-end 2 × 25 bp with an additional 8 cycles for each index. The data were aligned using Tophat version 2.0.10,¹⁵ and transcripts were quantified by the Broad Technology Labs computational pipeline using the Cufflinks suite version 2.2.1.¹⁶

Single-Cell RNA Reverse Transcriptase Quantitative Polymerase Chain Reaction (RT-qPCR) Analysis

Single cells were sorted into wells of a twin-tec PCR 96-well plate (Eppendorf) as described above and analyzed following the protocol and parameters in the Single Cell-to-CT qRT-PCR kit (Thermo Fisher Scientific). The following Taqman gene expression assays were used: COL1A1 (Hs00164004_m1), COL3A1 (Hs0094389_m1), B2M (Hs99999907_m1) (Thermo Fisher Scientific).

Flow Cytometry Analysis

Cells were trypsinized, washed with PBS, and stained with CD44-FITC (Southern Biotech, Birmingham, AL), CD90-FITC (Molecular Probes, Eugene, OR), HLA-ABC-APC (BD Biosciences, USA), HLA-DR-APC (eBioscience), CD105-FITC (Sigma-Aldrich), CD14-FITC (Sigma-Aldrich), CD73-PE (BD Biosciences), CD34-FITC (BD Biosciences), CD45-FITC (eBioscience), and CD19-APC (BD Biosciences) for 20 minutes at room temperature before being washed twice and resuspended in PBS. Mouse isotype-FITC/PE/APC antibodies (eBioscience) were used as controls. The cells were analyzed on the BD LSRFortessa cell analyzer (BD Biosciences) and FlowJo (V10, Flowjo LLC, Ashland, OR) for image acquisition.

Bioinformatics Analysis

Cufflinks provided 2 tables with quantified transcripts that were used for analysis: one table containing FPKM values and one table containing gene counts. Genes with average FPKM and gene counts below 10 was filtered away before analysis. PCA analysis was performed in Biojupies¹⁷ using the gene count table, and t-SNE analysis was performed using the SeqMonk software (Babraham Bioinformatics, UK). Cross R² values were calculated in R using fitting linear models function in stats package. Hierarchical clustering, using Pearson correlation, and average linkage clustering and differential expression analysis, using the LIMMA package, were performed using the Multiple Experiment Viewer (MEV 4.9.0). Graphpad Prism version 5 was used to plot FPKM values for chosen genes in Figure 1C .

Figure 1.

PCA and hierarchical analysis. (A) PCA analysis of ACs. (B) Hierarchical clustering of all genes detected in the ACs. The top of the heatmap was chosen for the main image, and again for the amplified image. The colored numbers above the color scale on the top represent log2 FPKM expression values. Cell numbers (1-86) are shown on the top of the dendogram while gene names are listed vertically on the right. (C) A plot of all 86 cells against FPKM expression values for the housekeeping genes GAPDH, B2M, and ACTB.

The single-cell RNA sequencing data have been uploaded under the NCBI SRA BioProject ID: PRJNA522964.

Results

Clustering Analysis Reveals Heterogeneous Gene Expression, but No Evidence of Subsets

In order to investigate the heterogeneity of the ACs, PCA and hierarchical clustering analysis were performed. PCA analysis showed that the cells clustered as a single population of cells, with no evidence of separate subpopulations ( Fig. 1A ). A t-SNE clustering analysis was also performed, but no distinct clusters of cells were observed (not shown). Figure 1B shows the upper part of a heatmap derived from hierarchical clustering analysis of all cells, with the very uppermost part amplified to facilitate visualization. The hierarchical cluster analysis is based on the expression of all genes. The cluster tree on the top (dendogram) indicated different clusters, but the graphical representation (heatmap) did not indicate any obvious clusters. It should be noted that there will always be cluster trees on top of heatmaps, even when there are no true subsets of cells. When differential gene expression analysis was performed between the 2 most separated clusters (cells below the orange and blue lines in Fig. 1B ), no genes were differently expressed (not shown). This suggests that there were no systematic differences in gene expression between the different clusters. Although these analyses pointed away from the existence of separate subsets of ACs, a close look at the heatmap revealed that the expression of many genes varied considerably between the cells ( Fig. 1B ). This variability was investigated in 2 ways. First, reasoning that frequently used housekeeping genes should be stably expressed, we investigated the expression of glyceraldehyde-3-phosphate dehydrogenase (GAPDH), beta-2-microglobulin (B2M), and actin beta (ACTB). The results, shown in Figure 1C , found this to be only partially true. For GAPDH, the cell with the highest level of expression showed 2.9-fold higher expression than the cell with the lowest level of expression. For B2M and ACTB the differences between the highest and lowest expressing cells were 4.9-fold and 30-fold, respectively. Second, to determine how different the cells were to each other based on the expression of all genes, a cross R² analysis was used. Here, all the cells are compared with each other and R² values = 1.0 indicates identical gene expression profiles, while gradually lower values indicate gradually greater differences in gene expression profiles. In Figure 2 , the right bar = 1 and represents each cell compared with itself. When all the cells were compared against each other the mean R² value was 0.25, showing that the individual cells had very different gene expression profiles ( Fig. 2A ). When the analysis was repeated, but restricted to highly expressed genes, the mean R² value increased to 0.6, but the spread of the R² values was much greater ( Fig. 2B ).

Figure 2.

Cross R² analysis. Gene expression comparison of all cells against each other represented by R² values. (A) All detected genes, (B) Genes with average FPKM > 500.

Recently, the circadian rhythm pathway has been shown to control cartilage homeostasis and integrity and also to control the expression of a large proportion of genes in ACs undergoing synchronized culture in vitro.^18,19 As lack of synchronization could potentially explain some of the variability in gene expression between the cells studied here, we examined the expression of some of the genes most essentially involved in the positive and negative loop regulating circadian rhythm. A hierarchical clustering analysis was performed ( Fig. 3 ). No distinct clusters of cells were found. This lack of rhythmic synchronization may therefore explain some of the variability in gene expression between the cells.

Figure 3.

Circadian rhythm. A heatmap of genes associated with circadian rhythm. The numbers above the color scale on the top represent log2 FPKM expression values.

Human ACs Fulfil the Criteria for Being MSCs

Stem or progenitor cells have been identified in articular cartilage in several studies based on expression of different stem cell markers and telomerase and proliferation assays.^2,20-24 To investigate the expression of stem cell markers in the present dataset, we first sought the genes coding for the cell surface markers defining the minimal criteria for MSCs ( Fig. 4A ).²⁵ Again there were no distinct clusters of cells. THY1 (CD90) was expressed at high levels in all the cells, and NT5E (CD73) and ENG (CD105) were also expressed in practically all the cells. No cells expressed PTPRC (CD45), 2 cells expressed low levels of CD34, some cells expressed CD14, 1 cell expressed CD19, and 2 cells expressed HLA-DRA. All these markers, in addition to HLA-ABC and CD44 that are also expressed by MSCs, were validated at the protein level using flow cytometry ( Fig. 4B ). Another criterion defining MSCs is plastic adherence, also fulfilled by all in vitro expanded ACs. Finally, MSCs must be able to differentiate to osteoblasts, adipocytes, and chondrocytes in vitro. For obvious reasons this cannot be tested at the single-cell level. However, polyclonal populations of the same cells as were examined by single-cell sequencing were differentiated into adipocytes and osteocytes for 4 weeks and stained with Oil Red O to visualize lipid droplets and Alizarin Red to visualize calcium deposition. The dedifferentiated ACs had both adipogenic and osteogenic differentiation properties ( Fig. 4C ). As the cells were inherently chondrocytes, this shows that the cultured ACs had trilineage differentiation capability. Thus, by the minimal criteria for defining MSCs, the vast majority of in vitro expanded ACs are MSCs.

Figure 4.

Stem cells markers and in vitro differentiation. (A) A heatmap of genes used to define mesenchymal stem cells. Numbers above the color scale on the top represent log2 FPKM expression values. (B) Flow cytometry of CD markers. (C) Osteogenic and adipogenic differentiation of ACs. (D) Heatmap of other stem cells markers. Numbers above the color scale on top indicate log2 FPKM expression values.

A number of other markers have been used to define MSCs. Binding of the STRO-1 antibody has been used as an important marker.²⁶ The receptor for this antibody is heat shock cognate 70, encoded by the gene HSPA8, which was expressed on all the ACs ( Fig. 4A ).²⁷ CD271 has been identified as perhaps the best marker for prospective isolation of MSCs from bone marrow aspirates. However, the prospectively isolated CD271⁺ cells lose CD271 expression on in vitro culture.²⁸ NGFR, the gene encoding CD271, was not expressed. In mice, Gremlin 1 has been shown to identify a skeletal stem cell with bone, cartilage, and reticular stromal potential.²⁹ GREM1 was expressed in almost half of the ACs. Also in mice, expression of leptin receptor has been shown to identify MSCs with trilineage differentiation capability in vivo.³⁰ A large proportion of the ACs expressed LEPR, and at least 15 of the cells co-expressed GREM1 and LEPR. The expression of a number of other genes encoding MSC markers can be seen in Figure 4A , and the expression of markers of other stem cell subsets is shown in Figure 4D .

Expression of Genes Encoding ECM Molecules and Soluble Factors

The expression of genes important for cartilage ECM synthesis and functionality is shown in Figure 5A , and those encoding relevant growth factors and cytokines in Figure 5B . Again the heatmap did not reveal any distinct group of cells. However, the ACs show gene expression typical of dedifferentiated chondrocytes, with practically all cells co-expressing COL1A1, COL1A2, and COL3A1 while COL2A1 was expressed in only one cell. The chondrogenic transcription factors SOX5, SOX6, and SOX9 were expressed in a variable number of cells, but were rarely co-expressed. The genes encoding the ECM degrading enzymes MMP1, MMP3, and MMP13 were expressed in a number of cells, sometimes at high levels, and these were not infrequently co-expressed. ADAMTS5 was expressed in almost all the cells. Interestingly, while the growth factor thought to be essential for the development of articular cartilage from interzone cells, GDF5, was expressed at moderate to high levels in the vast majority of the cells, its receptor, BMPR1B, was hardly expressed in any of the cells. Finally, while none of the cells expressed the genes encoding IL1β and TNFα, most of the cells expressed genes encoding their receptors. To validate some of these results, we analyzed expression of the highly expressed genes COL1A1 and COL3A1, in addition to the housekeeping gene B2M, using single-cell RT-qPCR. Similar to the sequencing data there were considerable variation in in gene expression between individual cells ( Fig. 5C ).

Figure 5.

Cartilage and OA associated genes. (A) A heatmap of genes encoding extracellular matrix proteins and transcription factors. (B) A heatmap of genes encoding cytokines and inflammatory molecules. Numbers above the color scale on top indicate log2 FPKM expression values. (C) Single-cell RT-qPCR analysis of COL1A1, COL3A1, and B2M.

Discussion

In order to obtain enough cells for a good clinical outcome from ACI, the ACs need to be expanded over many population doublings in vitro. It has long been recognized that proliferation and plastic adherence leads to dedifferentiation of the ACs.^6,7 The dedifferentiation process may be delayed by allowing the ACs to proliferate as rounded cells attached to its own pericellular matrix,³¹ but as soon as the cells become plastic adherent, dedifferentiation is inevitable. As dedifferentiated ACs will produce inferior repair cartilage, several research groups have searched for subpopulations among the culture expanded ACs, which might produce better articular cartilage on implantation. This has led to the claim that stem cells exist within articular cartilage and that these cells may be identified among the cultured cells when ACs are expanded in vitro.^2,20-23 To determine if global gene expression analysis at the single-cell level could help us identify subpopulations of chondrocytes, stem cells, or progenitor cells within the cultured AC population, we have performed RNA sequencing on 86 individual ACs collected at the stage where the ACs would be suitable for implantation. The data show that most of these cells qualify as MSCs, but no subset of cells expressing stem cell genes over and above the others could be identified. This does not rule out the presence of stem cells among the ACs—for that, a much larger population of cells would have to be examined—but our data show that if a stem cell population is present, it is very small.

Since the advent of single-cell RNA sequencing procedures, it has become clear that cells from the same tissue, belonging to the same population, vary considerably in their expression of a large number of genes.^8-10 We show here that in vitro expanded chondrocytes consist of one single population of cells with considerable variability between the cells with regard to expression of a large number of genes, including the 3 housekeeping genes GAPDH, B2M, and ACTB. A relatively small variation of GAPDH expression, 2.9-fold differences between the cell with the lowest and highest expression, was observed. For B2M and ACTB, the differences were 4.9-fold and 30-fold, respectively. This suggests that all the cells had the same requirement for energy consumption (GAPDH), while the 30-fold difference in ACTB expression may reflect differences in cell shape and motility during cell division as ACTB is involved in cell membrane architecture and cell mobility. The differences for B2M may be due to stochastic variation or other unknown factors.

The molecular basis for heterogeneity in gene expression between cells within single populations is not fully known. Some factors act at the level of gene transcription. Here, closely related cells may show biallelic or monoallelic expression, and the level of expression may be differentially regulated on the 2 alleles.³² Also, many genes are transcribed in bursts, adding to the likelihood of variability in gene expression between cells.³³ Then there are stochastic effects occurring during the processes of RNA processing and decay.³⁴ On an organismic level, circadian rhythms have been shown to lead to changes in the expression of genes in many organs.³⁵ This has been shown to be the case also for articular cartilage.¹⁹ More surprisingly, perhaps, this mechanism was shown to be of importance also in cultured human OA chondrocytes.¹⁸ The rhythmicity of core clock genes NR1D1 and BMAL1 could be synchronized in populations of ACs by addition of dexamethasone. For our 86 cells these core genes and other genes of importance for circadian rhythms were not synchronized, adding to the variability in gene expression between these cells. However, if we accept that these cultured cells, with such variability in gene expression actually belong to one population, we still have to account for the fact that these cells were initially isolated from different layers of native articular cartilage containing cells with different morphology, gene expression, and production of matrix molecules.³

The variability between chondrocytes in native cartilage was again documented in a recent study, where single-cell RNA sequencing was performed on uncultured ACs from the different cartilage layers. In this study, 7 subpopulations of cells were identified after sequencing of 1464 cells from 10 patients at different stages of OA.³⁶ At first sight, this is in stark contrast to the results presented here. The differences may be explained in many ways: by differences in sequencing strategy, by some differences in choice of bioinformatic analyses, and by the fact that the Ji et al. data were based on the cumulated results from cells from 10 patients, while we have sequenced cells from a single patient. However, the most important explanation is probably the fact that Ji et al. sequenced uncultured cells, while we sequenced cultured cells ready for implantation.

When ACs dedifferentiate they lose the genotype and phenotype typical of specialized hyaline ACs and obtain characteristics more typical for progenitor cells. In fact, we show here that practically all the cultured ACs satisfied all the criteria for MSCs: by marker expression, adherence to plastic, and trilineage differentiation. Others have shown that dedifferentiated ACs and cultured bone marrow MSCs (BM-MSCs) have quite similar secretory profiles, and that dedifferentiated ACs are basically identical to BM-MSC when comparing expression of 50 CD markers and differentiation potential.^37,38

MSCs and dedifferentiated ACs also behave similarly under chondrogenic differentiation conditions in hydrogel scaffolds.³⁹ In fact, one of the few important differences found in that study was for type X collagen, where mRNA levels were 50 to 100 times higher in MSCs, and practically all the differentiated MSCs produced the protein, while only a few percent of the redifferentiated ACs did.³⁹ This latter observation was supported in the present study, where 4 out of the 86 cells expressed COL10A1. This may suggest that, in ACs, the COL10A1 locus is epigenetically closed in all but a few of the cells, perhaps cells from the calcified zone known to synthesize type X collagen in vivo, while in MSCs the locus is opened during differentiation in all the cells, most likely by changes in histone modifications.⁴⁰ Thus, if uncultured ACs may be statistically subdivided into subpopulations,³⁶ it is likely that they lose their subpopulation characteristics in the course of cell culture, leaving a dedifferentiated, MSC-like population of cells without subdivisions. If these cells retain their subpopulation characteristics after implantation remains to be investigated, but it seems unlikely, as transplanted ACs do not tend to redifferentiate to hyaline matrix-producing cells following ACI.

The ACs used in this study were isolated from a biopsy containing cells from all cartilage zones. Except for perhaps COL10A1 we did not detect clear zonal differences between the cells. The reason for this may be explained by the fact that we only analyzed 86 cells. If cells from the superficial zone was in abundance and/or proliferated faster than cells from the other zones they would quickly outgrow cells from other zones. Thus, there is a chance that most of the cells originated from the superficial zone as they all expressed high levels of COL1A1, which is expressed by cells in the superficial zone. Another explanation could be the effect of dedifferentiation. This is unknown, but perhaps the dedifferentiation process forced all the cells to express the same gene expression profile independent of zonal origin.

In conclusion, single-cell RNA sequencing of cultured ACs at a population doubling number similar to cells used for ACI show that these cells all qualify as MSCs, with no evidence of separate cell clusters, but with considerable heterogeneity in gene expression between the cells. We argue that this heterogeneity is caused mainly by stochastic influences on factors regulating expression of genes, perhaps also by lack of synchronization of circadian rhythms. They all express very high levels of fibrocartilage collagen genes, which may be disconcerting for subset selection strategies to improve the histological and clinical outcome of ACI.

Footnotes

Author Contributions

Tommy A. Karlsen: Conception and design, analysis and interpretation of the data, drafting of the article, critical revision of the article for important intellectual content, final approval of the article, and collection and assembly of data.

Arvind Y. M. Sundaram: Analysis and interpretation of the data, drafting of the article, critical revision of the article for important intellectual content, final approval of the article, and statistical expertise.

Jan E. Brinchmann: Conception and design, analysis and interpretation of the data, drafting of the article, critical revision of the article for important intellectual content, final approval of the article, provision of study materials or patients, obtaining of funding, administrative, technical, or logistic support, and collection and assembly of data.

Acknowledgments and Funding

We would like to acknowledge Yan Zhang and Anna Bujko, Department of Pathology, Oslo University Hospital, for help with sorting of single cells into 96-well plates; Hege Brincker Fjerdingstad, The National Core Facility for Human Pluripotent Stem Cells at the Norwegian Center for Stem Cell Research, Oslo University Hospital, for technical assistance with flow cytometry; and Bjørn Ødegaard, Department of Orthopedic Surgery, Lovisenberg Diakonale Hospital, Oslo, for providing the cartilage tissue sample. The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This study was funded by a grant from the South-Eastern Norway Regional Health Authority, Grant Number: 2013077. The funding source (Southern and Eastern Norway Regional Health Authority) was not involved in the study design, collection and interpretation of data, writing of the manuscript, or in the decision to submit the manuscript.

Declaration of Conflicting Interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Ethical Approval

Ethical approval for this study was obtained from the Regional Commitee for Medical Research Ethics, Southern Norway, Section A. All methods and experiments were performed in accordance with the relvant guidelines and regulations.

Informed Consent

Written informed consent was obtained from all subjects before the study.

ORCID iD

Tommy A. Karlsen

References

Poole

Kojima

Yasuda

Mwale

Kobayashi

Laverty

Composition and structure of articular cartilage: a template for tissue repair. Clin Orthop Relat Res. 2001;(391 suppl):S26-S33.

Jiang

Tuan

RS.

Origin and function of cartilage stem/progenitor cells in osteoarthritis. Nat Rev Rheumatol. 2015;11:206-12. doi:10.1038/nrrheum.2014.200

Grogan

Duffy

Pauli

Koziol

D’Lima

, et al. Zone-specific gene expression patterns in articular cartilage. Arthritis Rheum. 2013;65:418-28. doi:10.1002/art.37760

Brittberg

Lindahl

Nilsson

Ohlsson

Isaksson

Peterson

Treatment of deep cartilage defects in the knee with autologous chondrocyte transplantation. N Engl J Med. 1994;331:889-95.

Knutsen

Engebretsen

Ludvigsen

Drogset

Grøntvedt

Solheim

, et al. Autologous chondrocyte implantation compared with microfracture in the knee. A randomized trial. J Bone Joint Surg Am. 2004;86-A:455-64.

Binette

McQuaid

Haudenschild

Yaeger

McPherson

Tubo

Expression of a stable articular cartilage phenotype without evidence of hypertrophy by adult human articular chondrocytes in vitro. J Orthop Res. 1998;16:207-16. doi:10.1002/jor.1100160208

Holtzer

Abbott

Lash

Holtzer

. The loss of phenotypic traits by differentiated cells in vitro, I. dedifferentiation of cartilage cells. Proc Natl Acad Sci U S A. 1960;46:1533-42.

Athanasiadis

Botthof

Andres

Ferreira

Lio

Cvejic

Single-cell RNA-sequencing uncovers transcriptional states and fate decisions in haematopoiesis. Nat Commun. 2017;8:2045. doi:10.1038/s41467-017-02305-6

Freeman

Jung

Ogle

BM.

Single-cell RNA-seq of bone marrow-derived mesenchymal stem cells reveals unique profiles of lineage priming. PLoS One. 2015;10:e0136199. doi:10.1371/journal.pone.0136199

10.

Glass

Calero-Nieto

Jawaid

Larraufie

Kay

Göttgens

, et al. Single-cell RNA-sequencing reveals a distinct population of proglucagon-expressing cells specific to the mouse upper small intestine. Mol Metab. 2017;6:1296-303. doi:10.1016/j.molmet.2017.07.014

11.

Karlsen

Jakobsen

Mikkelsen

Brinchmann

JE.

microRNA-140 targets RALA and regulates chondrogenic differentiation of human mesenchymal stem cells by translational enhancement of SOX9 and ACAN. Stem Cells Dev. 2014;23:290-304. doi:10.1089/scd.2013.0209

12.

Schallmoser

Strunk

Preparation of pooled human platelet lysate (pHPL) as an efficient supplement for animal serum-free human stem cell cultures. J Vis Exp. 2009;(30):1523. doi:10.3791/1523

13.

Picelli

Björklund

ÅK

Faridani

Sagasser

Winberg

Sandberg

Smart-seq2 for sensitive full-length transcriptome profiling in single cells. Nat Methods. 2013;10:1096-8. doi:10.1038/nmeth.2639

14.

Trombetta

Gennert

Satija

Shalek

Regev

Preparation of single-cell RNA-seq libraries for next generation sequencing. Curr Protoc Mol Biol. 2014;107:4.22.21-17. doi:10.1002/0471142727.mb0422s107

15.

Kim

Pertea

Trapnell

Pimentel

Kelley

Salzberg

SL.

TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14:R36. doi:10.1186/gb-2013-14-4-r36

16.

Trapnell

Roberts

Goff

Pertea

Kim

Kelley

, et al. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc. 2012;7:562-78. doi:10.1038/nprot.2012.016

17.

Torre

Lachmann

Ma’ayan

BioJupies: automated generation of interactive notebooks for RNA-seq data analysis in the cloud. Cell Syst. 2018;7:556-61. doi:10.1016/j.cels.2018.10.007

18.

Akagi

Akatsu

Fisch

Alvarez-Garcia

Teramura

Muramatsu

, et al. Dysregulated circadian rhythm pathway in human osteoarthritis: NR1D1 and BMAL1 suppression alters TGF-β signaling in chondrocytes. Osteoarthritis Cartilage. 2017;25:943-51. doi:10.1016/j.joca.2016.11.007

19.

Dudek

Gossan

Yang

Ruckshanthi

Yoshitane

, et al. The chondrocyte clock gene Bmal1 controls cartilage homeostasis and integrity. J Clin Invest. 2016;126:365-76. doi:10.1172/JCI82755

20.

Grogan

Miyaki

Asahara

D’Lima

Lotz

MK.

Mesenchymal progenitor cell markers in human articular cartilage: normal distribution and changes in osteoarthritis. Arthritis Res Ther. 2009;11:R85. doi:10.1186/ar2719

21.

Hiraoka

Grogan

Olee

Lotz

Mesenchymal progenitor cells in adult human articular cartilage. Biorheology. 2006;43:447-54.

22.

Khan

Bishop

Gilbert

Archer

CW.

Clonal chondroprogenitors maintain telomerase activity and Sox9 expression during extended monolayer culture and retain chondrogenic potential. Osteoarthritis Cartilage. 2009;17:518-28. doi:10.1016/j.joca.2008.08.002

23.

Williams

Khan

Richardson

Nelson

McCarthy

Analbelsi

, et al. Identification and clonal characterisation of a progenitor cell sub-population in normal human articular cartilage. PLoS One. 2010;5:e13246. doi:10.1371/journal.pone.0013246

24.

Koelling

Kruegel

Irmer

Path

Sadowski

Miro

, et al. Migratory chondrogenic progenitor cells from repair tissue during the later stages of human osteoarthritis. Cell Stem Cell. 2009;4:324-35. doi:10.1016/j.stem.2009.01.015

25.

Dominici

Le Blanc

Mueller

Slaper-Cortenbach

Marini

Krause

, et al. Minimal criteria for defining multipotent mesenchymal stromal cells. The International Society for Cellular Therapy Position Statement. Cytotherapy. 2006;8:315-7. doi:10.1080/14653240600855905

26.

Simmons

Torok-Storb

Identification of stromal cell precursors in human bone marrow by a novel monoclonal antibody, STRO-1. Blood. 1991;78:55-62.

27.

Fitter

Gronthos

Ooi

Zannettino

AC.

The mesenchymal precursor cell marker antibody STRO-1 binds to cell surface heat shock cognate 70. Stem Cells. 2017;35:940-51. doi:10.1002/stem.2560

28.

Ghazanfari

Zacharaki

Lim

Scheding

Human non-hematopoietic CD271^pos/CD140a^low/neg bone marrow stroma cells fulfill stringent stem cell criteria in serial transplantations. Stem Cells Dev. 2016;25:1652-58. doi:10.1089/scd.2016.0169

29.

Worthley

Churchill

Compton

Tailor

Rao

, et al. Gremlin 1 identifies a skeletal stem cell with bone, cartilage, and reticular stromal potential. Cell. 2015;160:269-84. doi:10.1016/j.cell.2014.11.042

30.

Zhou

Yue

Murphy

Peyer

Morrison

SJ.

Leptin-receptor-expressing mesenchymal stromal cells represent the main source of bone formed by adult bone marrow. Cell Stem Cell. 2014;15:154-68. doi:10.1016/j.stem.2014.06.008

31.

Shahdadfar

Løken

Dahl

Tunheim

Collas

Reinholt

, et al. Persistence of collagen type II synthesis and secretion in rapidly proliferating human articular chondrocytes in vitro. Tissue Eng Part A. 2008;14:1999-2007.

32.

Deng

Ramskold

Reinius

Sandberg

Single-cell RNA-seq reveals dynamic, random monoallelic gene expression in mammalian cells. Science. 2014;343:193-6. doi:10.1126/science.1245316

33.

Kumar

Singh

Kulkarni

RV.

Transcriptional bursting in gene expression: analytical results for general stochastic models. PLoS Comput Biol. 2015;11:e1004292. doi:10.1371/journal.pcbi.1004292

34.

Patange

Girvan

Larson

DR.

Single-cell systems biology: probing the basic unit of information flow. Curr Opin Syst Biol. 2018;8:7-15. doi:10.1016/j.coisb.2017.11.011

35.

Andreani

Itoh

Yildirim

Hwangbo

Allada

Genetics of circadian rhythms. Sleep Med Clin. 2015;10:413-21. doi:10.1016/j.jsmc.2015.08.007

36.

Zheng

Zhang

Fan

Hou

, et al. Single-cell RNA-seq analysis reveals the progression of human osteoarthritis. Ann Rheum Dis. 2019;78:100-10. doi:10.1136/annrheumdis-2017-212863

37.

de la Fuente

Abad

García-Castro

Fernández-Miguel

Petriz

Rubio

, et al. Dedifferentiated adult articular chondrocytes: a population of human multipotent primitive cells. Exp Cell Res. 2004;297:313-28. doi:10.1016/j.yexcr.2004.02.026

38.

Polacek

Bruun

Elvenes

Figenschau

Martinez

The secretory profiles of cultured human articular chondrocytes and mesenchymal stem cells: implications for autologous cell transplantation strategies. Cell Transplant. 2011;20:1381-93. doi:10.3727/096368910X550215

39.

Fernandes

Herlofsen

Karlsen

Küchler

Fløisand

Brinchmann

JE.

Similar properties of chondrocytes from osteoarthritis joints and mesenchymal stem cells from healthy donors for tissue engineering of articular cartilage. PLoS One. 2013;8:e62994. doi:10.1371/journal.pone.0062994

40.

Herlofsen

Bryne

Høiby

Wang

Issner

Zhang

, et al. Genome-wide map of quantified epigenetic changes during in vitro chondrogenic differentiation of primary human mesenchymal stem cells. BMC Genomics. 2013;14:105. doi:10.1186/1471-2164-14-105