【第4回】打者の「能力」を数値化してみる 【因子分析】
こんにちは。肉なし回鍋肉です。今回は「打者の「能力」を数値化してみる」の最終回です。前回まででタイトルにある能力を数値化する式が完成しました。今回は数値化した能力から何が分かるのかをご紹介します。
全体の目次
今回の目次
能力の成長と衰退
プロ野球選手は、厳しい練習の積み重ねで日々成長しています。一方で年齢を重ねるにつれて体が思うように動かなくなることもあるでしょう。選手が最も成績を残せる年齢を知ることはチーム作りにおいて最も関心のある話題でしょう。実際、日本プロ野球の大手データ分析サイト(?)である1.02さんのコラムでも能力の推移を取り扱っています。
選手のキャリアを予測する際、年齢による経年変化は欠かせない判断材料となる。今後NPBでもこうした研究を進めることが、より正確に選手の将来を予測することにつながっていくだろう。
このコラムでは総合指標(野手の場合は打撃走塁守備すべてを含む)のWARと攻撃WARの分析をしている。このコラムと同じような結果が出るのか。それとも違う結果になるのか。
まず、さきほども説明したように攻撃、守備の総合である野手WARは26歳から29歳にかけてピークが出現している。その構成要素である攻撃WAR（紫）は29歳前後にピークを迎え、そこから数値を下げていく。
分析対象
分析対象は、2年連続で規定打席に到達した打者です。ただ単純に規定打席到達者を年齢で分けて各年齢の選手ごとに平均を出すと偏りが生まれます。例えば「20歳で出場した選手の集団」と「21歳で出場した選手の集団」で分けるとそこに含まれる選手が異なるため、純粋に年齢による変化だけを見ることができません。
年齢による変化を観測するには、年齢以外の条件を揃える必要があります。なので、2年連続で規定打席に到達した選手のデータを使うことで「実際に20歳で出た選手の集団」と「21歳で実際に出場した選手の集団」は同一選手によって作られた集団になります。こうして求められた因子得点の平均を比較すれば、年齢による違いのみが反映されたデータとなるって仕組みです。ちなみに今回は全部Excelでやったのでコードはありません。
変換と計算
これは先ほどの条件で集計した年齢別因子得点のデータです。20歳の欄の因子得点は「20歳から21歳において連続出場した選手の集団」の20才時点での因子得点で、次年度の因子得点は同集団の21歳における因子得点です。人数の欄はその集団の人数を示してます。
各年齢における因子得点の平均は算出できましたが、因子得点のままだと負の数が混ざるため変化率を計算するのには都合が悪いです。そこで、因子得点に10を掛けて50を足す1次変換によって因子得点を偏差値に変換します。これなら変化率が計算できます。
このデータを使って年齢間の変化率を求めます。変化率は、「変化率=次年度の成績÷当該年度の成績」で求めます。
これが打者能力の変化率です。例えば25歳の欄のミート力は1.03となってますが、これは25歳から26歳になると平均的にはミート力が1.03倍になるという意味です。1.00を超えている場合は能力が成長し、下回る場合は能力が衰えていると言えます。
年齢変化による成績変化の倍率を求めることができたので、平均的な打者の各年齢における因子得点を求めることができます。「変換後の年齢別因子得点平均」より、最も集団の人数が多い27歳から28歳の集団における27歳の因子得点を基準とします。28歳以上の成績は当該年次の前年の成績に当該年次の変化率を掛け、26歳以下の成績は当該年次の翌年の成績を当該年次の変化率で割ることで打撃能力偏差値の経年変化を求めると以下のようになります。
これは偏差値なので、元々の因子得点に再変換すれば、年齢による因子得点の変化がわかります。
表だと分かりにくいので、図にしてみます。
先ほどのコラムは攻撃WARとして一括りにしていましたが、このグラフでは攻撃WARに関連する能力を4つに分解して推移を検証しています。
ミート力は27歳あたりまで上昇していき、29歳を境に大きく下降していきます。
走力は若いうちから年齢を重ねるごとに緩やかに劣化しています。
選球眼は、27歳あたりまで順調に上昇を続け、その後は小さいながらも緩やかに下降しています。
長打力は23歳まで急激に成長し、その後は能力を維持して29歳から劣化を始めています。
走力以外は29歳前後がピークとなっています。先述のコラムとおおよそ一致した結果であると言えるでしょう。このデータから、例えば長打に長けた選手であれば年齢を重ねても多少は通用するが、ミート力を武器とする選手は劣化のペースが早まる可能性が指摘できます。
各年齢の集団のサンプル数が十分に足りてないため確実性のあるデータではありませんが、ある程度の成績の変化は掴めたでしょう。また、この変化はあくまでも一般的な傾向であり、全ての選手に当てはまるものではないが、一般的な傾向でもそれなりに有意義だと考えています。実際に参考になるかは知らんけど。
今回のまとめ
- ミート力は29歳を境に衰退がスタートする。
- 走力はじわじわ衰えていく。
- 選球眼は27歳頃まで成長し、その後はほぼ横ばい。
- 長打力は29歳頃から緩やかに下降する。
総まとめと展望
因子分析によって、打者の持つ能力を「ミート力」、「走力」、「選球眼」、「長打力」の4つに分けることに成功しました。また、与えられた因子得点を目的変数として回帰分析を行い、簡易的に打者の能力を数値化することができるようになりました。さらに、年齢経過による能力の変化についても検証し、一般的な選手の成長や劣化を予測できるようになった。これはセイバーメトリクスの課題の一つである未来の成績予測の手助けになると思われます。
今回紹介した例以外にも、過去の選手の能力別ランキングや、各能力の相関関係などについて分析するのも面白いですね。また、パワプロの能力値も今回求めた能力値から計算できるかも？（特殊能力による調整等があるから難しいかも）
この研究が、少しでも皆様の野球ライフを盛り上げることができたら幸いです。最後までお読みくださりありがとうございました！
【第3回】打者の「能力」を数値化してみる 【因子分析】
こんにちは。肉なし回鍋肉です。今回は「打者の「能力」を数値化してみる」の3回目です。前回は繰り返し因子の抽出を行い、満足いく結果が出せました。今回は因子の命名と下位尺度の構成、因子得点の予測をします。いよいよ実際に数値化できるようになります！
全体の目次
今回の目次
因子の命名
まずは前回の結果をおさらいします。
第1因子は、打率、安打、BABIP、二塁打で因子負荷量が大きくなりました。これらの指標に共通するような能力が、第1因子の名前になります。
これらの項目はボールをフェアグラウンドに飛ばすことによって得られる出塁に関する項目であり、バットをボールに的確に当てる能力と考えられるので「ミート力」と命名します。
続いて第2因子です。先ほどの表は1, 4, 3, 2の順番なのでお間違いなく。ここでは盗塁、盗塁刺、三塁打で因子負荷量が大きいですね。これら３項目はいずれも走塁に関する項目であるので「走力」と命名しましょう。
盗塁は出塁後に発生するプレーであり、打席での打撃行為とは関係ないプレーです。しかし、足の早さによって内野ゴロが安打、安打が二塁打、二塁打が三塁打になることがあります。また、足の早い選手ならセーフティバントという選択肢も生まれます。このように、足の早さも打撃結果に影響を与えることから「走力」も打撃に影響を与える因子であると考えることにします。
ちなみに盗塁死が走力に関係していることになっていますが、足の遅い選手はそもそも盗塁をしないため盗塁死しないということでしょう。詳しくは後ほど。
第3因子はBB%、BB/K、出塁率で因子負荷量が大きくなりました。これらの項目はボール球を見送って四球を選ぶことで得られる出塁に関する項目であり、ストライクかボールかの判定を正しく行う能力と考えられるので「選球眼」と命名します。
最後に第4因子では、K%、ISO、HR%で因子負荷量が大きくなっています。ISOとHR%は長打に関する指標で、K%は三振に関する指標です。K%とISO、K%とHR%の相関係数はそれぞれ .522と .547であり、長打を打つにはある程度の三振のリスクを背負う必要があると言えるでしょう。ここではK%は三振を恐れずに長打を狙った結果喫した三振の率と捉えることにします。つまりK%も長打に関係する指標と考えて、この因子はボールをより遠くに飛ばす能力であるとして「長打力」と命名します。ただし、K%とISO、K%とHR%の相関は中程度であり、長打を打てるが三振の少ない打者、長打を打てないが三振の多い打者というのも一定数存在することは確認しておくべきでしょう。
下位尺度の構成
ここでは各因子を代表する観測変数を選出します。素直に因子負荷量が高い項目を選ぶと以下の通りになります。
- 第1因子（ミート力）：打率・安打・BABIP、二塁打
- 第2因子（走力）：盗塁・盗塁刺・三塁打
- 第3因子（選球眼）：BB%・BB/K・出塁率
- 第4因子（長打力）：K%・ISO・HR%
これらの項目に対する信頼性を、α係数というのを用いて検証します。これが1に近づくほど信頼性（内部一貫性）が高いと言えます。だいたい0.7を超えてると良いらしいですが、項目数によって値の出やすさが変わるため、一概には言えないです。
#下位尺度候補まとめ
l.f1 <- c("打率", "安打", "二塁打",
"BABIP")
l.f2 <- c("三塁打", "盗塁", "盗塁刺")
l.f3 <- c("出塁率", "BB.", "BB.K")
l.f4 <- c("K.", "HR.", "ISO")
d.f1 <- zdatas[l.f1]
d.f2 <- zdatas[l.f2]
d.f3 <- zdatas[l.f3]
d.f4 <- zdatas[l.f4]
#α係数
alpha(d.f1)
alpha(d.f2)
alpha(d.f3)
alpha(d.f4)
こちらが結果です。各因子が３項目または４項目と少ない項目から構成されているにも関わらず、全因子でα係数が .80を超えました。これは十分な内部一貫性を有していると言えるでしょう。また、下位尺度を構成するある項目を削除した時のα係数も算出しました。このα係数が通常のそれよりも高い時は、その項目がなければより内部一貫性が高くなったという意味になります。今回、長打力からK%を削除した時に内部一貫性が大きく上昇してますが、まあもともと0.87あるならいいでしょう。
因子得点の予測
因子得点とは
因子得点とは、各ケース（ここでは打者）が各因子に対してどれほどの重みを持っているのかを計算したものです。今回は回帰法で因子得点を求めてあります。
#因子得点 score <- f5$scores sdatas <- data.frame(datas, score) zsdatas <- data.frame(zdatas, score)
回帰分析
回帰分析とは、目的変数Yを説明変数Xで説明するような式を作る作業です。今回のケースでは、Yが因子得点、Xが観測変数になります。つまり、観測変数から因子得点を計算する式を作ろうってことです。
実際の因子得点は、全観測変数を使用して計算します。しかし、回帰分析で全観測変数を説明変数にすると式がめちゃくちゃになります。実用性皆無です。そこで、今回は説明変数とする観測変数は各下位尺度を構成する観測変数に限定して回帰分析をします。なお、実際の成績から計算できることを目的としているので、説明変数は標準化していないデータで計算します。
まずはミート力です。
#重回帰分析 meet <- lm(GLS1 ~ 打率 + 安打 + 二塁打 + BABIP, sdatas) summary(meet) vif(meet)
Call:
lm(formula = GLS1 ~ 打率 + 安打 + 二塁打 + BABIP, data = sdatas)
Residuals:
Min 1Q Median 3Q Max
-0.16935 -0.05575 -0.01354 0.04597 0.32813
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -9.1813662 0.0329127 -278.96 <2e-16 ***
打率 21.2402840 0.2225161 95.45 <2e-16 ***
安打 0.0088365 0.0002014 43.87 <2e-16 ***
二塁打 0.0410495 0.0005103 80.44 <2e-16 ***
BABIP 3.1694483 0.1621077 19.55 <2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.08058 on 751 degrees of freedom
Multiple R-squared: 0.9935, Adjusted R-squared: 0.9935
F-statistic: 2.864e+04 on 4 and 751 DF, p-value: < 2.2e-16出力結果そのままだと分かりにくいので、以降は綺麗にした結果のみの紹介にします。
表の「係数」は、各成績にかける係数です。「決定係数R2」は当てはまりの良さを表す値で、大きいほど良いです。「VIF」(分散拡大係数)は、多重共線性というのを判断する値で、10以下なら大丈夫です。「Pr(>|t|)」ってのは小さければ有意な説明変数ってことになります。
まあ結論としては、「ミート力=21.249×打率+0.009×安打+0.041×二塁打+3.171×"BABIP"-9.183」という式でミート力が求められるってことです。
次は走力です。Rのコードは先ほどのミート力と同じ感じなので省略。
結論の式は「走力=0.182×三塁打+0.030×盗塁+0.082×盗塁刺-1.024」となります。
盗塁死が多いほど走力が高くなることになっていますが、これは足の速い選手ほど盗塁試行回数が増えることが原因でしょう。しかし、この式では盗塁死しまくっても回数を積み重ねれば走力のある選手ということになってしまいます。盗塁成功率の高い選手が正当に評価されないことにな流ので、この点は改善の余地があると言えるでしょう。
次は選球眼です。
選球眼の結果は文句なしです。予測式は「選球眼=9.467×出塁率+12.924×"BB%"+1.349×"BB/K"-5.259」
最後に長打力なのですが・・・。因子の命名で、長打力とはボールをより遠くに飛ばす能力であること、長打を打てるが三振の少ない打者、長打を打てないが三振の多い打者というのも一定数存在することを説明しました。K%を説明変数として予測式を作成すると、三振が多ければ長打力が高くなる構造となっちゃいます。逆に、三振の少ないパワーヒッターが正当に評価されないことになります。α係数でも、K%を削除した方が内部一貫性は高くなることが示されており、このままK%を説明変数として回帰分析を行うのは不適当であると判断しました。よって、K%を除き、HR%、ISOを説明変数として重回帰分析を行った。
(盗塁刺は消すとα係数が結構落ちちゃうので消さなかったんだと思います)
今回も問題なし・・・ではありませんね。VIFがでかいので、どちらかの指標を削除して計算し直すことになります。HR%は本塁打限定の指標ですので、外野の頭を大きく超えた打球でもホームランにならなければ単打と同じ扱いになります。一方、ISOは二塁打や三塁打も評価するので長打力をより適切に表しているのはISOであると考えます。そこで、HR%を削除し、ISOのみを説明変数として今度は単回帰分析を行います。
これは問題ないですね。予測式は「長打力=12.778×"ISO"-1.881」となります。
すべての能力（因子）について回帰式をまとめると、
- ミート力=21.249×打率+0.009×安打+0.041×二塁打+3.171×"BABIP"-9.183
- 走力=0.182×三塁打+0.030×盗塁+0.082×盗塁刺-1.024
- 選球眼=9.467×出塁率+12.924×"BB%"+1.349×"BB/K"-5.259
- 長打力=12.778×"ISO"-1.881
となります。これで、プロ野球の打者の能力を数値化できるようになりました！
因子得点は平均0、標準偏差1の分布になっています。なので、求めた因子得点を10倍し、50を足せば偏差値に変換することができます。偏差値(というか正規分布)の性質を知れば、より楽しめるでしょう。
ただし、注意点として、これらの式は「規定打席に到達した打者」を分析して作った式であることがあります。この式を規定打席未満の選手に当てはめても正常な数値を出せるとは限りません。高校野球にも使えません。今後適応範囲を広げられると良いですね。
【第2回】打者の「能力」を数値化してみる 【因子分析】
こんにちは。肉なし回鍋肉です。今回は「打者の「能力」を数値化してみる」の2回目です。実際に因子分析ってやつをやります。ちなみに今回は基本統計量、ヒストグラム、相関行列の作り方はやりませんが、やっとくと良い（らしい）です。
全体の目次
今回の目次
因子分析の解説
前回お話ししたように、因子分析とは、多次元のデータを分析する手法（多変量解析）の一つです。因子分析は、観測された変数の相関関係を使い、それらの背後に因子と呼ばれる共通の原因があると仮定して探し出す手法です。最低限これだけ分かってれば良いので、結果だけ知りたい人は「分析2 因子分析」まで飛ばしてもらって構いません。
因子分析の詳細
因子分析は、(教育)心理学とかマーケティングとかでよく使われる手法です。5教科の成績から文系能力と理系能力を定めるみたいな感じです。
用語とかを簡単に解説します。筆者も詳しく理解してるわけではないので多少間違いがあるかもしれません。
これは前回も使用した因子分析の簡単なイメージです。この図の「成績」は、私たちが実際に目に見える形で得られたデータです。このようなデータを「観測変数」と言います。図の中で「能力」としているものは、実際には目に見えない存在です。このように、観測変数に影響を与えている、目に見えない潜在的な要因を「因子」と言います。
因子には2種類あります。一つはこの図の「能力」のように、複数の観測変数に影響する因子です。これは「共通因子」と呼ばれます。もう一つは、この図には書かれていないのですが1つの観測変数にのみ影響を与えている因子です。例としては、「打点」におけるランナーの存在みたいな感じです。ランナーがいても(だいたい)打率やホームラン数には影響しないですが、打点はランナーの存在で大きく変わりますよね。イメージはこんな感じです。
各因子が、それぞれの説明変数を同じ程度説明しているわけではなく、そこには強弱があります。その強弱の値を「因子負荷量」って言います。詳しい資料は記事の最後の「参考」に載せます。
因子分析の流れ
ここでは因子分析をどういう流れで行うか説明します。
1.観測変数の選出 : 前回紹介した29変数すべてが使えるわけではありません。分析すると不要な変数が出てきたりします。ということで、まずは変数選びからです。
2.共通因子数の決定 : そもそも因子が何個あるとして分析するかを決めます。少ない共通因子数で多くの観測変数を説明できると天晴れですので、多すぎず少なすぎずな因子数を決めることは大切なんです。適当に決めてはまずいのです。具体的な方法はいくつかありますが、今回はスクリー法によって決定します。スクリー法とは、固有値をプロットした図であるスクリープロットにおいて、固有値の減少がなだらかになる因子の直前の因子までを採用する方法です。固有値とは、各因子が観測変数をどれだけ説明しているかを表す数値です。例えばある因子の固有値が3.7ならば、その因子は観測変数3.7個分を説明していることになります。理論上、因子は最大で観測変数の数だけ存在して、全因子の固有値の和は観測変数の総数になります。（なんでなだらかになる直前の因子までを採用するとかは知りません）
3.因子の抽出 : 色々計算するパートです。計算方法にも色々あるのですが、今回は一般化最小二乗法を採用します。最尤法ってのもポピュラーなんですけど、これはデータの正規性を仮定した手法なんです。今回のデータは正規分布じゃないデータがいっぱいあるので今回は一般化最小二乗法です。
4.因子を回転させる : 突然回転とか言い出して横転。抽出したままの結果は、そのままだと意味不明です。そこで、因子軸を回転させても元の観測変数の関係は不変である性質を利用して、抽出された因子の解釈を簡単にするために因子軸を回転させます。ここの説明は詳しくできるほど理解してないので「結果を分かりやすくする」程度の認識です。
回転方法には大きく分けると直交回転と斜交回転があります。直交回転は、因子同士が無相関であると仮定した回転法で、斜交回転は因子間に相関があると仮定した回転法です。今回は因子間の相関を仮定して、斜交回転の一つであるプロマックス回転を採用します。
これを満足するまで繰り返します。
分析2 因子分析
ようやく、ここから本番です。
1回目
まずは因子数を決定します。
#因子分析1
#分析に使う観測変数を決める
statsNames1 <- c("打率", "試合",
"打席", "打数", "得点", "安打",
"二塁打", "三塁打", "本塁打",
"塁打", "打点", "盗塁", "盗塁刺",
"犠打", "犠飛", "四球", "故意四球",
"死球", "三振", "併殺打", "長打率",
"出塁率", "OPS", "BABIP", "BB.",
"K.", "HR.", "BB.K", "ISO")
zstats1 <- zdatas[statsNames1]
#スクリープロットの出力
png("screePlot1.png", width = 500, height = 500)
VSS.scree(zstats1)
dev.off()
第6因子からなだらかですね。ということで因子数を5にします。
f1 <- fa(zstats1, nfactors = 5, fm = "gls", rotate = "promax", scores = T) print(f1, sort = T, digit = 3)
結果は以下の通り。
> f1 <- fa(zstats1, nfactors = 5, fm = "gls", + rotate = "promax", scores = T) solve.default(S) でエラー: システムは数値的に特異です: 条件数の逆数 = 1.13542e-18 追加情報: 警告メッセージ: 1: cor.smooth(R) で: Matrix was not positive definite, smoothing was done 2: cor.smooth(R) で: Matrix was not positive definite, smoothing was done > print(f1, sort = T, digit = 3) print(f1, sort = T, digit = 3) でエラー: オブジェクト 'f1' がありません
エラー出ちゃった！原因としては、観測変数同士で相関係数が高い組み合わせがあったことです。計算の過程でおかしくなるらしい。ということで、相関係数の絶対値が0.9を超える組み合わせが以下の通り。
相関係数の高い組み合わせ：打席と打数、本塁打とHR%、本塁打とISO、四球とBB%、三振とK%、長打率とOPS、長打率とISO、HR%とISO
これらの組み合わせと解釈性を考慮して、試合、打席、打数、本塁打、四球、故意四球、三振、長打率、OPSの9変数を削除して2回目の因子分析を行います。
2回目
ネタバレ(?)ですが、ここからは最後まで因子数は4です。ということで、
- 実行コード
- スクリープロット
- 出力結果
- コメント
の流れで進めます。
#因子分析2
statsNames2 <- c("打率", "得点", "安打",
"二塁打", "三塁打",
"塁打", "打点", "盗塁", "盗塁刺",
"犠打", "犠飛",
"死球", "併殺打",
"出塁率", "BABIP", "BB.",
"K.", "HR.", "BB.K", "ISO")
zstats2 <- zdatas[statsNames2]
png("screePlot2.png", width = 500, height = 500)
VSS.scree(zstats2)
dev.off()
f2 <- fa(zstats2, nfactors = 4, fm = "gls",
rotate = "promax", scores = T)
print(f2, sort = T, digit = 3)
> f2 <- fa(zstats2, nfactors = 4, fm = "gls",
+ rotate = "promax", scores = T)
Loading required namespace: GPArotation
警告メッセージ:
fac(r = r, nfactors = nfactors, n.obs = n.obs, rotate = rotate, で:
A loading greater than abs(1) was detected. Examine the loadings carefully.
> print(f2, sort = T, digit = 3)
Factor Analysis using method = gls
Call: fa(r = zstats2, nfactors = 4, rotate = "promax", scores = T,
fm = "gls")
Warning: A Heywood case was detected.
Standardized loadings (pattern matrix) based upon correlation matrix
item GLS4 GLS1 GLS2 GLS3 h2 u2 com
安打 3 1.015 -0.111 0.054 -0.133 0.8817 0.1183 1.06
打率 1 0.929 -0.266 -0.031 0.127 0.8308 0.1692 1.21
BABIP 15 0.721 -0.099 0.342 0.036 0.6058 0.3942 1.48
二塁打 4 0.699 0.108 -0.222 -0.162 0.5735 0.4265 1.37
塁打 6 0.632 0.526 -0.135 -0.095 0.9516 0.0484 2.09
得点 2 0.545 0.440 0.368 0.085 0.7192 0.2808 2.78
K. 17 -0.450 1.005 0.211 -0.191 0.7786 0.2214 1.57
ISO 20 0.010 0.826 -0.235 0.043 0.9103 0.0897 1.17
HR. 18 -0.074 0.822 -0.300 0.013 0.9110 0.0890 1.28
打点 7 0.262 0.571 -0.409 -0.016 0.8574 0.1426 2.27
死球 12 0.047 0.210 0.063 0.024 0.0512 0.9488 1.32
盗塁 8 0.194 0.044 0.804 -0.022 0.6476 0.3524 1.12
盗塁刺 9 0.187 -0.037 0.751 -0.042 0.5996 0.4004 1.14
三塁打 5 0.202 0.058 0.709 -0.061 0.5072 0.4928 1.19
併殺打 13 0.084 0.007 -0.555 -0.077 0.3194 0.6806 1.09
犠飛 11 0.203 -0.012 -0.408 -0.014 0.2095 0.7905 1.47
犠打 10 -0.092 -0.302 0.379 -0.151 0.4159 0.5841 2.39
BB. 16 -0.324 0.377 0.052 0.927 0.9359 0.0641 1.60
BB.K 19 0.043 -0.374 -0.122 0.914 0.8315 0.1685 1.37
出塁率 14 0.423 0.080 0.070 0.698 0.9354 0.0646 1.70本当はもっと色々出てくるのですが、あまりにも長いので省略させていただきます。。簡単に解説すると、GLS1~4は第1因子から第4因子。各項目の数字は因子負荷量です。「h2」は共通性といって、その項目がどれだけ共通因子の影響を受けているかを表しています。小さいと因子とは関係ないってことです。「u2」は独自性です。1-共通性で求められ、どれだけ独自因子の影響が強いかです。「com」は複雑性と言います。これが1に近いほど、一つの因子からのみ影響を受けていると言えます。
h2の小さい観測変数は、各因子から大きな影響を受けていないとみなせます。これらを外すとさらに良くなるので、0.4以下を基準とし、死球、併殺打、犠飛を除いて３回目の因子分析を行います。
3回目
#因子分析3
statsNames3 <- c("打率", "得点", "安打",
"二塁打", "三塁打",
"塁打", "打点", "盗塁", "盗塁刺",
"犠打",
"出塁率", "BABIP", "BB.",
"K.", "HR.", "BB.K", "ISO")
zstats3 <- zdatas[statsNames3]
png("screePlot3.png", width = 500, height = 500)
VSS.scree(zstats3)
dev.off()
f3 <- fa(zstats3, nfactors = 4, fm = "gls",
rotate = "promax", scores = T)
print(f3, sort = T, digit = 3)
> print(f3, sort = T, digit = 3)
Factor Analysis using method = gls
Call: fa(r = zstats3, nfactors = 4, rotate = "promax", scores = T,
fm = "gls")
Standardized loadings (pattern matrix) based upon correlation matrix
item GLS4 GLS1 GLS3 GLS2 h2 u2 com
安打 3 0.968 -0.040 -0.134 0.119 0.869 0.1312 1.07
打率 1 0.951 -0.245 0.087 -0.071 0.864 0.1362 1.16
BABIP 12 0.724 -0.160 0.010 0.262 0.608 0.3916 1.37
二塁打 4 0.704 0.200 -0.174 -0.149 0.563 0.4371 1.39
塁打 6 0.619 0.616 -0.083 -0.009 0.951 0.0489 2.04
K. 14 -0.444 0.954 -0.156 0.268 0.743 0.2575 1.66
HR. 15 -0.030 0.880 0.023 -0.197 0.919 0.0806 1.10
ISO 17 0.044 0.879 0.055 -0.131 0.918 0.0817 1.06
打点 7 0.270 0.688 -0.002 -0.244 0.822 0.1784 1.58
犠打 10 -0.125 -0.376 -0.145 0.291 0.402 0.5979 2.48
BB. 13 -0.352 0.360 0.952 0.118 0.926 0.0742 1.61
BB.K 16 0.014 -0.353 0.920 -0.110 0.835 0.1650 1.32
出塁率 11 0.415 0.063 0.692 0.072 0.921 0.0788 1.68
盗塁 8 0.061 -0.049 0.035 0.805 0.685 0.3155 1.02
盗塁刺 9 0.064 -0.125 0.009 0.742 0.632 0.3679 1.07
三塁打 5 0.085 -0.016 -0.014 0.709 0.528 0.4720 1.03
得点 2 0.457 0.445 0.129 0.458 0.730 0.2697 3.15犠打が怪しいですけど、いずれの観測変数も共通性が0.4を超えましたね。しかし、これだとまだ各因子がどういう能力なのか分かりにくいですね。共通性の高い観測変数を除くことで分かりやすくなりますので、2以上を基準に犠打、塁打、得点を取り除いて再分析します。
4回目
#因子分析4
statsNames4 <- c("打率", "安打",
"二塁打", "三塁打", "打点", "盗塁", "盗塁刺",
"出塁率", "BABIP", "BB.",
"K.", "HR.", "BB.K", "ISO")
zstats4 <- zdatas[statsNames4]
png("screePlot4.png", width = 500, height = 500)
VSS.scree(zstats4)
dev.off()
f4 <- fa(zstats4, nfactors = 4, fm = "gls",
rotate = "promax", scores = T)
print(f4, sort = T, digit = 3)
> print(f4, sort = T, digit = 3)
Factor Analysis using method = gls
Call: fa(r = zstats4, nfactors = 4, rotate = "promax", scores = T,
fm = "gls")
Warning: A Heywood case was detected.
Standardized loadings (pattern matrix) based upon correlation matrix
item GLS2 GLS1 GLS3 GLS4 h2 u2 com
打率 1 0.927 -0.172 0.089 0.020 0.890 0.1097 1.09
安打 2 0.912 -0.072 -0.097 0.104 0.782 0.2183 1.06
BABIP 9 0.723 0.019 -0.002 0.427 0.720 0.2801 1.62
二塁打 3 0.718 0.154 -0.149 -0.167 0.561 0.4386 1.30
K. 11 -0.326 1.014 -0.139 0.323 0.841 0.1593 1.47
ISO 14 0.142 0.814 0.097 -0.165 0.914 0.0860 1.18
HR. 12 0.062 0.798 0.065 -0.241 0.907 0.0929 1.21
打点 5 0.330 0.594 0.041 -0.288 0.785 0.2146 2.08
BB. 10 -0.278 0.359 0.946 0.119 0.924 0.0757 1.51
BB.K 13 -0.006 -0.394 0.904 -0.143 0.863 0.1371 1.42
出塁率 8 0.447 0.119 0.689 0.134 0.946 0.0543 1.88
盗塁 6 0.057 -0.002 0.031 0.787 0.624 0.3764 1.01
盗塁刺 7 0.048 -0.072 0.006 0.743 0.605 0.3946 1.03
三塁打 4 0.083 0.051 -0.009 0.732 0.520 0.4800 1.04打点の複雑性が2を超えているので、取り除きます。
ちなみに、K%の第1因子の因子負荷量、1を超えていることに気がつきましたか？本来は1を超えてはいけないのですが、Heywood caseって問題が発生しているためこのような結果になってます(結果の最初の方に警告があります)。因子負荷量が1を超えたからって、K%を取り除けばいいわけではないので注意です。むしろ1に近いなら重要な観測変数と言えます。
5回目
#因子分析5
statsNames5 <- c("打率", "安打",
"二塁打", "三塁打", "盗塁", "盗塁刺",
"出塁率", "BABIP", "BB.",
"K.", "HR.", "BB.K", "ISO")
zstats5 <- zdatas[statsNames5]
png("screePlot5.png", width = 500, height = 500)
VSS.scree(zstats5)
dev.off()
f5 <- fa(zstats5, nfactors = 4, fm = "gls",
rotate = "promax", scores = T)
print(f5, sort = T, digit = 3)
> print(f5, sort = T, digit = 3)
Factor Analysis using method = gls
Call: fa(r = zstats5, nfactors = 4, rotate = "promax", scores = T,
fm = "gls")
Standardized loadings (pattern matrix) based upon correlation matrix
item GLS1 GLS4 GLS3 GLS2 h2 u2 com
打率 1 0.926 -0.172 0.083 -0.013 0.899 0.1008 1.09
安打 2 0.888 -0.092 -0.086 0.084 0.761 0.2391 1.06
BABIP 8 0.740 0.041 -0.016 0.384 0.734 0.2662 1.51
二塁打 3 0.729 0.133 -0.138 -0.207 0.552 0.4477 1.31
K. 10 -0.249 0.995 -0.129 0.259 0.858 0.1424 1.30
ISO 13 0.202 0.762 0.118 -0.227 0.896 0.1040 1.38
HR. 11 0.121 0.746 0.086 -0.299 0.885 0.1148 1.40
BB. 9 -0.242 0.328 0.948 0.104 0.924 0.0764 1.41
BB.K 12 -0.026 -0.407 0.898 -0.114 0.869 0.1308 1.43
出塁率 7 0.471 0.101 0.686 0.101 0.951 0.0491 1.88
盗塁 5 0.030 0.011 0.033 0.799 0.634 0.3661 1.01
盗塁刺 6 0.018 -0.056 0.006 0.757 0.611 0.3892 1.01
三塁打 4 0.057 0.059 -0.003 0.743 0.533 0.4671 1.02
GLS1 GLS4 GLS3 GLS2
SS loadings 3.132 2.480 2.258 2.236
Proportion Var 0.241 0.191 0.174 0.172
Cumulative Var 0.241 0.432 0.605 0.777
Proportion Explained 0.310 0.245 0.223 0.221
Cumulative Proportion 0.310 0.555 0.779 1.000
With factor correlations of
GLS1 GLS4 GLS3 GLS2
GLS1 1.000 0.127 0.344 0.086
GLS4 0.127 1.000 0.117 -0.391
GLS3 0.344 0.117 1.000 -0.087
GLS2 0.086 -0.391 -0.087 1.000
Mean item complexity = 1.3
Test of the hypothesis that 4 factors are sufficient.
The degrees of freedom for the null model are 78 and the objective function was 20.011 with Chi Square of 15004.7
The degrees of freedom for the model are 32 and the objective function was 9.761
The root mean square of the residuals (RMSR) is 0.066
The df corrected root mean square of the residuals is 0.104
The harmonic number of observations is 756 with the empirical chi square 518.525 with prob < 3.31e-89
The total number of observations was 756 with Likelihood Chi Square = 7293.031 with prob < 0
Tucker Lewis Index of factoring reliability = -0.19
RMSEA index = 0.5511 and the 90 % confidence intervals are 0.5377 0.5588
BIC = 7080.934
Fit based upon off diagonal values = 0.966
Measures of factor score adequacy
GLS1 GLS4
Correlation of (regression) scores with factors 0.996 0.999
Multiple R square of scores with factors 0.992 0.997
Minimum correlation of possible factor scores 0.983 0.995
GLS3 GLS2
Correlation of (regression) scores with factors 0.989 0.967
Multiple R square of scores with factors 0.978 0.934
Minimum correlation of possible factor scores 0.957 0.869今回の結果は、特にいちゃもんつける点がありません。ということで、これで分析パートは終了です。最後なので後に続くいろんな出力結果も載せておきます。
分かりやすくしたのはこちら↓
ここからは、各因子がどのような能力なのか分析するのですが・・・
すでに結構ボリュームあるので今回はここまでにします。
今回のまとめ
- 因子数決める→計算→回転
- 共通性や複雑性に気をつけながら、納得いくまで繰り返す
- Heywood caseに注意
今回はここまでです。次回は因子の命名と下位尺度の構成、因子得点の予測をします。
参考
当時参考にしてたページを貼っておきます。高校の図書室には参考になる本が無かったので（そりゃあそう）、ネットでいろんな大学のpdfとか読んで0から勉強しました。
授業関係
Rの練習帳
この二つは特に重宝してました。気が狂うほど読みました。
http://www.team1mile.com/asarin/hus/ip-hus/materials/factor.pdf
http://www.u.tsukuba.ac.jp/~hirai.akiyo.ft/forstudents/eigokyouikugaku7.files/M4_9(N.K.).pdf
http://cogpsy.educ.kyoto-u.ac.jp/personal/Kusumi/datasem13/masuda.pdf
https://koumurayama.com/koujapanese/factor.htm
http://psy.isc.chubu.ac.jp/~oshiolab/teaching_folder/datakaiseki_folder/08_folder/da08_01.html
この他にも色々あります。興味あったら探してみてください。
【第1回】打者の「能力」を数値化してみる 【因子分析】
こんにちは。夏休みを迎えた肉なし回鍋肉です。今回から4回(予定)に分けて、高校の卒業研究でやった内容を書きます。流石に時効でしょってことで・・・。内容は「打者の「能力」を数値化してみる 」です。統計学的に色々やってるのですが、統計をよく知らない人でも分かるように進めていくので、どうぞよろしくお願いします。
今回は第1回ということで、どういう研究なのか、どういう手法でやるのか等の導入・準備回です。実際の分析は統計ソフトのRを使用しましたので、その具体的なコードも併せて紹介します。なんかの参考になれば幸いです。
全体の目次
今回の目次
序論
野球におけるバッターには、いろんな特徴を持った選手がいます。例えば、イチローみたいなヒットを量産するバッターや、一撃必殺のホームランバッターみたいな感じです。このような選手の特徴は、打撃成績を見れば分かります。先ほどの例で言えば打率や本塁打ですね。
ところで、打撃成績について、次のように考えることができませんか？「打撃成績とは、打者がもともと持っている能力から得られた結果である」と。言い換えると、それぞれの打者の持っている「能力」が、打撃成績を決定しているのだと。図で表すとこんな感じです。
パワプロって各選手にパワーが86、ミートが14とか能力が割り振られて、オーペナで回すとその能力相応の成績を残しますよね。イメージとしてはまさにこの関係です。
この図に出てくる「能力」と「成績」のうち、我々が実際に観測できるのは「成績」だけです。与えられた成績から能力を数値化するってのが、今回の目標になります。
方法
分析対象
ここでは能力を導き出す「成績」を選びます。まずは対象の選手です。データは多いに越したことはないので、日本プロ野球の2005~2017年シーズンで規定打席に到達した選手756名とします。出典はNPBのホームページからです。
対象の選手を決めたら、次は対象となる成績です。NPBで発表されている成績は以下の22成績です。
打率、試合、打席、打数、得点、安打、二塁打、三塁打、本塁打、塁打、打点、盗塁、盗塁刺、犠打、犠飛、四球、故意四、死球、三振、併殺打、長打率、出塁率
ただ、これだけじゃあ物足りないですね。いわゆる「セイバーメトリクス」の指標は含まれていないのです。ということで、これらの22成績から計算できる以下の7成績を加えます。
- BB%=(四球+故意四球)÷打席
- K%=三振÷打席
- HR%=本塁打÷打席
- BB/K=(四球+故意四球)÷三振
- ISO=長打率-打率
今考えると四球に故意四球が含まれてるかもしれないなと思うのですが、これは考えないことにします。各指標の意味は省略します。
以上756名の29成績について、分析していきます。
因子分析とは
因子分析とは、多次元のデータを分析する手法（多変量解析）の一つです。因子分析は、観測された変数の相関関係を使い、それらの背後に因子と呼ばれる共通の原因があると仮定して探し出す手法です。ちょうど、序論で話したイメージと似た感じですよね。詳細は次回で。
統計ソフトR
今回の分析には、Rという統計ソフトを使いました。インストールや設定などの解説は行いませんので、必要であればご自身でお願いします。
分析1 データを取り込む・パッケージの読み込み
ここからはいよいよ実際の分析です！と言ってもまずはデータを読み込む作業なのですが。
自分の場合、実際に使用するデータはExcelに全部コピペして、csv形式で保存しました。ガチ勢はプログラムでデータすら回収してしまうらしいのですが、当時の私にはその発想すらなかったです。次のコードで、csvファイルを読み込みます。
#データ読み込み
datas <- read.csv("motodata.csv", header = T,
fileEncoding = "CP932")
zdatas <- read.csv("hyojyunkadata.csv", header = T,
fileEncoding = "CP932")read.csvで、作業ディレクトリ内のcsvファイルを読み込むことができます。header=Tはcsvファイルの1行目が列名であるという意味です。自分の環境(Mac)では、「fileEncoding = "CP932"」ってのが無いとうまく読み込めなかったです。これは人によっては不要かもしれません。datasには純粋な成績、zdatasには標準化された成績が入ってます。
これだけでデータは読み込めました。ただ、Rの都合上、列名の「K%」や「BB/K」は「K.」「BB.K」として読み込まれてしまいます。グラフなどを作るときにこれだと困るので、正しい成績名のデータフレーム(さっきのdatasとかのこと)を作ります。ついで、正しく無い成績名のデータフレームも作っちゃいます（あとで使うので）。
#指標名データフレーム
allStatsNames <- c("打率", "試合", "打席", "打数", "得点", "安打", "二塁打", "三塁打", "本塁打", "塁打", "打点", "盗塁", "盗塁刺", "犠打", "犠飛", "四球", "故意四球", "死球", "三振", "併殺打", "長打率", "出塁率", "OPS", "BABIP", "BB.", "K.", "HR.", "BB.K", "ISO")
allNames <- c("打率", "試合", "打席", "打数", "得点", "安打", "二塁打", "三塁打", "本塁打", "塁打", "打点", "盗塁", "盗塁刺", "犠打", "犠飛", "四球", "故意四球", "死球", "三振", "併殺打", "長打率", "出塁率", "OPS", "BABIP", "BB%", "K%", "HR%", "BB/K", "ISO")すごい見にくいけど勘弁してください。
あとは、サンプルサイズの値をnumbersに記録しておきます。
#サンプル数 numbers <- 756
こうすると、もし後からデータ数が変化した時も簡単に対応できます。
最後に、今回使用するパッケージを読み込みます。
#パッケージ読み込み library(psych) library(MASS) library(parallel) library(car)
以上で準備は完了となります。
今回のまとめ
- 打者の成績は打者の能力によって決定されていると考える！
- 「因子分析」で能力を数値化してみる！
次回は実際に因子分析を行います。続きもぜひご覧ください。それではさようなら。
NEW GAME!メドレーの簡単な解説
初めまして。肉なし回鍋肉です。
今回は、3月中旬にYouTube、ニコニコ動画に投稿した動画『「NEW GAME!」全曲メドレー【作業用】』が最近になって再生される機会が増えてきたので、簡単に曲順の意図など書き残しておきます。書いてたら意外と長くなっちゃったけど、簡単にです。
まずは聞いてくれ
そもそも聞いたことがないとこれ以降の文章を読んでも伝わらないので、どうぞ。
1時間の長いメドレーなので、ぜひ作業しながらでも聞いてください。自分はなかなか眠れない時に聞いたりもします。
曲を大して聞いてない奴が作ってると信用できない(?)ので、スクショに入る範囲で再生数載せておきます。最少は「サヨナラ飛行機雲」(212回)でした。全曲シャッフル再生してるので、最新の曲は少なめ。
作り方
知の共有目的で、一応作り方もざっと紹介します。
- 曲を集め、BPMを調べる
BPMとはBeats Per Minute、ざっくり言えばテンポです。テンポが近い曲を並べると綺麗につながりやすいらしいです。自分はMixMeister BPM Analyzerってのを使いました。 - 順番を決める
後述します - 録音する
Garage Band(以下GB)というMacに標準搭載されているアプリを使用しました。GBに曲をD&Dして、左の音量スライダーみたいなのをマウスで操作しながら、うまくつながるように並べます。（写真参照）
全部準備できたら、macの機能で録音します。
- 動画に編集する
編集はこちらもmacに標準搭載されてるiMovieで編集しました。適当な画像を準備して、字幕入れて、完成です。
曲順の意図
次に、これらをBPM順で並べ替えてみます。
基本はこの順番のままです。赤字になってるfourfolium楽曲は比較的知名度があると考えられるのですが、それが序盤に固まっています。また、「YOU-TOPIA!」と「Now Loading!!!!」のBPM差が、うまく繋ぎにくい程度には離れてしまっています。
うまく説明できないのですが、BPMは2倍差なら一周回って近いらしいので、この性質を生かして一部のfourfolium楽曲は終盤に回す方針にしました。
こんな感じで、前半14曲は素直にBPMに従って並べました。「春色Runway」から「Piece of Peace」は、fourfoliumの4人が並ぶように調整しました。
15曲目以降は、BPMよりも曲の雰囲気みたいなのを重視して順番を決めました。
「Blindly Sky」と「ふたり空」をくっつけたい、この目的のために順番をごちゃ混ぜにしました。「カラフル☆フュージョン」から「Rainbow Days!!」は、BPMが1.5倍の関係(今後多用する)なので、まあ許容範囲かなってことで。ここからBPMを下げていき「Blindly Sky」まで持っていきます。「Rainbow Days!!」→「イッチョクセン!」の青葉繋ぎが熱い。
「ふたり空」で曲調がしんみり系(?)になるので、ここに同系統の曲を固めました。これ以降コウの曲はないので、ここで「サヨナラ飛行機雲」(本質曲)を投入。「Green Leaves」→「CONTINUE...?」の青葉挿入歌繋ぎが熱い。マジで。さらに、そこからの「ユメイロコンパス」で______________
「ユメイロコンパス」以降はラストのfourfoliumラッシュに向けてテンションを上げて行きます。「そうじろうのうた」→「Dreaming Stars」のひふみ(+青葉)繋ぎから、「あめのちほしぞら」を経由し、「BUG! BUG! SURVIVAL!」、「YOU-TOPIA!」へ。この2曲でねね、うみこ(境目で一瞬)、しずくからfourfoliumへ。
本動画の目的は、「曲を知ってもらい買ってもらう」ことなので、基本的に1番で終わりに統一してたのですが、「STEP by STEP UP↑↑↑↑」と「Now Loading!!!!」はクライマックスってことで長めに流しました。(編集しやすかったのもある。)個人的には「ユメイロコンパス」もラスサビ流したかったのですが、技術不足で断念。
画像の意図
動画を見てる人が飽きないように、各曲に合うそれっぽい画像を動画に入れてあります。基本的に1期の曲は1期から、2期の曲は2期から選んでいます。歌ってる人が映るようにしてます（当然）が、「Blindly Sky」と「ふたり空」だけはあえて逆にしてます。ここはこだわりポイントです。暇があれば、どの画像が何話のものか探して見てください。にゅげ力のトレーニングになりますよ(???)。
終わりに
長々と書いてしまいました(照)。途中でも書いたように、本動画を作った目的は「曲を知ってもらい買ってもらう」ことです。また、先日のNEW GAME!展の開催を受け、もっと多くの人に「NEW GAME!」を知ってほしいという願いもあり、この動画を作りました。気に入った曲があれば是非買って、フルで聞いてほしいです。
最後まで読んでくださった方、ありがとうございました。ブログは不慣れなので、何かアドバイスありましたら是非コメントなどにお願いします。それでは。
