1. load libraries

2. Load Seurat Object


All_samples_Merged <- readRDS("../../../PHD_3rd_YEAR_Analysis/0-Seurat_RDS_OBJECT_FINAL/All_samples_Merged_with_Renamed_Clusters_final-26-10-2025.rds")

pbmc <- All_samples_Merged

rm(All_samples_Merged)

gc()
             used    (Mb) gc trigger    (Mb)   max used    (Mb)
Ncells   12854706   686.6   24322128  1299.0   21953026  1172.5
Vcells 3677004921 28053.4 6031079703 46013.5 4904709583 37420.0
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. 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. run UMAP and check clusters of all the cells in the object.

RunDimPlot(object = my.pbmc.pi)

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

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

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

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

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

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

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
)
setwd("~/cluster_home/0-Cluster_to_Computer_Transfer_files_folder/IdentifyTopMarkers_and_rename_clusters-26OCt-2025/RAGAS_Visualizarion")


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)

6. Marker analysis and matrix plot (Predictions)


my.pbmc.pi <- RunFindAllMarkers(
  my.pbmc.pi,
  ident = "renamed_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|renamed_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|renamed_clusters|AllMarkers|test.use=wilcox`$data <- markers_clusters_filtered

##. Matrix plot (Predictions)

p2


library(ComplexHeatmap)

# --- Save as PNG ---
png(filename = "Top5_Markers_Heatmap_p2.png", width = 5000, height = 3000, res = 300)
draw(p2)
dev.off()
png 
  2 
# --- Save as PDF ---
pdf(file = "Top5_Markers_Heatmap_p2.pdf", width = 50, height = 30)  # width/height in inches
draw(p2)
dev.off()
png 
  2 

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