【CELL】compass用單細(xì)胞轉(zhuǎn)錄組+流平衡分析預(yù)測(cè)細(xì)胞代謝狀態(tài)
細(xì)胞代謝調(diào)控正常細(xì)胞功能以及多種疾病狀態(tài)的病理生理。最近，免疫細(xì)胞代謝研究（immunometabolism）成為一個(gè)研究熱點(diǎn)，揭示了包括抗病毒免疫、自身免疫和抗腫瘤反應(yīng)在內(nèi)的炎癥代謝調(diào)節(jié)。然而，由于代謝網(wǎng)絡(luò)的規(guī)模和復(fù)雜性，某個(gè)代謝的擾動(dòng)（perturbation）能夠產(chǎn)生級(jí)聯(lián)效應(yīng)(cascade effects)并改變代謝網(wǎng)絡(luò)中看似較遠(yuǎn)的部分或者讓經(jīng)典通路“彎道超車”。因此，需要在系統(tǒng)水平上觀測(cè)特定反應(yīng)/酶，以理解代謝及其在疾病中的異常調(diào)控。在CELL文章《Metabolic modeling of single Th17 cells reveals regulators of autoimmunity》中作者提出了一種基于單細(xì)胞RNA測(cè)序和流平衡分析（flux balance analysis）描述細(xì)胞代謝狀態(tài)的算法 – compass（https://github.com/YosefLab/Compass）。
流平衡分析（flux balance analysis）將代謝網(wǎng)絡(luò)拓?fù)鋵W(xué)和化學(xué)計(jì)量學(xué)知識(shí)翻譯成數(shù)學(xué)對(duì)象，并用在代謝流預(yù)測(cè)上。提供了一種非常有用的場(chǎng)景：根據(jù)單細(xì)胞基因表達(dá)譜或者傳統(tǒng)RNA-seq基因表達(dá)譜研究細(xì)胞的代謝異質(zhì)性。
圖1. COMPASS算法
如圖1所示，使用單細(xì)胞轉(zhuǎn)錄組（Th17p細(xì)胞和Th17n細(xì)胞）和代謝模型作為COMPASS的輸入，經(jīng)過(guò)流平衡分析，找到了兩群細(xì)胞間異常的代謝通路（圖2）及相關(guān)的反應(yīng)/酶，并進(jìn)行了分子和功能驗(yàn)證，以及小鼠體內(nèi)驗(yàn)證。
圖2.4個(gè)差異通路
如圖2所示：左上：糖酵解通路（Glycolysis）在Th17p細(xì)胞中更加活躍（反應(yīng)/酶幾乎都分布在虛線右側(cè)）；右上：三羧酸循環(huán)（TCA cycle）中大部分反應(yīng)/酶在Th17p細(xì)胞中更活躍。左下：脂肪酸氧化的3個(gè)酶在Th17n細(xì)胞中更活躍；右下：氨基酸代謝中，絲氨酸羥甲基化轉(zhuǎn)移酶在Th17p細(xì)胞中更活躍，而精氨酸脫羧酶在Th17n細(xì)胞中更活躍。
看到這里，大家應(yīng)該對(duì)文章有了基本的理解。細(xì)胞代謝網(wǎng)絡(luò)太復(fù)雜，難以研究，然而我們手上不是有現(xiàn)成的一堆單細(xì)胞/bulk RNA-seq的表達(dá)數(shù)據(jù)么，我們可以通過(guò)大量樣品和嚴(yán)格的假設(shè)檢驗(yàn)以酶的mRNA表達(dá)為橋梁來(lái)研究細(xì)胞代謝，將基因表達(dá)-代謝反應(yīng)/酶活性-細(xì)胞代謝狀態(tài)聯(lián)系起來(lái)。例如我們知道了某反應(yīng)在癌細(xì)胞中活躍，就可以通過(guò)抑制劑阻斷（或者平衡）這個(gè)反應(yīng)，改變免疫細(xì)胞的代謝狀態(tài)，從而達(dá)到治療癌癥的目的。
1．安裝
Compass是用python編寫的，需要python 3環(huán)境，numpy，scipy，matplotlib等包。
由于計(jì)算量太大，還需要一個(gè)IBM CPLEX Optimization Studio優(yōu)化器（學(xué)術(shù)用戶可以注冊(cè)下載，注意版本）
conda 個(gè)python3環(huán)境
pip install numpy, scipy, matplotlib
下載IBM CPLEX Optimization Studio(需要學(xué)術(shù)版)
./xxx.bin 裝到software目錄（自己指定），一路enter或者yes
cd /home/software/CPLEX_Studio129/cplex/python/3.6/x86-64_linux
python setup.py install
python -m pip install git+https://github.com/yoseflab/Compass.git --upgrade
compass -h 測(cè)試安裝是否成功
2.運(yùn)行
輸入：表達(dá)矩陣（tpm，cpm等，一個(gè)基因一行，一個(gè)樣品一列），單細(xì)胞的話配上genes和sample文件。運(yùn)行速度很慢，可以多線程加速
圖3. Expression.csv輸入
compass --data expression.csv --num-processes 10 --species homo_sapiens （常規(guī)轉(zhuǎn)錄組）
compass --data-mtx expression.mtx genes.tsv sample_names.tsv --num-processes 10 --lambda 0.25 --species homo_sapiens （單細(xì)胞）
3．結(jié)果
運(yùn)行完上面的命令后，會(huì)在當(dāng)前目錄生成reaction.csv文件(圖4)。其中的值是反應(yīng)罰分（reaction penalties），高得分表示該反應(yīng)的可能性較低。
圖4. reaction.csv輸出
4. 后處理
步驟3獲得的是每個(gè)反應(yīng)在每個(gè)樣品中的反應(yīng)罰分，我們需要根據(jù)Recon2 meta信息（圖5）將A列的反應(yīng)編碼轉(zhuǎn)成人們可讀的代謝通路名字，將反應(yīng)罰分進(jìn)行轉(zhuǎn)化（轉(zhuǎn)成數(shù)值大的表示反應(yīng)活性高，加1，取-log）,然后使用wilcoxon或者Cohen’d進(jìn)行差異分析。
圖5. 后處理需要的3個(gè)文件
圖6. 最終結(jié)果表格，根據(jù)這個(gè)表格尋找更興趣的反應(yīng)/代謝，繪圖
各列說(shuō)明：
A：reaction編碼， “_pos” 表示化學(xué)反應(yīng)往右側(cè)，“_neg”表示反應(yīng)向左側(cè)
B/C/E：wilcoxon統(tǒng)計(jì)結(jié)果
D：cohens_d：cohen’s d用于計(jì)算兩組均值差異的效應(yīng)量， d = 0.2(效果小) ; d = 0.5(效果中) ; d = 0.8(效果大)。供參考
F：reaction編碼
G：reaction名字
H：reaction公式
I：相關(guān)基因
J：reaction所處的subsystem（例如酪氨酸代謝）
E：酶的EC號(hào)
L：置信度， 4 = most confident; 0 = unassigned confidence
注：1，由于單細(xì)胞計(jì)算量非常大，作者提供了一個(gè)micropooling腳本，將細(xì)胞劃分為簇，以簇的平均值代表這類細(xì)胞，然后進(jìn)行分析。
2，You may also apply Compass to bulk transcriptome data sets (e.g. bulk RNA-seq or microarray data sets) if there are enough observations (samples) to gain statistical power。可以用在常規(guī)轉(zhuǎn)錄組，甚至芯片數(shù)據(jù)集上。
有單細(xì)胞測(cè)序數(shù)據(jù)、常規(guī)轉(zhuǎn)錄組測(cè)序數(shù)據(jù)、芯片數(shù)據(jù)、代謝組數(shù)據(jù)的小伙伴可以試試這款軟件。又打開(kāi)了一個(gè)新的免疫代謝相關(guān)發(fā)文思路！
微生信助力高分文章，用戶63000+，引用830+
