Data Visualization (Ragas)-cell_state_visualization

1 1. load libraries

2 2. Load Seurat Object

All_samples_Merged <- readRDS("../../../PHD_3rd_YEAR_Analysis/0-Seurat_RDS_OBJECT_FINAL/All_samples_Merged_with_Renamed_Clusters_Cell_state-03-12-2025.rds.rds")

pbmc <- All_samples_Merged

rm(All_samples_Merged)

gc()

             used   (Mb) gc trigger    (Mb)   max used   (Mb)
Ncells   10768195  575.1   19539392  1043.6   13168465  703.3
Vcells 1236377650 9432.9 1535543907 11715.3 1236510487 9433.9

pbmc

An object of class Seurat 
62900 features across 49305 samples within 6 assays 
Active assay: RNA (36601 features, 0 variable features)
 2 layers present: data, counts
 5 other assays present: ADT, prediction.score.celltype.l1, prediction.score.celltype.l2, prediction.score.celltype.l3, SCT
 5 dimensional reductions calculated: integrated_dr, ref.umap, pca, umap, harmony

3 3. The post-integration (Pi) object

my.pbmc.pi <- CreatePostIntegrationObject(object = pbmc)

my.pbmc.pi

An object of class Pi 
6 fields in the object: seurat.obj, exp.freq, markers, ds, cell.prop, parent.meta.data.
The following field has been processed:
    seurat.obj: A Seurat object of 36601 features and 49305 cells.
        6 assays: RNA, ADT, prediction.score.celltype.l1, prediction.score.celltype.l2, prediction.score.celltype.l3, SCT, and 5 reductions: integrated_dr, ref.umap, pca, umap, harmony
Metadata from the parent object provided? No 
Subclusters integrated? No

4 4. run UMAP and check clusters of all the cells in the object.

RunDimPlot(object = my.pbmc.pi)

5 4. run UMAP and check Prediction of all the cells in the object.

RunDimPlot(object = my.pbmc.pi,
           group.by = "cell_state", label = F)

5.1 run UMAP and check Prediction of all the cells in the object.

RunDimPlot(object = my.pbmc.pi,
           group.by = "Prediction",label = F)

5.2 run UMAP and check Prediction of all the cells in the object.

RunDimPlot(object = my.pbmc.pi,
           group.by = "predicted.celltype.l2",label = F)

6 5. Marker analysis and matrix plot(Clusters)

filter_markers <- function(markers_df) {
  blacklist_patterns <- c(
    "^TRAV", "^TRBV", "^TRGV", "^TRDV", "^TRBC", "^TRAC", "^TRDC", "^TRGC", # TCR
    "^IGH", "^IGK", "^IGL", "^IGJ",                                         # Ig genes
    "^RPL", "^RPS",                                                         # ribosomal
    "^MT-",                                                                 # mitochondria
    "^HBA", "^HBB", "^HB[ABZ]",                                             # hemoglobins
    "^NEAT1$", "^MALAT1$",                                                  # optional lncRNAs
    "^XIST$"
  )
  blacklist_regex <- paste(blacklist_patterns, collapse = "|")
  
  markers_filtered <- markers_df %>%
    filter(!grepl(blacklist_regex, gene, ignore.case = TRUE))
  
  return(markers_filtered)
}



my.pbmc.pi <- RunFindAllMarkers(
  my.pbmc.pi,
  ident = "seurat_clusters",
  logfc.threshold = 0.25,
  min.pct = 0.25,
  min.pct.diff = 0.2,
  only.pos = TRUE,
  return.thresh = 0.05
)


markers_clusters <- my.pbmc.pi$markers$`Markers|seurat_clusters|AllMarkers|test.use=wilcox`

# Extract the data frame
markers_clusters_df <- markers_clusters$data



markers_clusters_filtered <- filter_markers(markers_clusters_df)


my.pbmc.pi$markers$`Markers|seurat_clusters|AllMarkers|test.use=wilcox`$data <- markers_clusters_filtered

##. Matrix plot(Clusters)

p1 <- RunMatrixPlot(my.pbmc.pi,
              markers.key = "Markers|seurat_clusters|AllMarkers|test.use=wilcox", 
              column.anno.name.rot = 45, 
              heatmap.height = 7)

Performing relative-counts-normalization

  |                                                                                     
  |                                                                               |   0%
  |                                                                                     
  |===============================================================================| 100%

p1

library(ComplexHeatmap)

# --- Save as PNG ---
png(filename = "Top5_Markers_Heatmap.png", width = 3000, height = 1000, res = 300)

draw(p1)
dev.off()

png 
  2 

# --- Save as PDF ---
pdf(file = "Top5_Markers_Heatmap.pdf", width = 30, height = 10)  # width/height in inches
draw(p1)

dev.off()

png 
  2 

message("✅ Heatmap saved successfully as both PNG and PDF.")

7 6. Marker analysis and matrix plot (Predictions)

my.pbmc.pi <- RunFindAllMarkers(
  my.pbmc.pi,
  ident = "cell_state",
  logfc.threshold = 0.25,
  min.pct = 0.25,
  min.pct.diff = 0.2,
  only.pos = TRUE,
  return.thresh = 0.05
)


markers_clusters <- my.pbmc.pi$markers$`Markers|renamed_clusters|AllMarkers|test.use=wilcox`

# Extract the data frame
markers_clusters_df <- markers_clusters$data



markers_clusters_filtered <- filter_markers(markers_clusters_df)

Error in UseMethod("filter") : 
  no applicable method for 'filter' applied to an object of class "NULL"

##. Matrix plot (Predictions)

p2 <- RunMatrixPlot(my.pbmc.pi,
              markers.key = "Markers|cell_state|AllMarkers|test.use=wilcox", 
              column.anno.name.rot = 90, 
              heatmap.height = 9)

Performing relative-counts-normalization

  |                                                                                     
  |                                                                               |   0%
  |                                                                                     
  |===============================================================================| 100%

p2


library(ComplexHeatmap)

# --- Save as PNG ---
png(filename = "Top5_Markers_Heatmap_p2_cell_state.png", width = 5000, height = 3000, res = 300)

draw(p2)
dev.off()

png 
  2 

# --- Save as PDF ---
pdf(file = "Top5_Markers_Heatmap_p2_cell_state.pdf", width = 50, height = 30)  # width/height in inches
draw(p2)

dev.off()

png 
  2 

message("✅ Heatmap saved successfully as both PNG and PDF.")