Single-Cell RNA Sequencing Reveals the Migration of Osteoclasts in Giant Cell Tumor of Bone

Giant cell tumor of bone (GCTB) is benign tumor that can cause significant osteolysis and bone destruction at the epiphysis of long bones. Osteoclasts are thought to be highly associated with osteolysis in GCTB. However, the migration of osteoclasts in GCTB remains unclear. A deeper understanding of the complex tumor microenvironment is required in order to delineate the migration of osteoclasts in GCTB. In this study, samples were isolated from one patient diagnosed with GCTB. Single-cell RNA sequencing (scRNA-seq) was used to detect the heterogeneity of GCTB. Multiplex immunofluorescence staining was used to evaluate the cell subtypes identified by scRNA-seq. A total of 8,033 cells were obtained from one patient diagnosed with GCTB, which were divided into eight major cell types as depicted by a single-cell transcriptional map. The osteoclasts were divided into three subsets, and their differentiation trajectory and migration status were further analyzed. Osteoclast migration may be regulated via a series of genes associated with cell migration. Furthermore, four signaling pathways (RANKL, PARs, CD137 and SMEA3 signaling pathway) were found to be highly associated with osteoclast migration. This comprehensive single-cell transcriptome analysis of GCTB identified a series of genes associated with cell migration as well as four major signaling pathways that were highly related to the migration of osteoclasts in GCTB. Our findings broaden the understanding of GCTB bionetworks and provides a theoretical basis for anti-osteolysis therapy against GCTB in the future.

Introduction

Giant cell tumor of bone (GCTB) is destructive osteolytic benign tumor that often affect the epiphysis of long bones and can lead to severe motor dysfunction (1). GCTB usually occurs in young adults between the ages of 30 and 40 years (2). The main pathologic feature of GCTB is severe osteolysis that is thought to be caused by osteoclasts, which are multinucleated cells derived from the monocyte/macrophage lineage (3). Surgical resection is the main treatment for GCTB. However, repetitive operations after the primary surgery due to local recurrence can lead to serious functional complications (4). The use of drugs designed to inhibit the function of osteoclasts, such as denosumab or zoledronic acid, can contribute to symptomatic relief but do not inhibit the development of GCTB (5). An in-depth understanding of the mechanism behind osteoclast-mediated osteolysis by GCTB is therefore of critical importance.

The combination of receptor activator of nuclear factor-κ B (RANK) and its ligand (RANKL) regulates osteoclast activity, which involves osteoclast attachment to the bone surface and digestion of the bone matrix via the secretion of acid (H⁺ and Cl^-) and proteinases (cathepsin K) (6). These phenomena are not only observed under normal physiologic conditions but also in the setting of pathologic bone resorption, such as by bone tumors such as GCTB (3). However, why osteoclasts migrate directly to the lesion and why osteoclast-mediated osteolysis occurs remains unknown. The answer may be related to the tumor microenvironment, in which cells directly communicate with each other through cellular contact and interaction with the extracellular matrix (7). Analysis of the unique types of GCTB tumor microenvironments can reveal key factors involved in osteoclast migration and osteoclast-mediated osteolysis that may serve as future targets of novel therapeutic strategies. We hypothesized that osteoclast migration is influenced by the complex and dynamic features of the GCTB tumor microenvironment.

Given that conventional bulk RNA-sequencing (RNA-seq) is based on the hypothesis that every gene is expressed equally by every cell, it is therefore unable to accurately characterize the heterogeneity of tumor microenvironments at the cell-type level (8). With advancement in the scRNA-seq technique, heterogeneous tissues can be delineated at the single-cell level (9). This technology permits the massively parallel characterization of thousands of cells at the transcriptome level and can better explain cell–cell interactions. However, relatively few studies have been published on the single-cell transcriptome of GCTB. Here we used scRNA-seq to investigate the intratumoral heterogeneity of GCTB. A relative comprehensive single-cell transcriptome profiling of GCTB was performed, and gene features and cellular dynamics associated with osteoclast migration were identified.

Materials and Methods

Patient and Sample Collection

This work was approved by the First Affiliated Hospital of Guangxi Medical University (No.2019KY-E-153) and complied with all relevant ethical regulations. The patient with GCTB whose clinical and cellular data was used in this study provided informed consent. The patient’s basic information, X-ray and MRI imaging, and pathology slides are shown in Supplementary Figure 1 and Supplementary Table 1. Fresh specimens acquired at the time of surgical resection were collected in cold Hank’s balanced salt solution (cat. no. 311-512-CL; Wisent Bio Products) containing 1% antibiotic-antimycotic (cat. no. 15240062; Thermo Fisher) and transported to the laboratory as soon as possible.

Preparation of the Single-Cell Suspension

Fresh tumor specimens were mechanically isolated and enzymatically digested to produce single-cell suspensions. The tumor was washed two times with Dulbecco’s phosphate-buffered saline (DPBS; cat. no. 14190250; Gibco) and minced on ice. The tumor was then subjected to enzymatic type II (cat. no. 17101015; Thermo Fisher) digestion in a water bath at 37 ° to generate GCTB cell suspensions. A 100-μm nylon cell strainer (cat. no. 352340; Falcon) was used to filter impurities. Red blood cells (RBC) were lysed with RBC lysis buffer (cat. no. 1966634; Invitrogen) containing DNase I (1 unit/ml) and removed via centrifugation. The cell suspensions were refiltered with a 40 μm nylon cell strainer (cat. no. 352340; Falcon) to capture the isolated cells. The dissociated cells were stained with trypan blue (0.4%; cat. no. 420301; Thermo Fisher) to calculate cell viability and diluted with DPBS containing 1% fetal bovine serum (FBS; cat. no. 10091148; Thermo Fisher) into appropriate concentrations for the next step.

Preparation of Single-Cell Suspensions for Library Construction and scRNA-Sequencing

Single-cell suspensions of GCTB were uploaded into an emulsion droplet using 10× Genomics Chromium Controller (version 3) to generate single-cell gel bead-in-emulsions (GEMs). mRNA in drops was subjected to reverse-transcription reactions, and cDNA amplification was performed according to the manufacturer’s instructions. The 10× libraries were sequenced with the HiSeq Xten (Illumina, San Diego, CA) sequencing platform.

Pre-Processing of scRNA-Sequencing Data

CellRanger (version 4.0.0) was used to convert the preliminary sequencing results (bcl files) generated from HiSeq Xten into fastq files. The fastq files were then aligned to the human genome reference sequence GRCh38. A Seurat package (version 3.1.1) in R software (version 3.6.3, R-Foundation, Vienna, Austria) was used to generate a gene-barcode matrix containing the barcoded cells and gene expression counts. Low-quality cells (gene numbers < 200 or > 4,000 and with a percentage of mitochondrial genes of > 10%) were directly filtered. A total of 8,033 cells were ultimately included for further bioinformatic analysis.

Cell Clustering Analysis, Visualization, and Annotation

Cell-clustering and sub-clustering analyses were performed with the FindClusters function of the Seurat package (resolution = 0.1). Identified cell clusters and sub-clusters were presented with uniform manifold approximation and projection (UMAP) analysis. The cell clusters were annotated based on highly expressed genes and some canonical cellular markers from the literature.

Pseudotime Trajectory Analysis

The evolutionary processes of osteoclasts and mononuclear macrophages were analyzed using the Monocle 3 package (version 0.2.3.0; https://coletrapnelllab.github.io/monocle3/). Cells were organized into discontinuous trajectories. The identified genes varied in their expression over these trajectories. In the selection of cell sequencing parameters, the starting point of the cell was selected based on the results of RNA velocity analysis. Six genes (CD74, HLA-DRA, CD14, ACP5, CTSK, and ATP6V0D2) that play an important role in the maturation of osteoclasts were also subjected to pseudotime trajectory analysis. The evolutionary process of CD8⁺T cells was analyzed using the Monocle 2 package (version 2.14.0; http://cole-trapnell-lab.github.io/monocle-release/) (10). The following parameters were set: Mean expression >0.1, num_cells_expressed >=10. GZMK and GZMA also underwent pseudotime trajectory analysis.

Cell–Cell Communication Analysis

CellPhoneDB analysis, which is based on the interaction between ligands and receptors (11, 12), was performed using a CellPhoneDB Python package (version 0.22; https://github.com/Teichlab/cellphonedb). We also identified some relevant cell type-specific interactions based on ligand–receptor pairs (P value < 0.05).

CellChat Analysis

CellChat, an intercellular interaction analysis tool that studies ligand receptor action in specific signaling pathways (13), was performed using the CellChat R package (version 0.0.1; https://github.com/sqjin/CellChat). A CellChat object was created using the R package process. Cell information was added into the meta slot of the object. The ligand–receptor interaction database was set, and the matching receptor inference calculation was performed. The graphic visualization parameter was nPatterns = 5.

Functional Enrichment Analysis

Gene set variation analysis (GSVA) was used to identify the gene set that was significantly enriched in each subset of mononuclear macrophages (14). We downloaded all gene sets from the Molecular Signatures Database MSigDB (https://www.gseamsigdb.org/gsea/downloads.jsp).

Single-Cell Regulatory Network Inference and Clustering (SCENIC) Analysis

SCENIC is a suitable method for using transcription factors (TFs) to guide the discovery of cellular states from scRNA-seq data (15). The R package for SCENIC (version 1.2.4) was used as previously described to perform a transcription factor network inference analysis (16). Briefly, a read-count matrix from the Seruat S4 object was inputted. The matrix was filtered using default parameters and used to establish a gene regulatory network. Differentially activated TFs in different cell types were identified using the Wilcoxon rank sum test. TFs with an adjusted P-value <0.05 and logFC >0.1 were considered significant.

RNA Velocity Analysis

RNA velocity predicts dynamic changes in transcription and the future state of individual cells (17). We obtained.loom files through velocyto.py (version 0.17.17). The.loom files were uploaded to the velocyto.R package (version 0.6.0). The following parameters were set: fit.quantile = 0.02; kCells = 25; and DeltaT = 1. Arrows indicating velocity vector were projected onto the UMAP plot.

Cellular Spatial Organization Mapper (CSOmap) Analysis

CSOmap (version 1.0) is an analytical method for calculating the spatial information of cells by using the interactions between ligand receptors (18). We generated an affinity matrix between cells by integrating thousands of ligand-receptor pairs. The resulting high-dimensional affinity matrix is embedded in three-dimensional space. The contribution of ligand receptor genes to spatial information was then calculated. Finally, genes highly expressed by C2_Mature osteoclasts were selected for computer overexpression and knockdown to observe changes in spatial structure.

Multiplex Immunohistochemistry (IHC) Staining

Multiplex immunohistochemistry staining was performed according to the manufacturer’s instructions using the PANO 7-plex IHC kit (Panovue, Beijing, China) (19, 20). Sections with a thickness of 3 μm were incubated overnight at 4°C with primary antibodies: anti-ACP5 (cat. no. 11594-1-AP; Rabbit; 1:300; Proteintech), anti-ATP6V0D2 (cat. no. bs-12548R; Rabbit; 1:300; Bioss), and anti-CKLF (cat. no. ab250213; 1:300; Rabbit; Abcam). Secondary antibodies were then used to incubate the sections at room temperature for 15 min, after which the tyramide signal amplification (TSA) plus working solution was applied. Other primary antibodies were applied to the slides, and the steps mentioned above were repeated until the last antibody was used. Finally, 4-6-diamidino-2-phenylindole (DAPI; SigmaAldrich) was used to stain the nuclei and multispectral images were collected with a confocal laser-scanning microscope (LSM880; Zeiss).

Statistical Analysis

All statistical analysis and figures were generated using R software (version 3.6.3). A p-value < 0.05 was considered statistically significant.

Results

GCTB Cellular Contribution

scRNA-seq analysis was performed on samples obtained from a patient diagnosed with GCTB during tumor resection (Figure 1A). Following initial quality control assessment, the single-cell transcriptomic data of 8,033 cells from the primary GCTB lesion were used for further analysis. Eight main cell clusters were identified (Figure 1B): mononuclear macrophages, osteoblasts, NK/T cells, osteoclasts, pericytes, proliferating cells, endothelial cells, and chondrocytes (Figure 1C). The canonical markers of each cell cluster were as follows (Figure 1D): [1] LYZ, CD163, CD14, MRC1, MSR1, C1QA, and C1QB for mononuclear macrophages (20); [2] RUNX2 and IBSP for osteoblasts (21, 22); [3] CD3D, CD3E, CD3G, and NKG7 for NK/T cells (23–25); [4] CTSK and ACP5 for osteoclasts (26); [5] RGS5 and ACTA2 for pericytes (24); [6] MKI67 and TOP2A for proliferative cells (24, 27); [7] EGFL7 and PLVAP for endothelial cells (28, 29); and [8] ACAN and COL10A1 for chondrocytes (24, 30) (Figure 1E).

FIGURE 1

Figure 1 scRNA-seq clustering analysis of GCTB. (A) Workflow depicting the collection and processing of GCTB tumor specimens for scRNA-seq. (B) UMAP plot of 8,033 cells demonstrating the eight main cell types in GCTB. (C) Cell number and percentage of the assigned cell types. (D) The heat-map shows and highlights the differentially expressed genes of each cluster. (E) UMAP plots show the expression of representative well-known markers across cell types in GCTB. scRNA-seq, single-cell RNA sequencing; GCTB, giant cell tumor of bone; UMAP, uniform manifold approximation and projection.

Heterogeneity of Mononuclear Macrophages

Previous studies found that mononuclear macrophages are critical abundant components of the tumor microenvironment and have been widely implicated in tumor stimulating and suppressing activities (31, 32). Two distinct subclusters comprised of C1_Mononuclear macrophages and C2_Mononuclear macrophages were identified in the UMAP analysis (Figure 2A). C1_Mononuclear macrophages were 95.26% of all mononuclear macrophages in the specimen, while 4.74% were C2_Mononuclear macrophages (Figure 2B). C1_Mononuclear macrophages had a high expression level of HLA-related genes, indicating that they were antigen presenting (Figure 2C). C1_Mononuclear macrophages were correlated with inflammation due to their higher expression of anti- and pro-inflammatory genes, while C2_Mononuclear macrophages were less correlated with inflammation (Figure 2D). C1_Mononuclear macrophages were enriched in gene sets associated with the inflammatory response, including the inflammatory response pathway and the inflammatory response to an antigen stimulus. C1_Mononuclear macrophages may be highly involved in the GCTB inflammatory microenvironment. C2_Mononuclear macrophages were enriched in gene sets associated with an immune response, including regulation of the innate immune response and the immune response to tumor cells, suggesting that they may be immune cells that have anti-GCTB properties (Figure 2E). The trajectory and RNA velocity analysis of mononuclear macrophages showed that the differentiation trajectory began with C1_Mononuclear macrophages, which can differentiate into C2_Mononuclear macrophages (Figures 2F–H).

FIGURE 2

Figure 2 Heterogeneity of macrophage populations in GCTB. (A) C1_Mononuclear macrophages and C2_Mononuclear macrophages were identified by UMAP analysis. (B) A pie chart depicting macrophage cell composition. (C) Heatmap of genes involved in antigen presenting by C1_Mononuclear macrophages and C2_Mononuclear macrophages. (D) Heatmap of genes associated with anti- and pro-inflammation properties of C1_Mononuclear macrophages and C2_Mononuclear macrophages. (E) GSVA analysis shows the function of C1_Mononuclear macrophages and C2_Mononuclear macrophages. (F, G) The differentiation and developmental trajectory of C1_Mononuclear macrophages and C2_Mononuclear macrophages in GCTB. (H) Velocity field projected onto the UMAP plot of C1_Mononuclear macrophages and C2_Mononuclear macrophages. Arrows indicate the direction of the differentiation and the average velocity. GCTB, giant cell tumor of bone; UMAP, uniform manifold approximation and projection for dimension reduction.

Heterogeneity of Osteoclasts in GCTB

GCTB, also known as “osteoclastoma,” is thought to be highly associated with osteoclasts, which play a critical role in the osteolysis of GCTB (33). Three distinct subtypes of osteoclasts were identified according to known marker genes (Figures 3A–C). These subclusters were described as follows: (1) C1_Progenitor osteoclasts, which accounted for most osteoclasts (52.25%) and expressed high levels of myeloid cell markers such as C1QA, C1QB, HLA-DRA, major histocompatibility complex II (MHC-II; CD74), and CD14; (2) C2_Mature osteoclasts, which accounted for 35.21% of osteoclasts and demonstrated high levels of the osteoclast markers tartrate acid phosphatase (ACP5), cathepsin K (CTSK), and the D2 isoform of vacuolar (H⁺) ATPase H⁺ Transporting V0 Subunit D2 (ATP6V0D2), which are necessary for the maturation of osteoclasts (34, 35); and (3) C3_Dysfunctional osteoclasts, which accounted for 12.54% of all osteoclasts and poorly expressed genes related to osteoclast function. Osteoclasts followed a differentiation trajectory that mainly started from the initial cluster of C1_Progenitor osteoclasts, which differentiated into C2_Mature osteoclasts and C3_Dysfunctional osteoclasts (Figures 3D–F). The cell trajectory analysis of marker genes confirmed these results. CD74, HLA-DRA, and CD14 were located at the initial position of the pseudotime and the expression levels of CTSK, ACP5, and ATP6V0D2 gradually increased along the pseudotime trajectory (Figure 3G). The regulatory network underlying each subset of osteoclasts was examined with SCENIC, and specific TF regulons for each osteoclast subset were identified. As shown in C2_Mature osteoclasts, four genes of TFs (JDP2, NFATC1, MLX, and ESRRA) were significantly up-regulated and four TFs were activated (Figures 3H, I). Jun dimerization protein 2 (JDP2) is a transcription factor of the AP-1 family that regulates osteoclast differentiation by RANKL (36, 37). The nuclear factor of activated T cells, cytoplasmic 1 (NFATc1), which belongs to the NFAT family, plays a vital role in osteoclast formation and function (38). MLX can promote myogenesis by inducing the expression of several myokines, including insulin-like growth factor 2 (IGF2) (39, 40). ESRRA is a transcription factor that is involved in tumorigenesis, such as oral squamous cell carcinoma (41).

FIGURE 3

Figure 3 Heterogeneity of osteoclasts in GCTB lesions. (A) An UMAP plot shows three main subtypes of osteoclasts: C1_Progenitor osteoclasts, C2_Mature osteoclasts, and C3_Dysfunctional osteoclasts. (B) A pie chart depicts the cell composition of the osteoclasts. (C) Violin plots showing relevant marker genes in the osteoclast subtypes (**P < 0.01, ***P <0.001, NS: No statistical significance). (D, E) The dynamics of the osteoclast subclusters showed with a monocle 3 trajectory plot. (F) Velocity field projected onto the UMAP plot of the osteoclasts. Arrows indicate the direction of the differentiation and its average velocity. (G) The expression levels of marker genes (CD74, HLA-DRA, CD14, ACP5, CTSK, and ATP6V0D2) related to the differentiation of osteoclasts. (H) Heatmap of the AUC scores of transcription factor expression regulation by SCENIC. (I) UMAP plot of osteoclasts color-coded for expression level (up) and the AUC of the estimated regulon activity of these transcription factors (down). GCTB, giant cell tumor of bone; UMAP, uniform manifold approximation and projection for dimension reduction; SCENIC, single-cell regulatory network inference and clustering; AUC, area under the curve.

Spatial information of osteoclasts and other cells (mononuclear macrophages, osteoblasts, NK/T cells, pericytes, proliferating cells, endothelial cells, and chondrocytes) was obtained with CSOmap analysis (Figure 4A). Four pairs of ligand receptors contributed more than 10% to the spatial reconstruction of the cell (Figure 4B). Based on the analysis of the expression levels of these ligand receptor pairs (MMP9-CD44, MMP9-LRP1, TIMP1-CD63, RPS19-C5AR1), it was found that CD63 and MMP9 were highly expressed in C2_mature osteoclasts (Figure 4C). Through computer simulation of CD63 overexpression osteoclasts were found to have a closer spatial structure, suggesting that they were in a state of aggregation at this time (Figure 4D). When CD63 was knocked down, osteoclasts had a looser spatial structure (Figure 4E), suggesting a state of dispersion. We quantified the distance between the osteoclast cells and the pseudo-space center to further confirm the visual features (Figure 4H). Similarly, through computer simulation the overexpression of MMP9 led to a closer osteoclast spatial structure (Figure 4F), while knockdown of MMP9 led to a looser spatial structure (Figure 4G). However, the change in the distance between osteoclast cells and the pseudo-space center after MMP9 knock-down was not statistically significant (Figure 4I).

FIGURE 4

Figure 4 CSOmap revealing osteoclast spatial information. (A) Global (left, lower right) and cross-sectional (upper right) views of the spatial information of osteoclasts and other cells. (B) Contributions of all ligand-receptor genes to the interactions between osteoclasts. (C) Heatmap of the 4 pairs of ligand receptors: (MMP9-CD44, MMP9-LRP1, TIMP1-CD63, RPS19-C5AR1). (D, E) Global (left, lower right) and cross-sectional (upper right) views of the spatial information on osteoclasts and other cells after the overexpression and knockdown of CD63. (F, G) Global (left, lower right) and cross-sectional (upper right) views of the spatial information of osteoclasts and other cells after the overexpression and knockdown of MMP9. (H) The distance of osteoclasts (normal control/CD63up/down) to the center (*P < 0.05, ***P < 0.001). (I) The distance of osteoclasts (normal control/CD63up/down) to the center (***P < 0.001; NS, No statistical significance). CSOmap, Cellular Spatial Organization mapper.

Heterogeneity of NK/T Cells

NK/T cells play an essential role in mediating response to chemotherapy and improving clinical outcomes in various cancers, including liver cancer (42), ovarian cancer (43), and breast cancer (44). Five subclusters of cells were identified from the GCTB lesion: C1_CD8⁺ T cells, C2_T cells, C3_T cells, C4_NK cells and C5_ CD8⁺ T cells (Figures 5A, B). C1_CD8⁺ T cells, and C2_CD8⁺ T cells shared specific T-cell markers genes, such as CD3D, CD3E, CD3G, CD8A, and CD8B (19, 45). C2_T cells and C3_T cells expressed CD3D, CD3E, and CD3G. C4_NK cells significantly expressed NKG7, GNLY, GZMA, GZMB, and KLRD1, which are thought to be the canonical marker genes of NK cells (23, 46) (Figure 5C). Natural killer cell granule protein 7 (NKG7) is a marker of NK cells that is critical for controlling cancer initiation, growth, and metastasis. NKG7 was highly expressed in C1_CD8⁺ T cells, suggesting the abundant infiltration of NK/T cells in GCTB.

FIGURE 5

Figure 5 Heterogeneity of NK/T cells in GCTB. (A) UMAP plot showing different subtypes of NK/T cells. (B) A pie chart depicting the subtypes composition of NK/T cells. (C) Violin plots displaying the expression of representative, well-known markers across the NK/T cells types identified in GCTB (*P <0.05, **P <0.01, ***P <0.001, NS: No statistical significance). (D) UMAP plot showing different subtypes of CD8⁺ T cells. (E) A pie chart indicating the subtypes composition of CD8⁺ T cells. (F) Dot-plot of marker genes in subtypes of CD8⁺ T cells. Shades of red represent the expression level, and the dot sizes represent the relative abundance of each gene. (G, H) The dynamics of CD8⁺ T cell subtype as shown by the monocle 2 trajectory plot. (I) An RNA velocity field projected onto the UMAP plot of the subtype of CD8⁺ T cells. Arrows indicate the direction of the differentiation and its average velocity. (J, K) The expression levels of the marker genes GZMK and GZMA, which are related to the differentiation of a subtype of CD8⁺ T cells. (L) Heatmap of the AUC scores of expression regulation by transcription factors as estimated by SCENIC. (M) UMAP plot of CD8⁺ T cell subtype, the expression level (up) and for the AUC of the estimated regulon activity of these transcription factors (down). GCTB, giant cell tumor of bone; UMAP, uniform manifold approximation and projection for dimension reduction; NK cells, natural killer cells; GZMK, Granzyme K; GZMA, Granzyme A.

Three distinct subtypes of CD8⁺ T cells were identified according to known marker genes (Figures 5D–F). These subclusters were as follows: (1) subC1_CD8⁺ T cells accounted for 44.14% of all CD8⁺ T cells and expressed low levels of markers related to T cell function such as PDCD1, CTLA4, LAG3, TIGIT, and HAVCR2. However, subC1_CD8⁺ T cells had a relatively high expression levels of CD69, GZMM, GZMK, GZMA, CST7, HLA-DPA1, HLA-DRB1 and HLA-DRA; (2) subC2_CD8⁺ T cells, which accounted for 44.14% of all CD8+ T cells and expressed the same markers as subC1_CD8⁺ T cells, but at lower levels; and (3) subC3_CD8⁺ T cells, which accounted for 11.72% of all CD8+ T cells and expressed high levels of FOSB, a marker gene of early-stage CD8⁺ T cells. CD8⁺ T cells followed a differentiation trajectory that mainly started from the initial cluster of partial subC3_CD8⁺ T cells, which differentiated into subC2_CD8⁺ T cells and subC1_CD8⁺ T cells (Figures 5G–I). The expression levels of GZMK and GZMA gradually increased along the pseudotime trajectory (Figurse 5J, K). Four TFs (FOSB, JUN, KLF6, and JUND) were slightly up-regulated and 4 TFs were activated in subC1_CD8⁺ T cells, while three TFs (POLR2A, ETS1 and JUNB) were slightly up-regulated and 3 TFs were activated in subC3_CD8⁺ T cells (Figures 5L, M).

FOSB, a member of the Fos transcription factor family and a component of the AP-1 complex, has been implicated in diverse biological processes including bone-forming tumors and a subset of vascular tumors (47–49). JUN, an extensively studied protein of the AP-1 complex, is involved in numerous cell activities such as proliferation, apoptosis, survival, tumorigenesis, and tissue morphogenesis (50). Kruppel like factor 6 (KLF6), a transcription factor of the zinc finger family, has been reported to play an important role in lipid homeostasis regulation in clear cell renal cell carcinoma (51). Jund proto-oncogene subunit (JunD), a prominent AP-1 component, is a bone formation inhibitor that contributes to low bone mass induced by estrogen depletion (52). RNA polymerase II subunit A (POLR2A), one of the subunits of the human RNA polymerase II complex, encodes the largest subunit that is indispensable for polymerase activity during messenger RNA synthesis. Prior work suggested that inhibiting POLR2A with a-amanitin could lead to extensive cell death (53). ETS1 belongs to a large family of the ETS domain family of transcription factors and is thought to be linked to poor survival during cancer progression (54) as well as transcription networks associated with CD8 T cell differentiation (55). JUNB, an AP-1 transcription factor, is expressed by eTreg cells and promotes an IRF4-dependent transcription program (56).

Complex Intercellular and Molecular Interaction Networks in GCTB

CellPhoneDB analysis was used to calculate the number of ligand–receptor pairs among all of the cell types in our data (Figure 6A). These ligand-receptor relationships are shown in Figures 6B, C. Thirty-three migration-associated genes were expressed by C2_Mature osteoclasts (Figure 6D, Supplementary Table 2). Cells of interest were selected, and their ligand–receptor relationships were identified (Figure 6E). All these ligand–receptor pairs were expressed individually. All cells expressed the ligand MIF, which is a macrophage migration inhibitory factor that regulates the function of mononuclear macrophages in host defense at sites of inflammation, while its receptor TNFRSF14 was expressed in C2_Mature osteoclasts. Osteoblasts significantly expressed TNFSF11 (RANKL), while the corresponding receptor TNFRSF11A (RANK) was expressed by C2_Mature osteoclasts. The TNFSF11/TNFRSF11A pathway regulates the activation of osteoclasts and induces the migration of tumor cells, notably in breast and lung cancer (57). CKLF (chemokine-like factor) is a novel cloned chemotactic cytokine that regulates the migration of immune cells (58) and was significantly expressed by C2_Mature osteoclasts while its receptor LPR6 was expressed in chondrocytes, endothelial cells, and osteoblasts. We also measured the expression levels of migration-related genes in the different cells (Supplementary Figure 2A, Supplementary Table 3). ACP5, CTSK, and ATP6V0D2, which are associated with the function of osteoclasts, were highly expressed in C2_Mature osteoclasts. Further, genes associated with cell migration were mostly expressed by C1_Progenitor osteoclasts (for instance, PTGER4, HBEGF, VEGFA, CD44, and VSIR) and C2_Mature osteoclasts (such as TNFRSF14, CD58, DPP4, EFNB1, SEMA7A, GRN, TNFRSF11A, APP, and CKLF). CKLF was more highly expressed by C2_Mature osteoclasts than by C1_Progenitor osteoclasts and C2_Mononuclear macrophages. Multiplex immunofluorescence experiments also revealed the presence of C2_Mature osteoclasts (Supplementary Figure 2B).

FIGURE 6

Figure 6 Cell–cell communication network of CellPhoneDB in the GCTB. (A) Heatmap showing the number of potential ligand–receptor pairs among the predicted cell types. (B, C) Interaction network constructed by CellPhoneDB. Thicker lines indicate more interaction with other types of cells. (D) Venn diagram of genes related to cells migration obtained from GSEA, and those that acted as receptors expressed in C2_Mature osteoclasts in our data. (E) Ligand–receptor pairs detected with CellPhoneDB are shown in a bubble plot. GCTB, giant cell tumor of bone; UMAP, uniform manifold approximation and projection for dimension reduction; GSEA, gene set enrichment analysis.

We algorithmically divided all of the cells into five patterns via CellChat analysis. C2_Mature osteoclasts mainly gathered in Pattern3, which expressed the RANKL, PARs, CD137 RANKL, and SEMA3 signaling pathways (Figures 7A, B). The bubble diagram shows the expression of these signaling pathways by C2_Mature osteoclasts (Figure 7C). Signaling pathways of interest were then selected, and their interaction networks are shown (Figure 7D). The RANKL signaling pathway was found to be highly associated with osteoclast formation and migration (57). Proteinase-activated receptors (PARs), which are key components of the PARs signaling pathway, are expressed in various epithelial tumors and can regulate progression and metastasis (59). The CD137 signaling pathway is critically involved in the promotion of breast cancer bone metastasis via enhancement of the migration and differentiation of monocytes/macrophages into osteoclast (60). SEMA3, a subfamily of signaling molecules in the SEMA3 signaling pathway, plays diverse roles in the regulation of cell proliferation, differentiation, and migration (61).

FIGURE 7

Figure 7 Cell–cell communication network of CellChat in GCTB. (A) Heatmap showing the enrichment of different cell types in different modules. The deeper red the color is, the higher its degree of enrichment. (B) Heatmap showing the contribution of different signaling pathways in different expression modules. The deeper the red color is, the higher the contribution of the signaling pathways. (C) A bubble plot showing the signaling pathways of interest in GCTB as detected by CellChat. (D) An interaction network shows the ligand–receptor pairs that have a known biological significance. GCTB, giant cell tumor of bone.

The interaction between NK/T cells and osteoclasts was analyzed (Figures 8A–C). Cells of interest were selected, and their ligand-receptor interactions were identified (Figure 8D). All immune cells expressed ligand HLA-C, a type of human leucocyte antigen that is correlated with the enhanced potent cytotoxic response of T cells immunity (62). Its receptor FAM3C, a member of the FAM3C gene family that is related to the inflammatory response and bone remodeling (63), was expressed in all subtypes of osteoclasts. C2_Mature osteoclasts expressed higher levels of TNFRSF14, while its receptor MIF, which is macrophage migration inhibitory factor that regulates the function of mononuclear macrophages in host defense at sites of inflammation, was expressed in C2_T cells, C3_T cells, C4_NK cells, subC2_CD8⁺T cells, and subC3_CD8⁺T cells. In addition to HLA-C_FAM3C, CD2_CD58 was also more highly expressed by these immune cells than osteoclasts. CXCR6 was significantly expressed by C2_T cells, C3_T cells and subC2_CD8⁺ T cells, while its receptor CXCL16 was expressed in all sub types of osteoclasts. CXCR6 has been shown to draw CD8⁺ T cells to the liver in graft-versus-host disease and is required for the maintenance of liver-resident CD8⁺ T cells following infection (64, 65). CXCL16 is a chemokine that binds to CXCR6 on Th1 and activated CD8 effector T cells and plays an important role in their recruitment to sites of inflammation (66, 67).

FIGURE 8

Figure 8 Cell-cell communication network between NK/T cells and osteoclasts in GCTB. (A) Heatmap showing the number of potential ligand-receptor pairs among the selected cell types. (B, C) Interaction network constructed by CellPhoneDB. Thicker lines mean more interaction with other types of cells. (D) Ligand-receptor pairs are shown in the bubble plot. GCTB, giant cell tumor of bone.

Discussion

The existing literature on GCTB is mostly comprised of genetic and genomic studies and clinical case reports (68, 69). These studies provide a basis for understanding the current treatment strategies for GCTB. However, accurate analysis of the cellular heterogeneity of GCTB remains challenging. To address the heterogeneity of GCTB, the present study evaluated a patient who was diagnosed with GCTB. To the best of our knowledge, this is the first comprehensive single-cell transcriptome study on a patient with GCTB. Khazaei et al. (70) characterized the transcriptomic of G34W GCTB and reported that a G34W mutation is necessary for GCTB tumorigenesis. However, the number of GCTB subtypes, their distinctive properties, and their heterogeneity remains unclear. In the present report, we identified eight major cell clusters in GCTB with UMAP clustering. The complexity of the GCTB cellular ecosystem and the intratumor heterogeneity of GCTB were revealed by scRNA-seq.

The main clinical manifestation of GCTB is osteolysis, which can lead to significant motor dysfunction. Osteoclasts, the only cells that are known to be involved in osteolysis in the human body, are believed to be intensively involved with GCTB (71). An increasingly large body of evidence has shown that key genes such as ACP5, CTSK, and ATP6V0D2 (72) are involved in the bone resorption function of osteoclasts. The directional migration of osteoclasts at the lesion site also attracted our interest. C1_Progenitor osteoclasts and C2_Mature osteoclasts have their own unique genes that regulate their migration and their highly migratory state. Interestingly, our data showed the CKLF, a newly cloned chemotactic cytokine, was relatively increased in subtypes of C2_Mature osteoclasts, indicating improved osteoclast mobility during maturation. In addition, the key genes that affected the aggregation of osteoclasts in GCTB were identified and may be key targets of future therapies, although more experiments are needed to further confirm this. As the high migration state of osteoclasts was confirmed, analysis of key ligand–receptor pairs that regulate the interactions of the cells within the tumor environment of GCTB is important. Osteoclasts were regulated by several cells in GCTB. We mainly focused on the key ligand-receptor pairs that interacted with and regulated C2_Mature osteoclasts and other cells in the GCTB tumor microenvironment. Previous reports showed that TNFRSF11A_TNFSF11, a classical ligand-receptor pair that regulates osteoclast maturation and function (73), was also included in our data. These key ligand–receptor pairs provide strong evidence of molecular crosstalk in the bone microenvironment (74) and can serve as guidance for follow-up studies on osteoclast migration and function. Moreover, four key signaling pathways involved in the formation and migration of C2_Mature osteoclasts were identified: the RANKL signaling pathway network, the CD137 signaling pathway network, the PARs signaling pathway network, and the SMEA3 signaling pathway network. The results provide a theoretical basis for future studies on the osteolytic effects of GCTB.

There are several limitations of this study: 1. scRNA-seq only contains about 0.1 pg of mRNA per cell on average, so there are too few materials available for sequencing (75); 2. the single-cell suspension digestion scheme directly affects the composition ratio of cells, resulting in the loss of rare cells, increasing cell stress and affecting cell gene expression (9); 3. restrictions on the depth of the sequencing, as at present, most scRNA-seq can only detect 10%-20% of mRNA (76, 77); 4. the small number of GCTB cases used in this work can only reflect the heterogeneity of the cells rather than between patients, and the results may be biased due to the small number of included cases.

Overall, this study characterized the heterogeneity of GCTB, provided a valuable single-cell transcription atlas and novel insights into GCTB, and identified a new mechanism and target for clinical GCTB treatment.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/, GSE168664.

Ethics Statement

The studies involving human participants were reviewed and approved by The First Affiliated Hospital of Guangxi Medical University (No.2019KY-E-153). The patients/participants provided their written informed consent to participate in this study.

Author Contributions

QW, YL, and JM were responsible for the conception of the study. WF, ZQ, QH, and MH harvested the GCTB tissue and performed the scRNA-seq experiments. WF, MH, XJ, and CL analyzed the data. MH, TX, HL, and SL performed the immunofluorescence staining. WF, MH, and JH wrote the paper. QW, YL, JM, and JX contributed to the revision and discussion of the paper. All authors contributed to the article and approved the submitted version.

Funding

This project was supported in part by the Guangxi Key Laboratory of Genomics and Personalized Medicine. This project was also supported by the Nanning Scientific Research and Technology Development Plan (grant no. 20183012), Guangxi Key R&D Program (grant no.AB17292073), the National Natural Science Foundation of China (grant no. 81960768), Guangxi Key Laboratory of Genomics and Personalized Medicine open project (grant no. GXGPMC201903), Natural Science Foundation of Guangxi Province (grant no. 2020GXNSFAA259088), Guangxi Natural Science Foundation of Youth Science Foundation project (grant no. 2017GXNSFBA198098), the “Medical Excellence Award” Funded by the Creative Research Development Grant from the First Affiliated Hospital of Guangxi Medical University, and the Youth Science and Technology Project of the First Affiliated Hospital of Guangxi Medical University (grant no. 201903038). YL and SL were visiting scholars to the Molecule Biology Lab at the University of Western Australia.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We are grateful for the participation and cooperation of the GCTB patient. We are also grateful for the support of the Guangxi Key Laboratory of Genomics and Personalized Medicine.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2021.715552/full#supplementary-material

Supplementary Figure 1 | (A) 99mTc-MDP Bone Scan shows a clumpy reflexive-concentrated focus in the lower right distal femur. (B, C) Plain radiographs shows an osteolytic lesion in the distal femur. (D) A sagittal T1-weighted image shows the tumor with heterogeneous signal intensity, and (E) cystic changes are visible on a fat-suppressed T2-weighted image. (F) The pathologic images show local tumor cell infiltration (100×). GCTB, giant cell tumor of bone.

Supplementary Figure 2 | (A) Expression levels of migration-related genes between the different cells (*P < 0.05, **P < 0.01, ***P < 0.001; NS, No statistical significance). (B) Multiplex IHC staining of GCTB tissue. Arrows represent the co-expression positions of the relevant indicators (DAPI, ACP5, ATP6V0D2, CKLF) in Multiplex IHC (scale bar=200 μm). GCTB, giant cell tumor of bone; IHC, immunohistochemistry

Supplementary Table 1 | Characteristics of the patient with GCTB included in this study.

Supplementary Table 2 | Migration-associated genes.

Supplementary Table 3 | Expression levels of migration-related genes.

References

1. Ueda T, Morioka H, Nishida Y, Kakunaga S, Tsuchiya H, Matsumoto Y, et al. Objective Tumor Response to Denosumab in Patients With Giant Cell Tumor of Bone: A Multicenter Phase II Trial. Ann Oncol (2015) 26:2149–54. doi: 10.1093/annonc/mdv307

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Lipplaa A, Dijkstra S, Gelderblom H. Challenges of Denosumab in Giant Cell Tumor of Bone, and Other Giant Cell-Rich Tumors of Bone. Curr Opin Oncol (2019) 31:329–35. doi: 10.1097/CCO.0000000000000529

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Shibuya I, Takami M, Miyamoto A, Karakawa A, Dezawa A, Nakamura S, et al. In Vitro Study of the Effects of Denosumab on Giant Cell Tumor of Bone: Comparison With Zoledronic Acid. Pathol Oncol Res (2019) 25:409–19. doi: 10.1007/s12253-017-0362-8

PubMed Abstract | CrossRef Full Text | Google Scholar

4. van der Heijden L, Dijkstra S, van de Sande M, Gelderblom H. Current Concepts in the Treatment of Giant Cell Tumour of Bone. Curr Opin Oncol (2020) 32:332–8. doi: 10.1097/CCO.0000000000000645

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Li H, Gao J, Gao Y, Lin N, Zheng M, Ye Z. Denosumab in Giant Cell Tumor of Bone: Current Status and Pitfalls. Front Oncol (2020) 10:580605. doi: 10.3389/fonc.2020.580605

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Branstetter DG, Nelson SD, Manivel JC, Blay JY, Chawla S, Thomas DM, et al. Denosumab Induces Tumor Reduction and Bone Formation in Patients With Giant-Cell Tumor of Bone. Clin Cancer Res (2012) 18:4415–24. doi: 10.1158/1078-0432.CCR-12-0578

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Hui L, Chen Y. Tumor Microenvironment: Sanctuary of the Devil. Cancer Lett (2015) 368:7–13. doi: 10.1016/j.canlet.2015.07.039

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Balko JM, Cook RS, Vaught DB, Kuba MG, Miller TW, Bhola NE, et al. Profiling of Residual Breast Cancers After Neoadjuvant Chemotherapy Identifies DUSP4 Deficiency as a Mechanism of Drug Resistance. Nat Med (2012) 18:1052–9. doi: 10.1038/nm.2795

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Potter SS. Single-Cell RNA Sequencing for the Study of Development, Physiology and Disease. Nat Rev Nephrol (2018) 14:479–92. doi: 10.1038/s41581-018-0021-7

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Qiu X, Mao Q, Tang Y, Wang L, Chawla R, Pliner HA, et al. Reversed Graph Embedding Resolves Complex Single-Cell Trajectories. Nat Methods (2017) 14:979–82. doi: 10.1038/nmeth.4402

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Efremova M, Vento-Tormo M, Teichmann SA, Vento-Tormo R. CellPhoneDB: Inferring Cell-Cell Communication From Combined Expression of Multi-Subunit Ligand-Receptor Complexes. Nat Protoc (2020) 15:1484–506. doi: 10.1038/s41596-020-0292-x

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Vento-Tormo R, Efremova M, Botting RA, Turco MY, Vento-Tormo M, Meyer KB, et al. Single-Cell Reconstruction of the Early Maternal-Fetal Interface in Humans. Nature (2018) 563:347–53. doi: 10.1038/s41586-018-0698-6

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Jin S, Guerrero-Juarez CF, Zhang L, Chang I, Ramos R, Kuan CH, et al. Inference and Analysis of Cell-Cell Communication Using CellChat. Nat Commun (2021) 12:1088. doi: 10.1038/s41467-021-21246-9

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Hanzelmann S, Castelo R, Guinney J. GSVA: Gene Set Variation Analysis for Microarray and RNA-Seq Data. BMC Bioinf (2013) 14:7. doi: 10.1186/1471-2105-14-7

CrossRef Full Text | Google Scholar

15. Aibar S, Gonzalez-Blas CB, Moerman T, Huynh-Thu VA, Imrichova H, Hulselmans G, et al. SCENIC: Single-Cell Regulatory Network Inference and Clustering. Nat Methods (2017) 14:1083–6. doi: 10.1038/nmeth.4463

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Li Y, Zhang QY, Sun BF, Ma Y, Zhang Y, Wang M, et al. Single-Cell Transcriptome Profiling of the Vaginal Wall in Women With Severe Anterior Vaginal Prolapse. Nat Commun (2021) 12:87. doi: 10.1038/s41467-020-20358-y

PubMed Abstract | CrossRef Full Text | Google Scholar

17. La Manno G, Soldatov R, Zeisel A, Braun E, Hochgerner H, Petukhov V, et al. RNA Velocity of Single Cells. Nature (2018) 560:494–8. doi: 10.1038/s41586-018-0414-6

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Ren X, Zhong G, Zhang Q, Zhang L, Sun Y, Zhang Z. Reconstruction of Cell Spatial Organization From Single-Cell RNA Sequencing Data Based on Ligand-Receptor Mediated Self-Assembly. Cell Res (2020) 30:763–78. doi: 10.1038/s41422-020-0353-2

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Liu Y, He S, Wang XL, Peng W, Chen QY, Chi DM, et al. Tumour Heterogeneity and Intercellular Networks of Nasopharyngeal Carcinoma at Single Cell Resolution. Nat Commun (2021) 12:741. doi: 10.1038/s41467-021-21043-4

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Chen YP, Yin JH, Li WF, Li HJ, Chen DP, Zhang CJ, et al. Single-Cell Transcriptomics Reveals Regulators Underlying Immune Cell Diversity and Immune Subtypes Associated With Prognosis in Nasopharyngeal Carcinoma. Cell Res (2020) 30:1024–42. doi: 10.1038/s41422-020-0374-x

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Clark BS, Stein-O’Brien GL, Shiau F, Cannon GH, Davis-Marcisak E, Sherman T, et al. Single-Cell RNA-Seq Analysis of Retinal Development Identifies NFI Factors as Regulating Mitotic Exit and Late-Born Cell Specification. Neuron (2019) 102:1111–26.e5. doi: 10.1016/j.neuron.2019.04.010

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Baird A, Lindsay T, Everett A, Iyemere V, Paterson YZ, McClellan A, et al. Osteoblast Differentiation of Equine Induced Pluripotent Stem Cells. Biol Open (2018) 7:bio033514. doi: 10.1242/bio.033514

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Azizi E, Carr AJ, Plitas G, Cornish AE, Konopacki C, Prabhakaran S, et al. Single-Cell Map of Diverse Immune Phenotypes in the Breast Tumor Microenvironment. Cell (2018) 174:1293–308.e36. doi: 10.1016/j.cell.2018.05.060

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Zhou Y, Yang D, Yang Q, Lv X, Huang W, Zhou Z, et al. Single-Cell RNA Landscape of Intratumoral Heterogeneity and Immunosuppressive Microenvironment in Advanced Osteosarcoma. Nat Commun (2020) 11:6322. doi: 10.1038/s41467-020-20059-6

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Peng J, Sun BF, Chen CY, Zhou JY, Chen YS, Chen H, et al. Single-Cell RNA-Seq Highlights Intra-Tumoral Heterogeneity and Malignant Progression in Pancreatic Ductal Adenocarcinoma. Cell Res (2019) 29:725–38. doi: 10.1038/s41422-019-0195-y

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Nordstrand A, Bovinder Ylitalo E, Thysell E, Jernberg E, Crnalic S, Widmark A, et al. Bone Cell Activity in Clinical Prostate Cancer Bone Metastasis and Its Inverse Relation to Tumor Cell Androgen Receptor Activity. Int J Mol Sci (2018) 19:1223. doi: 10.3390/ijms19041223

CrossRef Full Text | Google Scholar

27. Liu T, Zhang H, Yi S, Gu L, Zhou M. Mutual Regulation of MDM4 and TOP2A in Cancer Cell Proliferation. Mol Oncol (2019) 13:1047–58. doi: 10.1002/1878-0261.12457

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Chim SM, Qin A, Tickner J, Pavlos N, Davey T, Wang H, et al. EGFL6 Promotes Endothelial Cell Migration and Angiogenesis Through the Activation of Extracellular Signal-Regulated Kinase. J Biol Chem (2011) 286:22035–46. doi: 10.1074/jbc.M110.187633

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Auvinen K, Lokka E, Mokkala E, Jappinen N, Tyystjarvi S, Saine H, et al. Fenestral Diaphragms and PLVAP Associations in Liver Sinusoidal Endothelial Cells Are Developmentally Regulated. Sci Rep (2019) 9:15698. doi: 10.1038/s41598-019-52068-x

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Chen J, Chen F, Bian H, Wang Q, Zhang X, Sun L, et al. Hypertrophic Chondrocyte-Specific Col10a1 Controlling Elements in Cre Recombinase Transgenic Studies. Am J Trans Res (2019) 11:6672–9.

Google Scholar

31. Dazzi F. Cancer Makes New Friends With Old Tricks. Blood (2013) 122:1093–4. doi: 10.1182/blood-2013-06-509620

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Arneth B. Tumor Microenvironment. Med (Kaunas) (2019) 56:15. doi: 10.3390/medicina56010015

CrossRef Full Text | Google Scholar

33. Hoch B, Inwards C, Sundaram M, Rosenberg AE. Multicentric Giant Cell Tumor of Bone. Clinicopathologic Analysis of Thirty Cases. J Bone Joint Surg Am (2006) 88:1998–2008. doi: 10.2106/JBJS.E.01111

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Lotinun S, Ishihara Y, Nagano K, Kiviranta R, Carpentier VT, Neff L, et al. Cathepsin K-Deficient Osteocytes Prevent Lactation-Induced Bone Loss and Parathyroid Hormone Suppression. J Clin Invest (2019) 129:3058–71. doi: 10.1172/JCI122936

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Kim T, Ha HI, Kim N, Yi O, Lee SH, Choi Y. Adrm1 Interacts With Atp6v0d2 and Regulates Osteoclast Differentiation. Biochem Biophys Res Commun (2009) 390:585–90. doi: 10.1016/j.bbrc.2009.10.010

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Kawaida R, Ohtsuka T, Okutsu J, Takahashi T, Kadono Y, Oda H, et al. Jun Dimerization Protein 2 (JDP2), a Member of the AP-1 Family of Transcription Factor, Mediates Osteoclast Differentiation Induced by RANKL. J Exp Med (2003) 197:1029–35. doi: 10.1084/jem.20021321

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Maruyama K, Fukasaka M, Vandenbon A, Saitoh T, Kawasaki T, Kondo T, et al. The Transcription Factor Jdp2 Controls Bone Homeostasis and Antibacterial Immunity by Regulating Osteoclast and Neutrophil Differentiation. Immunity (2012) 37:1024–36. doi: 10.1016/j.immuni.2012.08.022

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Yang Y, Chung MR, Zhou S, Gong X, Xu H, Hong Y, et al. STAT3 Controls Osteoclast Differentiation and Bone Homeostasis by Regulating NFATc1 Transcription. J Biol Chem (2019) 294:15395–407. doi: 10.1074/jbc.RA119.010139

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Hunt LC, Xu B, Finkelstein D, Fan Y, Carroll PA, Cheng PF, et al. The Glucose-Sensing Transcription Factor MLX Promotes Myogenesis via Myokine Signaling. Genes Dev (2015) 29:2475–89. doi: 10.1101/gad.267419.115

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Mejhert N, Kuruvilla L, Gabriel KR, Elliott SD, Guie MA, Wang H, et al. Partitioning of MLX-Family Transcription Factors to Lipid Droplets Regulates Metabolic Gene Expression. Mol Cell (2020) 77:1251–64 e9. doi: 10.1016/j.molcel.2020.01.014

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Tiwari A, Swamy S, Gopinath KS, Kumar A. Genomic Amplification Upregulates Estrogen-Related Receptor Alpha and Its Depletion Inhibits Oral Squamous Cell Carcinoma Tumors In Vivo. Sci Rep (2015) 5:17621. doi: 10.1038/srep17621

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Zheng C, Zheng L, Yoo JK, Guo H, Zhang Y, Guo X, et al. Landscape of Infiltrating T Cells in Liver Cancer Revealed by Single-Cell Sequencing. Cell (2017) 169:1342–56 e16. doi: 10.1016/j.cell.2017.05.035

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Santoiemma PP, Powell DJ Jr. Tumor Infiltrating Lymphocytes in Ovarian Cancer. Cancer Biol Ther (2015) 16:807–20. doi: 10.1080/15384047.2015.1040960

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Stanton SE, Disis ML. Clinical Significance of Tumor-Infiltrating Lymphocytes in Breast Cancer. J Immunother Cancer (2016) 4:59. doi: 10.1186/s40425-016-0165-6

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Zhang S, Liu Z, Wu D, Chen L, Xie L. Single-Cell RNA-Seq Analysis Reveals Microenvironmental Infiltration of Plasma Cells and Hepatocytic Prognostic Markers in HCC With Cirrhosis. Front Oncol (2020) 10:596318. doi: 10.3389/fonc.2020.596318

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Kim N, Kim HK, Lee K, Hong Y, Cho JH, Choi JW, et al. Single-Cell RNA Sequencing Demonstrates the Molecular and Cellular Reprogramming of Metastatic Lung Adenocarcinoma. Nat Commun (2020) 11:2285. doi: 10.1038/s41467-020-16164-1

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Hung YP, Fletcher CD, Hornick JL. FOSB is a Useful Diagnostic Marker for Pseudomyogenic Hemangioendothelioma. Am J Surg Pathol (2017) 41:596–606. doi: 10.1097/PAS.0000000000000795

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Fittall MW, Mifsud W, Pillay N, Ye H, Strobl AC, Verfaillie A, et al. Recurrent Rearrangements of FOS and FOSB Define Osteoblastoma. Nat Commun (2018) 9:2150. doi: 10.1038/s41467-018-04530-z

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Antonescu CR, Chen HW, Zhang L, Sung YS, Panicek D, Agaram NP, et al. ZFP36-FOSB Fusion Defines a Subset of Epithelioid Hemangioma With Atypical Features. Genes Chromosomes Cancer (2014) 53:951–9. doi: 10.1002/gcc.22206

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Meng Q, Xia Y. C-Jun, at the Crossroad of the Signaling Network. Protein Cell (2011) 2:889–98. doi: 10.1007/s13238-011-1113-3

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Syafruddin SE, Rodrigues P, Vojtasova E, Patel SA, Zaini MN, Burge J, et al. A KLF6-Driven Transcriptional Network Links Lipid Homeostasis and Tumour Growth in Renal Carcinoma. Nat Commun (2019) 10:1152. doi: 10.1038/s41467-019-09116-x

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Kawamata A, Izu Y, Yokoyama H, Amagasa T, Wagner EF, Nakashima K, et al. JunD Suppresses Bone Formation and Contributes to Low Bone Mass Induced by Estrogen Depletion. J Cell Biochem (2008) 103:1037–45. doi: 10.1002/jcb.21660

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Lindell TJ, Weinberg F, Morris PW, Roeder RG, Rutter WJ. Specific Inhibition of Nuclear RNA Polymerase II by Alpha-Amanitin. Science (1970) 170:447–9. doi: 10.1126/science.170.3956.447

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Dittmer J. The Role of the Transcription Factor Ets1 in Carcinoma. Semin Cancer Biol (2015) 35:20–38. doi: 10.1016/j.semcancer.2015.09.010

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Moskowitz DM, Zhang DW, Hu B, Le Saux S, Yanes RE, Ye Z, et al. Epigenomics of Human CD8 T Cell Differentiation and Aging. Sci Immunol (2017) 2:eaag0192. doi: 10.1126/sciimmunol.aag0192

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Koizumi SI, Sasaki D, Hsieh TH, Taira N, Arakaki N, Yamasaki S, et al. JunB Regulates Homeostasis and Suppressive Functions of Effector Regulatory T Cells. Nat Commun (2018) 9:5344. doi: 10.1038/s41467-018-07735-4

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Wang Y, Song Y, Che X, Zhang L, Wang Q, Zhang X, et al. Caveolin1 Enhances RANKLinduced Gastric Cancer Cell Migration. Oncol Rep (2018) 40:1287–96. doi: 10.3892/or.2018.6550

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Wang ZZ, Yuan YH, Zhang Y, Wang XF, Chu SF, Han N, et al. Chemokine-Like Factor 1 Promotes the Migration of Rat Primary Cortical Neurons by the Induction of Actin Polymerization. Neuroreport (2014) 25:1221–6. doi: 10.1097/WNR.0000000000000252

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Kaufmann R, Hollenberg MD. Proteinase-Activated Receptors (PARs) and Calcium Signaling in Cancer. Adv Exp Med Biol (2012) 740:979–1000. doi: 10.1007/978-94-007-2888-2_45

PubMed Abstract | CrossRef Full Text | Google Scholar

60. Jiang P, Gao W, Ma T, Wang R, Piao Y, Dong X, et al. CD137 Promotes Bone Metastasis of Breast Cancer by Enhancing the Migration and Osteoclast Differentiation of Monocytes/Macrophages. Theranostics (2019) 9:2950–66. doi: 10.7150/thno.29617

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Fu X, Xu J, Chaturvedi P, Liu H, Jiang R, Lan Y. Identification of Osr2 Transcriptional Target Genes in Palate Development. J Dent Res (2017) 96:1451–8. doi: 10.1177/0022034517719749

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Sips M, Liu Q, Draghi M, Ghebremichael M, Berger CT, Suscovich TJ, et al. HLA-C Levels Impact Natural Killer Cell Subset Distribution and Function. Hum Immunol (2016) 77:1147–53. doi: 10.1016/j.humimm.2016.08.004

PubMed Abstract | CrossRef Full Text | Google Scholar

63. Montalvany-Antonucci CC, Zicker MC, Ferreira AVM, Macari S, Ramos-Junior ES, Gomez RS, et al. High-Fat Diet Disrupts Bone Remodeling by Inducing Local and Systemic Alterations. J Nutr Biochem (2018) 59:93–103. doi: 10.1016/j.jnutbio.2018.06.006

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Tse SW, Radtke AJ, Espinosa DA, Cockburn IA, Zavala F. The Chemokine Receptor CXCR6 is Required for the Maintenance of Liver Memory CD8(+) T Cells Specific for Infectious Pathogens. J Infect Dis (2014) 210:1508–16. doi: 10.1093/infdis/jiu281

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Wein AN, McMaster SR, Takamura S, Dunbar PR, Cartwright EK, Hayward SL, et al. CXCR6 Regulates Localization of Tissue-Resident Memory CD8 T Cells to the Airways. J Exp Med (2019) 216:2748–62. doi: 10.1084/jem.20181308

PubMed Abstract | CrossRef Full Text | Google Scholar

66. Matsumura S, Wang B, Kawashima N, Braunstein S, Badura M, Cameron TO, et al. Radiation-Induced CXCL16 Release by Breast Cancer Cells Attracts Effector T Cells. J Immunol (Baltimore Md: 1950) (2008) 181:3099–107. doi: 10.4049/jimmunol.181.5.3099

CrossRef Full Text | Google Scholar

67. Jiang X, Sun W, Zhu L, Guo D, Jiang H, Ma D, et al. Expression of CXCR6 on CD8(+) T Cells was Up-Regulated in Allograft Rejection. Transpl Immunol (2010) 22:179–83. doi: 10.1016/j.trim.2009.12.001

PubMed Abstract | CrossRef Full Text | Google Scholar

68. Yoshida KI, Nakano Y, Honda-Kitahara M, Wakai S, Motoi T, Ogura K, et al. Absence of H3F3A Mutation in a Subset of Malignant Giant Cell Tumor of Bone. Mod Pathol (2019) 32:1751–61. doi: 10.1038/s41379-019-0318-5

PubMed Abstract | CrossRef Full Text | Google Scholar

69. Tang H, Moro A, Feng W, Lai Y, Xiao Z, Liu Y, et al. Giant Cell Tumors Combined With Secondary Aneurysmal Bone Cysts Are More Likely to Develop Postoperative Recurrence: A Retrospective Study of 256 Cases. J Surg Oncol (2019) 120:359–65. doi: 10.1002/jso.25588

PubMed Abstract | CrossRef Full Text | Google Scholar

70. Khazaei S, De Jay N, Deshmukh S, Hendrikse LD, Jawhar W, Chen CCL, et al. H3.3 G34W Promotes Growth and Impedes Differentiation of Osteoblast-Like Mesenchymal Progenitors in Giant Cell Tumor of Bone. Cancer Discovery (2020) 10:1968–87. doi: 10.1158/2159-8290.CD-20-0461

PubMed Abstract | CrossRef Full Text | Google Scholar

71. Zheng MH, Robbins P, Xu J, Huang L, Wood DJ, Papadimitriou JM. The Histogenesis of Giant Cell Tumour of Bone: A Model of Interaction Between Neoplastic Cells and Osteoclasts. Histol Histopathol (2001) 16:297–307. doi: 10.14670/HH-16.297

PubMed Abstract | CrossRef Full Text | Google Scholar

72. Feng H, Cheng T, Steer JH, Joyce DA, Pavlos NJ, Leong C, et al. Myocyte Enhancer Factor 2 and Microphthalmia-Associated Transcription Factor Cooperate With NFATc1 to Transactivate the V-ATPase D2 Promoter During RANKL-Induced Osteoclastogenesis. J Biol Chem (2009) 284:14667–76. doi: 10.1074/jbc.M901670200

PubMed Abstract | CrossRef Full Text | Google Scholar

73. Wada T, Nakashima T, Hiroshi N, Penninger JM. RANKL-RANK Signaling in Osteoclastogenesis and Bone Disease. Trends Mol Med (2006) 12:17–25. doi: 10.1016/j.molmed.2005.11.007

PubMed Abstract | CrossRef Full Text | Google Scholar

74. Kular J, Tickner J, Chim SM, Xu J. An Overview of the Regulation of Bone Remodelling at the Cellular Level. Clin Biochem (2012) 45:863–73. doi: 10.1016/j.clinbiochem.2012.03.021

PubMed Abstract | CrossRef Full Text | Google Scholar

75. Wang Y, Navin NE. Advances and Applications of Single-Cell Sequencing Technologies. Mol Cell (2015) 58:598–609. doi: 10.1016/j.molcel.2015.05.005

PubMed Abstract | CrossRef Full Text | Google Scholar

76. Islam S, Zeisel A, Joost S, La Manno G, Zajac P, Kasper M, et al. Quantitative Single-Cell RNA-Seq With Unique Molecular Identifiers. Nat Methods (2014) 11:163–6. doi: 10.1038/nmeth.2772

PubMed Abstract | CrossRef Full Text | Google Scholar

77. Svensson V, Natarajan KN, Ly LH, Miragaia RJ, Labalette C, Macaulay IC, et al. Power Analysis of Single-Cell RNA-Sequencing Experiments. Nat Methods (2017) 14:381–7. doi: 10.1038/nmeth.4220

PubMed Abstract | CrossRef Full Text | Google Scholar