2019年12月11日水曜日

Rによる主成分分析

主成分分析(Principal Component Analysis)は、多数のデータを俯瞰的に見るためのツールです。今回はフリーで使える統計解析ソフトのRを使って、主成分分析を行う方法を説明します。

データの準備

図1
ここでは例として有名な アヤメという植物のデータで説明します。まずExcelなどを用いて図1のようにデータをまとめます。図のとおり1行に1個体、各列にそれぞれの形態データという形でまとめてください。このデータには、sepal.length(がく片の長さ)、sepal.width(がく片の幅)、petal.length(花弁の長さ)、petal.width(花弁の幅)の4つの量的変数と、Species(種、setosa, versicolor, virginicaの3種類)という質的変数が含まれています。まとめ終わったら、「ファイル」>「名前を付けて保存」を選択し、csv形式で保存してください。 ちなみに、ファイル名は半角英数字と「_」だけで構成することを強くおすすめします。
図2

 csvとは、Comma-Separated Valuesの略で、Excelのそれぞれの列を「,」で区切った形となります。保存したファイルをメモ帳などのテキストエディタで開くと、図2のように表示されます。この形式はRにデータを読み込ませる上で重要なので覚えておいてください。

データの読み込み

データが準備できたら、Rを起動します。起動できたら以下のコマンドでデータを読み込ませます。


x <- read.table("C:/Users/shinka/Documents/iris.csv", header=T, sep=",", na.strings="NA", 
dec=".", strip.white=T)
それぞれのコマンドの意味は、以下の通りです。
  • x <- read.table(): ()内の表を、x(任意の文字)に代入します。
  • "C:/Users/shinka/Documents/iris.csv": 読み込むファイルの場所と名前を指定します。ダブルコーテーション(")で囲ってください。
  • header=T: ファイルの先頭行をラベルとして扱います。今回の例では先頭に変数名(sepal.lengthなど)を書いているので、T(=True)です。
  • sep=",": 列の区切りに「,」を使っていることを示します。ダブルコーテーションで囲ってください。
  • na.strings="NA": 空欄がある場合に「NA」という文字で埋めます。これがないと空欄に対して、システムが勝手な解釈で値を振ってしまうそうです。
    ダブルコーテーションで囲ってください。
  • dec=".": 小数点記号を指定します。ダブルコーテーションで囲ってください。
  • strip.white=T: データに含まれる空白文字を自動で削除します。
上のコマンドを打ち終えたらEnterを押してください。もしエラーがでた場合はコマンドのどこかが間違っています。確認してください。 読み込みができたら、先に指定した変数名(ここでは「x」)を入力してEnterを押してみましょう。読み込ませた表が表示されるはずです。さらに、以下のコマンドで表の要約を見ることができます(図3)。

summary(x)


図3

主成分分析の実行

メインの主成分分析に移ります。まず、今回主成分分析に使うprcompというコマンドは質的変数(今回のデータではSpeciesが該当する)に対しては使用できません。そこで、「Species」の列を削ります。具体的には、「Species」以外の列を新しい変数に代入します。

y <- x[,-5] 

上のコマンドはxに入っている表の5列目以外をyという新しい変数に代入することを意味しています。[]内のコンマを忘れないでください。もし3列目を削りたければ、「-3」とします。
 実際にprcompを使います。ここでは、主成分分析の結果をpcaという変数に代入します。summary(pca)で結果の要約が見れます(図4)。

pca <- prcomp(y, scale=T)
summary(pca)

図4

図の「Proportion of Variance」は寄与率と呼ばれるものです。簡単にいえば、「一つの主成分軸で、データ全体の何割を説明できるか」を示しています。図の例では、PC1という主成分軸で全体の70%以上が説明可能です。また、下の「Cumulative Proportion」は累積寄与率で、寄与率を第一主成分から順に足していったものです。PC2の寄与率は22.85%なので、第二主成分までで、データの95.81%が説明できます。
 最後に主成分分析の結果を図示します。図示にはggbiplotを用います。これはラボのパソコンによってはインストールする必要があります。インストールは以下のコマンドで可能ですが、わからなければ先生か誰かに聞いてください。

install.packages("devtools")
devtools::install_github("vqv/ggbiplot")

ggbiplotを以下のコマンドで呼び出してから実行します。ここで変数xに入れているSpeciesを利用することで、種ごとに色分けすることができて便利です(groups=x$Speciesで指定する。x$Speciesは変数xのSpeciesというラベルを意味する)。

library(ggbiplot)
ggbiplot(pca, obs.scale=1, var.scale=1, groups=x$Species)


図5
結果は右のようになります(図5)。アヤメの種ごとにはっきりと分かれていることがみてとれます。矢印は各形態の方向性を示しています。これを見るとpetal.lengthとpetal.widthはともに第一主成分軸に沿っており、両者は「ほぼ同じ意味をもつ形態である」といえます。また、sepal.lengthも少し逸れてはいるものの、petalとよく似た方向性であり、主にこの3つが第一主成分軸を構成していることがわかります。一方で、sepal.widthはこの3つとは全く別の方向へ伸びており、これにより第二主成分軸が形成されています。以上より、ざっくりとですが、アヤメの分類では、sepal.widthとその他の形態というまとまりでデータを見ることができそうであることがわかりました。


以上が主成分分析の一連の方法となります。主成分分析について勉強したい方は以下のブログがわかりやすいと思います。

https://logics-of-blue.com/principal-components-analysis/

0 件のコメント:

コメントを投稿

MaxEnt解析手順メモ

前回までの記事 で琵琶湖南湖の湖沼図から地形に関するラスターデータを作成しました。今回はこれを利用してMaxEnt解析を試みます。MaxEntはオープンソースのソフトで ここ からダウンロードできます。なお、ここではバージョン3.4.4を使用しています。サンプルとして、ブルーギル...