隠れ変数の探索 を、Rで進める時のレシピです。
因果推論 として 相関関係による仮説の探索 をやろうと思い、 Rによる変数の類似度の分析 をしたものの、注目している変数と他の変数の間には、何も相関が見つからないことがあります。
このページの方法は、そのような時に、あきらめずにさらにデータを探索するための方法になります。
主成分分析 と 独立成分分析 の2つの方法がありますが、同じような結果になることもありますし、 どちらかの方が、良いこともあります。
主成分分析を使う方法の方が、パラメタの設定がないですし、計算時間も短いです。 主成分分析では、隠れ変数が見つからない時に、独立成分分析を試すようにするのが良いようです。
下記のコードが注目している変数が「Y」ということにしています。 他の変数については、名前の付け方に決まりはありません。
library(dummies) # ライブラリを読み込み
library(ggplot2) # ライブラリを読み込み
library(tidyr) # ライブラリを読み込み
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
Data1 <- Data # 加工用のデータを作る
Y <- Data1$Y # Yの列を別に作っておく
Data1$Y <- NULL # Yの列を消す
Data2 <- dummy.data.frame(Data1)# 質的変数はダミー変換
pc <- prcomp(Data2, scale=TRUE)# 主成分分析をして主成分を抽出
Data4 <- cbind(Y, Data2, pc$x)# 加工したデータとYの列を合わせる
Data_long <- tidyr::gather(Data4, key="X", value = Xs, -Y) # 縦型に変換(Yの列以外を積み上げる)
ggplot(Data_long, aes(x=Xs,y=Y)) + geom_point() + facet_wrap(~X,scales="free")# 散布図を大量に描く
この結果を見ると、まず、Yと元のデータの他の変数(X01、X02、X03)には何も相関が見られないことがわかります。 次に、Yと主成分の関係を見ると、YとPC3には相関があることがわかります。 隠れ変数のようです。
この後ですが、隠れ変数が元の変数のどれに影響しているのかを調べます。
pc2 <- sweep(pc$rotation, MARGIN=2, pc$sdev, FUN="*") # 因子負荷量を計算
pc2# 結果の数値の出力
library(heatmaply) # ライブラリを読み込み
heatmaply(pc2) # クラスター分析をしてからヒートマップを作成
この結果を見ると、X03とPC3は、因子負荷量の絶対値が0に近いことから、隠れ変数のPC3は、X03には影響していないらしいことがわかりました。 PC3は、X01とX02に少し影響しているような因子のようです。
library(fastICA) # ライブラリを読み込み
library(dummies) # ライブラリを読み込み
library(ggplot2) # ライブラリを読み込み
library(tidyr) # ライブラリを読み込み
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
Data1 <- Data # 加工用のデータを作る
Y <- Data1$Y # Yの列を別に作っておく
Data1$Y <- NULL # Yの列を消す
Data2 <- dummy.data.frame(Data1)# 質的変数はダミー変換
summary(prcomp(Data2, scale=TRUE))# 主成分分析
Data3 <- fastICA(Data2, 3)$S# 独立成分分析。成分は3つ抽出の場合
Data4 <- cbind(Y, Data2, Data3)# 加工したデータとYの列を合わせる
Data_long <- tidyr::gather(Data4, key="X", value = Xs, -Y) # 縦型に変換(Yの列以外を積み上げる)
ggplot(Data_long, aes(x=Xs,y=Y)) + geom_point() + facet_wrap(~X,scales="free")# 散布図を大量に描く
Yと独立成分の関係を見ると、Yと1には相関があるようなので、隠れ変数があるようです。
この後ですが、隠れ変数が元の変数のどれに影響しているのかを調べたいのですが、 主成分分析の時のようなやり方を筆者がわからないです。 わかったらここに追記するつもりです。
主成分分析や、独立成分分析では、隠れ変数がわかったとしても、それが見えている変数にどのように影響しているかの考察が難しいです。 そのデータがどのようなデータなのかを考察する必要があります。
ここで 因子分析が役に立ちます。 Rによる主成分回帰分析 でも紹介している方法です。
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
factanal(Data,2)$loadings # 潜在変数を2個で仮定して、因子分析
すると、Yと共通の潜在変数が含まれているのは、X03とX04の2つであることがわかります。
psychというライブラリでも同じ結論を出せます。psych
だと、変数が少な過ぎるということが起きにくいです。
library(psych) # ライブラリの読み込み
library(GPArotation) # ライブラリの読み込み
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
fa(Data, nfactors = 2, fm = "ml", rotate = "varimax")$loadings # 潜在変数を2と仮定して、因子負荷量を出力
因子分析に、Yも入れるのがポイントです。 こうすることで、「変数が少なくて分析できない」ということが少し起きにくくなります。 また、Yが入っていると、Yを潜在変数のひとつの基準として、他の潜在変数を見つけるような計算になりますので、隠れ変数がうまくあぶり出せます。 Yを入れないと、主成分分析の時とあまり変わらない計算になり、考察が進めにくくなります。
library(dummies) # ライブラリを読み込み
library(ggplot2) # ライブラリを読み込み
library(tidyr) # ライブラリを読み込み
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
Data1 <- Data # 加工用のデータを作る
Y <- Data1$Y # Yの列を別に作っておく
Data1$Y <- NULL # Yの列を消す
Data2 <- dummy.data.frame(Data1)# 質的変数はダミー変換
pc <- prcomp(Data2, scale=TRUE)# 主成分分析をして主成分を抽出
Data4 <- cbind(Y, Data2, pc$x)# 加工したデータとYの列を合わせる
Data_long <- tidyr::gather(Data4, key="X", value = Xs, -Y) # 縦型に変換(Yの列以外を積み上げる)
ggplot(Data_long, aes(x=Y,y=Xs)) + geom_boxplot() + facet_wrap(~X,scales="free")# 箱ひげ図
次に、Yと主成分の関係を見ると、PC3ではYがかなり明確に分離されます。 隠れ変数のようです。
隠れ変数が元の変数のどれに影響しているのかの調べ方は、注目している変数が量的変数の場合で主成分分析を使った時と同じです。
library(fastICA) # ライブラリを読み込み
library(dummies) # ライブラリを読み込み
library(ggplot2) # ライブラリを読み込み
library(tidyr) # ライブラリを読み込み
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
Data1 <- Data # 加工用のデータを作る
Y <- Data1$Y # Yの列を別に作っておく
Data1$Y <- NULL # Yの列を消す
Data2 <- dummy.data.frame(Data1)# 質的変数はダミー変換
summary(prcomp(Data2, scale=TRUE))# 主成分分析
Data3 <- fastICA(Data2, 3)$S# 独立成分分析。成分は3つ抽出の場合
Data4 <- cbind(Y, Data2, Data3)# 加工したデータとYの列を合わせる
Data_long <- tidyr::gather(Data4, key="X", value = Xs, -Y) # 縦型に変換(Yの列以外を積み上げる)
ggplot(Data_long, aes(x=Y,y=Xs)) + geom_boxplot() + facet_wrap(~X,scales="free")# 箱ひげ図
Yと独立成分の関係を見ると、Yと1には相関があるようなので、隠れ変数があるようです。
注目している変数が質的変数の場合、量的変数の時のように、すべての変数を使って、因子分析をする方法ができないです。
このため、うまく行きません。