细胞因子-热图

rm(list = ls())

library(limma) 

exp=副本ATP.细胞因子.all
exp=t(exp)
dat=normalizeBetweenArrays(exp)
boxplot(dat,outline=FALSE, notch=T,col=9, las=2)

group_list=c(rep('NC',4),rep('ATP',4))

library(pheatmap)
# 默认绘图

batch <- c(rep("1",3),rep("2",1),rep("1",3),rep("2",1))
expr <- removeBatchEffect(dat, batch)

expr=normalizeBetweenArrays(expr)
boxplot(expr,outline=FALSE, notch=T,col=9, las=2)
expr <- expr[-49,]

df = expr[apply(expr, 1, function(x) sd(x)!=0),]
df = df[,apply(df, 2, function(x) sd(x)!=0)]

pheatmap(df)
# scale = "row"参数对行进行归一化
pheatmap(df, scale = "row")

# clustering_method参数设定不同聚类方法，默认为"complete",可以设定为'ward', 'ward.D', 'ward.D2', 'single', 'complete', 'average', 'mcquitty', 'median' or 'centroid'
pheatmap(df,scale = "row", clustering_method = "average")

# clustering_distance_rows = "correlation"参数设定行聚类距离方法为Pearson corralation，默认为欧氏距离"euclidean"
pheatmap(df, scale = "row", clustering_distance_rows = "correlation")
# cluster_row = FALSE参数设定不对行进行聚类
pheatmap(df, cluster_row = FALSE)

# legend_breaks参数设定图例显示范围，legend_labels参数添加图例标签
pheatmap(df, legend_breaks = c(1:5), legend_labels = c("1.0","2.0","3.0","4.0","5.0"))

# display_numbers = TRUE参数设定在每个热图格子中显示相应的数值，number_color参数设置数值字体的颜色
pheatmap(df, display_numbers = TRUE,number_color = "blue")

library(corrplot)

M <- cor(df)
M <-df
M=t(M)

library(pheatmap)
anno <- data.frame(sampleType=group_list)
rownames(anno) <- rownames(M)#满足条件
pheatmap(M)
pheatmap(M,scale = "row", clustering_method = "ward.D2",annotation = anno,treeheight_row=0, treeheight_col=0,border_color=NA)

install.packages("ggcorrplot")
install.packages("ggthemes")
library(ggcorrplot)
library(ggthemes)
#计算这几个样本的相关性
sample_cor <- round(cor(M1),3)
#计算P值，看显著性
sample_P <- round(cor_pmat(M1),3）
ggcorrplot(sample_cor,#计算的相关系数矩阵
method = 'circle',#热图显示为圆圈
hc.order = T,#聚类
hc.method = "ward.D",#聚类方法
outline.color = "white",
ggtheme = theme_bw(),
type = "upper",#显示上三角
colors = c("#6D9EC1","white","#E46726"),#热图颜色
lab = T,#热图上显示相关系数
lab_size = 2.5)

计算差异基因

exprSet = M
ex <- exprSet
qx <- as.numeric(quantile(ex, c(0., 0.25, 0.5, 0.75, 0.99, 1.0), na.rm=T))
LogC <- (qx[5] > 100) ||
(qx[6]-qx[1] > 50 && qx[2] > 0) ||
(qx[2] > 0 && qx[2] < 1 && qx[4] > 1 && qx[4] < 2)

if (LogC) { ex[which(ex <= 0)] <- NaN
exprSet <- log2(ex)
print("log2 transform finished")}else{print("log2 transform not needed")}

group_list=factor(group_list)
design= model.matrix(~0 + group_list)
colnames(design) <- levels(group_list)
design

library(limma)
contrast.matrix=makeContrasts("ATP-NC",
levels=design)
contrast.matrix
fit <- lmFit(exprSet, design)
fit2 <- contrasts.fit(fit, contrast.matrix)
fit2 <- eBayes(fit2)

allDiff=topTable(fit2,adjust='fdr',coef=1,number=Inf)
write.csv(allDiff,file="Diff_ATP_cytokins2.csv")

library(ggplot2)#加载ggplot2包
library(ggrepel)#用来添加名字

Cytokin=t(副本ATP.细胞因子.all)
write.csv(Cytokin,file="Cytokin.csv")
write.csv(exprSet,file="Cytokin-removebatch-logtrans.csv")

原文地址：https://www.cnblogs.com/jolinstory/p/16001267.html