Dataset DOI: 10.5061/dryad.s4mw6m9jh

Description of the data and file structure

This repository contains a Seurat object (named “seurat_object”) generated from CosMx SMI using the manufacturer's 1000 plex panel. The object enables researchers to reproduce or extend the analyses performed in our study. All cell metadata has been de-identified according to our IRB-approved protocol, and no personally identifiable information (PII) is included. The data comes from 2 tissue microarrays of formalin fixed paraffin embedded tissue of human cutaneous SCC specimens taken before, during, or after immunotherapy. 

Files and variables

Seurat_Object_cSCC.RDS

 – A serialized R object in Seurat format, which includes: 
   – Filtered and normalized gene expression data. 
   – Cell-level metadata (excluding direct identifiers). 
   – Cluster assignments and other annotations used in the study.

Data De-identification

• In compliance with applicable regulations and IRB guidelines, all patient‐specific identifiers have been removed or replaced with hashed or random codes.

• Only variables relevant to the analyses are retained within the Seurat object.

Using the data

To load the Seurat object in R:

library(Seurat)
# Replace "path/to/" with the local path where you placed the .rds file:
seurat_object <- readRDS("path/to/seurat_object.rds")`

After loading, you can explore gene expression levels, run Seurat’s built-in visualization tools, or incorporate the object into your own pipelines.

Requirements

• R (version 4.0 or higher recommended)

• Seurat (version 4 or higher recommended)

Code/software

Creating_object_and_clustering_cells.R

This code demonstrates how to preprocess, transform, and cluster single‐cell (or spatial) transcriptomic data using Seurat. It includes filtering low‐count cells, performing SCTransform, detecting/removing doublets, and using PCA/UMAP for dimension reduction. The pipeline also applies InSituType for immunofluorescence‐based classification and identifies top differentially expressed markers in each cluster for downstream biological insights.

GSEA.R

This code demonstrates a pipeline for performing spatially resolved GSEA, starting from normalized gene expression data and metadata about cell locations. It uses GSVA to compute per‐cell and summarized enrichment scores for hallmark gene sets and then compares them across samples, clusters, and time points. Finally, the code generates heatmaps, plots, and statistical outputs to help interpret how these gene sets vary within and across patients and biopsy timepoints.

InSituCor_spatial_gene_modules.R

This code uses InSituCor to identify gene modules that are spatially co‐expressed across the tissue. After normalizing and preparing the data (including single‐cell metadata and neighboring information), it fits correlation modules and calculates both per‐cell and environment scores. Finally, it visualizes the modules via heatmaps and correlation networks, providing a functional look at how gene expression groups are spatially organized.

Moran_s_I.R

This code calculates bivariate Moran’s I statistics to measure spatial autocorrelation between ligand–receptor pairs (e.g., IFN-γ pathways) across multiple fields of view (FOVs). It subsets metadata and gene expression, identifies nearest neighbors for each cell, and then computes Moran’s I along with bootstrap confidence intervals. Finally, the script merges these results with patient metadata for comparison and visualization of spatial ligand–receptor interaction patterns by sample or clinical grouping.

ROC_analysis.R

This code performs ROC (receiver operating characteristic) analysis to assess the ability of various features—namely the proportion of “favorable niches” (spatial cluster composition) and PD‐L1 expression levels—to distinguish between responders and nonresponders. It merges and filters relevant metadata, calculates AUC values via the pROC package, identifies optimal cutoffs (sensitivity/specificity), and uses statistical tests (e.g., Wilcoxon) to compare distributions of those features by response category, ultimately creating ROC plots to visualize classifier performance.

Spatial_Clustering.R

This code demonstrates an approach for clustering cells based on their local neighborhood gene‐expression profiles. It first computes a nearest‐neighbor graph, extracts neighborhood‐level expression values (e.g., average per gene in each cell’s local neighborhood), and then performs PCA on these neighborhood features. Next, it uses Mclust to define six spatial clusters, which are integrated into the Seurat metadata. Finally, the script generates bar plots and summary tables examining how these clusters vary across biopsy timepoints, patients, and clinical response groups, and identifies spatial‐cluster‐specific marker genes.

Transcriptome_Correlation.R

This code computes an inter-patient Pearson correlation matrix of average transcriptome profiles. Specifically, it calculates gene‐expression averages for each patient and then generates a correlation heatmap, splitting rows and columns by response type. This helps visualize which patients have similar global expression patterns and whether these clusters align with clinical outcomes.

Human subjects data

All patient data used in this study were de-identified in accordance with institutional review board (IRB) requirements and relevant data privacy regulations (e.g., HIPAA). Specifically, the research team removed or replaced all direct and indirect identifiers, including names, dates of birth, addresses, and medical record numbers, with coded pseudonyms or aggregated values. No link to identifiable personal information was retained, ensuring that the data cannot be traced back to individual participants. All de-identified data were stored in secure, access-controlled repositories with strict oversight, and were only used for the purposes of this study.