SeekSpace小鼠脑demo数据分析¶
背景信息¶
本⽂档是关于寻因SeekSpace技术所得到的数据及其seurat分析⽅法的说明。
寻因SeekSpace技术可以检测到单细胞精度的基因表达数据,同时也能定位出每个细胞在组织中的空间坐标。
SeekSpace的数据分析简便易学,能够很方便地兼容常见单细胞转录组的分析软件,例如seurat和scanpy。
SeekSpace技术中一个像素点的大小是约为 0.2653 微米,将像素点的坐标乘以 0.2653 即可转换计算细胞在真实空间上的距离。
本文档的输入数据为小鼠脑数据的分析结果。首先使用SeekSpace tools,在默认参数下进行分析;然后,使用seurat进行聚类;进一步,将具有高比例线粒体cluster的细胞过滤掉,过滤后的结果作为本文档输入数据。
⽂件说明¶
SeekSpace技术的配套的基础数据处理的软件是SeekSpace tools,该软件可以从测序文库中识别细胞的表达信息,同时可以识别出每个细胞在空间上的位置。
经过SeekSpace tools软件处理之后所得到的结果文件格式如下:
├── WTH1092_filtered_feature_bc_matrix # 表达矩阵⽬录,可以使⽤seurat的Read10X命令进⾏读取
│ ├── barcodes.tsv.gz
│ ├── features.tsv.gz
│ ├── matrix.mtx.gz
│ └── cell_locations.tsv.gz 为⼩⿏脑测序数据中细胞的空间坐标⽂件。第1列是barcode,顺序与matrix中的filtered_feature_bc_matrix/barcode中细胞的顺序⼀致;第2列和第3列分别为该barcode所代表的细胞的空间位置(即空间芯⽚上的像素坐标)。
├── WTH1092_aligned_DAPI.png # 为⼩⿏脑组织切⽚的DAPI染⾊图⽚
├── WTH1092_aligned_HE_TIMG.png # 为⼩⿏脑组织切⽚的DAPI染⾊图⽚
├── seekspace_of_seurat.ipynb # 为使⽤seurat分析该⼩⿏脑空间数据的jupyter示例⽂件
└── seekspace_of_scanpy.ipynb # 为使⽤scanpy分析该⼩⿏脑空间数据的jupyter示例⽂件
SeekSpace技术中一个像素点的大小是约为 0.2653 微米,将像素点的坐标乘以 0.2653 即可转换计算细胞在真实空间上的距离。
library(Seurat)
library(tidyverse)
library(ggplot2)
Attaching SeuratObject Seurat v4 was just loaded with SeuratObject v5; disabling v5 assays and validation routines, and ensuring assays work in strict v3/v4 compatibility mode ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ── ✔ dplyr 1.1.4 ✔ readr 2.1.5 ✔ forcats 1.0.0 ✔ stringr 1.5.1 ✔ ggplot2 3.5.1 ✔ tibble 3.2.1 ✔ lubridate 1.9.3 ✔ tidyr 1.3.1 ✔ purrr 1.0.2 ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ── ✖ dplyr::filter() masks stats::filter() ✖ dplyr::lag() masks stats::lag() ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
创建seurat对象¶
读取SeekSpace的矩阵⽂件,这里直接使用seurat常用的读入和创建矩阵的方法即可。
mouse_brain.data <- Read10X('./Outs/WTH1092_filtered_feature_bc_matrix')
mouse_brain <- CreateSeuratObject(counts=mouse_brain.data,project='mouse_brain')
降维聚类¶
这里使用常见的默认参数进行降维聚类,实际分析的时候,可以按照样品的特征进行相应的修改。
mouse_brain <- NormalizeData(mouse_brain, normalization.method = "LogNormalize", scale.factor = 10000) %>%
FindVariableFeatures(selection.method = "vst", nfeatures = 2000) %>%
ScaleData() %>%
RunPCA() %>%
FindNeighbors(dims = 1:15) %>%
FindClusters(resolution = seq(0.2, 1.4, 0.3)) %>%
RunUMAP(dims = 1:15)
Idents(mouse_brain) <- 'RNA_snn_res.0.8'
Centering and scaling data matrix PC_ 1 Positive: Rbfox1, Kcnip4, Grin2b, Snap25, Tenm2, Grin2a, Nrxn3, Atp1a3, Opcml, Dlgap2 Ppp3r1, Erc2, Atp1b1, Lrrc7, Grm5, Olfm1, Ncdn, Nrgn, Aldoa, Cacna1e Cntnap2, Atp2b1, Cnksr2, Enc1, Ppp3cb, Ppp3ca, Kctd16, Nell2, Frmpd4, Slit3 Negative: Plp1, Trf, Epas1, Glul, Mag, Atp1a2, Cldn11, Mbp, Cnp, Selenop Mal, Mobp, Ptgds, Flt1, Prr5l, Slc1a3, Gab1, St18, Car2, Fa2h Apod, Myrf, Cryab, Plpp3, Slco1c1, Mog, Npas3, Kcnj10, Ptprb, Ugt8a PC_ 2 Positive: Nrgn, 2010300C02Rik, Enc1, Ptk2b, Cnksr2, Wipf3, Cacna1e, Slc8a2, Neurod2, Hpca Ppp3ca, Tafa1, Slit3, Nptx1, Camk2a, Cpne6, Chrm1, Phactr1, Zbtb20, Egr3 Iqgap2, Nptxr, Slc17a7, Ano3, Rprml, C1ql3, Neurod6, Pde1a, Slit1, Dlgap2 Negative: Rbms3, Tcf7l2, Ntng1, Synpo2, Kcnc1, Zfhx3, Kcnc2, Btbd11, Nxph1, Prkcd Grik1, Grip1, Robo1, Gad1, Pcp4l1, Nap1l5, Gad2, Cit, Inpp4b, Ank1 Slc32a1, Slc17a6, Amotl1, 6330411D24Rik, Alk, Rgs16, Adarb2, Adarb1, Rcan2, Gria4 PC_ 3 Positive: Plp1, Tmeff2, Mag, Pde4b, Trf, Cldn11, St18, Cnp, Fa2h, Mal Myrf, Enpp2, Mbp, Prr5l, Slc24a2, Mog, Dock10, Cryab, Mobp, Tspan2 Ugt8a, Gpr37, Rnf220, Plekhh1, Pex5l, Ptprd, D7Ertd443e, Phldb1, Gjc3, Gm16168 Negative: Flt1, Epas1, Ptprb, Rgs5, Slc1a2, Atp1a2, Utrn, Slco1a4, Fn1, Pltp Slco1c1, Eng, Ebf1, Tns1, Vtn, Ahnak, Myl9, Igfbp7, Podxl, Adgrl4 Cxcl12, Itga1, Cdh5, Igf2, Adgrf5, Ltbp4, Pitpnc1, Mecom, Pdgfrb, Col4a1 PC_ 4 Positive: Dpp10, Egr1, Ngef, Mef2c, Nrg1, Satb2, Pdzrn3, Kcnh7, Bc1, Car10 Stx1a, Tbr1, Sgcz, Atp1a1, Camk2n1, Homer1, Tmem132d, Lamp5, Efna5, Unc5d Rgs6, Ipcef1, Slc1a2, Kcnq5, Kcnh5, Atp2b4, Etl4, Garnl3, Snap25, Nrgn Negative: Shisa6, Zbtb20, Nr3c2, Adarb2, Gria1, Tanc1, Galnt17, Prox1, Epha6, Gm20754 St6galnac5, Nkain3, Cacnb2, Slc44a5, Wipf3, Epha7, Sema5a, Lsamp, Iqgap2, Sipa1l3 Grm1, C1ql2, Dscaml1, Synpr, Slit1, Glis3, Ryr3, Pip5k1b, Dpp6, Dgkg PC_ 5 Positive: Lsamp, Gpc5, Npas3, Gm3764, Gad1, Trpm3, Sox6, Luzp2, Slc6a1, Gad2 Nxph1, Ptprt, Ptprz1, Slc32a1, Slc4a4, Rgs20, Gli2, Erbb4, Plpp3, Grip1 Mertk, Maml2, Slc1a3, Gli3, Egfr, Slc6a11, Nwd1, Reln, Htra1, Kcnmb2 Negative: Lef1, Flt1, Ptprb, Rgs5, Synpo2, Epas1, Rnf220, Heg1, Slco1a4, Fn1 Plekhg1, Prkcd, Plp1, Adgrf5, Igfbp7, Kitl, Eng, Pakap.1, Vtn, Ebf1 Ltbp4, Pltp, Myl9, Cxcl12, Podxl, Jcad, Tns1, Adgrl4, Cdh5, Amotl1 Computing nearest neighbor graph Computing SNN
Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck Number of nodes: 22870 Number of edges: 834752 Running Louvain algorithm... Maximum modularity in 10 random starts: 0.9677 Number of communities: 15 Elapsed time: 3 seconds Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck Number of nodes: 22870 Number of edges: 834752 Running Louvain algorithm... Maximum modularity in 10 random starts: 0.9429 Number of communities: 17 Elapsed time: 3 seconds Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck Number of nodes: 22870 Number of edges: 834752 Running Louvain algorithm... Maximum modularity in 10 random starts: 0.9234 Number of communities: 22 Elapsed time: 3 seconds Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck Number of nodes: 22870 Number of edges: 834752 Running Louvain algorithm... Maximum modularity in 10 random starts: 0.9076 Number of communities: 26 Elapsed time: 3 seconds Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck Number of nodes: 22870 Number of edges: 834752 Running Louvain algorithm... Maximum modularity in 10 random starts: 0.8952 Number of communities: 32 Elapsed time: 3 seconds
Warning message: “The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation' This message will be shown once per session” 10:46:03 UMAP embedding parameters a = 0.9922 b = 1.112 10:46:03 Read 22870 rows and found 15 numeric columns 10:46:03 Using Annoy for neighbor search, n_neighbors = 30 10:46:03 Building Annoy index with metric = cosine, n_trees = 50 0% 10 20 30 40 50 60 70 80 90 100% [----|----|----|----|----|----|----|----|----|----| * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * | 10:46:06 Writing NN index file to temp file /tmp/RtmpCReRpJ/file10268882ade 10:46:06 Searching Annoy index using 1 thread, search_k = 3000 10:46:15 Annoy recall = 100% 10:46:15 Commencing smooth kNN distance calibration using 1 thread with target n_neighbors = 30 10:46:16 Initializing from normalized Laplacian + noise (using irlba) 10:46:19 Commencing optimization for 200 epochs, with 1000636 positive edges 10:46:33 Optimization finished
添加空间坐标¶
接下来我们使用seurat的CreateDimReducObject函数,给seurat对象中的每个细胞添加一个空间的坐标。
注意!!
CreateDimReducObject函数读入空间坐标矩阵的细胞的顺序必须与mouse_brain里面的meta.data里面细胞的顺序相一致,否则在多样品整合分析的时候会出现细胞顺序错乱的情况。
这里我们使用了一个简单的一行代码对空间坐标矩阵进行了排序。
spatial_df <- read.table('./Outs/WTH1092_filtered_feature_bc_matrix/cell_locations.tsv.gz', row.names = 1, sep = '\t',header = T)
colnames(spatial_df) <- c("spatial_1","spatial_2")
spatial_matrix <- as.matrix(spatial_df)
spatial_matrix_sorted <- spatial_matrix[match(row.names(mouse_brain@meta.data),row.names(spatial_matrix)), ]
mouse_brain@reductions$spatial <- CreateDimReducObject(embeddings = spatial_matrix_sorted, key='spatial_', assay='RNA')
head(mouse_brain@reductions$spatial@cell.embeddings)
spatial_1 | spatial_2 | |
---|---|---|
AAGTCGTCGACTCTAGAGCTCTGCTTC | 34480 | 17330 |
ACAGGGATACTTCCAGTCCTAATAACG | 31634 | 12314 |
ACCACAAATGGCGTTTCGCTCTGCTTC | 13407 | 16626 |
ACTATGGTACAACTCCTCTAAGTACGA | 29719 | 4091 |
AGACTCAATGTTCTGGCTAGTCATACG | 32670 | 4484 |
AGCCAATCGAACGGGATCCTAATAACG | 13874 | 15676 |
经过这步处理之后,我们可以在seurat对象中看到一个新的坐标系”spatial“,这个坐标系即为每个细胞在空间位置上的坐标。
str(mouse_brain@reductions$spatial)
Formal class 'DimReduc' [package "SeuratObject"] with 9 slots ..@ cell.embeddings : int [1:22870, 1:2] 34480 31634 13407 29719 32670 13874 19186 18630 22301 37179 ... .. ..- attr(*, "dimnames")=List of 2 .. .. ..$ : chr [1:22870] "AAGTCGTCGACTCTAGAGCTCTGCTTC" "ACAGGGATACTTCCAGTCCTAATAACG" "ACCACAAATGGCGTTTCGCTCTGCTTC" "ACTATGGTACAACTCCTCTAAGTACGA" ... .. .. ..$ : chr [1:2] "spatial_1" "spatial_2" ..@ feature.loadings : num[0 , 0 ] ..@ feature.loadings.projected: num[0 , 0 ] ..@ assay.used : chr "RNA" ..@ global : logi FALSE ..@ stdev : num(0) ..@ jackstraw :Formal class 'JackStrawData' [package "SeuratObject"] with 4 slots .. .. ..@ empirical.p.values : num[0 , 0 ] .. .. ..@ fake.reduction.scores : num[0 , 0 ] .. .. ..@ empirical.p.values.full: num[0 , 0 ] .. .. ..@ overall.p.values : num[0 , 0 ] ..@ misc : list() ..@ key : chr "spatial_"
添加组织图像信息¶
下面我们展示一下在空间数据上添加背景图片。
samplename = 'WTH1092'
size_x = 55128
size_y = 19906
#如果芯片编号的倒数第二个字母是A,则size_x=55128,size_y=19906;
#如果芯片编号的倒数第二个字母是B,则size_x=55050,size_y=19906;
#如果芯片编号的倒数第二个字母是C,则size_x=55050,size_y=19906.
#(可以在samplename_report.html质控报告的Summary Table中查看芯片编号)
mouse_brain@misc$info[[`samplename`]]$size_x = as.integer(size_x)
mouse_brain@misc$info[[`samplename`]]$size_y = as.integer(size_y)
添加DAPI图片¶
img = './Outs/WTH1092_aligned_DAPI.png'
#base64格式
img_64 = base64enc::dataURI(file = img)
mouse_brain@misc$info[[`samplename`]]$img = img_64
#png格式
img_gg <- png::readPNG(img)
img_grob <- grid::rasterGrob(img_gg, interpolate = FALSE, width = grid::unit(1,"npc"), height = grid::unit(1, "npc"))
mouse_brain@misc$info[[`samplename`]]$img_gg = img_grob
添加HE图片¶
img = './Outs/WTH1092_aligned_HE_TIMG.png'
#base64格式
img_64 = base64enc::dataURI(file = img)
mouse_brain@misc$info[[`samplename`]]$img_he = img_64
#png格式
img_gg <- png::readPNG(img)
img_grob <- grid::rasterGrob(img_gg, interpolate = FALSE, width = grid::unit(1,"npc"), height = grid::unit(1, "npc"))
mouse_brain@misc$info[[`samplename`]]$img_he_gg = img_grob
str(mouse_brain@misc$info)
List of 1 $ WTH1092:List of 6 ..$ size_x : int 55128 ..$ size_y : int 19906 ..$ img : chr "data:;base64,iVBORw0KGgoAAAANSUhEUgAAD6AAAAWnCAAAAAByYdMoAAAgAElEQVR4AezBaa9l6WGe5/t537X22vOZ6lSdququrq5ukt2kRI"| __truncated__ ..$ img_gg :List of 12 .. ..$ raster : 'raster' chr [1:1447, 1:4000] "#000000" "#000000" "#000000" "#000000" ... .. ..$ x : 'simpleUnit' num 0.5npc .. .. ..- attr(*, "unit")= int 0 .. ..$ y : 'simpleUnit' num 0.5npc .. .. ..- attr(*, "unit")= int 0 .. ..$ width : 'simpleUnit' num 1npc .. .. ..- attr(*, "unit")= int 0 .. ..$ height : 'simpleUnit' num 1npc .. .. ..- attr(*, "unit")= int 0 .. ..$ just : chr "centre" .. ..$ hjust : NULL .. ..$ vjust : NULL .. ..$ interpolate: logi FALSE .. ..$ name : chr "GRID.rastergrob.11" .. ..$ gp : list() .. .. ..- attr(*, "class")= chr "gpar" .. ..$ vp : NULL .. ..- attr(*, "class")= chr [1:3] "rastergrob" "grob" "gDesc" ..$ img_he : chr "data:;base64,iVBORw0KGgoAAAANSUhEUgAAD6AAAAWnCAIAAADYaBujAAAgAElEQVR4AZThUZYlx4JlR8rW4Jy5+NvzJBcHkq6HZtcdQACJV5"| __truncated__ ..$ img_he_gg:List of 12 .. ..$ raster : 'raster' chr [1:1447, 1:4000] "#DEDDE0" "#777982" "#74757F" "#727179" ... .. ..$ x : 'simpleUnit' num 0.5npc .. .. ..- attr(*, "unit")= int 0 .. ..$ y : 'simpleUnit' num 0.5npc .. .. ..- attr(*, "unit")= int 0 .. ..$ width : 'simpleUnit' num 1npc .. .. ..- attr(*, "unit")= int 0 .. ..$ height : 'simpleUnit' num 1npc .. .. ..- attr(*, "unit")= int 0 .. ..$ just : chr "centre" .. ..$ hjust : NULL .. ..$ vjust : NULL .. ..$ interpolate: logi FALSE .. ..$ name : chr "GRID.rastergrob.12" .. ..$ gp : list() .. .. ..- attr(*, "class")= chr "gpar" .. ..$ vp : NULL .. ..- attr(*, "class")= chr [1:3] "rastergrob" "grob" "gDesc"
细胞注释结果¶
接下来我们再把细胞的注释结果添加到每个细胞上,这样我们就可以看不同的细胞类型在空间上的分布情况了。
SeekSpace细胞注释的过程,跟普通的单细胞注释的过程完全一致。
我们之前已经对demo数据进行了大群和亚群的注释,注释的结果文件放在了anntation.csv文件中。
接下来我们简单的读入处理一下。
anno <- read.csv("./annotation.csv",row.names = 1,header = TRUE)
mouse_brain <- AddMetaData(mouse_brain, anno)
head(mouse_brain)
orig.ident | nCount_RNA | nFeature_RNA | RNA_snn_res.0.2 | RNA_snn_res.0.5 | RNA_snn_res.0.8 | RNA_snn_res.1.1 | RNA_snn_res.1.4 | seurat_clusters | sc_barcode | Main_CellType | Sub_CellType | |
---|---|---|---|---|---|---|---|---|---|---|---|---|
<fct> | <dbl> | <int> | <fct> | <fct> | <fct> | <fct> | <fct> | <fct> | <chr> | <chr> | <chr> | |
AAGTCGTCGACTCTAGAGCTCTGCTTC | mouse_brain | 200 | 156 | 11 | 12 | 15 | 15 | 17 | 17 | AAGTCGTCGACTCTAGAGCTCTGCTTC | OPC | OPC |
ACAGGGATACTTCCAGTCCTAATAACG | mouse_brain | 200 | 146 | 3 | 3 | 1 | 0 | 2 | 2 | ACAGGGATACTTCCAGTCCTAATAACG | Astro | Astro |
ACCACAAATGGCGTTTCGCTCTGCTTC | mouse_brain | 200 | 145 | 14 | 16 | 21 | 24 | 29 | 29 | ACCACAAATGGCGTTTCGCTCTGCTTC | Inh | Ependymal |
ACTATGGTACAACTCCTCTAAGTACGA | mouse_brain | 200 | 157 | 1 | 1 | 11 | 7 | 8 | 8 | ACTATGGTACAACTCCTCTAAGTACGA | Oligo | Oligo |
AGACTCAATGTTCTGGCTAGTCATACG | mouse_brain | 200 | 120 | 4 | 15 | 18 | 16 | 23 | 23 | AGACTCAATGTTCTGGCTAGTCATACG | Ext | nb |
AGCCAATCGAACGGGATCCTAATAACG | mouse_brain | 200 | 135 | 1 | 1 | 3 | 7 | 5 | 5 | AGCCAATCGAACGGGATCCTAATAACG | Oligo | Oligo |
AGGTTGTCGAGCGTATTGCTCTGCTTC | mouse_brain | 200 | 150 | 3 | 3 | 1 | 0 | 2 | 2 | AGGTTGTCGAGCGTATTGCTCTGCTTC | Astro | Astro |
AGTAGCTCGACATCGCGGGACTGTGGA | mouse_brain | 200 | 140 | 4 | 5 | 2 | 2 | 7 | 7 | AGTAGCTCGACATCGCGGGACTGTGGA | Endo | Endo |
ATGGATCTACTACGGTAAGCGACGGAT | mouse_brain | 200 | 126 | 1 | 1 | 3 | 7 | 8 | 8 | ATGGATCTACTACGGTAAGCGACGGAT | Oligo | Oligo |
CCTAAGAGCTACCGACCAAGGCACTAT | mouse_brain | 200 | 120 | 4 | 5 | 2 | 2 | 7 | 7 | CCTAAGAGCTACCGACCAAGGCACTAT | Endo | Endo |
saveRDS(mouse_brain,"WTH1092_demo_mouse_brain.rds")
#dimplot umap
options(repr.plot.height=7, repr.plot.width=7)
DimPlot(mouse_brain, reduction = 'umap')
options(repr.plot.height=10, repr.plot.width=10)
FeaturePlot(mouse_brain, reduction = 'umap', features=c('Mbp','Mobp', 'Olig1','Plp1'), pt.size = 1)
当我们把 DimPlot 函数中的 reduction 参数设置为 "spatial" 的时候,我们即可将细胞绘制在空间水平上。
options(repr.plot.height=7, repr.plot.width=15)
DimPlot(mouse_brain, reduction = 'spatial',pt.size = 1)
options(repr.plot.height=7, repr.plot.width=15)
DimPlot(mouse_brain, reduction = 'spatial',pt.size = 1,group.by = "Sub_CellType")
同样,我们可以在使用FeaturePlot 函数时,将reduction参数设置为"spatial",来查看基因在空间位置上的表达情况。
options(repr.plot.height=10, repr.plot.width=20)
FeaturePlot(mouse_brain, reduction = 'spatial', features = c('Mbp','Mobp', 'Olig1','Plp1'),pt.size = 0.7)
SeekSpace数据对Seurat的兼容性比较好,取子集等等函数均可运行。
options(repr.plot.height=7, repr.plot.width=15)
sub_mouse_brain <- subset(mouse_brain, Main_CellType == 'Ext')
DimPlot(sub_mouse_brain, reduction = 'spatial',pt.size = 1.5, group.by = "Sub_CellType")
带有DAPI/HE的图片绘制¶
自定义了ImageSpacePlot和FeatureSpacePlot两个函数绘制细胞在空间上的分群信息和细胞连续变量的指标
######################绘图细胞分群函数
ImageSpacePlot = function(obj, group_by, type="DAPI", sample=names(obj@misc$info)[1], size=1, alpha=1,color=MYCOLOR){
MYCOLOR=c(
"#6394ce", "#2a4c87", "#eed500", "#ed5858",
"#f6cbc2", "#f5a2a2", "#3ca676", "#6cc9d8",
"#ef4db0", "#992269", "#bcb34a", "#74acf3",
"#3e275b", "#fbec7e", "#ec4d3d", "#ee807e",
"#f7bdb5", "#dbdde6", "#f591e1", "#51678c",
"#2fbcd3", "#80cfc3", "#fbefd1", "#edb8b5",
"#5678a8", "#2fb290", "#a6b5cd", "#90d1c1",
"#a4e0ea", "#837fd3", "#5dce8b", "#c5cdd9",
"#f9e2d6", "#c64ea4", "#b2dfd6", "#dbdfe7",
"#dff2ec", "#cce8f3", "#e74d51", "#f7c9c4",
"#f29c81", "#c9e6e0", "#c1c5de", "#750000"
)
raster_type <- switch(type,
HE = "img_he_gg",
DAPI = "img_gg",
stop("Invalid type. Must be 'HE' or 'DAPI'.")
)
spatial_coord1 <- as.data.frame(obj[[group_by]])
colnames(spatial_coord1) <- group_by
spatial_coord2 <- as.data.frame(obj@reductions$spatial@cell.embeddings)
spatial_coord <-cbind(spatial_coord2,spatial_coord1)
ImageSpacePlot <- ggplot2::ggplot() + ggplot2::annotation_custom(grob = obj@misc$info[[sample]][[raster_type]],
xmin = 0, xmax = obj@misc$info[[sample]]$size_x,
ymin = 0, ymax = obj@misc$info[[sample]]$size_y) +
ggplot2::geom_point(data = spatial_coord, ggplot2::aes(x = spatial_1,y = spatial_2, color = !!sym(group_by),
fill = !!sym(group_by)), size=size, alpha=alpha)+
labs(size = group_by) + guides(alpha = "none")+
ggplot2::theme_classic()+
scale_color_manual(values = color)+ coord_fixed()
return(ImageSpacePlot)
}
###################绘图基因表达函数
FeatureSpacePlot = function(obj, feature, type="DAPI", sample=names(obj@misc$info)[1], size=1, alpha=c(1,1),color=c("lightgrey","blue")){
raster_type <- switch(type,
HE = "img_he_gg",
DAPI = "img_gg",
stop("Invalid type. Must be 'HE' or 'DAPI'.")
)
spatial_coord1 <- as.data.frame(obj@reductions$spatial@cell.embeddings)
spatial_coord2 <- FetchData(obj,feature)
colnames(spatial_coord2) <- feature
spatial_coord <-cbind(spatial_coord1,spatial_coord2)
FeatureSpacePlot <-ggplot2::ggplot() + ggplot2::annotation_custom(grob = obj@misc$info[[sample]][[raster_type]],
xmin = 0, xmax = obj@misc$info[[sample]]$size_x, ymin = 0, ymax = obj@misc$info[[sample]]$size_y) +
ggplot2::geom_point(data = spatial_coord, ggplot2::aes(x = spatial_1, y = spatial_2,color = !!sym(feature),alpha = !!sym(feature)), size=size)+
labs(color = feature)+
guides(alpha = "none")+
ggplot2::theme_classic()+
ggplot2::scale_alpha_continuous(range=alpha)+
scale_color_gradient(low=color[1],high = color[2])+ coord_fixed()
return(FeatureSpacePlot)
}
# DAPI背景的细胞分群图
options(repr.plot.height=7, repr.plot.width=15)
ImageSpacePlot(obj=mouse_brain, group_by = "Sub_CellType",type="DAPI",size=0.7)
# HE背景的细胞分群图
options(repr.plot.height=7, repr.plot.width=15)
ImageSpacePlot(obj=mouse_brain, group_by = "Sub_CellType",type="HE")
# DAPI背景的基因表达图
options(repr.plot.height=7, repr.plot.width=15)
FeatureSpacePlot(obj=mouse_brain, feature="Hpca",type="DAPI")
# HE背景的基因表达图
options(repr.plot.height=7, repr.plot.width=15)
FeatureSpacePlot(obj=mouse_brain, feature="Hpca",type="HE")