FeaturePlot important Genes

Prerequisites: The analysis assumes a Seurat object named All_samples_Merged is loaded, containing both RNA (SCT assay) and ADT (ADT assay) data.

# Replace with the actual path to your 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") 

# Verify the object structure
print(All_samples_Merged)
# Set default assay to SCT for initial RNA-based context
DefaultAssay(All_samples_Merged) <- "SCT"

DimPlot(All_samples_Merged, group.by = "orig.ident", label = T)

1 Featureplots

DefaultAssay(All_samples_Merged) <- "SCT"

Idents(All_samples_Merged) <- "seurat_clusters"

all_markers_SS <- c(

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 1. CD4 T-CELL LINEAGE / MEMORY / ACTIVATION
  #–––––––––––––––––––––––––––––––––––––––––––––
  "CD4","CD3D","CD3E","CD3G","TRAC","TRBC1","TRBC2",
  "IL7R","CCR7","SELL","TCF7","LTB","CD28","CD40LG",
  "TNFRSF4","TNFRSF9","ICOS","TOX2",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 2. EXHAUSTION / CHRONIC ACTIVATION
  #–––––––––––––––––––––––––––––––––––––––––––––
  "PDCD1","CTLA4","LAG3","TIGIT","HAVCR2","ENTPD1",
  "CD244","CXCL13","TOX","BATF",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 3. SÉZARY-SPECIFIC MALIGNANT SIGNATURES
  #–––––––––––––––––––––––––––––––––––––––––––––
  "GATA3","IKZF2","TOX","TOX2","DUSP4","CCR4",
  "KIR3DL2","KLRB1","JUN","JUNB","FOS","FOSB","DUSP2",
  "CD7","DPP4","CD26","CD3D","CD3E","CD3G",
  "NAMPT","LAMP3","HSPA1A","HSPA1B","TNFRSF1B",
  "FCRL3","HLA-DRA","HLA-DRB1","B2M",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 4. MALIGNANT DOWNREGULATED / TUMOR SUPPRESSION LOSS
  #–––––––––––––––––––––––––––––––––––––––––––––
  "SELL","TCF7","CD7","DPP4","STAT4","BACH2","CD26",
  "CD45RA","IL7R","CCR7",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 5. Th1 / Th2 / Th17 / Tfh SUBSET MARKERS
  #–––––––––––––––––––––––––––––––––––––––––––––
  # Th1
  "TBX21","IFNG","CXCR3",
  # Th2
  "GATA3","IL4","IL5","IL13",
  # Th17
  "RORC","IL17A","IL17F","CCR6",
  # Tfh
  "BCL6","CXCL13","PDCD1","ICOS",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 6. TREG PANEL
  #–––––––––––––––––––––––––––––––––––––––––––––
  "FOXP3","IL2RA","CTLA4","TIGIT","IKZF2","CCR8","TNFRSF18",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 7. CELL CYCLE / PROLIFERATION
  #–––––––––––––––––––––––––––––––––––––––––––––
  "MKI67","TOP2A","CENPF","PCNA","MCM4","MCM6",
  "UBE2C","CDK1","AURKA","AURKB","CCNB1","CCNB2",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 8. METABOLIC / STRESS RESPONSE (RELEVANT IN SS)
  #–––––––––––––––––––––––––––––––––––––––––––––
  "NAMPT","LGALS3","HSP90AA1","HSPB1","SOD2","TXN",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 9. CYTOKINES / CHEMOKINES FREQUENTLY ACTIVE IN SS
  #–––––––––––––––––––––––––––––––––––––––––––––
  "CCL5","CCL4","CCL22","CCL17","IL32","IL2","IL7",

  #–––––––––––––––––––––––––––––––––––––––––––––
  # 10. CONTAMINATION CONTROL — MYELOID / NK / B
  #–––––––––––––––––––––––––––––––––––––––––––––
  # Myeloid
  "LYZ","S100A8","S100A9","CD14","FCGR3A","MS4A7","CD68","CD163","MRC1",
  # NK
  "NKG7","GNLY","KLRD1","PRF1","GZMB","TRDC","TRGC2","CD160",
  # B cells / plasma
  "MS4A1","CD79A","CD79B","MZB1","JCHAIN"
)

FeatPlot_canonical <- lapply(all_markers_SS, function(x){
    
    name = x
    print(name)
    jpeg(paste0("immune_featPlot_", name, ".jpeg"), width=10, height=10, units = 'in', res = 300)
    print(FeaturePlot(All_samples_Merged, features = x, reduction = "umap", label = TRUE, min.cutoff = 0, cols = c("navyblue", "darkseagreen","yellow","firebrick1")))
    dev.off()
    
}) 

2 Featureplots

DefaultAssay(All_samples_Merged) <- "SCT"
Idents(All_samples_Merged) <- "orig.ident"

library(viridis)

FeatPlot_canonical <- lapply(all_markers_SS, function(gene){

    message("Plotting: ", gene)

    # Extract gene expression
    expr <- FetchData(All_samples_Merged, vars = gene)[, 1]
    expr[!is.finite(expr)] <- NA   # remove Inf, NaN

    # Case 1: gene not expressed
    if (all(is.na(expr)) || all(expr == 0, na.rm = TRUE)) {

        warning(paste0("Gene ", gene, " has zero/NA expression. Plot will be blank."))

        p <- FeaturePlot(
            All_samples_Merged,
            features = gene,
            reduction = "umap",
            min.cutoff = NA,
            max.cutoff = NA,
            label = TRUE
        ) +
        scale_colour_viridis(option = "magma", na.value = "grey90") +
        ggtitle(paste0(gene, " (no expression)")) +
        theme_void()

        print(p)
        return(NULL)
    }

    # Case 2: normal gene → NO quantile cutoffs (avoids your error)
    p <- FeaturePlot(
        All_samples_Merged,
        features = gene,
        reduction = "umap",
        min.cutoff = NA,
        max.cutoff = NA,
        label = TRUE
    ) +
    scale_colour_viridis(option = "magma", na.value = "grey90") +
    ggtitle(gene) +
    theme_void() +
    theme(
        plot.title = element_text(size = 16, face = "bold"),
        legend.position = "right"
    )

    print(p)
})

3 Featureplots-Clusters

DefaultAssay(All_samples_Merged) <- "SCT"
Idents(All_samples_Merged) <- "seurat_clusters"

library(viridis)

FeatPlot_canonical <- lapply(all_markers_SS, function(gene){

    message("Plotting: ", gene)

    # Extract gene expression
    expr <- FetchData(All_samples_Merged, vars = gene)[, 1]
    expr[!is.finite(expr)] <- NA   # remove Inf, NaN

    # Case 1: gene not expressed
    if (all(is.na(expr)) || all(expr == 0, na.rm = TRUE)) {

        warning(paste0("Gene ", gene, " has zero/NA expression. Plot will be blank."))

        p2 <- FeaturePlot(
            All_samples_Merged,
            features = gene,
            reduction = "umap",
            min.cutoff = NA,
            max.cutoff = NA,
            label = TRUE
        ) +
        scale_colour_viridis(option = "magma", na.value = "grey90") +
        ggtitle(paste0(gene, " (no expression)")) +
        theme_void()

        print(p2)
        return(NULL)
    }

    # Case 2: normal gene → NO quantile cutoffs (avoids your error)
    p2 <- FeaturePlot(
        All_samples_Merged,
        features = gene,
        reduction = "umap",
        min.cutoff = NA,
        max.cutoff = NA,
        label = TRUE
    ) +
    scale_colour_viridis(option = "magma", na.value = "grey90") +
    ggtitle(gene) +
    theme_void() +
    theme(
        plot.title = element_text(size = 16, face = "bold"),
        legend.position = "right"
    )

    print(p2)
})