Julia で VolcanoPlot 完結編

Julia での VolcanoPlot を描くプログラムが見つからなかったので，以下を参考にして，プログラムを書いた。

1) は非常に参考になった。ただ近接するデータ点のアノテーションは Julia では難しいので，描画結果を時後に部分選択・拡大を行うことで代替することとした。

参考サイト

1) 19.11 Volcano plots
   https://biocorecrg.github.io/CRG_RIntroduction/volcano-plots.html

2) Rでvolcano plotを描いてみる
   https://qiita.com/tomo_metal/items/619e08f6dc18100884ee

3) RNA-seq解析 - R言語でvolcano plotを描く
   https://qiita.com/shibanosuke/items/b4aa286993fd23a5c2ea

4) Pythonで作るVolcano plot
   https://qiita.com/insilicomab/items/c382f10208e8bd2482b5

プログラム

function volcano_plot( # 最初の3つの引数さえ与えればよい。描画色の指定をするなら4番目の引数まで指定する。あとは，お好み次第。
            log2FoldChange, # log2(変化倍率)
            pvalue, # p 値
            name; # 名前
            color=[:black, :blue, :red], # 描画点の色, "NO", "DOWN", "UP" の順。':' に続いて一般的な色の名称を指定
            markersize=3, # 描画点の大きさ
            alpha = 0.5, # 描画点のα値
            grid=false, # グリッドを書くか？ false/true
            tickfontsize=6, # 目盛りの文字サイズ
            xlabel="log2FoldChange", # 横軸の名前
            ylabel="-log10(pvalue)", # 縦軸の名前
            guidefontsize=8, # 軸名の文字サイズ
            vlinecolor=:red, # 追加の縦線の色の指定
            vlinewidth=1, # 追加の縦線の線の幅の指定
            vlinestyle=:dot, # :solid, :dash, :dot, :dashdot, :dashdotdot # 線の種類
            hlinecolor=:red, # vline* に準ずる
            hlinewidth=1,
            hlinestyle=:dot, # :solid, :dash, :dot, :dashdot, :dashdotdot
            annotatefontsize=8, # 名前のフォントサイズ
            width=500, height=600, # 出力画像ファイルの横幅と高さをピクセル単位で指定
            threshold_log2FoldChange = 0.6, # log2FoldChange の閾値 ± threshold_log2FoldChange
            threshold_pvalue = 0.05, # 水平線の描画位置
            leftright=15, # 名前 name を右揃えにするか揃えにするかの閾値
            ypos = 999 # "DOWN"/"UP" を描く y 座標（デフォルトでは自動算出）
            )
      assigncolor(x) = x == "NO" ? color[1] : x == "DOWN" ? color[2] : color[3]
      pyplot(grid=grid, size=(width, height), label="")
      n = length(log2FoldChange)
      diffexpressed = fill("NO", n);
      diffexpressed[(log2FoldChange .>  threshold_log2FoldChange) .& (pvalue .< threshold_pvalue)] .= "UP";
      diffexpressed[(log2FoldChange .< -threshold_log2FoldChange) .& (pvalue .< threshold_pvalue)] .= "DOWN";
      plt = scatter(log2FoldChange, -log10.(pvalue),
            color=assigncolor.(diffexpressed),
            markersize=markersize, markerstrokewidth=0, alpha=alpha,
            tick_direction=:out, tickfontsize=tickfontsize,
            xlabel=xlabel, ylabel=ylabel,
            guidefontsize=guidefontsize)
      vline!([-0.6, 0.6], color=vlinecolor, linewidth=vlinewidth, linestyle=vlinestyle)
      hline!([-log10(0.05)], color=hlinecolor, linewidth=hlinewidth, linestyle=hlinestyle)
      ypos = ypos == 999 ? (maximum(-log10.(de.pvalue))-log10.(0.05)) / 1.5 : ypos
      annotate!(-threshold_log2FoldChange, ypos, text(" DOWN", annotatefontsize, color[2], :left))
      annotate!( threshold_log2FoldChange, ypos, text("UP ",   annotatefontsize, color[3],  :right))
      for i = 1:size(de, 1)
            x = log2FoldChange[i]
            y = -log10(pvalue[i])
            d = diffexpressed[i]
            if d != "NO"
                  annotate!(x, y,
                        text(" $(name[i]) ", 6, assigncolor(d),
                        y > leftright ? (d == "DOWN" ? :left : :right) : (d == "DOWN" ? :right : :left)))
            end
      end
      plt
end

使用例

using RData, DataFrames, Plots

データ読み込み

# Download the data we will use for plotting
# download.file("https://raw.githubusercontent.com/biocorecrg/CRG_RIntroduction/master/de_df_for_volcano.rds", "de_df_for_volcano.rds", method="curl"

tmp = load("de_df_for_volcano.rds");de = dropmissing(tmp);

引数に与えるベクトルデータの設定

log2FoldChange = de[!, 2];
pvalue = de[!, 3];
name = de[!, 1];

関数呼び出し

plt = volcano_plot(log2FoldChange, pvalue, name)
savefig("fig4.png")

データ点が密になる箇所がある場合は，関数の戻り値を変数に代入し（例では plt），x 範囲，y 範囲を指定して再描画する。
この操作（範囲の再設定・再描画）は何回でも行うことができるので，描画結果を見ながら最適な描画範囲を設定できる。

plot(plt, xlims=(-0.9, -0.5), ylims=(6, 18), size=(400, 400))
savefig("fig5.png")