Abstract

Background

   Papillary thyroid carcinoma (PTC) is the main pathological type of
   thyroid carcinoma (TC). Gender is a prominent background parameter for
   patients with PTC. Here, we aimed to delineate the differences in cell
   clusters and immune microenvironment in relation to gender in PTC.

Methods

   We generated 6720, 14,666, and 33,373 single-cell transcriptomes that
   were pooled from the tissues of four male patients with PTC, seven
   female patients with PTC, and three patients with nodular goiter,
   respectively. We performed single-cell RNA-sequencing (scRNA-seq) based
   on BD Rhapsody and characterized the first single-cell transcriptomic
   landscape of PTC involving gender. The differential cell clusters and
   their gene profiles were identified and analyzed via a multi-resolution
   network in male and female patients. The interactions of fibroblasts
   and endothelial cells with malignant epithelial cells and the
   difference in the immune infiltration of B and T lymphocytes according
   to gender were assessed.

Results

   Malignant epithelial cells were divided into two distinct subsets in
   male and female patients with PTC. Moreover, significant differences
   involving inferred copy-number variations (CNVs), gene profiles, and
   cell differentiation were detected between male and female patients.
   Regarding the interactions of fibroblasts and endothelial cells with
   malignant epithelial cells, members of the human leukocyte antigen
   (HLA) family and their receptors were considered as typical in female
   patients with PTC, while transforming growth factor beta 1 (TGFB1) and
   its receptors were typical of male patients with PTC. The
   characteristics of B cells, including cell clusters, cell
   differentiation, and dominant gene sets, were significantly different
   between genders.

Conclusions

   Our data revealed the detailed differences in cell clusters and immune
   microenvironment in PTC according to gender at the single-cell level,
   which provided new insights into the understanding of the impact of
   gender on PTC.

   Keywords: Papillary thyroid carcinoma, Single‐cell RNA-sequencing, Cell
   differentiation, Immune microenvironment, Gender‐driven

Introduction

   Thyroid cancer is the most common endocrine-system malignancy worldwide
   [[45]1, [46]2]. PTC is the most common histopathological type of TC,
   accounting for about 75–85% [[47]3]. Intriguingly, obvious differences
   were observed between male and female patients [[48]1, [49]4–[50]6].
   The incidence of PTC in women is significantly higher than that in men
   and the female-to-male ratio is about three or more [[51]1, [52]4].
   Moreover, epidemiological data indicate that male patients with PTC
   have a more aggressive phenotype and worse clinical outcomes [[53]6].
   Although females exhibit a higher incidence of PTC, male gender has
   been shown to be associated with higher rates of advanced disease stage
   and a worse disease prognosis. Kurozumi et al. showed that
   extrathyroidal extension was closely related to male gender [[54]7].
   Machens and colleagues performed a retrospective cohort analysis and
   the results showed that the primary tumor size of PTC in male patients
   was significantly larger than that in female patients and male patients
   had a higher incidence of lymph node metastases and distant metastases
   [[55]8]. Nevertheless, the differences induced by gender in PTC are
   little understood at the cellular level.

   The tumor immune microenvironment (TIME) is a collection of immune
   cells and related immune factors (including innate and adaptive immune
   factors) located in and around tumor tissues [[56]9]. The interaction
   between malignant cells and the proximal immune cells leads to an
   environment that promotes tumor growth and metastasis [[57]10, [58]11].
   The therapeutic strategies targeting the TIME are expected to be highly
   promising in view of the powerful regulatory ability of the TIME on
   tumor cells [[59]12–[60]15]. Single-cell transcriptomic sequencing is a
   technology that developed rapidly in the last decade that can
   accurately distinguish and explore cell clusters and related gene
   profiles in complex mixtures [[61]16, [62]17].

   Here, the single-cell transcriptomic landscape of PTC was characterized
   for the first time and the differences in cell clusters between male
   and female patients with PTC were identified. We also analyzed the
   interactions of fibroblasts and endothelial cells with malignant
   epithelial cells, as well as the differences in the immune infiltration
   of B and T lymphocytes.

Methods

Patients and demographic data

   Human samples of either nodular goiter or PTC were collected from
   patients who received surgery at the First Affiliated Hospital of Sun
   Yat-sen University (Guangzhou, China) from March 2020 to April 2020.
   All patients were diagnosed with the corresponding disease types, and
   informed consent was obtained before sample collection. Patients with
   autoimmune thyroid disease (Hashimoto’s thyroiditis or lymphocytic
   thyroiditis) were excluded. In total, fourteen patients (three females
   with nodular goiter, seven females with PTC, and four males with PTC)
   were enrolled in the study. Age of males with PTC averaged 38.0 years
   with a range of 30–50 years. Age of females with PTC averaged
   36.6 years with a range of 20–70 years. Age of patients with nodule
   averaged 39.0 years with a range of 29–56 years. Tumor size of male
   averaged 1.53 cm with a range of 0.5–3.0 cm and female averaged 1.04 cm
   with a range of 0.3–1.2 cm. Nodule size of averaged 2.43 cm with a
   range of 0.4–4.6 cm. There were no multicentricity, extrathyroidal
   invasion and distant metastasis in PTC. Almost all the patients with
   PTC at TNM stage I except one male with PTC at stage II. All the above
   operations were in accordance with the relevant ethical norms and were
   approved by the Institutional Research Ethics Committee.

Preparation of single‐cell samples

   Human tissues were rinsed twice with cold saline and then minced using
   scissors. Tissue fragments were collected for enzymatic digestion in
   Dulbecco’s Modified Eagle’s Medium (DMEM, Gibco, Gaithersburg, MD, USA)
   containing 40 U/ml DNase type I (Gibco, Carlsbad, CA, USA) and
   1.0 ml/ml collagenase type IV (Sigma-Aldrich, St. Louis, Missouri,
   USA). The digestive system was incubated at 37 °C for 1 h, and cells
   were then passed through a 70 µm nylon cell strainer (Millipore,
   Billerica, MA, USA). The cells were centrifuged and resuspended with
   fetal bovine serum (FBS, Gibco, Grand Island, NY, USA). The prior
   viability and cell density were detected by trypan blue, and the
   average viability rate of cells before loading was no less than 90%. We
   loaded PTC cells from male patients into one channel, cells from female
   patients into a second channel and cells from patients with goiter into
   the third channel. Each of these three pooled samples was prepared for
   the next step of single-cell library preparation.

Single‐cell RNA-sequencing

   We performed single-cell RNA-sequencing (scRNA-seq) based on BD
   Rhapsody and characterized the first single-cell transcriptomic
   landscape of PTC involving gender. Single-cell libraries were prepared
   using the BD Rhapsody Single-Cell Analysis System (BD, USA), according
   to the manufacturer’s protocol. Libraries were sequenced using multiple
   runs on an Illumina NextSeq platform. The differential cell clusters
   and their gene profiles were identified and analyzed via a
   multi-resolution network in male and female patients. The interactions
   of fibroblasts and endothelial cells with malignant epithelial cells
   and the difference in the immune infiltration of B as well as T
   lymphocytes according to gender were assessed.

Unsupervised clustering of cells and t-distributed Stochastic Neighbor
Embedding (tSNE) visualization

   The Seurat package (version 3.1.2) implemented in R was used to perform
   the unsupervised clustering of cells. Genes that were expressed in less
   than two cells were filtered out. Cells with expression of > 200 genes
   and < 10% mitochondrial genes were further processed. The variation
   coefficient of genes was calculated using the Seurat arithmetic.
   Dimensionality reduction of all data was performed by principle
   component analysis (PCA) according to the first 1500 highest alterable
   genes. A k-nearest neighbor graph was executed with Euclidean distances
   in the space of the first 10 significant principal components. The
   cells were clustered by the Louvain Modularity optimization algorithm
   and visualized by tSNE visualization.

Marker gene identification and cell annotation

   To ascertain marker genes, the single-cell gene expression of each
   cluster was calculated using the “bimod” test implemented in the Seurat
   Find Markers function. Only genes with an expression with a log2
   average difference of 0.585 and P < 0.05 were classified as marker
   genes. ImmGen and Enrichr were used to represent cell clusters, which
   were annotated by canonical markers of known cell types.

Infer CNV analysis

   According to the mRNA copy detected by single-cell sequencing, the gene
   was located on the genome, and the heat map of inferred CNV was drawn.
   Taking the nodular goiter group as the reference, the differences were
   analyzed according to the copy number of mRNA and the degree of
   sequence variation, and the tumor cells were identified.

Hallmark signal pathway analysis

   The hallmark gene sets were obtained from the public Molecular
   Signature Database, and the signal pathway activity of each cell was
   scored using gene set variation analysis methods, according to a
   standardized strategy. Student’s t-test was used to calculate the
   differences in pathways between two groups. Significance was set at
   P < 0.05.

KEGG enrichment analysis

   KEGG analysis of selected genes was performed using the KOBAS software,
   with Benjamini-Hochberg testing adjustments.

Statistical analysis of population shifts

   The significant alteration of the percent composition of cell clusters
   in samples was evaluated using the Benjamini-Hochberg correction and
   the nonparametric Kruskal-Wallis test with Dunn’s multiple comparisons
   test.

Data statistics

   Data are presented as the mean ± SEM, and differential analysis was
   performed via either variance (three groups) or Student’s t-test (two
   groups) using GraphPad Prism 7 Significance was set at P < 0.05.

Results

Determination of the human PTC transcriptomic landscape

   Cells were isolated respectively from three patients with nodular
   goiter, four male patients with PTC, and seven female patients with
   PTC, followed by single-cell transcriptomic sequencing based on BD
   Rhapsody (Fig. [63]1a). All cell clusters were differentiated according
   to the three disease statuses (nodular, male PTC, and female PTC),
   which showed that the dominant clusters were different between male and
   female patients (Fig. [64]1b). The proportion of endothelial cells,
   fibroblasts, and B cells in female PTC patients were greatly higher
   than these in male PTC patients (Fig. [65]1c). Next, we focused on
   malignant epithelial cells, fibroblasts, endothelial cells, T cells and
   B cells by analyzing the differences in cell differentiation and gene
   expression profiles, as well as their interactions according to gender
   (Fig. [66]1d).

Fig. 1.

   [67]Fig. 1
   [68]Open in a new tab

   Transcriptome profiling of the three groups at the single-cell level. a
   Process and principle of the experiment. Tissues were collected and
   isolated for the preparation of single-cell libraries, which were then
   sequenced and all data were analyzed. b tSNE plot visualization showed
   cell fractions in three groups. The subjects of the three groups are
   displayed using different colors. c The frequency of representative
   cell clusters in the three groups is shown in the bar graph. d
   Subsequent data analysis strategy for inferCNV, including gene mapping
   in opposite-sex malignant cells and differentiation in opposite-sex T
   and B cells

Differences in cell differentiation and gene profiles of malignant epithelial
cells between male and female patients with PTC

   CNV is the premise of normal cell cancerization. Epithelial cells from
   nodular goiter tissues were used as references, and the CNVs of