Objective: The objective of this analysis is to explore the immune landscape and tumor microenvironment (TME) in cervical cancer through single-cell RNA sequencing (scRNA-seq) data. By focusing on the expression of key immune markers, I aimed to characterize the different immune cell populations present within the tumor and evaluate their potential roles in tumor biology and immune response. —— Loading the important Library

library(hdf5r)
library(Seurat)
## Loading required package: SeuratObject
## Loading required package: sp
## 
## Attaching package: 'SeuratObject'
## The following objects are masked from 'package:base':
## 
##     intersect, t
library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

Loading and processing our Data (scRNA-seq data from cervical cancer tissue from 10x genomic dataset)

cervcal <- Read10X_h5("Data/9k_Cervical_Cancer_scFFPE_count_filtered_feature_bc_matrix.h5")
file
## function (description = "", open = "", blocking = TRUE, encoding = getOption("encoding"), 
##     raw = FALSE, method = getOption("url.method", "default")) 
## {
##     .Internal(file(description, open, blocking, encoding, method, 
##         raw))
## }
## <bytecode: 0x0000022eac07ad10>
## <environment: namespace:base>
 ####creating a seurat object
cervical_object <- CreateSeuratObject(counts = cervcal, min.cells=3, min.feature= 200)

Quality Control

###addMetadata: percentage of mitochondrial gene
cervical_object[["percent.mt"]] <- PercentageFeatureSet(cervical_object, pattern = "^MT-")

#####violine plot viz of QC 
VlnPlot(cervical_object, feature = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol=3)
## Warning: Default search for "data" layer in "RNA" assay yielded no results;
## utilizing "counts" layer instead.

####
FeatureScatter(cervical_object, feature1 = "nCount_RNA", feature2 = "nFeature_RNA")+
  geom_smooth(method ="lm")
## `geom_smooth()` using formula = 'y ~ x'

filtring cell base on QC metrics

cervical_object <- subset(cervical_object, subset = nFeature_RNA > 200 & nFeature_RNA < 2500 & percent.mt < 5)

Normalize The Data

cervical_object <- NormalizeData(cervical_object)
## Normalizing layer: counts
####variable feature identification
cervical_object <- FindVariableFeatures(cervical_object, selection.method = "vst", nfeatures = 2000)
## Finding variable features for layer counts
# viewing top10 most variable gene
top10 <- head(VariableFeatures(cervical_object), 10)
variable_plot <- VariableFeaturePlot(cervical_object)
###gene labeling
LabelPoints(plot = variable_plot, points = top10, repel = TRUE)
## When using repel, set xnudge and ynudge to 0 for optimal results
## Warning in scale_x_log10(): log-10 transformation introduced infinite values.

Dimensionality reduction

cervical_object <- ScaleData(cervical_object)
## Centering and scaling data matrix
cervical_object <- RunPCA(cervical_object, features= VariableFeatures(object = cervical_object))
## PC_ 1 
## Positive:  LTB, TRBC1, PTGDS, COL4A2, PIM2, JAML, RGS1, VCAN, ALOX5, FBN1 
##     PGR, ITGA9, KLRK1, BATF, SERPINB9, PDGFRB, CCL5, ITGB2, STAB1, CTLA4 
##     THBS2, ADAMTS9, BLK, MCOLN2, IRF8, ADAMTS4, DUSP2, MS4A1, TWIST2, CD163 
## Negative:  TACSTD2, PRSS8, ELF3, PPL, PRSS22, JUP, SPINT1, KRT80, CLDN4, TMPRSS11E 
##     PITX1, TMPRSS2, LCN2, SDC1, KLF5, NCCRP1, PERP, CDKN2B, EPHA2, TRIM29 
##     DSP, KRT17, SPRR1B, MYH14, MUC21, CLDN1, KRT16, MUC4, TMPRSS11D, NECTIN4 
## PC_ 2 
## Positive:  LYZ, IFI30, SPI1, CD163, ALOX5, ITGAX, C1QA, CD14, TYROBP, C1QC 
##     CSF1R, C1QB, MS4A7, CYBB, CSF2RA, STAB1, MS4A6A, LAIR1, SLCO2B1, NCF2 
##     SLC11A1, OLR1, CLEC7A, MSR1, MRC1, HSPA6, FCER1G, FCGR2A, TLR2, CD68 
## Negative:  BGN, AEBP1, COL4A2, CALD1, SPARC, THBS2, TPM2, TNC, C1S, MMP2 
##     IGFBP7, LAMB1, LUM, MGP, SULF1, COL5A2, COL7A1, SFRP2, C11orf96, PRRX1 
##     FN1, MYL9, COL4A1, CTHRC1, SERPINH1, DCN, FBN1, TAGLN, COL12A1, INHBA 
## PC_ 3 
## Positive:  FN1, VCAN, SPARC, AEBP1, THBS2, BGN, MMP2, LUM, FTH1, SULF1 
##     SFRP2, MGP, COL4A2, THBS1, COL5A2, C1S, PRRX1, C11orf96, CTSB, CCDC71L 
##     CHI3L1, C3, DCN, IGFBP7, TAGLN, PSAP, LAMB1, MYL9, POSTN, RARRES2 
## Negative:  FAT2, SFN, TNS4, TP63, ITGB4, DSC3, COL17A1, S100A2, FXYD3, ITGB6 
##     TPX2, LAMB3, FABP5, FGFR3, RBP1, E2F1, RNF152, NDUFA4L2, KRT14, PLP2 
##     SLC2A1, CRABP2, DSP, KRT8, PTHLH, MCM4, TRIM29, HIST1H1B, CLCA2, CALML3 
## PC_ 4 
## Positive:  IFI30, LYZ, CD163, C1QA, CD14, CTSB, C1QC, GLUL, CTSZ, C1QB 
##     PSAP, APOE, SPI1, TYROBP, CSF1R, MS4A7, CAPG, ITGAX, OLR1, TGFBI 
##     ALOX5, STAB1, PLTP, SLC11A1, SLCO2B1, MRC1, FTL, GPNMB, FCER1G, MSR1 
## Negative:  TMPRSS2, NCCRP1, MUC21, KRT78, PRSS27, EPS8L1, CEACAM5, TMPRSS11B, SCEL, ADGRF1 
##     CEACAM6, SPRR2A, LCN2, PSCA, S100P, KRT80, FUT3, A2ML1, SPRR2D, LYPD2 
##     TMPRSS11E, SPNS2, SLC6A14, PRSS22, SPRR3, MUC4, IL36RN, MAB21L4, CLCA4, SPINK5 
## PC_ 5 
## Positive:  EGFL7, SLCO2A1, MCF2L, AQP1, RAMP3, PTPRB, TMEM255B, VWF, PLVAP, PECAM1 
##     CDH5, CLDN5, ACKR1, NPDC1, FLT1, SPRY1, HOXD4, LIMS2, ADAMTS9, SHANK3 
##     HOXD9, KANK3, MMRN2, ECSCR, CALCRL, CLEC14A, ADCY4, ADGRL4, PODXL, PIK3R3 
## Negative:  THBS2, LUM, C1S, PRRX1, BGN, SFRP4, DCN, AEBP1, PDGFRA, SPON2 
##     GREM1, VCAN, COL12A1, RARRES2, SFRP2, CCDC71L, CHI3L1, MMP11, MXRA8, CCN4 
##     COL5A2, SULF1, TNC, LTB, COL7A1, SRPX2, PDGFRB, C3, NKD2, SERPINF1
####vizulaization
ElbowPlot(cervical_object)

Clustring Analysis to have an idea of gene expression patterns and identify the types of cells present in the tumor microenvironment. This will help us explore the cellular composition of cervical cancer tissue.

####Finding Nearest Neigbour
cervical_object <- FindNeighbors(cervical_object, dims =1:10)
## Computing nearest neighbor graph
## Computing SNN
cervical_object <- FindClusters(cervical_object, resolution = 0.5)
## Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
## 
## Number of nodes: 8404
## Number of edges: 276287
## 
## Running Louvain algorithm...
## Maximum modularity in 10 random starts: 0.9120
## Number of communities: 12
## Elapsed time: 1 seconds
####Vizulizing using  UMAP
cervical_object <- RunUMAP(cervical_object,dims= 1:10)
## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per session
## 12:15:38 UMAP embedding parameters a = 0.9922 b = 1.112
## 12:15:38 Read 8404 rows and found 10 numeric columns
## 12:15:38 Using Annoy for neighbor search, n_neighbors = 30
## 12:15:38 Building Annoy index with metric = cosine, n_trees = 50
## 0%   10   20   30   40   50   60   70   80   90   100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 12:15:39 Writing NN index file to temp file C:\Users\iyand\AppData\Local\Temp\RtmpiY5JNf\file6280290a60ae
## 12:15:39 Searching Annoy index using 1 thread, search_k = 3000
## 12:15:42 Annoy recall = 100%
## 12:15:42 Commencing smooth kNN distance calibration using 1 thread with target n_neighbors = 30
## 12:15:43 Initializing from normalized Laplacian + noise (using RSpectra)
## 12:15:44 Commencing optimization for 500 epochs, with 363192 positive edges
## 12:16:03 Optimization finished
DimPlot(cervical_object, reduction = "umap", label = TRUE)

we need to identify marker gene by annotation Annotating the clusters helps us understand what types of cells are present in the tumor, such as cancerous cells, immune cells, or stromal cells. This gives insights into the tumor microenvironment and cancer biology.

marker <- FindAllMarkers(cervical_object, only.pos= TRUE, min.pct = 0.25, logfc.threshold = 0.25)
## Calculating cluster 0
## For a (much!) faster implementation of the Wilcoxon Rank Sum Test,
## (default method for FindMarkers) please install the presto package
## --------------------------------------------
## install.packages('devtools')
## devtools::install_github('immunogenomics/presto')
## --------------------------------------------
## After installation of presto, Seurat will automatically use the more 
## efficient implementation (no further action necessary).
## This message will be shown once per session
## Calculating cluster 1
## Calculating cluster 2
## Calculating cluster 3
## Calculating cluster 4
## Calculating cluster 5
## Calculating cluster 6
## Calculating cluster 7
## Calculating cluster 8
## Calculating cluster 9
## Calculating cluster 10
## Calculating cluster 11
###display top 10 markers
top_markers <- marker %>% group_by(cluster) %>% top_n(n= 10, wt = avg_log2FC)

FeaturePlot(cervical_object, features = c("EPCAM", "CD3D", "PTPRC", "CD68", "CD163", "MS4A1", "KRT19"))

FeaturePlot(cervical_object, features = "EPCAM")

#FeaturePlot(cervical_object, features = "KRT18")
FeaturePlot(cervical_object, features = "CD3D")

FeaturePlot(cervical_object, features = "PTPRC")

FeaturePlot(cervical_object, features = "CD68")

FeaturePlot(cervical_object, features = "CD163")

FeaturePlot(cervical_object, features = "MS4A1")

FeaturePlot(cervical_object, features = "KRT19")

FeaturePlot(cervical_object, features = "VIM")

####
DoHeatmap(cervical_object, features = top_markers$gene, size = 20)+NoLegend() 
## Warning in DoHeatmap(cervical_object, features = top_markers$gene, size = 20):
## The following features were omitted as they were not found in the scale.data
## slot for the RNA assay: TPM1, PLPP3, ADAM33, COL27A1, CRISPLD2, CLEC2D, TMC8,
## ETS1, SPOCK2, TRBC2, TRAC, ZAP70, IKZF1, IL7R

Looking into the immunological landscape and tumour micro envrionment

library(ggplot2)
library(ggrepel)
immune_markers <- c("CD3D", "CD4", "CD8A", "FOXP3", "CD19", "MS4A1", 
                    "CD68", "CD163", "LYZ", "CD163", "IFI30", "C1QA", "TRBC1", "ITGAX")

marker_names <- c("T-cell receptor CD3 delta", "T-cell co-receptor CD4", 
           "T-cell co-receptor CD8A", "Forkhead box P3", 
           "B-lymphocyte marker CD19", "B-lymphocyte marker MS4A1", 
           "Macrophage marker CD68", "Macrophage marker CD163", 
           "Lysozyme", "Interferon gamma-inducible protein 30", 
           "Complement component C1q A chain", 
           "T-cell receptor beta chain", 
           "Integrin alpha X" )
for (i in seq_along(immune_markers)) {
  plot <- FeaturePlot(cervical_object, features = immune_markers[i], pt.size = 1) +
    ggtitle(marker_names[i]) +  
    theme_minimal() +  
    theme(plot.title = element_text(hjust = 0.5)) 
  print(plot)
}

Quantitative Analysis (Average expression of selected immune marker)

immune_markers <- c("CD3D", "CD4", "CD8A", "FOXP3", "CD19", "MS4A1", 
                    "CD68", "CD163", "LYZ", "CD163", "IFI30", "C1QA", "TRBC1", "ITGAX")
avg_expression <- AverageExpression(cervical_object, features = immune_markers)
## As of Seurat v5, we recommend using AggregateExpression to perform pseudo-bulk analysis.
## First group.by variable `ident` starts with a number, appending `g` to ensure valid variable names
## This message is displayed once per session.
# Convert to a data frame for easier visualization
avg_expression_df <- as.data.frame(avg_expression$RNA)
print(avg_expression_df)
##               g0         g1         g2          g3          g4        g5
## C1QA   0.3996149 0.35619481 0.44643158 0.447712111  8.56844783 0.2830056
## CD8A   5.6470090 0.32634211 0.40406524 0.159537399  0.42612160 0.3094180
## TRBC1 12.2354256 0.84746057 0.89134201 0.292503841  0.68784562 0.3447863
## MS4A1  1.4387382 0.51210539 0.56019587 0.075757379  0.75480796 0.2001858
## CD3D   3.6737802 0.27935001 0.30312693 0.173255185  0.16804421 0.1548049
## CD4    3.1475932 0.26684849 0.30995636 0.257535107  5.50659722 0.3445585
## CD163  0.2526186 0.56025744 0.16882076 0.615739552 13.67868361 0.2251964
## LYZ    0.9224197 0.87096206 0.72021423 1.266798108 29.32077828 1.2352009
## CD19   0.1746619 0.07724593 0.03438756 0.008303717  0.07520176 0.0000000
## ITGAX  1.2334088 0.68304987 0.46638363 0.628477714 35.93755228 0.4507245
## CD68   0.2254856 0.35030460 0.69904925 0.502906109  3.19412913 0.1785921
## IFI30  1.0683809 1.11655199 1.17984637 1.440126733 22.37539611 0.8083363
## FOXP3  2.4217344 0.11480146 0.22098862 0.058838883  0.12847649 0.0244810
##                g6         g7         g8          g9       g10        g11
## C1QA  25.83731260  0.2178339 0.14001310  0.32080730  4.744000 0.00000000
## CD8A   0.10345553  0.7049577 0.08847937  0.07418926  1.831371 0.19778481
## TRBC1  0.28052935  2.0411570 0.10517468  0.81337060  3.525386 0.94651792
## MS4A1  0.02522787 20.9053813 0.00000000  0.36350307  1.716631 0.00000000
## CD3D   0.01764677  0.5969689 0.04725959  0.09058221  1.367831 0.00000000
## CD4    8.13400719  0.6438944 0.09072342 11.22677464  5.287958 0.04876821
## CD163 31.16546445  0.1191825 0.06272044  0.52766648  4.174224 0.23162027
## LYZ   64.60123233  0.8516187 0.75223869  1.77921094 20.266275 0.63669663
## CD19   0.02166796  2.2921282 0.00000000  0.06641517  0.382196 0.00000000
## ITGAX 27.91317370  3.0092579 0.48750027  0.48008222  7.456404 4.41746253
## CD68   9.54198894  0.5302786 3.60192494  1.65635746  1.644860 0.51458018
## IFI30 43.18270773  2.9614044 0.62113940  1.10452243  8.703795 0.34818855
## FOXP3  0.00000000  0.5094456 0.14380659  0.00000000  1.561381 0.15229045

Printed average expression value organized into New DataFrame

### data matrix of the average value
data_matrix <- matrix(c(
  0.3996149, 0.35619481, 0.44643158, 0.447712111, 8.56844783, 0.2830056,
  5.6470090, 0.32634211, 0.40406524, 0.159537399, 0.42612160, 0.3094180,
  12.2354256, 0.84746057, 0.89134201, 0.292503841, 0.68784562, 0.3447863,
  1.4387382, 0.51210539, 0.56019587, 0.075757379, 0.75480796, 0.2001858,
  3.6737802, 0.27935001, 0.30312693, 0.173255185, 0.16804421, 0.1548049,
  3.1475932, 0.26684849, 0.30995636, 0.257535107, 5.50659722, 0.3445585,
  0.2526186, 0.56025744, 0.16882076, 0.615739552, 13.67868361, 0.2251964,
  0.9224197, 0.87096206, 0.72021423, 1.266798108, 29.32077828, 1.2352009,
  0.1746619, 0.07724593, 0.03438756, 0.008303717, 0.07520176, 0.0000000,
  1.2334088, 0.68304987, 0.46638363, 0.628477714, 35.93755228, 0.4507245,
  0.2254856, 0.35030460, 0.69904925, 0.502906109, 3.19412913, 0.1785921,
  1.0683809, 1.11655199, 1.17984637, 1.440126733, 22.37539611, 0.8083363,
  2.4217344, 0.11480146, 0.22098862, 0.058838883, 0.12847649, 0.0244810,
  
  25.83731260, 0.2178339, 0.14001310, 0.32080730, 4.744000, 0.00000000,
  0.10345553, 0.7049577, 0.08847937, 0.07418926, 1.831371, 0.19778481,
  0.28052935, 2.0411570, 0.10517468, 0.81337060, 3.525386, 0.94651792,
  0.02522787, 20.9053813, 0.00000000, 0.36350307, 1.716631, 0.00000000,
  0.01764677, 0.5969689, 0.04725959, 0.09058221, 1.367831, 0.00000000,
  8.13400719, 0.6438944, 0.09072342, 11.22677464, 5.287958, 0.04876821,
  31.16546445, 0.1191825, 0.06272044, 0.52766648, 4.174224, 0.23162027,
  64.60123233, 0.8516187, 0.75223869, 1.77921094, 20.266275, 0.63669663,
  0.02166796, 2.2921282, 0.00000000, 0.06641517, 0.382196, 0.00000000,
  27.91317370, 3.0092579, 0.48750027, 0.48008222, 7.456404, 4.41746253,
  9.54198894, 0.5302786, 3.60192494, 1.65635746, 1.644860, 0.51458018,
  43.18270773, 2.9614044, 0.62113940, 1.10452243, 8.703795, 0.34818855,
  0.00000000, 0.5094456, 0.14380659, 0.00000000, 1.561381, 0.15229045
), nrow = 13, byrow = TRUE)

###Assigning row and column names
rownames(data_matrix) <- c("C1QA", "CD8A", "TRBC1", "MS4A1", "CD3D", "CD4", "CD163", "LYZ", "CD19", "ITGAX", "CD68", "IFI30", "FOXP3")
colnames(data_matrix) <- c("g0", "g1", "g2", "g3", "g4", "g5", "g6", "g7", "g8", "g9", "g10", "g11")

#### creating the dataFrame
expression_df <- as.data.frame(data_matrix)

####Print the table
print(expression_df)
##                g0          g1         g2           g3          g4         g5
## C1QA   0.39961490  0.35619481 0.44643158  0.447712111  8.56844783 0.28300560
## CD8A  12.23542560  0.84746057 0.89134201  0.292503841  0.68784562 0.34478630
## TRBC1  3.67378020  0.27935001 0.30312693  0.173255185  0.16804421 0.15480490
## MS4A1  0.25261860  0.56025744 0.16882076  0.615739552 13.67868361 0.22519640
## CD3D   0.17466190  0.07724593 0.03438756  0.008303717  0.07520176 0.00000000
## CD4    0.22548560  0.35030460 0.69904925  0.502906109  3.19412913 0.17859210
## CD163  2.42173440  0.11480146 0.22098862  0.058838883  0.12847649 0.02448100
## LYZ    0.10345553  0.70495770 0.08847937  0.074189260  1.83137100 0.19778481
## CD19   0.02522787 20.90538130 0.00000000  0.363503070  1.71663100 0.00000000
## ITGAX  8.13400719  0.64389440 0.09072342 11.226774640  5.28795800 0.04876821
## CD68  64.60123233  0.85161870 0.75223869  1.779210940 20.26627500 0.63669663
## IFI30 27.91317370  3.00925790 0.48750027  0.480082220  7.45640400 4.41746253
## FOXP3 43.18270773  2.96140440 0.62113940  1.104522430  8.70379500 0.34818855
##                g6        g7         g8         g9        g10       g11
## C1QA   5.64700900 0.3263421 0.40406524 0.15953740  0.4261216 0.3094180
## CD8A   1.43873820 0.5121054 0.56019587 0.07575738  0.7548080 0.2001858
## TRBC1  3.14759320 0.2668485 0.30995636 0.25753511  5.5065972 0.3445585
## MS4A1  0.92241970 0.8709621 0.72021423 1.26679811 29.3207783 1.2352009
## CD3D   1.23340880 0.6830499 0.46638363 0.62847771 35.9375523 0.4507245
## CD4    1.06838090 1.1165520 1.17984637 1.44012673 22.3753961 0.8083363
## CD163 25.83731260 0.2178339 0.14001310 0.32080730  4.7440000 0.0000000
## LYZ    0.28052935 2.0411570 0.10517468 0.81337060  3.5253860 0.9465179
## CD19   0.01764677 0.5969689 0.04725959 0.09058221  1.3678310 0.0000000
## ITGAX 31.16546445 0.1191825 0.06272044 0.52766648  4.1742240 0.2316203
## CD68   0.02166796 2.2921282 0.00000000 0.06641517  0.3821960 0.0000000
## IFI30  9.54198894 0.5302786 3.60192494 1.65635746  1.6448600 0.5145802
## FOXP3  0.00000000 0.5094456 0.14380659 0.00000000  1.5613810 0.1522904

Barplot

library(reshape2)
expression_df$Gene <- rownames(expression_df)
expression_long <- melt(expression_df, id.vars = "Gene")
mean_expression <- expression_long %>%
  group_by(Gene, variable) %>%
  summarize(Mean = mean(value), .groups = 'drop')

ggplot(mean_expression, aes(x = Gene, y = Mean, fill = variable)) +
  geom_bar(stat = "identity", position = "dodge") +
  labs(title = "Average Gene Expression by Group", x = "Gene", y = "Average Expression Value") +
  theme_minimal() +
  scale_fill_brewer(palette = "Paired") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

BoxPlot

colnames(expression_long) <- c("Gene", "Group", "Expression")
dim(expression_long)
## [1] 156   3
head(expression_long)
ggplot(expression_long, aes(x = Gene, y = Expression)) +
  geom_boxplot() +
  labs(title = "Gene Expression Boxplot", x = "Gene", y = "Expression Value") +
  theme_minimal()

Conclusion Overall Immune Landscape Discussion The expression patterns observed across clusters illustrate a dynamic immune landscape characterized by distinct populations of T cells, macrophages, and potentially B cells. CD68 emerged as the most highly expressed marker, particularly in Clusters 4 and 6, suggesting that macrophages play a central role in the immune response. The presence of both effector T cells (CD3D, CD8A) and regulatory T cells (FOXP3) within Clusters 0 and 4 indicates an active interplay between immune activation and regulation.

FOXP3 and IFI30 expression levels, particularly in Cluster 0, underscore the potential for regulatory mechanisms in modulating immune responses. This duality suggests that while the immune system is active, there are also regulatory controls in place to prevent overactivation and potential tissue damage.

more opinion Cluster 0: Characterized by strong T cell markers (CD3D, CD8A) and moderate FOXP3 expression, suggesting active immune responses and potential cytotoxic activity.

Cluster 1: Displays lower expression levels across markers, possibly representing a resting or less activated state of immune cells.

Cluster 4: Rich in various T cell markers (CD4, TRBC1) and macrophage markers (CD68, CD163), indicating an active immune environment.

Cluster 6: Characterized by high expression of T cell and macrophage markers, suggesting a significant role in orchestrating immune responses.

Cluster 7: Notable for MS4A1 expression, suggesting a specialized immune cell subset, potentially related to B cell activity or alternative macrophage activation.

Clusters 8 and 9: Intermediate expression levels for various markers, suggesting potential transitional states or mixed cell populations.