上次的視頻中已花費(fèi)大量時(shí)間講解過單樣本分析的基本流程，所以這節(jié)課的學(xué)習(xí)需要有上節(jié)課的基礎(chǔ)饮六，希望大家按順序觀看其垄。此次的內(nèi)容較簡單、篇幅也較小卤橄，代碼與視頻請(qǐng)看下文绿满，測試數(shù)據(jù)集與代碼存于文末鏈接之中。由于測試數(shù)據(jù)比較特殊窟扑，并沒有展示出去批次的精妙之處喇颁，留一個(gè)懸念給大家吧，可以用自己的數(shù)據(jù)集測試一下嚎货。
（B站同步播出橘霎，先看一遍視頻再跟著代碼一起操作，建議每個(gè)視頻至少看三遍）
###########單純的merge#################
library(Seurat)
library(multtest)
library(dplyr)
library(ggplot2)
library(patchwork)
##########準(zhǔn)備用于拆分的數(shù)據(jù)集##########
#pbmc <- subset(pbmc, downsample = 50)
ifnb <- readRDS('pbmcrenamed.rds')
ifnb.list <- SplitObject(ifnb, split.by = "group")
C57 <- ifnb.list$C57
AS1 <- ifnb.list$AS1
######簡單merge########
#不具有去批次效應(yīng)功能
pbmc <- merge(C57, y = c(AS1), add.cell.ids = c("C57", "AS1"), project = "ALL")
pbmc
head(colnames(pbmc))
unique(sapply(X = strsplit(colnames(pbmc), split = "_"), FUN = "[", 1))
table(pbmc$orig.ident)
##############anchor###############
library(Seurat)
library(tidyverse)
### testA ----
myfunction1 <- function(testA.seu){
testA.seu <- NormalizeData(testA.seu, normalization.method = "LogNormalize", scale.factor = 10000)
testA.seu <- FindVariableFeatures(testA.seu, selection.method = "vst", nfeatures = 2000)
return(testA.seu)
}
C57 <- myfunction1(C57)
AS1 <- myfunction1(AS1)
### Integration ----
testAB.anchors <- FindIntegrationAnchors(object.list = list(C57,AS1), dims = 1:20)
testAB.integrated <- IntegrateData(anchorset = testAB.anchors, dims = 1:20)
#需要注意的是：上面的整合步驟相對(duì)于harmony整合方法殖属，對(duì)于較大的數(shù)據(jù)集（幾萬個(gè)細(xì)胞）
#非常消耗內(nèi)存和時(shí)間姐叁，大約9G的數(shù)據(jù)32G的內(nèi)存就已經(jīng)無法運(yùn)行；
#當(dāng)存在某一個(gè)Seurat對(duì)象細(xì)胞數(shù)很少（印象中200以下這樣子），
#會(huì)報(bào)錯(cuò)外潜，這時(shí)建議用第二種整合方法
DefaultAssay(testAB.integrated) <- "integrated"
# # Run the standard workflow for visualization and clustering
testAB.integrated <- ScaleData(testAB.integrated, features = rownames(testAB.integrated))
testAB.integrated <- RunPCA(testAB.integrated, npcs = 50, verbose = FALSE)
testAB.integrated <- FindNeighbors(testAB.integrated, dims = 1:30)
testAB.integrated <- FindClusters(testAB.integrated, resolution = 0.5)
testAB.integrated <- RunUMAP(testAB.integrated, dims = 1:30)
testAB.integrated <- RunTSNE(testAB.integrated, dims = 1:30)
p1<- DimPlot(testAB.integrated,label = T,split.by = 'group')#integrated
DefaultAssay(testAB.integrated) <- "RNA"
testAB.integrated <- ScaleData(testAB.integrated, features = rownames(testAB.integrated))
testAB.integrated <- RunPCA(testAB.integrated, npcs = 50, verbose = FALSE)
testAB.integrated <- FindNeighbors(testAB.integrated, dims = 1:30)
testAB.integrated <- FindClusters(testAB.integrated, resolution = 0.5)
testAB.integrated <- RunUMAP(testAB.integrated, dims = 1:30)
testAB.integrated <- RunTSNE(testAB.integrated, dims = 1:30)
p2 <- DimPlot(testAB.integrated,label = T,split.by = 'group')
p1|p2
###########harmony 速度快原环、內(nèi)存少################
if(!require(harmony))devtools::install_github("immunogenomics/harmony")
test.seu <- pbmc
test.seu <- test.seu%>%
Seurat::NormalizeData() %>%
FindVariableFeatures(selection.method = "vst", nfeatures = 2000) %>%
ScaleData()
test.seu <- RunPCA(test.seu, npcs = 50, verbose = FALSE)
#####run 到PCA再進(jìn)行harmony，相當(dāng)于降維########
test.seu=test.seu %>% RunHarmony("group", plot_convergence = TRUE)
test.seu <- test.seu %>%
RunUMAP(reduction = "harmony", dims = 1:30) %>%
FindNeighbors(reduction = "harmony", dims = 1:30) %>%
FindClusters(resolution = 0.5) %>%
identity()
test.seu <- test.seu %>%
RunTSNE(reduction = "harmony", dims = 1:30)
p3 <- DimPlot(test.seu, reduction = "tsne", group.by = "group", pt.size=0.5)+theme(
axis.line = element_blank(),
axis.ticks = element_blank(),axis.text = element_blank()
)
p4 <- DimPlot(test.seu, reduction = "tsne", group.by = "ident", pt.size=0.5, label = TRUE,repel = TRUE)+theme(
axis.line = element_blank(),
axis.ticks = element_blank(),axis.text = element_blank()
)
p3|p4
本系列其他課程
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（一）——緒論
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（二）——各類輸入文件讀取
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（三）——單樣本分析
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（四）——多樣本整合
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（五）——細(xì)胞類型注釋
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（六）——組間差異分析及可視化
手把手教你做單細(xì)胞測序數(shù)據(jù)分析（七）——基因集富集分析
歡迎關(guān)注同名公眾號(hào)~
