ggscatter (ggpubr package) ~Rを用いた図の描画~

##ggpubr##
ggpubr は ggplot2 をベースに描かれた R パッケージで、比較的容易に様々な作図が可能です。（使ってみて、めちゃ便利でした。）ggpubr package に関しては以下参照。

・[ggpubr 論文記事１](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/)
・[ggpubr 論文記事２](https://rpkgs.datanovia.com/ggpubr/)

今回は、ggpubr package の **ggscatter 関数** を使って、様々な形式の散布図を作成しようと思います。
まずは、パッケージのインストールから行います。
R のインストールは [前の記事 (RainCloud plots)](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/) を参考にどうぞ。

```R:
# 必要なツール（パッケージ）のインストール (この作業は１回きりで大丈夫です)
install.packages("dplyr")
install.packages("tidyverse")
install.packages("ggpubr")
# 必要なツールを使えるようにする (この作業は R を立ち上げるたびに行いましょう)
library("dplyr")
library("tidyverse")
library("ggpubr")

# 乱数発生のための関数
cor_get <- function(r, n, family = 'normal', rep = 1000, seed = 1234, 
                    mean = 0, sd = 1, lambda = 1, size = 10, prob = 0.5, 
                    min = 0, max = 1){
  ...　以下省略 (https://www.rpubs.com/yyoneda/random_correlation) ...
}
```
乱数発生の関数は [この記事](https://www.rpubs.com/yyoneda/random_correlation) から引用しました。


次にインプットデータの作成です。この記事では、野生型 (WT)、変異体A (MutantA)、変異体B (MutantB) の背丈と重さのデータを **cor_get 関数**を用いて作成します。

```R:
# 野生型 (WT) のデータ作成
WT_HW <- cor_get(0.9, 100, mean=20, sd=0.2)
WT_HW_1 <- data.frame(WT_HW[[1]]$x, WT_HW[[1]]$y, rep("WT", 100))
colnames(WT_HW_1) <- c("Height", "Weight", "Class")

# 変異体A (MutantA) のデータ作成
MutantA_HW <- cor_get(0.5, 100, mean=20, sd=0.2)
MutantA_HW_1 <- data.frame(MutantA_HW[[1]]$x, MutantA_HW[[1]]$y, rep("MutantA", 100))
colnames(MutantA_HW_1) <- c("Height", "Weight", "Class")

# 変異体B (MutantB) のデータ作成
MutantB_HW <- cor_get(0.3, 100, mean=20, sd=0.2)
MutantB_HW_1 <- data.frame(MutantB_HW[[1]]$x, MutantB_HW[[1]]$y, rep("MutantB", 100))
colnames(MutantB_HW_1) <- c("Height", "Weight", "Class")

# ３つのデータの結合
AllClass_HW <- rbind(WT_HW_1, MutantA_HW_1, MutantB_HW_1)
head(AllClass_HW)

    Height   Weight Class
1 19.75859 26.55619    WT
2 20.05549 26.92860    WT
3 20.21689 26.88711    WT
4 19.53086 26.41052    WT
5 20.08582 26.85866    WT
6 20.10121 26.81220    WT
```

このように、１列目に背丈、２列目に重さ、３列目に分類群（野生型か変異体か）があるデータフレームを作成しました。

さて、本題の散布図作成に移ります。
まずは、デフォルト（？）の散布図を作成します。
使用する関数は、散布図を作成する [ggscatter](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/78-perfect-scatter-plots-with-correlation-and-marginal-histograms/) と図を保存する [ggexport](https://www.rdocumentation.org/packages/ggpubr/versions/0.2.5/topics/ggexport) です。

```R:
# はじめ : 作成したデータフレーム, x : 作成したデータフレームのx軸に列 (Height)
# y : 作成したデータフレームのy軸に列 (Weight)
# color : 色の塗り分け（WT, MutantA, MutantB で塗り分けたい）
# pallet : 塗り分ける色のパレット (JCO pallet に関して https://nanx.me/ggsci/reference/pal_jco.html) 
# size : ドットサイズ, alpha : ドットの透明度
# ggtheme = theme_bw() : 軸やラベル位置、グリッドなどの包括設定

g <- ggscatter(AllClass_HW, x = "Height", y = "Weight",
                color = "Class", palette = "jco",
                size = 2, alpha = 0.6, ggtheme = theme_bw()) 

# width & height : 出力画像の幅と高さ (デフォルト各480)
# PDF で出力する場合は、filename="保存したいディレクトリ/プロットの名前.pdf" (width & height の記入はいらない)
g %>% ggexport(filename="保存したいディレクトリ/プロットの名前.png", width = 480, height = 480)
```

結果、このような図が出力されます。
![scatter1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/595527/67909ec1-e7c0-920e-62a6-fd8883684b72.png)

次は、この図のそれぞれの分類群に Pearson 相関解析を行い、図示してみましょう。相関解析には、[stat_cor 関数](https://www.rdocumentation.org/packages/ggpubr/versions/0.2.5/topics/stat_cor) を使います。

```R:
# stat_cor(aes(color = XXX)) : 相関解析するグループ (今回は色分けした分類群ごとに行う)
# label.x & label.y : 相関係数などをプロットする位置 (今回は x 軸は３つとも 20.25 の位置に、y 軸はそれぞれ26.2, 26.0, 25.8 の位置にプロット)
# 他にも method 引数があり、Pearson (デフォルト), Spearman, Kendall が指定できる

g <- ggscatter(AllClass_HW, x = "Height", y = "Weight",
                color = "Class", palette = "jco",
                size = 2, alpha = 0.6, ggtheme = theme_bw()) +
  stat_cor(aes(color = Class), label.x = 20.25, label.y = c(26.2, 26.0, 25.8)) 
g %>% ggexport(filename="保存したいディレクトリ/プロットの名前.png", width = 480, height = 480)
```

![scatter2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/595527/a15f53c5-ae6c-e4a5-ab58-fcc96e6d01ca.png)

次は、散布図に回帰直線と信頼区間の描写を行います。

```R:
# add  = "reg.line" : 散布図に加える図 (今回は reg.line: regression line を追記)
# conf.int = TRUE : 信頼区間の描写をオン

g <- ggscatter(AllClass_HW, x = "Height", y = "Weight",
                color = "Class", palette = "jco",
                size = 2, alpha = 0.6, add = "reg.line", conf.int = TRUE, ggtheme = theme_bw()) +
  stat_cor(aes(color = Class), label.x = 20.25, label.y = c(26.2, 26.0, 25.8)) 
g %>% ggexport(filename="保存したいディレクトリ/プロットの名前.png", width = 480, height = 480)
```
![scatter4.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/595527/98d86167-2af7-9e79-41d6-78d210d528c6.png)

次は、このグラフをそれぞれ単独に描写し、ひとつの図にまとめます。

```R:
# facet.by = "Class" : 分類群ごとに分けて作図

g <- ggscatter(AllClass_HW, x = "Height", y = "Weight",
                color = "Class", palette = "jco",
                size = 2, alpha = 0.6, add = "reg.line", conf.int = TRUE, ggtheme = theme_bw(),
          facet.by = "Class") +
  stat_cor(aes(color = Class), label.x = 20.25, label.y = c(26.2, 27.0, 25.8)) 

# nrow = 1, ncol = 3 : １行３列で作図

g %>% ggexport(filename="保存したいディレクトリ/プロットの名前.png", width = 480*3, height = 480, nrow=1, ncol=3)
```
![scatter5.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/595527/e77b9bde-16d7-e9c1-31d7-8529a7b80c64.png)

最後に、散布図にラグ（密度情報）の追加、ドット等のカラー変更を行います。

```R:
# palette = c("magenta", "darkgray", "skyblue") : 色の変更
# rug=TRUE : ラグの追記

g <- ggscatter(AllClass_HW, x = "Height", y = "Weight",
                color = "Class", palette = c("magenta", "darkgray", "skyblue"), rug=TRUE, 
                size = 2, alpha = 0.8, add = "reg.line", conf.int = TRUE, ggtheme = theme_bw()) +
  stat_cor(aes(color = Class), label.x = 20.25, label.y = c(26.2, 26.0, 25.8)) 
g %>% ggexport(filename="保存したいディレクトリ/プロットの名前.png", width = 480, height = 480)
```
![scatter7.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/595527/17397257-cae2-162a-c240-c65a0fe553da.png)

以上です。奥が深いです。

以下参照
・[ggpubr パッケージの論文記事１](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/)
・[ggpubr パッケージの論文記事２](https://rpkgs.datanovia.com/ggpubr/)
・[前の記事 (R のインストール等)](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/)
・[乱数発生の関数の引用](https://www.rpubs.com/yyoneda/random_correlation) 
・[ggscatter 関数](http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/78-perfect-scatter-plots-with-correlation-and-marginal-histograms/) 
・[ggexport 関数](https://www.rdocumentation.org/packages/ggpubr/versions/0.2.5/topics/ggexport) 
・[stat_cor 関数](https://www.rdocumentation.org/packages/ggpubr/versions/0.2.5/topics/stat_cor)