Rで主成分分析をやってみた

統計の勉強で主成分分析をやってみました。

分析概要

NBA選手のstats（各成績のデータ）を使用して、選手の総合力を測る指標を作成する
データは2016-2017年シーズンのものを使用
結果の大まかな検証として年収のデータと比較する

データ

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

年収のデータは以下のリンクから入手したものをcsvファイルにして読み込んで使用します。
https://hoopshype.com/salaries/players/2016-2017/

データの読み込み

まずは成績のデータと年収のデータを読み込みます。
成績のデータは項目数が多いのですが、今回はtargetsで指定した4変数のみ使用します。

library(formattable)
library(ggbiplot)
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
  "PTS"  # points
)

# change data type
temp1 <- stats %>% select(-targets)
temp2 <- stats %>% select(targets) %>% sapply(as.integer)
stats <- cbind(temp1, temp2)

# salary data
# https://hoopshype.com/salaries/players/2016-2017/
salaries <- read_csv("./salary_2016_2017.csv")

成績データと年収データのマージ

次に読み込んだ成績データと年収データをマージします。

stats <- stats %>% 
  mutate(flg_st = "1",
         name = tolower(Player)) 

salaries <- salaries %>% 
  mutate(flg_sa = "1",
         name = tolower(name))

stats_salaries <- left_join(stats, salaries, by = c("name"))

マージ結果を確認します。

stats_salaries %>% 
  group_by(flg_st, flg_sa) %>% 
  dplyr::summarise(n = n())

# A tibble: 2 x 3
# Groups:   flg_st [1]
  flg_st flg_sa     n
  <chr>  <chr>  <int>
1 1      NA        44
2 1      1        551

一部くっついていないデータがありますが、選手の母集団が違う可能性があるため、今回はこれでよしとしようと思います。
年収が取得できなかったレコードは除外します。

stats_salaries <- stats_salaries %>% 
  filter(!is.na(flg_sa))

少しデータを見てみようと思います。成績の相関係数と年収のTOP10です。

### check the data
#correlation matrix
cor <- data.frame(cor(stats_salaries[targets]))
formattable(cor)

# salary top5 bottom5
top5_bottom5 <- stats_salaries %>% 
  arrange(desc(salary)) %>% 
  dplyr::mutate(n = row_number(),
                salary = comma(salary)) %>% 
  filter(n <= 5 | n >= (nrow(stats_salaries) - 4)) %>% 
  select(n, Player, Pos, Age, Tm, salary)
  
formattable(top5_bottom5)

f:id:hadadada00:20190616164820p:plain — 相関係数行列

相関係数行列を見ると4つの変数間に相関は概ねありそうですが、BLK（ブロック数）だけ、やや相関が低いといえそうです。

f:id:hadadada00:20190616165241p:plain — 年収top5とbottom5

年収（USドル）のランキングを見ると、1位のLebron Jamesは3,000万ドル（約30億円）もらっているのに対して、一番低いDahntay Jonesは24,000ドル（約240万円）です。

主成分分析の実施

以下で主成分分析を実行します。

### principal component analysis
# pick up target variables
pca_vars <- stats_salaries[, targets]

# pca
pca_results <- prcomp(pca_vars, scale = T)

結果を確認します。

summary(pca_results)
pca_results

> summary(pca_results)
Importance of components:
                          PC1    PC2     PC3    PC4
Standard deviation     1.6865 0.8834 0.44539 0.4205
Proportion of Variance 0.7111 0.1951 0.04959 0.0442
Cumulative Proportion  0.7111 0.9062 0.95580 1.0000

> pca_results
Standard deviations (1, .., p=4):
[1] 1.6865263 0.8834317 0.4453887 0.4204833

Rotation (n x k) = (4 x 4):
          PC1         PC2         PC3        PC4
AST 0.5096261 -0.46438991  0.09538605 -0.7180005
STL 0.5480113 -0.13136200 -0.73544490  0.3762293
BLK 0.3675924  0.87550923 -0.04149933 -0.3108652
PTS 0.5521238 -0.02386692  0.66955217  0.4962757

累積寄与率（Cumulatrive Proportion）を見ると、第1主成分で71%ほど、第2主成分までで全体の91%の変動が説明できています。
結果をグラフでも確認してみます。

ggbiplot(pca_results, 
         obs.scale = 1, 
         var.scale = 1, 
         alpha=0.5)

f:id:hadadada00:20190616165621p:plain — 散布図（PC1とPC2）

グラフを見ると第1主成分は概ねPTS（総得点数）、STL（スティール数）、AST（アシスト数）と相関が高く、第2主成分はBLK（ブロック数）と相関が高いです。
第1主成分はいずれの成績（PTS, STL, AST, BLK）とも正の相関のため、選手の能力の総合指標として使えそうです。

主成分スコア（PC1）と年収の比較

作成した主成分スコアが選手の能力を表わす総合指標として使えそうか、主観的な評価と年収との比較で確認します。
まずは、PC1のトップ10選手を確認します。

### compare PC1 scores with Salaries
stats_salaries_scores <- cbind(stats_salaries, pca_results$x)

# make rank variables
stats_salaries_scores <- stats_salaries_scores  %>% 
  mutate(salary_rank = rank(desc(salary)),
         score1_rank = rank(desc(PC1)))

# top10 score1 players
stats_salaries_scores %>% 
  arrange(score1_rank) %>% 
  filter(score1_rank <= 10) %>% 
  select(Player, Pos, Age, Tm, score1_rank, salary_rank) %>% 
  formattable()

f:id:hadadada00:20190616171539p:plain — PCScore1 top10

ランクインしている選手を見ると、salary_rankでも上位のスター選手がいる一方、salary_rankがそれほど高くない選手もいます。
以下に挙げたように年報と選手の評価やパフォーマンスは必ずしも一致しないため注意が必要ですが、PC1は選手のパフォーマンスを表わす指標として的外れではなさそうです。

2016-2017年にブレイクした選手は、その評価が年俸に反映されていない。
若い選手は契約の制約でそもそも年俸が上がりづらい

主観的な感想としては、ここにあがっている選手は現役のスター選手ばかりのためPC1がパフォーマンス指標としてある程度使えそうだと思われます。

次にPC1とsalaryの相関を見ます。

# correlation matrix  
stats_salaries_scores %>% 
  select(score1_rank, salary_rank) %>% 
  cor()

> 
score1_rank   1.0000000   0.5784745
salary_rank   0.5784745   1.0000000

# scatter plot of score1_rank and salary_rank
stats_salaries_scores %>% 
  ggplot(aes(score1_rank, salary_rank)) +
  geom_point()

f:id:hadadada00:20190616172235p:plain — PC1Score rankとsalary_rank

相関係数は0.58程度のため、そこまで高くはありません。散布図をみてもかなりバラけています。

総合的に見ると、主観的にはPC1は選手の総合的な能力を表わす指標として使えそうです。しかし、年俸といった選手の現時点での金銭的な価値とはかならずしも相関が高くないといえそうです。

補足

今回は主成分分析を用いて、主成分スコアを選手の総合能力を表わす指標として使えないか考えて見ました。しかし、わざわざ主成分スコアを計算しなくとも、各成績の平均値ではだめなのか気になりました。
ためしに、標準化してスケールを調整した上での成績の平均値を計算し、主成分スコアと比較してみました。

### supplement
stats_salaries_scores %>% 
  mutate(AST_std = scale(AST),
         STL_std = scale(STL),
         BLK_std = scale(BLK),
         PTS_std = scale(PTS),
         std_sum = (AST_std + STL_std + BLK_std + PTS_std),
         std_sum_rank = rank(desc(std_sum))) %>% 
  ggplot(aes(score1_rank, std_sum_rank)) +
  geom_point()