2018-07-21

R で主成分分析 (2)

PCA, Principal Component Analysis(主成分分析)の解析例として、前回は iris のデータセットを使いました [1]

この解析では、「がく」と「花弁」の長さと幅に4種類の変数について、「種」別に描画したボックスプロットを示しました。また、4つの主成分についても同様にボックスプロットで示しました。ボックスプロットを使ったのはどちらの場合も4つの変数しかなかったからです。

しかし、よくよく考えてみると4つしか変数がないからこそ、全ての組み合わせを散布図にして眺めることができるとも言えます。ggplot2 でも pairs でやるようにいっぺんにプロットできればいいのにと思って調べたら、同等な機能を実現できるパッケージがありました。参考サイト [2] で紹介されているサンプルを参考にしてプロットを作成してみました。

R スクリプト

本記事で使用したスクリプトを以下に示しました。

require(ggplot2)
require(dplyr)
require(tidyr)

gatherpairs <- function(data, ...,
                        xkey = '.xkey', xvalue = '.xvalue', ykey = '.ykey', yvalue = '.yvalue',
                        na.rm = FALSE, convert = FALSE, factor_key = FALSE) {
  vars <- quos(...)
  xkey <- enquo(xkey)
  xvalue <- enquo(xvalue)
  ykey <- enquo(ykey)
  yvalue <- enquo(yvalue)
  
  data %>% {
    cbind(gather(., key = !!xkey,
                 value = !!xvalue, !!!vars, na.rm = na.rm,
                 convert = convert, factor_key = factor_key), select(.,!!!vars))
  } %>% gather(., key = !!ykey,
               value = !!yvalue, !!!vars, na.rm = na.rm,
               convert = convert, factor_key = factor_key)
}

# -----------------------------------------------------------------------------
#  MAIN
# -----------------------------------------------------------------------------
iris %>% gatherpairs(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width) %>% {
  p <- ggplot(., aes(x = .xvalue, y = .yvalue, color = Species))
  p <- p + geom_point()
  p <- p + facet_wrap(.xkey ~ .ykey, ncol = length(unique(.$.ykey)),
                         scales = 'free', labeller = label_both)
  plot(p)
}

df <- as.data.frame(iris)
row.names(df) <- paste(df$Species, row.names(df), sep="_") 
df$Species <- NULL
head(df)
#
df_pca <- prcomp(df)
df_out <- as.data.frame(df_pca$x)
df_out$Species <- sapply(strsplit(as.character(row.names(df)), "_"), "[[", 1)
row.names(df_out) <-NULL
head(df_out)

df_out %>% gatherpairs(PC1, PC2, PC3, PC4) %>% {
  p <- ggplot(., aes(x = .xvalue, y = .yvalue, color = Species))
  p <- p + geom_point()
  p <- p + facet_wrap(.xkey ~ .ykey, ncol = length(unique(.$.ykey)),
                         scales = 'free', labeller = label_both)
  plot(p)
}

生データの散布図

setosa のその他の「種」において、「がく」と「花弁」の長さと幅にそれぞれ集団が別れているのが見て取れます。

生データのプロット

主成分の散布図

主要因 PC1 の軸から見ると setosa の集団とその他の集団とではっきり分離されていますが、他の主要因はそれほどでもありません。生データを使った「種」別における「がく」と「花弁」の長さと幅の傾向は、データを主成分分析すると、変換した主成分 PC1 だけで概ね説明できることがわかります。

主成分のプロット

参考サイト

  1. bitWalk's: R で主成分分析 (1) [2018-07-10]
  2. r - Create a matrix of scatterplots (pairs() equivalent) in ggplot2 - Stack Overflow

 

ブログランキング・にほんブログ村へ
にほんブログ村