close
導語

批次效應一直都是一個很麻煩的事情,不管是bulk RNA-seq還是單細胞數據,我們沒辦法徹底消除批次效應,只能儘量的減少批次效應,例如在實驗設計之初,儘可能做到平衡的實驗設計:

形狀表示不同的樣品類型(如組織或患者)和顏色處理批次。在混雜設計中,由於處理批次的原因,不可能從生物變異中分離出變異。在平衡設計中,通過使用組織複製和混合它們的批次,有可能區分生物和批次相關的變異。

再者,我們得到下游數據之後,可以利用一些算法降低批次效應,例如單細胞數據主流的一些算法CCA, Harmony, LIGER, fastMNN, Conos, Scanorama。但是麻煩的是,有些時候算法整合的力度不夠,又有一些時候算法整合的力度太強了導致過擬合。這也就是為什麼從單細胞誕生之初到現在,大量的整合算法仍是層出不窮。

那麼哪個算法更好呢?其實主流的幾個算法都是各有其優點的,2021年Nat Biotechnol有一篇測評文章《A multicenter study benchmarking single-cell RNA sequencing technologies using reference samples》對目前主流的單細胞整合方法做了測評和總結,大致結論如下:

image-20220228155146635

為了更好的理解單細胞整合去批次,幾大主流的方法是需要大量的進行實踐和練習的。本文就Seurat R包提供的CCA算法做一介紹,並討論整合後分析所用插槽的使用選擇問題。

一. Seurat CCA整合去批次

參考教程:

https://satijalab.org/seurat/articles/integration_introduction.html

https://www.singlecellcourse.org/scrna-seq-dataset-integration.html

https://hbctraining.github.io/scRNA-seq_online/lessons/06_integration.html

三個建議:

整合前先check一下未整合前的批次情況,再做決定是否需要做整合分析;

如果在liunx環境下運行:

linux環境下開啟library(future)多線程運行,但是注意資源的合理分配;
熟練者可以把寫好的R腳本用Rscript批量運行,如下所示,這樣後台運行節省時間。
#cat>Rscript_seurat_CCA.Rrm(list=ls())library(Seurat)library(dplyr)library(patchwork)library(readr)library(ggplot2)library(RColorBrewer)#library(future)library(clustree)library(cowplot)library(stringr)library(SeuratData)#changethecurrentplantoaccessparallelization#plan("multisession",workers=2)#plan()#設置可用的內存#options(future.globals.maxSize=10*1024^3)#installdatasetInstallData("ifnb")#loaddatasetLoadData("ifnb")###CCA整合去批次##splitthedatasetintoalistoftwoseuratobjects(stimandCTRL)seurat_list<-SplitObject(ifnb.list,split.by="stim")##normalizeandidentifyvariablefeaturesforeachdatasetindependentlyseurat_list<-lapply(X=seurat_list,FUN=function(x){x<-NormalizeData(x)x<-FindVariableFeatures(x,selection.method="vst",nfeatures=2000)})##Performintegration#selectfeaturesthatarerepeatedlyvariableacrossdatasetsforintegrationfeatures<-SelectIntegrationFeatures(object.list=seurat_list)#找到整合的錨點anchorsseurat_anchors<-FindIntegrationAnchors(object.list=seurat_list,anchor.features=features,dims=1:30)#利用anchors進行整合seurat_int<-IntegrateData(anchorset=seurat_anchors,dims=1:30)names(seurat_int@assays)seurat_int@active.assay#移除不用的對象rm(seurat_list)rm(seurat_anchors)###5.降維#DefaultAssay(seurat_int)<-"integrated"features=row.names(seurat_int)seurat_int<-ScaleData(seurat_int,features=features,verbose=FALSE,vars.to.regress=c("percent.mt"))seurat_int<-RunPCA(seurat_int,npcs=30,verbose=FALSE)seurat_int<-RunUMAP(seurat_int,reduction="pca",dims=1:30)seurat_int<-FindNeighbors(seurat_int,reduction="pca",dims=1:30)for(resinc(0.05,0.1,0.2,0.3,0.5,0.8,1,1.2)){print(res)seurat_int<-FindClusters(seurat_int,graph.name="integrated_snn",resolution=res,algorithm=1)}#umap可視化cluster_umap<-wrap_plots(ncol=4,DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.0.05",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.0.1",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.0.2",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.0.3",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.0.5",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.0.8",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.1",label=T)&NoAxes(),DimPlot(seurat.qc,reduction="umap",group.by="integrated_snn_res.1.2",label=T)&NoAxes())#cluster_umapggsave(cluster_umap,filename="./Step4.After_inter.cluster_umap.pdf",width=20,height=16)write_rds(seurat_int,file="./Step4.After_int_Dimension_reduction_clustering.rds")

這是整合後的結果:

整合前長這樣:

img

如果是在Linux環境下,寫好腳本後,我習慣用Rscript後台運行:

##運行腳本nohupRscriptRscript_seurat_CCA.R1>seurat.CCA.log2>&1&二. 後續分析的插槽選擇

Seurat CCA整合分析後什麼情況下用RNAassay,什麼時候用intergratedassay?

為了回答這個問題,我google查了很多帖子:

https://github.com/satijalab/seurat/issues/2906
https://github.com/satijalab/seurat/issues/1659
https://github.com/satijalab/seurat/issues/1256
https://github.com/satijalab/seurat/issues/1091#issuecomment-457285643
https://github.com/satijalab/seurat/issues/1717
https://github.com/satijalab/seurat/issues/4274
https://github.com/satijalab/seurat/discussions/4032
https://github.com/satijalab/seurat/issues/1674

總結來說,

conserved cell type markers部分,用整合前的「RNA」或整合後的歸一化」integrated「分析不會造成差異;

作為一般規則,我們總是建議對原始的「RNA「執行差異分析,而不是對批次校正等值執行差異分析;

後續的分析建議建立在「RNA」上,因此在CCA之後,可以設置:

DefaultAssay(sce) <- "RNA"

https://github.com/satijalab/seurat/issues/1717 的一個總結:

To keep this simple: You should use the integratedassay when trying to 'align' cell states that are shared across datasets (i.e. for clustering, visualization, learning pseudotime, etc.)。
You should use the RNA assay when exploring the genes that change either across clusters, trajectories, or conditions. SingleR也是用RNA插槽。

https://github.com/satijalab/seurat/issues/1674

擬時序分析可以用整合後的,也可以用原始的Count插槽,具體視分析目的而定:

- END -
arrow
arrow
    全站熱搜
    創作者介紹
    創作者 鑽石舞台 的頭像
    鑽石舞台

    鑽石舞台

    鑽石舞台 發表在 痞客邦 留言(0) 人氣()