scWGCNA是一個(gè)生物信息學(xué)工作流程糖儡,是R包WGCNA的附加組件仅政,用于在單細(xì)胞或單核RNA-seq數(shù)據(jù)集中進(jìn)行加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析。WGCNA最初是為分析大量基因表達(dá)數(shù)據(jù)集而構(gòu)建的,由于scRNA-seq數(shù)據(jù)的固有稀疏性,vanilla WGCNA在單細(xì)胞數(shù)據(jù)上的性能受到限制。為了解決這個(gè)問題津坑,scWGCNA在運(yùn)行WGCNA分析之前將轉(zhuǎn)錄表達(dá)譜相似的細(xì)胞聚集成pseudo-bulk元細(xì)胞。
scWGCNA代碼比較簡短傲霸,對WGCNA的原有算法的一種改進(jìn),以適應(yīng)單細(xì)胞稀疏矩陣(很多表達(dá)為0);
該包只包含兩個(gè)函數(shù)metacells_by_groups和construct_metacells昙啄,construct_metacells()函數(shù)將 Seurat 對象穆役,基于knn算法,將相鄰細(xì)胞構(gòu)建平均的“元細(xì)胞”(averaged metacells)梳凛。
github網(wǎng)站:https://github.com/smorabit/scWGCNA
在線教程:https://smorabit.github.io/tutorials/9_scWGCNA_tutorial
作者在2021年nature genetics發(fā)表的文章Single-nucleus chromatin accessibility and transcriptomic characterization of Alzheimer's disease中應(yīng)用了此算法耿币。
1. 先決條件
運(yùn)行 scWGCNA,需要有一個(gè) Seurat 格式的單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)集韧拒,且已經(jīng)進(jìn)行聚類和降維淹接。
2. scWGCNA安裝
# 要運(yùn)行scWGCNA,需要安裝其他一些R包
install.packages('WGCNA')
install.packages('igraph')
install.packages('devtools')
# install Seurat, check their website for the most up-to-date instructions
install.packages('Seurat')
devtools::install_github('smorabit/scWGCNA')
library(scWGCNA)
3. scWGCNA實(shí)例教程
官方教程介紹從操作Seurat 對象到構(gòu)建metacells逐步展開分析叛溢。 如果你已經(jīng)有一個(gè)已聚類的數(shù)據(jù)集塑悼,你可以跳過接下來的前兩部分,從“構(gòu)建metacells”部分開始分析楷掉。
3.1 聚類和降維
首先厢蒜,從已發(fā)表的人類大腦健康和疾病的 snRNA-seq 研究中收集了幾個(gè)數(shù)據(jù)集,總計(jì)超過 50 萬個(gè)單核轉(zhuǎn)錄組細(xì)胞烹植。這些數(shù)據(jù)集都滿足scWGCNA 的先決條件斑鸦,并且可以通過許多不同的方式進(jìn)行聚類/降維。由于處理來自多個(gè)來源的數(shù)據(jù)草雕,因此使用online iNMF來構(gòu)建不受批次或數(shù)據(jù)來源影響的降維巷屿,以反映潛在的細(xì)胞異質(zhì)性。
代碼中總共包含9個(gè)文獻(xiàn)數(shù)據(jù)集墩虹,數(shù)據(jù)非常不容易下載攒庵,先放棄實(shí)操。
step1: 將9個(gè)seurat對象merge在一起败晴,生成seurat_obj對象浓冒,進(jìn)行Normalize歸一化和scale處理;
step2: 提取9個(gè)數(shù)據(jù)集的基因表達(dá)矩陣尖坤,生成expression_list稳懒,構(gòu)建Liger對象,Liger是一種去批次慢味,多數(shù)據(jù)集整合方法场梆,跟SCTransform-CCA,Harmony方法類似纯路。對Liger對象進(jìn)行如下操作:1)normalize或油,scale等預(yù)處理;2)執(zhí)行online iNMF驰唬;3) quantile歸一化顶岸。將liger_obj保存為9_datasets_online_liger.rds腔彰。
step3:將iNMF坐標(biāo)轉(zhuǎn)移至seurat_obj,我們之前的seurat_obj降維方式通常只有pca辖佣,umap霹抛,tsne,現(xiàn)在將增加iNMF降維卷谈,另外杯拐,用liger_obj的umap坐標(biāo)替換掉seurat_obj的坐標(biāo)。將修改后的seurat_obj保存為9_datasets_online_seurat.rds世蔗。
library(Seurat)
library(Matrix)
library(tidyverse)
library(rliger)
data_dir <- "data/"
# load individual seurat objects:
seurat_trem2 <- readRDS(file=paste0(data_dir, 'seurat_objects/trem2_nmed_unprocessed.rds'))
seurat_allen <- readRDS(file=paste0(data_dir, 'seurat_objects/allen_brain_map_2020_unprocessed.rds'))
seurat_jakel <- readRDS(file=paste0(data_dir, 'seurat_objects/jakel_2019_unprocessed.rds'))
seurat_grubman <- readRDS(file=paste0(data_dir, 'seurat_objects/grubman_2019_unprocessed.rds'))
seurat_schirmer <- readRDS(file=paste0(data_dir, 'seurat_objects/schirmer_2019_unprocessed.rds'))
seurat_velmeshev <- readRDS(file=paste0(data_dir, 'seurat_objects/velmeshev_2019_unprocessed.rds'))
seurat_swarup_splitseq <- readRDS(file=paste0(data_dir, 'seurat_objects/swarup_AD_splitseq_2019_unprocessed.rds'))
seurat_swarup <- readRDS(file=paste0(data_dir, 'seurat_objects/swarup_AD_2019_unprocessed.rds'))
seurat_tsai <- readRDS(file=paste0(data_dir, 'seurat_objects/tsai_AD_2019_unprocessed.rds'))
# make one big merged Seurat object:
seurat_obj <- merge(c(seurat_trem2, seurat_allen, seurat_jakel, seurat_grubman, seurat_schirmer, seurat_velmeshev, seurat_swarup_splitseq, seurat_swarup, seurat_tsai))
# Normalize and scale data:
seurat_obj <- NormalizeData(seurat_obj)
seurat_obj <- ScaleData(seurat_obj)
# make a list of expression matrices:
expression_list <- list(
'trem2' = GetAssayData(seurat_trem2, slot='counts'),
'allen'= GetAssayData(seurat_allen, slot='counts'),
'jakel'= GetAssayData(seurat_jakel, slot='counts'),
'grubman'= GetAssayData(seurat_grubman, slot='counts'),
'schirmer'= GetAssayData(seurat_schirmer, slot='counts'),
'velmeshev'= GetAssayData(seurat_velmeshev, slot='counts'),
'swarup_splitseq'= GetAssayData(seurat_swarup_splitseq, slot='counts'),
'swarup_AD'= GetAssayData(seurat_swarup, slot='counts'),
'mathys'= GetAssayData(seurat_tsai, slot='counts')
)
# append dataset names to barcodes:
for(dataset in names(expression_list)){
colnames(expression_list[[dataset]]) <- paste0(dataset, '_', colnames(expression_list[[dataset]]))
}
# combined metadata:
meta_columns <- c('orig.ident', 'nCount_RNA', 'nFeature_RNA', 'Age', 'Sex', 'Condition', 'Condition_specific', 'Batch', 'SampleID', 'DonorID', 'Original_cluster', 'Dataset', 'Tissue', 'Technology' )
seurat_meta <- rbind(
seurat_trem2@meta.data[,meta_columns],
seurat_allen@meta.data[,meta_columns],
seurat_jakel@meta.data[,meta_columns],
seurat_grubman@meta.data[,meta_columns],
seurat_schirmer@meta.data[,meta_columns],
seurat_velmeshev@meta.data[,meta_columns],
seurat_swarup_splitseq@meta.data[,meta_columns],
seurat_swarup@meta.data[,meta_columns],
seurat_tsai@meta.data[,meta_columns]
)
rownames(seurat_meta) <- paste0(seurat_meta$Dataset, '_', rownames(seurat_meta))
# remove individual seurat objects to save space:
rm(seurat_trem2, seurat_allen, seurat_jakel, seurat_grubman, seurat_schirmer, seurat_velmeshev, seurat_swarup, seurat_swarup_splitseq, seurat_tsai); gc();
# create liger object:
liger_obj <- createLiger(expression_list)
# pre-processing
liger_obj <- normalize(liger_obj)
liger_obj <- selectGenes(liger_obj,)
liger_obj <- scaleNotCenter(liger_obj)
# perform online iNMF
liger_obj = online_iNMF(liger_obj, k=50, max.epochs=5)
# quantile normalization
liger_obj = quantile_norm(liger_obj)
liger_obj = runUMAP(liger_obj)
saveRDS(liger_obj, file=paste0(data_dir, '9_datasets_online_liger.rds'))
# transfer iNMF matrix to seurat obj:
seurat_obj@reductions$iNMF <- CreateDimReducObject(
loadings=t(liger_obj@W),
embeddings=liger_obj@H.norm[colnames(seurat_obj),],
key="iNMF_",
assay="RNA"
)
# add UMAP:
seurat_obj@reductions$UMAP <- CreateDimReducObject(
embeddings = liger_obj@tsne.coords[colnames(seurat_obj),],
key='umap',
assay='RNA'
)
# add clusters from liger
seurat_obj@meta.data$liger_clusters <- liger_obj@clusters
saveRDS(seurat_obj, file=paste0(data_dir, '9_datasets_online_seurat.rds'))
問題1:liger_obj進(jìn)行runUMAP降維端逼,但是在后面的代碼中,出現(xiàn)liger_obj@tsne.coords污淋,有點(diǎn)看不懂顶滩,是不是誤寫?
seurat_obj@reductions$UMAP <- CreateDimReducObject(
embeddings = liger_obj@umap.coords[colnames(seurat_obj),],
key='umap',
assay='RNA'
)
找同行確認(rèn)了下芙沥,這應(yīng)該是作者誤寫诲祸。
問題2:seurat_meta整合9個(gè)樣本的meta.data信息(orig.ident,nCount_RNA等)而昨,但是代碼中沒有將seurat_meta賦值給liger_obj@meta.data或seurat_obj@meta.data救氯, 也沒有保存?
liger_obj@meta.data <- seurat_meta
seurat_obj@meta.data <- seurat_meta
下圖是用online iNMF整合方法得到的UMAP聚類圖
3.2 獲取感興趣的特定細(xì)胞類型子集
為了識別給定細(xì)胞類型內(nèi)的共表達(dá)模塊歌憨,作者僅使用數(shù)據(jù)的一個(gè)子集重復(fù)上述分析着憨。在這里,我們對感興趣的少突膠質(zhì)細(xì)胞譜系細(xì)胞進(jìn)行這種分析务嫡,包括少突膠質(zhì)祖細(xì)胞和成熟的少突膠質(zhì)細(xì)胞甲抖。
step1:從seurat_obj提取細(xì)胞類型為ODC,OPC的子集seurat_odc心铃;同上述操作准谚,對seurat_odc數(shù)據(jù)集進(jìn)行iNMF整合,降維去扣,將liger_obj保存為ODC_9_datasets_online_liger.rds柱衔。
step2:將iNMF坐標(biāo)轉(zhuǎn)移至seurat_odc,命名為ctiNMF愉棱,這將新增ctiNMF降維方式(類似于pca)唆铐,seurat_odc對象scale處理后,基于ctiNMF進(jìn)行umap降維奔滑,將生成的seurat_odc保存為9_datasets_ODC_seurat.rds艾岂。
# subset by neuronal clusters
seurat_odc <- subset(seurat_obj, celltype %in% c('ODC', 'OPC'))
# iNMF for just ODCs:
expr_matrix <- GetAssayData(seurat_odc, slot='counts')
expression_list <- list(
'zhou' = expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'zhou']],
'allen'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'allen']],
'jakel'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'jakel']],
'grubman'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'grubman']],
'schirmer'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'schirmer']],
'velmeshev'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'velmeshev']],
'swarup_splitseq'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'swarup_splitseq']],
'swarup_AD'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'swarup_AD_2020']],
'mathys'= expr_matrix[,colnames(seurat_odc)[seurat_odc$Dataset == 'mathys']]
)
seurat_meta <- rbind(
subset(seurat_odc@meta.data, Dataset=='zhou'),
subset(seurat_odc@meta.data, Dataset=='allen'),
subset(seurat_odc@meta.data, Dataset=='jakel'),
subset(seurat_odc@meta.data, Dataset=='grubman'),
subset(seurat_odc@meta.data, Dataset=='schirmer'),
subset(seurat_odc@meta.data, Dataset=='velmeshev'),
subset(seurat_odc@meta.data, Dataset=='swarup_splitseq'),
subset(seurat_odc@meta.data, Dataset=='swarup_AD_2020'),
subset(seurat_odc@meta.data, Dataset=='mathys')
)
# create liger object:
liger_obj <- createLiger(expression_list)
liger_obj <- normalize(liger_obj)
pdf("liger_variable_genes.pdf", width=8, height=8)
liger_obj <- selectGenes(
liger_obj,
var.thresh =c(0.05, 0.15, 0.025, 0.025, 0.075, 0.075, 0.05, 0.20, 0.025),
do.plot=T
)
dev.off()
liger_obj@var.genes %>% length
liger_obj <- scaleNotCenter(liger_obj)
# perform online iNMF
liger_obj = online_iNMF(liger_obj, k=20, max.epochs=5)
# quantile normalization
liger_obj = quantile_norm(liger_obj)
liger_obj = runUMAP(liger_obj)
saveRDS(liger_obj, file='ODC_9_datasets_online_liger.rds')
# transfer iNMF matrix to seurat obj:
seurat_odc@reductions$ctiNMF <- CreateDimReducObject(
loadings=t(liger_obj@W),
embeddings=liger_obj@H.norm[colnames(seurat_odc),],
key="ctiNMF_",
assay="RNA"
)
VariableFeatures(seurat_odc) <- liger_obj@var.genes
# scale expression data:
seurat_odc <- ScaleData(seurat_odc, features=VariableFeatures(seurat_odc))
# UMAP + clustering:
seurat_odc <- RunUMAP(seurat_odc, reduction='ctiNMF', dims=1:20)
# save the results
saveRDS(seurat_odc, file ='9_datasets_ODC_seurat.rds')
# optional: format for use in Scanpy
library(SeuratDisk)
SaveH5Seurat(cur_save, filename = "9_datasets_ODC_seurat.h5Seurat", overwrite=TRUE)
Convert("9_datasets_ODC_seurat.h5Seurat", dest='h5ad', overwrite=TRUE)
下面可以看到少突膠質(zhì)細(xì)胞譜系 細(xì)胞的UMAP聚類圖,作者根據(jù)已知標(biāo)記基因?qū)⑵浞譃?4 個(gè)主要cluster朋其。
3.3 構(gòu)建metacells
現(xiàn)在我們準(zhǔn)備構(gòu)建metacells王浴。為了節(jié)省內(nèi)存脆炎,我們只使用高可變基因來構(gòu)建metacells并執(zhí)行下游 WGCNA分析。 可以使用 HVG基因叼耙,也可以用所有基因腕窥,下游結(jié)果將基因數(shù)目的差異而有所不同粒没。因此筛婉,你可以重新運(yùn)行此步驟以嘗試不同的參數(shù)。
step1:這個(gè)數(shù)據(jù)集有幾種不同的疾病處理?xiàng)l件癞松,但是爽撒,我們只想聚合相同條件下的細(xì)胞,即odc_group和Condition屬性一致的細(xì)胞群响蓉。 因此硕勿,我們在 Seurat 對象中創(chuàng)建了一個(gè)新的metadata列metacell_group ,它將 ODC 亞群(odc_group)與處理?xiàng)l件(Condition)組合在一起枫甲。
step2:接下來源武,我們將來自同一 ODC 亞群和相同疾病處理組的細(xì)胞群構(gòu)建metacells。construct_metacells()函數(shù)有一個(gè)關(guān)鍵參數(shù)是 k想幻,有多少個(gè)細(xì)胞合并成一個(gè)metacells粱栖。 這里我們使用了一個(gè)較大數(shù)值 100,這僅僅是因?yàn)槲覀冞@里的數(shù)據(jù)集非常大脏毯。 你可以根據(jù)數(shù)據(jù)集中的細(xì)胞數(shù)來調(diào)整 k 值大小闹究。我只是通過 for 循環(huán)單線程運(yùn)行它,當(dāng)然食店,您可以更快速渣淤,例如,如果您要將每個(gè)調(diào)用construct_metacells函數(shù)的程序投入到HPC計(jì)算節(jié)點(diǎn)吉嫩,實(shí)現(xiàn)并行運(yùn)算价认。該construct_metacells() 函數(shù)輸出一個(gè)新的Seurat 對象,包含聚合后的基因表達(dá)譜自娩。
step3:將每個(gè)metacell_group構(gòu)建好的metacells(Seurat 對象)進(jìn)行merge用踩,scale縮放,Harmony整合椒功,umap聚類捶箱。
# load scWGCNA package
library(scWGCNA)
# construct metacells by ODC group, condition
seurat_odc$metacell_group <- paste0(
as.character(seurat_odc$odc_group), '_',
as.character(seurat_odc$Condition)
)
genes.keep <- VariableFeatures(seurat_odc)
# loop through each group and construct metacells
seurat_list <- list()
for(group in unique(seurat_odc$metacell_group)){
print(group)
cur_seurat <- subset(seurat_odc, metacell_group == group)
cur_seurat <- cur_seurat[genes.keep,]
cur_metacell_seurat <- scWGCNA::construct_metacells(
cur_seurat, name=group,
k=100, reduction='umap',
assay='RNA', slot='data'
)
cur_metacell_seurat$Condition <- as.character(unique(cur_seurat$Condition))
cur_metacell_seurat$odc_group <- as.character(unique(cur_seurat$odc_group))
seurat_list[[group]] <- cur_metacell_seurat
}
# merge all of the metacells objects
metacell_seurat <- merge(seurat_list[[1]], seurat_list[2:length(seurat_list)])
saveRDS(metacell_seurat, file='data/metacell_seurat.rds')
Optionally, we can run a dimensionality reduction on the metacell Seurat object to check if some cellular heterogeneity has been retained. The differences in the transcriptional profiles of the different disease groups are quite prominent in the UMAP space, so I am harmonizing the data on the basis of condition merely for visualization purposes.
metacell_subset <- ScaleData(metacell_subset, features = rownames(metacell_subset))
metacell_subset <- RunPCA(metacell_subset, features=rownames(metacell_subset))
metacell_subset <- RunHarmony(metacell_subset, group.by='Condition', dims=1:15)
metacell_subset <- RunUMAP(metacell_subset, reduction='harmony', dims=1:15)
pdf(paste0(fig_dir, "metacell_umap_group.pdf"), width=6, height=6)
DimPlot(metacell_subset, group.by='odc_group', reduction='umap', label=TRUE) + umap_theme
dev.off()
由下圖可以看出,通過construct_metacells處理后的seurat對象在umap圖中仍有很好地分離动漾。
前面代碼最關(guān)鍵的一步是construct_metacells丁屎,將cur_seurat轉(zhuǎn)成cur_metacell_seurat ,我們看下源碼旱眯,這步做了什么處理晨川。
construct_metacells <- function(seurat_obj, name='agg', k=50, reduction='umap', assay='RNA', slot='data'){
reduced_coordinates <- as.data.frame(seurat_obj@reductions[[reduction]]@cell.embeddings)
nn_map <- FNN::knn.index(reduced_coordinates, k = (k - 1))
row.names(nn_map) <- row.names(reduced_coordinates)
nn_map <- cbind(nn_map, seq_len(nrow(nn_map)))
good_choices <- seq_len(nrow(nn_map))
choice <- sample(seq_len(length(good_choices)), size = 1,
replace = FALSE)
chosen <- good_choices[choice]
good_choices <- good_choices[good_choices != good_choices[choice]]
it <- 0
k2 <- k * 2
get_shared <- function(other, this_choice) {
k2 - length(union(cell_sample[other, ], this_choice))
}
while (length(good_choices) > 0 & it < 5000) {
it <- it + 1
choice <- sample(seq_len(length(good_choices)), size = 1,
replace = FALSE)
new_chosen <- c(chosen, good_choices[choice])
good_choices <- good_choices[good_choices != good_choices[choice]]
cell_sample <- nn_map[new_chosen, ]
others <- seq_len(nrow(cell_sample) - 1)
this_choice <- cell_sample[nrow(cell_sample), ]
shared <- sapply(others, get_shared, this_choice = this_choice)
if (max(shared) < 0.9 * k) {
chosen <- new_chosen
}
}
cell_sample <- nn_map[chosen, ]
combs <- combn(nrow(cell_sample), 2)
shared <- apply(combs, 2, function(x) {
k2 - length(unique(as.vector(cell_sample[x, ])))
})
message(paste0("Overlap QC metrics:\nCells per bin: ",
k, "\nMaximum shared cells bin-bin: ", max(shared),
"\nMean shared cells bin-bin: ", mean(shared), "\nMedian shared cells bin-bin: ",
median(shared)))
if (mean(shared)/k > 0.1)
warning("On average, more than 10% of cells are shared between paired bins.")
exprs_old <- GetAssayData(seurat_obj, assay=assay, slot=slot)
mask <- sapply(seq_len(nrow(cell_sample)), function(x) seq_len(ncol(exprs_old)) %in%
cell_sample[x, , drop = FALSE])
mask <- Matrix::Matrix(mask)
new_exprs <- (exprs_old %*% mask) / k
colnames(new_exprs) <- paste0(name, '_', 1:ncol(new_exprs))
rownames(cell_sample) <- paste0(name, '_', 1:ncol(new_exprs))
colnames(cell_sample) <- paste0('knn_', 1:ncol(cell_sample))
# make seurat obj:
seurat_aggr <- CreateSeuratObject(
counts = new_exprs
)
seurat_aggr
}
代碼解析:
1.假設(shè)seurat_obj對象有5000個(gè)細(xì)胞证九,10000個(gè)基因,即5000*10000的基因表達(dá)矩陣共虑。提取seurat_obj對象的umap坐標(biāo)信息(umap1,umap2)愧怜,生成reduced_coordinates(5000*2矩陣)。
2.利用FNN快速最近鄰搜索算法妈拌,為每個(gè)細(xì)胞尋找99個(gè)最近鄰細(xì)胞拥坛,nn_map為5000*99矩陣;第1-99列均為最臨近細(xì)胞的次序號尘分;
添加第100列猜惋,為當(dāng)前細(xì)胞自身的次序號;nn_map為5000*100矩陣培愁;99個(gè)鄰近細(xì)胞加上自身細(xì)胞著摔,共100個(gè)。
3.隨機(jī)從5000個(gè)細(xì)胞中抽取一個(gè)細(xì)胞定续,如編號為2000的細(xì)胞N2000谍咆,然后從剩下的4999個(gè)細(xì)胞中再隨機(jī)抽取一個(gè)編號N100細(xì)胞,計(jì)算N2000和N100的最臨近細(xì)胞群數(shù)目(200個(gè)細(xì)胞編號)私股,對著200個(gè)細(xì)胞編號摹察,取其并集,如滿足max(shared) < 0.9 * k條件庇茫,也就是說港粱,如果這兩細(xì)胞越相似,他們的共享鄰居應(yīng)該越多旦签,把N2000和N100細(xì)胞合并成chosen組查坪,反之亦然。接著從剩下的細(xì)胞群抽取宁炫,與上一步的chosen 組比較偿曙,看相似程度,迭代5000次羔巢,把相似的細(xì)胞編號放入chosen組望忆。
4.通過上一步的找鄰居策略,chosen組為2000個(gè)竿秆,5000個(gè)總細(xì)胞群启摄,有2000個(gè)最相似的最近鄰居群, cell_sample為2000*100矩陣幽钢。
5.2000個(gè)細(xì)胞兩兩比較歉备,查看共有鄰居數(shù),統(tǒng)計(jì)max(shared)匪燕,median(shared)蕾羊,mean(shared)喧笔。
6.提取5000*10000的基因表達(dá)矩陣exprs_old ,(exprs_old %*% mask) / k的計(jì)算是2000個(gè)細(xì)胞鄰近細(xì)胞群的平均表達(dá)量龟再,轉(zhuǎn)換成2000*10000书闸,即5000個(gè)細(xì)胞由2000個(gè)細(xì)胞所代表,數(shù)據(jù)集也實(shí)現(xiàn)了降維利凑。
construct_metacells算法小結(jié):
1.seurat_obj對象有5000個(gè)細(xì)胞浆劲,10000個(gè)基因,是5000*10000的基因表達(dá)矩陣截碴;設(shè)置參數(shù)k=100梳侨,通過KNN算法查找100個(gè)最臨近鄰居群(包括自身細(xì)胞)蛉威,如A細(xì)胞日丹,找99個(gè)鄰居,加上自己蚯嫌,共100個(gè)哲虾。
2.指定一套篩選相似細(xì)胞策略,將細(xì)胞的鄰居群相互比較择示,保留相似的細(xì)胞束凑,5000個(gè)細(xì)胞保留了最相似的2000個(gè)細(xì)胞,用這2000個(gè)細(xì)胞代表整體細(xì)胞栅盲。
3.之前是5000*10000的基因表達(dá)矩陣汪诉,新的基因表達(dá)矩陣是:2000個(gè)細(xì)胞的鄰居群(100個(gè)細(xì)胞)的平均基因表達(dá)量,如A細(xì)胞的gene1表達(dá)量為100個(gè)鄰居群在gene1的表達(dá)平均值谈秫,因此新矩陣為2000*10000扒寄,剔除了噪音數(shù)據(jù),也實(shí)現(xiàn)了數(shù)據(jù)降維拟烫。這應(yīng)該也是處理單細(xì)胞表達(dá)譜的一種數(shù)據(jù)處理策略该编。