Rで因子分析をやってみた

Ｒによるデータ分析の練習として因子分析をやってみました。

概要

NBAの各選手の成績データを使用する
PT（ポイント数）やBLK（ブロック数）などの成績データから選手のパフォーマンスを説明する共通因子がないか確認する

データ

以下の記事で作成したデータを使用します。
hadadada00.hatenablog.com

データ準備

分析対象の変数はtargetsで指定した5つとします。

library(formattable)
library(nbastats)
library(tidyverse)

### data
# stats data
stats <- seasons_stats

# filter 2017 records
stats <- seasons_stats %>%
  filter(Year == 2017)

# targets for analysis
targets <- c(
  "AST", # assists
  "STL", # steals
  "BLK", # blocks
  "TRB", # total rebounds
  "PTS"  # points
)

# extract variables and change data type
vars <- stats %>% 
  select(targets) %>% 
  sapply(as.integer)

scree plot

因子数を決めるためにスクリープロットを作成します。

### scree plot
# correlation matrix
cor <- cor(vars)

# eigen value
ev <- eigen(cor)$values

# scree plot
plot(ev, type = "b")

f:id:hadadada00:20190624002651p:plain — スクリープロット

3つ目以降の固有値はほとんど説明力をもっていないため、因子数は2とします。

因子分析の実施

factanal関数を使用して因子分析を実施します。以下の記事を参考に、因子の抽出法は最尤法、回転法はプロマックス法を使用します。また、因子得点は回帰法で算出します。
因子分析の因子抽出方法について | Sunny side up!
因子分析における因子軸の回転法について | Sunny side up!

### factor analysis
# factor extraction : maximum likelihood
# rotation : promax
# factor score : regression 
fa <- factanal(x = vars, factors = 2, rotation = "promax", scores = "regression")
fa

結果は以下です。

>fa
Call:
factanal(x = vars, factors = 2, scores = "regression", rotation = "promax")

Uniquenesses:
  AST   STL   BLK   TRB   PTS 
0.137 0.225 0.300 0.047 0.192 

Loadings:
    Factor1 Factor2
AST  1.036  -0.202 
STL  0.753   0.190 
BLK -0.151   0.919 
TRB          0.923 
PTS  0.683   0.304 

               Factor1 Factor2
SS loadings      2.138   1.866
Proportion Var   0.428   0.373
Cumulative Var   0.428   0.801

Factor Correlations:
        Factor1 Factor2
Factor1   1.000  -0.601
Factor2  -0.601   1.000

Test of the hypothesis that 2 factors are sufficient.
The chi square statistic is 0.21 on 1 degree of freedom.
The p-value is 0.643

Cumulative Varを見ると、2つの因子で全体の80%を説明できています。

Factor1とFacto2の内訳を見ていると、Factor1はAST（アシスト数）、STL（スティール数）、PTS（ポイント数）が高い因子負荷量を示しており、Factor2はBLK（ブロック数）とTRB（リバウンド数）が高い因子負荷量を示しています。

因子の名前としては以下のようにつけることにしました。

因子1 ： non-power-play ability　
因子2 ： power-play ability

因子得点の確認

因子得点のトップ10を求めてどのような選手がそれぞれの因子得点が高いか確認します。まずはFactor1（non-power-play ability）の上位10選手を見てみます。

### check factor scores
stats <- cbind(stats, fa$scores)

# top10 non-power-play ability
tbl1 <- stats %>%
  arrange(desc(Factor1)) %>%
  mutate(n = row_number()) %>%
  filter(n <= 10) %>%
  select(Player, Pos, Factor1)

formattable(tbl1)

f:id:hadadada00:20190624003827p:plain — Factor1 top10 players

ランキングを見ると、ポイントガードのスター選手がランクインしています。Facto1はポイントガード能力といってもいいかもしれません。
次にFactor2（power-play ability）の上位10選手を見てみます。

# top10 power-play ability
tbl2 <- stats %>%
  arrange(desc(Factor2)) %>%
  mutate(n = row_number()) %>%
  filter(n <= 10) %>%
  select(Player, Pos, Factor2)

formattable(tbl2)