【ggplot】在R#语言之中绘制火山图的代码实现

https://github.com/rsharp-lang/ggplot

在生物信息学中的组学数据分析领域内，有一个非常常见的数据可视化图表：应用于可视化两两组别比对结果的火山图。在火山图之中，X坐标轴一般是log2FC，纵坐标Y轴，则一般是t检验的pvalue的-log10转换之后的值。由于fold change有大于1的值，A/B大于1，表示A的表达量高于B的表达量，反之小于一表示A的表达量低于B的表达量。这样子fold change经过log2转换之后，就会出现负数，散点一般呈轴对称分布在X=0的位置周围。这样子绘制出来的散点图就有点类似于火山喷发的样子了。

R#语言中绘制散点图

因为自己最近在进行R#语言的ggplot数据可视化程序包的编写工作。在这里通过火山图的绘制来为大家展示R#语言之中的ggplot程序包目前的开发进度。因为火山图本质上就是一个散点图，所以在这里我们可以很轻松的在R#语言之中通过geom_point图层绘制出来：

require(ggplot);

const volcano = read.csv(`${@dir}/log2FC.csv`);

volcano[, "p.value"] = -log10(volcano[, "p.value"]);

print("peeks of the raw data:");
print(head(volcano));

bitmap(file = `${@dir}/volcano.png`, size = [3000, 3600]) {
    ggplot(volcano, aes(x = "log2FC", y = "p.value"), padding = "padding:350px 100px 250px 350px;")
       + geom_point(color = "black", shape = "Circle", size = 21)
       + labs(x = "log2(FoldChange)", y = "-log10(P-value)")
       + ggtitle("Volcano Plot")
       + scale_x_continuous(labels = "F2")
       + scale_y_continuous(labels = "F2")
    ;
}

执行上面的代码，可以得到下面的一个基本的散点图。可以看得出来，有火山喷发的样式了：

到这里，我们已经完成了三分之一的可视化工作。那，上面的demo图看起来有点朴素，我们还需要为上面的散点图添加上下调颜色来增加信息量。

ggplot中通过scale_colour_manual进行颜色管理

如果我们需要对散点图进行上色，我们可以直接设置散点图图层函数geom_point的color参数。这样子会直接为散点设置了统一的颜色。但是这个效果很明显不是我们所需要的。除了通过color参数进行简单粗暴的设置颜色，在ggplot程序包之中，还可以通过aes数据映射函数中的color参数更加灵活的设置颜色。

在讲解aes函数进行颜色映射之前，我们需要先计算出原始数据之中的颜色因素，例如上下调关系以及是否差异基因等：

volcano[, "factor"]  = ifelse(volcano[, "log2FC"] >  log2(1.25), "Up", "Not Sig");
volcano[, "factor"]  = ifelse(volcano[, "log2FC"] < -log2(1.25), "Down", volcano[, "factor"]);
volcano[, "factor"]  = ifelse(volcano[, "p.value"] < 0.05, volcano[, "factor"], "Not Sig");
volcano[, "p.value"] = -log10(volcano[, "p.value"]);

print("peeks of the raw data:");
print(head(volcano));
print("count of the factors:");
print(`Up:      ${sum("Up"      == volcano[, "factor"])}`);
print(`Not Sig: ${sum("Not Sig" == volcano[, "factor"])}`);
print(`Down:    ${sum("Down"    == volcano[, "factor"])}`);

# [1] "peeks of the raw data:"
#         ID        p.value   log2FC    factor
# <mode>  <string>  <double>  <double>  <string>
# [1, ]   "Q5XJ10"   9.81     -0.205    "Not Sig"
# [2, ]   "A8WG05"   5.21      0.0472   "Not Sig"
# [3, ]   "Q8JH71"   11.6     -0.38     "Down"
# [4, ]   "Q7T3L3"   8.85      0.165    "Not Sig"
# [5, ]   "Q567C8"   21        0.837    "Up"
# [6, ]   "Q92005"   0.091    -0.00514  "Not Sig"
#
# [1] "count of the factors:"
# [1]     "Up:      275"
# [1]     "Not Sig: 2267"
# [1]     "Down:    284"

接着，我们就可以根据颜色因素，进行上下调颜色的设置了：

geom_point(aes(color = "factor"), color = "black", shape = "circle", size = 25)
+ scale_colour_manual(values = list(
     Up        = "red",
     "Not Sig" = "gray",
     Down      = "skyblue"
  ))

在这里，我们通过aes函数对散点图设定颜色映射所使用的颜色因素来自于原始数据中的factor列，并且通过scale_colour_manual手动设定了每一个factor所对应的颜色值。

scale_colour_manual功能的实现

我们首先来看一下这个映射函数在底层代码中的实现：函数的代码非常的简单，就是创建一个颜色映射对象用于作为调整ggplot绘图样式的选项参数对象。

<ExportAPI("scale_colour_manual")>
Public Function scale_colour_manual(<RRawVectorArgument> values As Object, Optional env As Environment = Nothing) As ggplotOption
    Return New ggplotColorProfile With {.profile = ggplotColorMap.CreateColorMap(values, env)}
End Function

Public Shared Function CreateColorMap(map As Object, env As Environment) As ggplotColorMap
    If TypeOf map Is String Then
        Return stringMap(DirectCast(map, String))
    ElseIf map.GetType.IsArray Then
        Dim strArray As String() = REnv.asVector(Of String)(map)

        If strArray.Length = 1 Then
            Return stringMap(strArray(Scan0))
        Else
            Return directMap(strArray)
        End If
    ElseIf TypeOf map Is list Then
        Return New ggplotColorFactorMap With {
            .colorMap = DirectCast(map, list).AsGeneric(Of String)(env)
        }
    Else
        Throw New NotImplementedException(map.GetType.FullName)
    End If
End Function

我们进行颜色映射设置，一般而言会有一下的几种情况：

• 单个字符串：可能是将原始数据的值直接作为颜色值进行映射，也可能是一个颜色谱的名称
• 字符串向量：对所有的原始数据点进行一对一的直接映射
• list列表：对原始数据中的factor分类数据进行一对一的映射，例如Up分类映射为红色，Down分类映射为蓝色，Not Sig分类映射为灰色

在这里我们传递进入了一个list对象，将对应的factor映射为某一种特定的颜色。在ggplot程序包的底层代码之中，实现这样子的factor颜色映射操作过程是通过一个映射对象来提供一个映射函数：

''' <summary>
''' used the field value as color factor
''' </summary>
Public Class ggplotColorFactorMap : Inherits ggplotColorMap

    Public Overrides Function ColorHandler(ggplot As ggplot) As Func(Of Object, String)
        Dim colorMap As Dictionary(Of String, String) = Me.colorMap
        Return Function(keyObj) colorMap.TryGetValue(any.ToString(keyObj), [default]:="black")
    End Function
End Class

这样子在底层代码中进行散点图的绘制的时候，就可以自动的通过这个映射函数来创建处一个颜色向量用于数据可视化了：

If useCustomColorMaps Then
    Dim factors As String() = REnv.asVector(Of String)(DirectCast(ggplot.data, dataframe).getColumnVector(reader.color))
    Dim maps = DirectCast(colorMap, ggplotColorFactorMap).ColorHandler(ggplot)

    colors = factors.Select(Function(factor) maps(factor)).ToArray
    legends = New legendGroupElement With {
        .legends = DirectCast(colorMap, ggplotColorFactorMap) _
            .GetLegends(shape, ggplot.ggplotTheme.legendLabelCSS) _
            .ToArray
    }

火山图的最终DEMO

那，基于上面所提到的scale_colour_manual函数进行颜色映射管理，我们可以得到最终的火山图可视化DEMO代码：

require(ggplot);

const volcano = read.csv(`${@dir}/log2FC.csv`);

volcano[, "factor"]  = ifelse(volcano[, "log2FC"] >  log2(1.5), "Up", "Not Sig");
volcano[, "factor"]  = ifelse(volcano[, "log2FC"] < -log2(1.5), "Down", volcano[, "factor"]);
volcano[, "factor"]  = ifelse(volcano[, "p.value"] < 0.05, volcano[, "factor"], "Not Sig");
volcano[, "p.value"] = -log10(volcano[, "p.value"]);

print("peeks of the raw data:");
print(head(volcano));
print("count of the factors:");
print(`Up:      ${sum("Up"      == volcano[, "factor"])}`);
print(`Not Sig: ${sum("Not Sig" == volcano[, "factor"])}`);
print(`Down:    ${sum("Down"    == volcano[, "factor"])}`);

bitmap(file = `${@dir}/volcano.png`, size = [3000, 3000]) {
    ggplot(volcano, aes(x = "log2FC", y = "p.value"), padding = "padding:250px 500px 250px 300px;")
        + geom_point(aes(color = "factor"), color = "black", shape = "circle", size = 25)
       + scale_colour_manual(values = list(
          Up        = "red",
          "Not Sig" = "gray",
          Down      = "steelblue"
       ))
       + labs(x = "log2(FoldChange)", y = "-log10(P-value)")
       + ggtitle("Volcano Plot (A vs B)")
       + scale_x_continuous(labels = "F2")
       + scale_y_continuous(labels = "F2")
    ;
}

运行一下这一段代码，输出的结果图效果挺不错的：

• Author
• Recent Posts

谢桂纲

高级数据科学家 at 苏州帕诺米克

Working on Engineered bacteria CAD design on its genome from scratch. Writing scientific computing software for Tianhe & Sunway TaihuLight supercomputer. Do scientific computing programming in R/R# language, he is also the programming language designer of the R# language on the .NET runtime.

Latest posts by 谢桂纲 (see all)

• 【MZKit】简单自动化组织分区 - 2023年11月5日
• 【MZKit教程】质谱成像原始数据文件查看 - 2023年6月29日
• 生物序列图嵌入算法 - 2023年6月29日