RNA-Seq分析中（二）——用R画带基因名标签的MA图

MA图

MA plot即M-versus-A plot，在芯片数据处理出现之前也称为Bland-Altman plot，是由发明者名字命名的，而MA plot是对M与A作图而得名，M是minus的缩写，代表两个值之差，A是add的缩写，代表两个值之和。有研究者也把MA plot称为Ratio-Intensity (RI) plots，同时MA也正好是micro-array的简写。

一、MA图简介

MA图主要应用在基因组数据可视化方面，实现数据分布情况的展示。早期主要应用于DNA芯片数据，现在常用于高通量测序数据中基因差异表达分析结果的展示。

其计算公式如下：

M一般做Y轴，A一般做X轴。

M常对应差异表达分析获得的差异对比组之间基因表达变化log2FC。

A可以利用差异对比组的FPKM进行计算，以R和G来表示差异对比组的话，可以取R组基因的平均FPKM和G组基因的平均FPKM进行计算。

二、MA图绘制

1. 数据准备（IL-1B_1_IL-1B_2_vs_Cu5_1_Cu5_2.all.xls）：

第一列为基因ID，*_FPKM为样品FPKM值（列名必须包含"FPKM"），倒数第二列列为FDR，最后一列列为log2FC。

2. 我们使用R语言来绘制MA图，代码如下：

# 安装ggpubr和ggthemes 包
#install.packages("ggpubr")
#install.packages("ggthemes")
# 加载ggpubr和ggthemes包
library(ggpubr)
library(ggthemes)

# 读取数据
#deg.data <- read.csv("IL-1B_1_IL-1B_2_vs_Cu5_1_Cu5_2.all.xls", header = T, sep = "\t")
deg.data <- read.table("IL-1B_1_IL-1B_2_vs_Cu5_1_Cu5_2.all.xls",header=TRUE,check.names =F, sep = "\t")
# 显示文件
head(deg.data)

colnames(deg.data)<-read.delim("IL-1B_1_IL-1B_2_vs_Cu5_1_Cu5_2.all.xls", header=F,check.names =F,stringsAsFactors=F,nrows=1)
f <- grepl("FPKM",colnames(deg.data)) #提取所有匹配FPKM的列，下面对所有该列值取平均数
fpkm <- deg.data[ , f] 
head(fpkm)

# 新加一列log10fpkm
deg.data$log10fpkm = log10(rowMeans(fpkm)) #对FPKM值取平均数后取log10()

head(deg.data)

# 对差异FDR值进行log10转换
deg.data$logQ <- -log10(deg.data$FDR)

# 绘制基本热图
ggscatter(deg.data, x = "log10fpkm", y = "log2FC") + theme_base()

# 新加一列Group
deg.data$Group = "normal"
# 将FDR小于0.05，logFC小于2的基因设为显著下调基
deg.data$Group[which( (deg.data$FDR < 0.01) & (deg.data$log2FC > 1) )] = "up"
deg.data$Group[which( (deg.data$FDR < 0.01) & (deg.data$log2FC < -1) )] = "down"
# 查看上调和下调基因数
table(deg.data$Group)

# 绘制新的MA图
ggscatter(deg.data, x = "log10fpkm", y = "log2FC",
          color = "Group") + theme_base()

# 改变MA图颜色(palette)和点的大小(size)
ggscatter(deg.data, x = "log10fpkm", y = "log2FC",
          color = "Group", 
          palette = c("green", "gray", "red"),
          size = 1) + theme_base()

# 新加一列Label，为了能在MA图中标记指定基因名
deg.data$Label = ""
# 对差异表达基因的log2FC值进行从大到小排序取top 20
deg.data <- deg.data[order(abs(deg.data$log2FC),decreasing = T), ]
log2FC.genes <- head(deg.data$ID, 20)
#head(deg.data)
# 对差异表达基因的log10(fpkm) 值进行从大到小排序取top 20
deg.data <- deg.data[order(abs(deg.data$log10fpkm),decreasing = T), ]
fpkm.genes <- head(deg.data$ID, 20)

# 将log2FC.genes 和fpkm.genes合并，并加入到Label
deg.top20.genes <- c(as.character(log2FC.genes), as.character(fpkm.genes))
deg.data$Label[match(deg.top20.genes, deg.data$ID)] <- deg.top20.genes
#print (deg.data$Label)


# 改变MA图点的颜色和坐标轴标注，使图片更美观
ggscatter(deg.data, x = "log10fpkm", y = "log2FC",
          color = "Group", 
          palette = c("#00BA38", "#BBBBBB", "#F8766D"),
          size = 2,
          label = deg.data$Label, 
          font.label = 8, 
          repel = T,
          xlab = "log10(fpkm)", 
          ylab = "log2(FC)") + 
  theme_base() 

三、MA图结果

MA图.png