作者，追風(fēng)少年i~
國(guó)慶前的最后一彈衣迷，分享一個(gè)簡(jiǎn)單的內(nèi)容畏鼓，空間軌跡向量場(chǎng)。
其中關(guān)于空間軌跡壶谒，我也寫了很多云矫，文章放在下面，供大家參考
圖片.png
首先我們來(lái)解讀以下這個(gè)圖片汗菜，這個(gè)地方類似于基因让禀、細(xì)胞類型或者通路的區(qū)域轉(zhuǎn)換（細(xì)胞遷移）挑社。為了探索代謝改變區(qū)域中遷移基因表達(dá)特征的富集，確定了特定基因表達(dá)特征的低富集和高富集之間的定向梯度的空間方向巡揍。 簡(jiǎn)化后痛阻，每個(gè)點(diǎn)的方向向量是基于其局部鄰域中所研究的基因表達(dá)特征的分級(jí)富集。這些向量場(chǎng)計(jì)算使我們能夠近似空間基因表達(dá)軌跡腮敌，從而能夠識(shí)別空間上相反的轉(zhuǎn)錄途徑阱当。基于這些矢量場(chǎng)計(jì)算糜工，報(bào)告缺氧響應(yīng)和遷移特征顯示反向空間軌跡（上圖C弊添、D）。 總之捌木，研究結(jié)果為代謝變化和氧化應(yīng)激是基因組多樣性的潛在互惠驅(qū)動(dòng)因素提供了證據(jù)油坝，從而導(dǎo)致 GBM 中的克隆進(jìn)化。
其中我們要實(shí)現(xiàn)的部分在
空間軌跡向量場(chǎng)
話不多說刨裆，我們直接來(lái)
library(ggplot2)
library(Seurat)
library(SPATA2)
library(dplyr)
source('runVectorFields.R')
source('plotVectorFields.R')
####腳本放在最后
x = readRDS(Seurat.spatial.rds)
y = transformSeuratToSpata(seurat_object = x, sample_name = 'FT')
####這個(gè)時(shí)候我們把需要分析的軌跡基因澈圈、通路或者細(xì)胞類型進(jìn)行分析，這里以CD3D為例
data = runVectorFields(y,'CD3D')
head(data)
data
有了data崔拥，就可以進(jìn)行向量場(chǎng)繪圖
p = plotVectorFields(data,'CD3D')
pdf('eg.pdf')
print(p)
dev.off()
就會(huì)得到向量場(chǎng)圖极舔。
圖片.png
其中的顏色，點(diǎn)的大小都可以更改链瓦，選擇自己喜歡的搭配拆魏，當(dāng)然了，我這里是拿一個(gè)基因作為展示慈俯，更為有生物學(xué)意義的是細(xì)胞類型和信號(hào)通路渤刃，照貓畫虎就可以了（就把對(duì)應(yīng)一個(gè)的基因值替換成你想要的細(xì)胞類型分?jǐn)?shù)或者通路得分）。
附上source的代碼
- runVectorFields.R
runVectorFields <- function(object,
features,
cut_off=NULL,
normalize=T,
smooth=T,
smooth_span=NULL,
dist.spot=10,
run.mcor=T,
workers=8,
ram=50){
# Get the data:
df <- SPATA2::hlpr_join_with_aes(object,
df = SPATA2::getCoordsDf(object),
color_by = features,
normalize = normalize,
smooth = smooth,
smooth_span=smooth_span)
if(!is.null(cut_off)){df[df[,features]<cut_off, features]=0}
# Prepare data ------------------------------------------------------------
NN.file <- SPATAwrappers::getSurroundedSpots(object)
if(run.mcor==T){
base::options(future.fork.enable = TRUE)
future::plan("multisession", workers = workers)
future::supportsMulticore()
base::options(future.globals.maxSize = ram *100* 1024^2)
message("... Run multicore ... ")
}
if(dist.spot<min(NN.file %>% filter(distance!=0) %>% pull(distance))){
dist.spot <-min(NN.file %>% filter(distance!=0) %>% pull(distance))
message(paste0("The distance was adopted for the minimal distance: ",dist.spot, "px"))}
VF <- furrr::future_map(.x=1:nrow(df), .f=function(i){
#Spot Def.
bc <- df[i, c("barcodes")]
cc <- df[i, c("x", "y")]
#Neighbour Spots
NN <-
NN.file %>%
dplyr::filter(xo < cc$x+dist.spot & xo > cc$x-dist.spot) %>%
dplyr::filter(yo < cc$y+dist.spot & yo > cc$y-dist.spot) %>%
dplyr::pull(bc_destination)
#Filter input DF
NN.df <- df %>% dplyr::filter(barcodes %in% NN) %>% as.data.frame()
parameter <- features
#Create Vector
V <- -c(as.numeric(cc) - c(NN.df$x[which.max(NN.df[,parameter])], NN.df$y[which.max(NN.df[,parameter])]))
if(length(V)==0){out <- data.frame(barcodes=bc, t.x=0, t.y=0)}else{out <- data.frame(barcodes=bc, t.x=V[1], t.y=V[2])}
return(out)
}, .progress = T) %>%
do.call(rbind, .) %>%
as.data.frame()
out <- cbind(df, VF)
out[is.na(out)] <- 0
return(out)
}
- plotVectorFields.R
plotVectorFields <- function(VF, parameter, pt.size=6,pt.alpha=0.8,color.extern=NULL,skip=1){
VF <-
VF %>%
dplyr::select(x,y,{{parameter}}, t.x, t.y) %>%
dplyr::rename("parameter":=!!sym(parameter))
color.points <- VF$parameter
if(!is.null(color.extern)){color.points <- color.extern}
if(color.points %>% class()=="factor"){
p <-
ggplot2::ggplot(data=VF, aes(x,y))+
ggplot2::geom_point(data=VF , mapping=aes(x,y, color=color.points), size=pt.size, alpha=pt.alpha)+
metR::geom_vector(aes(dx = t.x, dy = t.y),skip=skip) +
metR::scale_mag()+
ggplot2::theme_void()+
Seurat::NoLegend()
}else{
p <-
ggplot2::ggplot(data=VF, aes(x,y))+
ggplot2::geom_point(data=VF , mapping=aes(x,y, color=color.points), size=pt.size, alpha=pt.alpha)+
ggplot2::scale_color_viridis_c(guide = "none")+
metR::geom_vector(aes(dx = t.x, dy = t.y),skip=skip) +
ggplot2::theme_void()+
Seurat::NoLegend()+
metR::scale_mag()
}
return(p)
}
關(guān)于軌跡分析贴膘，推薦大家采用SPATA2卖子。
生活很好，有你更好,祝打賞的老板們國(guó)慶快樂
