作者，Evil Genius

關于RCTD，本來都不打算更新了，R版本的聯(lián)合分析我覺得大家自己寫寫就完了，現(xiàn)在看來，還是需要整理一下。

文章在Robust decomposition of cell type mixtures in spatial transcriptomics | Nature Biotechnology

參考鏈接在GitHub - dmcable/spacexr: Spatial-eXpression-R: Cell type identification (including cell type mixtures) and cell type-specific differential expression for spatial transcriptomics

至于引用該方法的文章，那就很多了。

RCTD有三種模式:

（1）doublet mode:每個spot分配1-2種細胞類型，推薦用于具有高空間分辨率的技術，如Stereo-seq、HD等（注意bin的大?。?

（2）full mode:每個spot分配任意數(shù)量的細胞類型，推薦用于具有低空間分辨率的技術，如Visium;

（3）multi mode : doublet mode的擴展，可以每個spot發(fā)現(xiàn)兩個以上的細胞類型，作為全模式的替代選項。

很多博主、公司推文寫了非常多的介紹和代碼示例，大家可以借鑒一下。

不過對于這些網絡寫手，我更喜歡網絡俠客這個稱呼。

官網的教程寫了很多，針對不同精度的平臺都有，列舉一下：

• Doublet mode: spatial transcriptomics vignette. Also, most other vignettes use doublet mode.
• Doublet mode on MERFISH: MERFISH nonparametric vignette.
• Full mode: full mode on Visium hippocampus
• Multi mode: multi mode on Visium hippocampus

不過官網分析的結果是真的丑

Seurat官網也借鑒了RCTD的方法，在Analysis, visualization, and integration of Visium HD spatial datasets with Seurat ? Seurat

示例代碼如下，采用了Doublet mode模式

library(spacexr)

# set up reference
ref <- readRDS("../data/mouse_hippocampus_reference.rds")
ref <- UpdateSeuratObject(ref)
Idents(ref) <- "celltype"

# extract information to pass to the RCTD Reference function
counts <- ref[["RNA"]]$counts
cluster <- as.factor(ref$celltype)
names(cluster) <- colnames(ref)
nUMI <- ref$nCount_RNA
names(nUMI) <- colnames(ref)
reference <- Reference(counts, cluster, nUMI)

# set up query with the RCTD function SpatialRNA
slide.seq <- SeuratData::LoadData("ssHippo")
counts <- slide.seq[["Spatial"]]$counts
coords <- GetTissueCoordinates(slide.seq)
colnames(coords) <- c("x", "y")
coords[is.na(colnames(coords))] <- NULL
query <- SpatialRNA(coords, counts, colSums(counts))


RCTD <- create.RCTD(query, reference, max_cores = 8)
RCTD <- run.RCTD(RCTD, doublet_mode = "doublet")
slide.seq <- AddMetaData(slide.seq, metadata = RCTD@results$results_df)

RCTD目前大多數(shù)用在高精度的平臺，比如Stereo-seq的bin20、bin30的情況，HD的8um情況等。

但是對于真正的文章，一般都要根據自己的課題進行修改。

下方是一個文章分析代碼示例

# Load required libraries
library(spacexr)
library(Seurat)
library(ggplot2)
library(dplyr)
library(tidyr)
library(pheatmap)
library(progeny)

# Data Preprocessing

# Load the Spatial data
pdac <- readRDS("/path/to/pdac_most_updated.rds")

# Load the Single cell data
sc <- readRDS("/path/to/sc_rctd.rds")

# Prepare data for RCTD
counts <- sc@assays$RNA@counts
Idents(sc) <- "SCT"
cluster <- as.factor(sc$celltype_nicheDE)
names(cluster) <- colnames(sc)
nUMI <- sc$nCount_RNA
names(nUMI) <- colnames(sc)

# Create the reference object
reference <- Reference(counts, cluster, nUMI)

# Prepare spatial transcriptomics data
counts <- pdac@assays$Spatial@counts

coordinates_list <- lapply(image_names, function(image_name) {
  pos <- GetTissueCoordinates(pdac, image = image_name)
  colnames(pos) <- c('x','y')
  return(pos)
})

coords <- do.call(rbind, coordinates_list)
rownames(coords) <- gsub("^.+\\.", "", rownames(coords))

# Create SpatialRNA object
query <- SpatialRNA(coords, counts, colSums(counts))

# Run RCTD
RCTD <- create.RCTD(query, reference, max_cores = 8)
RCTD.full <- run.RCTD(RCTD, doublet_mode = "full")

# Add RCTD results to Seurat object
pdac <- AddMetaData(pdac, metadata = RCTD.full@results$results_df)

# Normalize weights
weights <- RCTD.full@results$weights
norm_weights <- normalize_weights(weights)

# Add RCTD results as a new assay
pdac[["rctd_full"]] <- CreateAssayObject(data = t(as.matrix(norm_weights)))
if (length(pdac@assays$rctd_full@key) == 0) {
    pdac@assays$rctd_full@key <- "rctd_full_"
}

這樣就獲取了RCTD的解卷積空間細胞矩陣，當然了，通常分析到這里還沒結束，我們需要繼續(xù)分析空間細胞聚類，共定位等分析內容。

我們把RCTD、空間細胞聚類、共定位分析一起封裝起來，注意要適用各種空間平臺，寫好參數(shù)的設定。跟cell2location一樣的，最好匹配的樣本聯(lián)合分析。

#! usr/bin/R
### zhaoyunfei
### 20241111

suppressMessages({
library(Seurat)
library(compositions)
library(tidyverse)
library(clustree)
library(patchwork)
library(argparse)
library(robustbase)
library(ISCHIA)
library(factoextra)
library(dplyr)
library(scran)
library(ggplot2)
library(spacexr)
library(cluster)
library(showtext)
library(gridExtra)
library(pdftools)
})

parser = ArgumentParser()
parser$add_argument("--sc_rds", help="the sc data",required = T)
parser$add_argument("--spatial_rds", help="the sp data",required = T)
parser$add_argument("--sample", help="the sample name",required = T)
parser$add_argument("--outdir", help="the outdir",default = './')
parser$add_argument("--celltype", help="the annotation for celltype")
parser$add_argument("--normalization_method",default = 'SCT')
parser$add_argument("--reduction",help='Dimensional reduction to perform when finding anchors',choices = c('pcaproject','cca'),default = 'pcaproject')
parser$add_argument("--mode"，help='mode for RCTD',default = 'doublet',choices = c('doublet','full','multi'))
parser$add_argument("--cell_interest", help="the interest of cell type,eg : 'T,B:Fobro,epi'")
args <- parser$parse_args()
str(args)

sc_rds = args$sc_rds
spatial_rds = args$spatial_rds
outdir = args$outdir
celltype = args$celltype
normalization_method = args$normalization_method
reduction = args$reduction
sample = args$sample
mode = args$mode
cell_interest = args$cell_interest

接下來就是腳本主體