ICA降維簡介

ICA又稱為獨立成分分析，在單細胞轉錄組學里面作為一種線性降維參與到單細胞轉錄組的分析中。首先我們介紹一下觀測變量的線性組合：

簡而言之就是一組隨機觀測變量x，和另外一組隨機變量y，將它們用線性方式來表示
而ICA的目的是把一組數(shù)據(jù)拆分成不同的模塊，并用這些模塊線性表示這組數(shù)據(jù)

上式中 y 表示某一組隨機變量；而 x 表示為線性組合成 y 的模塊
而降維的一般原則就是盡可能在保證隨機變量 y 之間有差異的條件下，盡可能的減少模塊 x的個數(shù)
ICA的目的是將一組數(shù)據(jù)用其中最“獨立”的幾個模塊線性來表示

單細胞中的ICA

顯然，對應于單細胞轉錄組數(shù)據(jù)來說：

想利用ICA進行降維，那么對應于每一個cell（列向量）來說，我們可以依據(jù)每個cell的基因表達模式進行模塊的拆分
這一部分主要選擇cell間基因表達差異較大的基因作為模塊（類似于kmeans分基因表達模塊），那么我們將這些模塊線性表示為每個"cell"

每一個mod里面是一組基因表達量

舉個比較簡單的例子：

由此可見，每個細胞的數(shù)據(jù)特征便可以利用基因表達量及設計矩陣進行線性表示
那么，這里有五個基因，相當于個5維空間，假設我要降維成二維，意味著我僅僅保留兩個基因，卻能很好的將這5個cell區(qū)分開，那么我們勢必要得到其中各個細胞之間差異最大的兩個基因。比方說基因A和基因B在各個細胞中差異最大

那么此時保留住基因A和基因B即可，在設計矩陣中把其他基因的系數(shù)均改為0即可，那么所得式子如上，那么在二維坐標軸上即可可視化了，也就達到了降維的目的

在真正的單細胞表達矩陣里，我們往往根據(jù)基因的表達模塊（每個模塊含若干基因，類似于bulk-seq的kmeans分表達模塊）

其中m1表示的是cell_1對應mod_A的一組基因表達量向量

那么線性表示為上圖，假設mod_A和mod_B能較好的區(qū)分細胞，那么降維后：

將設計矩陣其他mod的系數(shù)均改為1
而這些模塊在生物學意義上表示為某些具有特定功能的一類基因，比方說mod_A控制信號通路，mod_B控制細胞增殖，等等
ICA降維的好處是，它能夠保留細胞間差異最大的模塊（基因集），有了這種差異我們就可以知道哪些細胞中的哪些基因先表達，哪一些后表達，根據(jù)這些特征追溯軌跡；而使用ICA降維，需要滿足數(shù)據(jù)呈現(xiàn)非高斯分布，而在單細胞數(shù)據(jù)里很難說明數(shù)據(jù)是否符合非高斯分布

細胞擬時分析

這一塊分析又稱為細胞軌跡追溯，為什么要做這一塊分析呢，是因為在同一個組織中，有不同類型的細胞；也就是說，哪怕是統(tǒng)一組織，細胞的分化程度也是不一樣的，有的細胞分化的快，成為了成熟的細胞，有的細胞分化慢，還處于初級狀態(tài)或者中間狀態(tài)，正是因為有這樣的差異，我們根據(jù)基因表達情況和biomarker可以判斷出一條細胞“演化”的軌跡

那么常規(guī)的追溯細胞規(guī)矩的方法依賴于降維，降維后根據(jù)距離來判斷軌跡，也就是，距離相近的點（細胞）算作為一類細胞

這里有兩種方法計算，一種是利益ICA或者PCA進行線性降維；另外一種是利用diffusion maps算法轉換為概率進行非線性降維

如果是采用線性降維，則需要用反向圖嵌入：

那么對于單細胞表達矩陣來說，先進性PCA或者ICA進行降維然后根據(jù)降維后的cluster對cell進行聚類，根據(jù)聚類情況“擬合”出一條細胞軌跡，有些較為遠離軌跡的點就近歸于最靠近的那條軌跡上，并且再次將數(shù)據(jù)映射到原有的高維空間，進行矯正
最終根據(jù)生物學知識定根，這樣的話，軌跡就計算出了
總結一些，采用線性降維畫軌跡圖，先聚類，用每一類細胞的均值作為節(jié)點來定軌跡

再有另外一種采用的是非線性降維方法diffusion maps
該方法將點（細胞）間的距離化為概率，類似于t-SNE：

圖片來自：here
然后計算每一條路徑的概率：

如上圖，從1號點到6號點的路徑有：一，1->2->6；二,1->4->5->6；三,1->4->7->6等
由于任意兩點（細胞）間的距離我們都用概率表示了，那么比如1->2->6的分數(shù)為 p(1->2)×p(2->6)，圖中p(1|2)表示1->2兩個點（細胞）的概率，即p(1->2)
又譬如1->4->5->6，1->4->7->6這兩條路徑，都是從 1 這個點（細胞）經(jīng)過 4 這個點（細胞）后，分成兩條路徑，即經(jīng)過 5 這個點（細胞）或 7 這個點（細胞）{相當于分類成兩條路徑}，最終回到 6 這個點（細胞）；那么其概率為 p(1->4)×p(4->5)×p(5->6) + p(1->4)×p(4->7)×p(7->6)，圖中p(1|4)表示1->4兩個點（細胞）的概率，即p(1->4)

接下來，采用最小生成樹來畫軌跡，即點與點（細胞與細胞間連線）概率之和最小的路徑，并以此作為軌跡，可參考：Single cell RNA-seq data analysis with R視頻學習筆記（八）

RNA velocity模型

這一個模型也是一個經(jīng)典的細胞軌跡模型，其基本原理是對于某個基因來說，其轉錄本可以分為兩類，一類是unspliced mRNA，另一類是spliced mRNA。
而根據(jù)生物學過程，unspliced mRNA一定發(fā)生在spliced mRNA之前，那么在某一時刻的細胞，有實驗室認為，可以從mRNA剪切前后的比例來估計mRNA豐度隨時間的變化來估算細胞發(fā)育軌跡
細胞為了保持活力，永遠不會停止不前。在任何給定時刻，細胞都會上調或下調數(shù)百或數(shù)千個基因，這個過程的特征是不同的mRNA表達水平，在這種細胞分化過程中，數(shù)百個基因開啟或關閉，這取決于這個祖細胞試圖變成何種類型的成體細胞。比方說A細胞分化為B細胞的過程中，ab這個基因會大量表達（ab基因在A細胞中表達低，在B細胞中表達高），那么ab基因增加的趨勢暗示著細胞分化的方向（若ab基因表達量增加，那么此時unspiced mRNA的比例會增加）

這是不同基因在不同細胞時期unspliced mRNA和spliced mRNA含量變化圖

根據(jù)文獻：
RNA velocity of single cells來看，unspliced mRNA和spliced mRNA滿足的動力學方程為：

簡化模型，即α(t)為常數(shù)，β(t)為 1，γ(t)為常數(shù)

目前有兩個假設來估計γ
第一種，Model I
假設：

那么在短時間內忽略mRNA的降解，那么有：

因此，測量不同時間點的s(t)與時間t做線性回歸可以得到γ
但是這種方法只能在極短時間內有效，因為忽略降解

第二種，Model II
我們還可以假設未剪切的分子數(shù)保持不變，即u(t)=u0，則：

同樣測量不同時間點的s(t)與時間t做擬合可以得到γ，因此帶入微分方程的解，即可求出u(t)和s(t)了

那么，每個細胞的發(fā)育方向該如何確定呢？
PCA降維
我們看一下單細胞的表達矩陣：

那么表達矩陣每一個細胞（每一列）可以看作為一個向量，對于PCA降維來說，由于我們已知 u(t)和s(t) 這兩個unspliced mRNA和splice mRNA表達量隨時間 t 變化的函數(shù) ，那么利用現(xiàn)在時刻的表達矩陣，預估未來一小段時間的各個細胞的基因表達量，即求出未來一段時間 t的u(t)和s(t)，從而形成一個新的細胞基因表達向量，那么發(fā)育箭頭指向就是未來時刻的細胞基因表達向量減去現(xiàn)在時刻的細胞基因表達向量，由現(xiàn)在時刻的細胞基因表達向量指向未來時刻的細胞基因表達向量

PCA降維

t-SNE降維
對于t-SNE降維來說，由于是非線性降維，定義其方向不能像線性降維那么簡單，作者首先定義：

即求出ri,j與di的相關性
其中：

而 di，在Model I 的條件下：

其中t表示未來一小段時間，也就是說 di 表示的是未來時刻的細胞基因表達向量減去現(xiàn)在時刻的細胞基因表達向量
即設現(xiàn)在時刻某細胞的各基因表達向量為 s0 ，未來某段時間該細胞的各基因表達向量為 s(t) = s0 + t·ds/dt ，那么未來時刻與現(xiàn)在時刻的差值為 t·ds/dt

作者定義單細胞表達矩陣的細胞基因表達向量為：