綜合簡(jiǎn)書作者們的文章以及自己操作，綜合整理如下宦棺，以備后面使用契讲。
有問題歡迎交流
Robust enumeration of cell subsets from tissue expression profiles | Nature Methods
以上為Cibersort的文章
使用Cibersort工具需要三個(gè)文件：
1、Cibersort.R
2冠场、LM22.txt
3秆剪、genes_exp.txt
1赊淑、Cibersort.R
此文件為源代碼，在使用之前請(qǐng)閱讀一下代碼中的注釋段仅讽，安裝一下前置包
在R中創(chuàng)建script陶缺，復(fù)制以下代碼保存為 “Cibersort.R”。
提示：不需要去理解代碼洁灵，直接復(fù)制粘貼饱岸，運(yùn)行就ok了。對(duì)代碼感興趣的話當(dāng)我沒說徽千。
#' CIBERSORT R script v1.03 (last updated 07-10-2015)
#' Note: Signature matrix construction is not currently available; use java version for full functionality.
#' Author: Aaron M. Newman, Stanford University (amnewman@stanford.edu)
#' Requirements:
#' R v3.0 or later. (dependencies below might not work properly with earlier versions)
#' install.packages('e1071')
#' install.pacakges('parallel')
#' install.packages('preprocessCore')
#' if preprocessCore is not available in the repositories you have selected, run the following:
#' source("http://bioconductor.org/biocLite.R")
#' biocLite("preprocessCore")
#' Windows users using the R GUI may need to Run as Administrator to install or update packages.
#' This script uses 3 parallel processes. Since Windows does not support forking, this script will run
#' single-threaded in Windows.
#'
#' Usage:
#' Navigate to directory containing R script
#'
#' In R:
#' source('CIBERSORT.R')
#' results <- CIBERSORT('sig_matrix_file.txt','mixture_file.txt', perm, QN)
#'
#' Options:
#' i) perm = No. permutations; set to >=100 to calculate p-values (default = 0)
#' ii) QN = Quantile normalization of input mixture (default = TRUE)
#'
#' Input: signature matrix and mixture file, formatted as specified at http://cibersort.stanford.edu/tutorial.php
#' Output: matrix object containing all results and tabular data written to disk 'CIBERSORT-Results.txt'
#' License: http://cibersort.stanford.edu/CIBERSORT_License.txt
#' Core algorithm
#' @param X cell-specific gene expression
#' @param y mixed expression per sample
#' @export
CoreAlg <- function(X, y){
#try different values of nu
svn_itor <- 3
res <- function(i){
if(i==1){nus <- 0.25}
if(i==2){nus <- 0.5}
if(i==3){nus <- 0.75}
model<-e1071::svm(X,y,type="nu-regression",kernel="linear",nu=nus,scale=F)
model
}
if(Sys.info()['sysname'] == 'Windows') out <- parallel::mclapply(1:svn_itor, res, mc.cores=1) else
out <- parallel::mclapply(1:svn_itor, res, mc.cores=svn_itor)
nusvm <- rep(0,svn_itor)
corrv <- rep(0,svn_itor)
#do cibersort
t <- 1
while(t <= svn_itor) {
weights = t(out[[t]]$coefs) %*% out[[t]]$SV
weights[which(weights<0)]<-0
w<-weights/sum(weights)
u <- sweep(X,MARGIN=2,w,'*')
k <- apply(u, 1, sum)
nusvm[t] <- sqrt((mean((k - y)^2)))
corrv[t] <- cor(k, y)
t <- t + 1
}
#pick best model
rmses <- nusvm
mn <- which.min(rmses)
model <- out[[mn]]
#get and normalize coefficients
q <- t(model$coefs) %*% model$SV
q[which(q<0)]<-0
w <- (q/sum(q))
mix_rmse <- rmses[mn]
mix_r <- corrv[mn]
newList <- list("w" = w, "mix_rmse" = mix_rmse, "mix_r" = mix_r)
}
#' do permutations
#' @param perm Number of permutations
#' @param X cell-specific gene expression
#' @param y mixed expression per sample
#' @export
doPerm <- function(perm, X, Y){
itor <- 1
Ylist <- as.list(data.matrix(Y))
dist <- matrix()
while(itor <= perm){
#print(itor)
#random mixture
yr <- as.numeric(Ylist[sample(length(Ylist),dim(X)[1])])
#standardize mixture
yr <- (yr - mean(yr)) / sd(yr)
#run CIBERSORT core algorithm
result <- CoreAlg(X, yr)
mix_r <- result$mix_r
#store correlation
if(itor == 1) {dist <- mix_r}
else {dist <- rbind(dist, mix_r)}
itor <- itor + 1
}
newList <- list("dist" = dist)
}
#' Main functions
#' @param sig_matrix file path to gene expression from isolated cells
#' @param mixture_file heterogenous mixed expression
#' @param perm Number of permutations
#' @param QN Perform quantile normalization or not (TRUE/FALSE)
#' @export
CIBERSORT <- function(sig_matrix, mixture_file, perm=0, QN=TRUE){
#read in data
X <- read.table(sig_matrix,header=T,sep="\t",row.names=1,check.names=F)
Y <- read.table(mixture_file, header=T, sep="\t", row.names=1,check.names=F)
X <- data.matrix(X)
Y <- data.matrix(Y)
#order
X <- X[order(rownames(X)),]
Y <- Y[order(rownames(Y)),]
P <- perm #number of permutations
#anti-log if max < 50 in mixture file
if(max(Y) < 50) {Y <- 2^Y}
#quantile normalization of mixture file
if(QN == TRUE){
tmpc <- colnames(Y)
tmpr <- rownames(Y)
Y <- preprocessCore::normalize.quantiles(Y)
colnames(Y) <- tmpc
rownames(Y) <- tmpr
}
#intersect genes
Xgns <- row.names(X)
Ygns <- row.names(Y)
YintX <- Ygns %in% Xgns
Y <- Y[YintX,]
XintY <- Xgns %in% row.names(Y)
X <- X[XintY,]
#standardize sig matrix
X <- (X - mean(X)) / sd(as.vector(X))
#empirical null distribution of correlation coefficients
if(P > 0) {nulldist <- sort(doPerm(P, X, Y)$dist)}
#print(nulldist)
header <- c('Mixture',colnames(X),"P-value","Correlation","RMSE")
#print(header)
output <- matrix()
itor <- 1
mixtures <- dim(Y)[2]
pval <- 9999
#iterate through mixtures
while(itor <= mixtures){
y <- Y[,itor]
#standardize mixture
y <- (y - mean(y)) / sd(y)
#run SVR core algorithm
result <- CoreAlg(X, y)
#get results
w <- result$w
mix_r <- result$mix_r
mix_rmse <- result$mix_rmse
#calculate p-value
if(P > 0) {pval <- 1 - (which.min(abs(nulldist - mix_r)) / length(nulldist))}
#print output
out <- c(colnames(Y)[itor],w,pval,mix_r,mix_rmse)
if(itor == 1) {output <- out}
else {output <- rbind(output, out)}
itor <- itor + 1
}
#save results
write.table(rbind(header,output), file="CIBERSORT-Results.txt", sep="\t", row.names=F, col.names=F, quote=F)
#return matrix object containing all results
obj <- rbind(header,output)
obj <- obj[,-1]
obj <- obj[-1,]
obj <- matrix(as.numeric(unlist(obj)),nrow=nrow(obj))
rownames(obj) <- colnames(Y)
colnames(obj) <- c(colnames(X),"P-value","Correlation","RMSE")
obj
}
2苫费、LM22.txt
此文件為22種免疫細(xì)胞的標(biāo)志基因表達(dá)量，是衡量細(xì)胞含量的標(biāo)準(zhǔn)双抽。
去Cibersort的文章里下載Supplementry table 1百框，下載后打開如下：
Supplementry table 1
只選取如下含有數(shù)據(jù)的部分（其他部分自行探索），如下：
所需數(shù)據(jù)
復(fù)制粘貼為txt文件牍汹，注意籃圈標(biāo)記部分铐维，后面自己文檔的基因列名要與此保持一致柬泽。如下圖：
txt文件
3、gene_exp.txt
此文件是自己的數(shù)據(jù)嫁蛇，在R中處理時(shí)锨并，導(dǎo)出為“sep=\t”的“.txt”文件，需要注意的地方主要有幾點(diǎn)：
1.基因名不能有重復(fù)
2.整個(gè)矩陣不能有空值
3.基因的列名和LM22文件保持一致
4.數(shù)據(jù)格式要和LM22保持一致睬棚，fpkm/tpm不要log處理
我自己的數(shù)據(jù)格式第煮，如下：
我的數(shù)據(jù)
如果有報(bào)錯(cuò)，就把這兩個(gè)表復(fù)制到excel上抑党，去檢查一下數(shù)據(jù)與我這個(gè)Excel文件有什么區(qū)別包警，還有LM22是不是有問題。
轉(zhuǎn)成Excel后
我犯過的問題就有：兩個(gè)表基因列的列名不一致新荤；LM22文件范圍沒選對(duì)揽趾；基因名有重復(fù)和空值出現(xiàn)
這個(gè)腳本報(bào)錯(cuò)信息不詳細(xì)台汇，遇到問題來這里看看苛骨，自己核對(duì)一下
4、最終步驟
將三個(gè)文件放到一個(gè)文件夾苟呐，然后將R當(dāng)前工作目錄轉(zhuǎn)到那個(gè)文件夾（setwd函數(shù)）之后直接輸入以下代碼痒芝，運(yùn)行Cibersort.R，然后等待一段不短的時(shí)間牵素，會(huì)自動(dòng)生成結(jié)果文件”"CIBERSORT-Results.txt"“严衬。如果需要處理多組數(shù)據(jù)，要及時(shí)對(duì)結(jié)果文件重命名笆呆，否則會(huì)重寫為新的分析結(jié)果请琳。
setwd("")
source("Cibersort.R")
result1 <- CIBERSORT("LM22.txt", "genes_exp.txt", perm = 1000, QN = T)
# perm置換次數(shù)=1000，QN分位數(shù)歸一化=TRUE
# 文件名可以自定義
關(guān)于數(shù)據(jù)格式：
fragments per kilobase per million (FPKM) and transcripts per kilobase million (TPM), are suitable for use with CIBERSORT—《Profiling Tumor Infiltrating Immune Cells with CIBERSORT》
