読者です 読者をやめる 読者になる 読者になる

廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

打者の調子の波のモデル化(幾何分布編)

仮説検定

Albert (2008) "Streaky Hitting in Baseball" ではベータ二項分布を用いて野球選手の調子の波を評価した。

Albert (2008) 打者の調子の波のモデル化 - 廿TT

下記はカルロス・ギーエンという選手の2005年の打撃成績のデータで、ヒットを 1、アウトを 0 とコード化してある。

GuillenC <- c(0, 1, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0,
1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0,
1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0,
0, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 0, 0, 1, 0, 1,
0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 1, 0,
0, 0, 0, 1, 1, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 1, 0,
0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0,
0, 1, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
1, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0,
0, 1, 1, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 0,
0, 0, 0, 1, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 1)

調子の波が存在しない選手は、常にコンスタントな打率でヒットを出すから、その打撃成績を上記のように 0 か 1 かに符号化すると、それはベルヌーイ過程になる。

帰無仮説:
ギーエンの0-1のプロセスがベルヌーイ過程である。”

対立仮説:
ギーエンの0-1のプロセスがベルヌーイ過程でない。”

として、仮説検定ができないだろうか。

ベルヌーイ試行で 1 が出るまでの待ち時間の分布は幾何分布になる。

ギーエンの場合、1 が出るまでの待ち時間(打席数)の最大値は 19 だった。

f:id:abrahamcow:20160416104815p:plain

spacings <-diff(which(c(1,GuillenC)==1))-1 #1が出るまでの待ち打席数
plot(table(spacings))

幾何分布の最大値の分布を帰無分布として、ギーエンの 19 という数字が得られる確率が十分に小さければ、ギーエンには調子の波が存在すると言えそうである。

帰無仮説の下で最大値が 19 以上(18を超える)になる確率が p-値になる。

n 個の標本の最大値が x 以下である確率は、

P\{x_{1}\leq x \mbox{ and } x_{2}\leq x … \mbox{ and } x_{n}\leq x\} \\
= {F(x)} ^n

で与えられる。

幾何分布のパラメータの最尤推定量は、

p(x) = p (1-p)^x\\
\sum_{i=1}^{n} log p + x_i log(1-p)\\
\frac{n}{p} - \frac{\sum_{i=1}^{n}  x}{1-p} =0\\
n(1-p) - p \sum_{i=1}^{n}  x =0\\
n-p(n +  \sum_{i=1}^{n}  x) =0\\
\hat p = n/(n +  \sum_{i=1}^{n}  x)
=1/(1 + \bar X)

である。

phat <- 1/(1+mean(spacings))
n=length(spacings)
1-pgeom(18,phat)^n

p-値は 0.067 で、通例使われる有意水準5%では、帰無仮説は棄却されず、ギーエンに調子の波が存在すると考える必要性はあまりなさそうだ。

ベータ幾何分布

ふつうの幾何分布では、ベルヌーイ試行の成功確率 p は一定だが、ベータ幾何分布は各試行ごとに p が変化すると解釈できる。

ベータ幾何分布の確率関数は、

P(x) = B(a+1,x+b)/B(a,b)

で与えられる。ここで B はベータ関数

これを、

P(x) = B(K\eta+1,x+K(1-\eta))/B(K\eta,K(1-\eta))
(K>0, 0 < \eta < 1)

と改めてパラメタライズすることで、η は打率の中心を決めるパラメータ、K は打率の精度を決めるパラメータと解釈できる。

K が大きいほどばらつきが小さくなる。

最尤推定でベータ幾何分布のパラメータを推定し、ふつうの幾何分布と当てはまりを比較する。

library(VGAM)
dbetageom2 <- function(x,K,eta,log=FALSE){
  dbetageom(x,K*eta,K*(1-eta),log = log)
}
points(0:19,dgeom(0:19,phat)*n,type="b",pch=4,col="red")
LL <- function(par){sum(dbetageom2(spacings,par[1],par[2],log = TRUE))}
fitbetageom <- optim(c(1,0.1),LL,control = list(fnscale=-1))
points(0:19,dbetageom2(0:19,fitbetageom$par[1],fitbetageom$par[2])*n,type="b",pch=4,col="blue")
legend("topright",legend=c("beta","beta-geometric"),pch=c(4,4),col=c("red","blue"))

f:id:abrahamcow:20160416111412p:plain

パラメータの推定値は、それぞれ、

 \hat K =27.96, \hat p =0.34

だった。

AIC は、

  • 幾何分布:420.93
  • ベータ幾何分布:421.67

だった。

-2*sum(dgeom(spacings,phat,log = TRUE))+2
-2*fitbetageom$value+2*2

わざわざベータ幾何分布を使ってモデルを複雑にしなくても、ふつうの幾何分布で間に合う。

やはり、ギーエンに調子の波が存在すると考える必要性はあまりなさそうだ。

ただし Albert (2008) "Streaky Hitting in Baseball" ではベータ二項分布を用いてベイズファクターを求め、このエントリとは逆の結論を導いている。

https://www.stat.berkeley.edu/~aldous/157/Papers/albert_streaky.pdf