今天我們來學習一個簡單的數(shù)學內(nèi)容票顾，NMF
NMF算法簡介
NMF是什么挽拂？
圖片
非負矩陣分解(Non-negative Matrix Factorization, NMF)本質(zhì)上說是一種矩陣分解的方法蝇更，對于任意給定的一個非負矩陣V余掖，NMF算法能夠?qū)ふ业揭粋€非負矩陣W和一個非負矩陣H祭陷，使得 V≈W*H 成立 峦椰，從而將一個非負的矩陣分解為左右兩個非負矩陣的乘積。
NMF的特點壁拉？
NMF最重要的特點是非負性約束谬俄。
矩陣分解的方法有很多種，如奇異值分解 (singular value decomposition, SVD) 扇商、獨立成分分析 (independent component analysis, ICA) 凤瘦、主成分分析 (principal component analysis, PCA) 等。這些方法的共同特點是案铺，即使初始矩陣 V 元素是非負的蔬芥，分解出來的因子 W 和 H 中的元素往往含有負值元素。從計算科學的角度來看控汉，分解出來的因子 W 和 H 中的元素含有負值元素并沒有問題笔诵， 但負值元素通常是無法解釋的。NMF約束了原始矩陣V和分解矩陣W姑子、H的非負性乎婿，這就意味著只能通過特征的相加來實現(xiàn)原始矩陣V的還原，最終導致的結(jié)果是：
- 非負性會引發(fā)稀疏
- 非負性會使計算過程進入部分分解
給大家對比一下PCA與NMF分解圖像的效果：
圖片.png
PCA分解之后街佑，每個主成分（PC）都會保留與其他PC不正交的全局特征谢翎，并且PC保留的特征是遞減的。
圖片.png
NMF分解之后沐旨，每個因子保留的都是局部特征森逮，它們的權重是基本平等的。通過這張圖可以看出磁携，很多因子能與面部特征一一對應起來褒侧，例如鼻子、眼睛谊迄、嘴巴都能找到相應的因子闷供。
NMF在單細胞研究中的優(yōu)勢
單細胞研究避免不了要回答兩個問題：組織中有哪些細胞類型，每個細胞類型又有哪些表達模式统诺？NMF解決這類問題具有天然的優(yōu)勢歪脏，因為它分解的因子很容易與細胞類型或表達模式對應起來。Github上有很多基于NMF和其變種算法的單細胞分析工具粮呢，我比較喜歡的有單細胞整合分析工具liger和空間轉(zhuǎn)錄組去卷積工具SPOTlight唾糯。應用NMF分析方法發(fā)表的高分文章也有很多怠硼，我給大家介紹一篇，更多的文章請自己搜索移怯。
Chen, YP., Yin, JH., Li, WF. et al. Single-cell transcriptomics reveals regulators underlying immune cell diversity and immune subtypes associated with prognosis in nasopharyngeal carcinoma. Cell Res 30, 1024–1042 (2020). https://doi.org/10.1038/s41422-020-0374-x
基礎NMF包的安裝與用法簡介
安裝NMF基礎包
BiocManager::install('Biobase')
nmf函數(shù)簡介
NMF包通過nmf()函數(shù)實現(xiàn)矩陣分解，它的用法及重要參數(shù)如下：
nmf(x, rank, method, seed, nrun, ...)
rank值一般是要通過測試評估后確定的这难，但是分析單細胞數(shù)據(jù)這是一個很難完成的工作舟误，5000個細胞的測試時間可能超過10個小時。替代辦法是使用經(jīng)驗或先驗知識指定姻乓，可以嘗試略多于細胞類型或細胞狀態(tài)（細胞亞群再聚類時）的一個數(shù)值嵌溢，例如我在本帖的PBMC數(shù)據(jù)分解中就指定為rank=10。因為NMF一般是從隨機數(shù)開始蹋岩，通過迭代算法收斂誤差的方法求出最優(yōu)W和H矩陣赖草，所以seed不同最后的結(jié)果也不同。為了減少seed的影響求得最優(yōu)解剪个，常規(guī)的辦法是通過nrun參數(shù)設置運行100-200次矩陣分解選取最優(yōu)值秧骑，也可以使用特殊的算法選擇一個最佳的seed（設置seed='nndsvd'或seed='ica'），這樣運行一次也能得到最優(yōu)解扣囊。下面我們測試一下不同方法的運行時間：
library(Biobase)
10X PBMC數(shù)據(jù)實測
測試數(shù)據(jù)
測試數(shù)據(jù)來源于10x genomics官網(wǎng)的示例數(shù)據(jù)集乎折。
圖片
數(shù)據(jù)下載鏈接：https://cf.10xgenomics.com/samples/cell-exp/3.0.2/5k_pbmc_v3_nextgem/5k_pbmc_v3_nextgem_filtered_feature_bc_matrix.h5
保存文件名：pbmc.h5
基于PCA分解的降維聚類
library(Seurat)
圖片
圖片
基于NMF分解的降維聚類
## 高變基因表達矩陣的分解
# pbmc大體可分成T，B侵歇，NK骂澄，CD14+Mono，CD16+Mono惕虑，DC坟冲，Platelet等類型，考慮冗余后設置rank=10
vm <- pbmc@assays$RNA@scale.data
res <- nmf(vm, 10, method = "snmf/r", seed = 'nndsvd')
runtime(res)
# 用戶 系統(tǒng) 流逝
#1063.147 78.019 1139.831
## 分解結(jié)果返回suerat對象
pbmc@reductions$nmf <- pbmc@reductions$pca
pbmc@reductions$nmf@cell.embeddings <- t(coef(res))
pbmc@reductions$nmf@feature.loadings <- basis(res)
## 使用nmf的分解結(jié)果降維聚類
set.seed(219)
pbmc.nmf <- RunUMAP(pbmc, reduction = 'nmf', dims = 1:10) %>%
FindNeighbors(reduction = 'nmf', dims = 1:10) %>% FindClusters()
## 結(jié)果可視化
p <- DimPlot(pbmc.nmf, label = T) + ggsci::scale_color_igv()
ggsave("pbmc_nmf.png", p, width = 9, height = 6)
p <- FeaturePlot(pbmc.nmf, features = c('CD3D', 'CD3E', 'MS4A1', 'CD79A', 'GNLY', 'NKG7', 'CD14',
'FCGR3A', 'PPBP', 'FCER1A', 'CD4', 'CD8A'), ncol = 4)
ggsave("pbmc_nmf_markers.png", p, width = 12, height = 8)
圖片.png
圖片
NMF因子可解釋性探索
對比PCA分析的結(jié)果溃蔫，NMF雖然毫不遜色健提，但是它的運行時間更長，我們?yōu)槭裁匆肗MF呢酒唉？一個很重要的原因是NMF的因子可解釋性更強矩桂，每個因子貢獻度最大的基因基本代表了某種或某個狀態(tài)細胞的表達模式，相比差異分析得到marker基因更有代表性痪伦。
NMF因子與細胞類型的關系
## 人工定義細胞類型
pbmc.nmf$celltype <- pbmc.nmf$seurat_clusters
pbmc.nmf$celltype <- recode(pbmc.nmf$celltype,
'1' = "B cells",
'5' = "NKs",
'10' = "CD8+ T",
'6' = "CD8+ T",
'4' = "CD4+ T",
'9' = "CD4+ T",
'0' = "CD4+ T",
'3' = "CD4+ T",
'7' = "CD4+ T",
'2' = "CD14+ Mono",
'8' = "CD14+ Mono",
'11' = "CD16+ Mono",
'12' = "DCs",
'13' = "Platelet",
'14' = "Unknown")
p <- DimPlot(pbmc.nmf, group.by = 'celltype', label = T, label.size = 3) + ggsci::scale_color_npg(alpha = 0.6)
ggsave("pbmc_nmf_celltype.png", p, width = 9, height = 6)
## 查看細胞因子上的荷載
tmp <- data.frame(t(coef(res)), check.names = F)
colnames(tmp) <- paste0("factor", 1:10)
pbmc.nmf <- AddMetaData(pbmc.nmf, metadata = tmp)
p <- FeaturePlot(pbmc.nmf, features = paste0("factor", 1:10), ncol = 4)
ggsave("pbmc_nmf_factors.png", p, width = 12, height = 8)
## 查看細胞主成分上的荷載
p <- FeaturePlot(pbmc.nmf, features = paste0("PC_", 1:12), ncol = 4)
ggsave("pbmc_nmf_PCs.png", p, width = 12, height = 8)
人工鑒定的細胞類型
圖片
細胞在因子上的值
圖片
細胞在PC軸上的值
圖片.png
對比上下兩張圖侄榴，很容易發(fā)現(xiàn)NMF的因子比PCA的PC軸解釋性更強。
提取celltype的signatures
## 提取每個因子貢獻度最大的20個基因
f <- extractFeatures(res, 20L)
f <- lapply(f, function(x) rownames(res)[x])
f <- do.call("rbind", f)
DT::datatable(t(f))
圖片
NMF的10個因子中网沾，很容易發(fā)現(xiàn)我們對細胞進行分類的marker基因癞蚕，是不是很神奇？
生活很好辉哥，等你超越
