パッケージ
まず使用するパッケージ群をロードします。
library(tidyverse) # 言わずもがな
library(psych) # 相関係数の総当たり計算に使用
library(tidygraph) # ネットワーク構築の計算
library(ggraph) # ネットワーク作図
データ用意
前提1
すでに相関係数を算出したい数値データの列のみになっている状態を仮定します。処理区等の相関係数算出に関係のないカテゴリデータの列は事前に落としておいてください。
前提2
欠損値がない状態を仮定します。欠損値は事前にna.omit()なり補完するなりしておいてください。
閾値設定&乱数設定
ネットワークに反映する相関係数およびFDRの有意水準の閾値を設定しておきます。
threshold_fdr <- 0.05 # 有意水準0.05
threshold_r <- 0.2 # 相関の有無の閾値を0.2とする
set.seed(100) # 乱数設定
相関係数を算出し、ネットワークのノードとエッジの情報を取得する
ノード:ネットワークの点
エッジ:ネットワークの線
以下解析の流れ
- まず、相関係数を算出しデータフレームにまとめます。前提に当てはまるデータフレームが変数「raw_data」に格納されていると仮定します。
- 変数「df_cor」 に総当たり相関の結果を格納します。(行数は 変数の数(列数)C2)
- 変数「df_edge」にエッジの情報を格納します。(行数はデータにより異なる)
- 変数「df_node」にノードの情報を格納します。(行数は変数の数(列数))
- 変数「d_graph」にtbl_graph形式(ggraphパッケージに使用)のデータを格納します。
「df_cor」「df_edge」「d_graph」はパイプ処理の途中経過を永続代入<<-
を用いてGlobal環境に保存しています。
※相関係数を閾値に基づき抽出する際、正のみ(0.2 < r < 1)負のみ (-1 < r < -0.2) 抽出する場合と正負両方 (-1 < r < -0.2, 0.2 < r < 1) 抽出する場合では途中のfilter関数の引数が若干異なります。今回は正負両方抽出を仮定し、正のみのケースはコメントアウトしています。
前提
今回はネットワークの中心性には媒介中心性、コミュニティ検出にはスピングラス法を用います。ここの方法はいろいろありますので他の方法を試したければ調べてください。Tips: コミュニティ検出
# 最終出力がdf_nodeとなる
df_node <-
raw_data %>%
# 総当たりの相関行列算出
corr.test(method = "pearson") %>%
# データフレームにする
print(short = F) %>%
mutate(
# 相関係数の正負を判別する列を追加
posi_nega = case_when(
raw.r > 0 ~ "positive",
raw.r < 0 ~ "negative",
),
# p値をFDRで補正
fdr = p.adjust(raw.p, method = "BH"),
# 行名(どの変数間の相関関係かの情報)を列として追加
name = rownames(.),
) %>%
# name列を"-"で分割し2列(A,B)に分ける
separate(
col = name,
into = c("A", "B"),
sep = "-"
) %>%
# 必要な列のみ選択
select(A, B, raw.r, raw.p, posi_nega, fdr) %>%
# ここまでの結果をdf_corに格納
{
df_cor <<- (.)
} %>%
# 設定した閾値に基づきフィルタリング
filter(
# 有意水準でフィルタリング
fdr < threshold_fdr,
# 相関係数でフィルタリング
# 正のみ抽出する場合
# raw.r > threshold_r
# 正負両方抽出する場合
raw.r < -(threshold_r) | raw.r > threshold_r
) %>%
# ここまでの結果をdf_edgeに格納
{
df_edge <<- (.)
} %>%
# tbl_graphオブジェクトに変換 (無向グラフなのでdirected = FALSE)
as_tbl_graph(directed = FALSE) %>%
mutate(
# 媒介中心性を追加
centrality = centrality_betweenness(),
# コミュニティ検出 (スピングラス法) 結果を追加
community = as.factor(group_spinglass())
) %>%
# ここまでの結果をd_graphに格納
{
d_graph <<- (.)
} %>%
# データフレームに変換
as.data.frame()
ここまでで、df_cor, df_edge, df_node, d_graphが出力されます。
補足情報
各パッケージに含まれる関数(tidyverseを除く)
psychパッケージ:corr.test()
tidygraphパッケージ:as_tbl_graph(), centrality_betweenness(), group_spinglass()
ネットワーク可視化
変数「d_graph」に格納したtbl_graphオブジェクトを用いてggraphパッケージによるネットワーク可視化を行います。
相関係数の正負でエッジを色分けしたい場合と全部同色で良い場合でgeom_edge_link()の引数が異なり、色指定にはscale_edge_color_manual()を使用します。今回は正負で色分けを仮定し、全部同色のケースはコメントアウトしています。
最後の方のテーマ設定はお好みで。今回は背景なしにしたかったのでggplot2の関数であるtheme_void()を使用し、同じくggplot2の関数であるtheme()で凡例を非表示にしました。このようにggraphではggplot2のテーマが適用可能です。
plot <-
d_graph %>%
# レイアウトの指定
ggraph(layout = "kk") +
# エッジの追加 (正負で色分けする場合)
geom_edge_link(aes(width = raw.r, color = posi_nega)) +
scale_edge_color_manual(values = c("red", "gray")) +
# エッジの追加 (全てgrayの場合)
# geom_edge_link(aes(width = raw.r), color = "gray") +
# エッジ太さ調節
scale_edge_width(range = c(0.1,1.2)) +
# ノードの追加 (コミュニティで色分け、媒介中心性を大きさで示す)
geom_node_point(aes(colour = community, size = centrality)) +
# ラベルの追加
geom_node_label(aes(label = name, colour = community, size = 100), repel = TRUE) +
# ノードの色分けを指定(コミュニティ数により数は異なる)
scale_color_manual(values = c("darkorange", "dodgerblue","darkgreen", "violet")) +
# テーマ設定 (ここ以下はご自由に)
theme_void() +
# 凡例非表示の場合
theme(legend.position = "none")
# プロット出力
print(plot)
ここまででネットワーク図が出力できるかと思います。途中で出力されるデータフレーム等はwrite.csv()等でよしなに保存してください。
Tips: コミュニティ検出
今回の例ではスピングラス法を適用していますが、コミュニティ検出法は様々なものがあります。
以下のgroup_spingrass()
部分に他の関数を適用することで様々なコミュニティ検出法を実装可能です。
# スピングラス法
community = as.factor(group_spinglass())
# 媒介中心性
group_edge_betweenness()
# 固有ベクトル
group_leading_eigen()
# ランダムウォーク
group_walktrap()
# 情報中心性
group_walktrap()
# グリーディアルゴリズム
group_fast_greedy()
参考