さあ、準備はととのった

MikuHatsune2011-02-20

バグ取りはたぶんうまくできたと思う。
 
\begin{matrix}&Simulation&Calculation\\(1)&hitting\hspace{3}probability&hitting\hspace{3}probability\\(2)&\downarrow&\downarrow\\(3)&one\hspace{3}inning\hspace{3}scores&one\hspace{3}inning\hspace{3}scores\\(4)&\downarrow&\downarrow\\(5)&N\hspace{3}innings\hspace{3}scores&N\hspace{3}innings\hspace{3}scores\\(6)&\downarrow&\downarrow\\(7)&2\hspace{3}teams\hspace{3}battle&2\hspace{3}teams\hspace{3}battle\\(8)&\downarrow&\downarrow\\(9)&pennant&pennant\end{matrix}
 
少ないアウト、少ないイニングで確かめたら、(1)〜(5)はうまくできている。
(5)以降はうまくできているので、(おそらく)大丈夫だろうと。
ただ、確率で(1)〜(5)を計算すると時間がかかるが、(5)以降は確率計算のほうが圧倒的に早いので、データ作りには(1)〜(5)をシミュレーション、(5)以降を確率計算に任せようと思う。
ということで今までやったやつをきちんと理由をつけてやり直してみる。
上のリンクから、#必要な関数、を読み込ませておいてスタート。
チームはアウト0.8、ヒット各0.05で、実際の野球に準ずる。

game<- 10000            # シミュレーション試合数。
battle.cmb<- combinations(people,2)
wdldata<- matrix(0,nr=nrow(battle.cmb),nc=3);dimnames(wdldata)<- list(1:nrow(wdldata),c("win","draw","lose"))
battledata<- matrix(0,nr=people,nc=3);dimnames(battledata)<- list(1:people,c("win","draw","lose"))
winprobdata<- matrix(0,nr=people,nc=game);dimnames(winprobdata)<- list(1:people,1:game)
names(dimnames(winprobdata))<- c("team.No","Game cycles")
# game回の得点分布を出しておいて、9チームで比較しあう。
WDL<- matrix(0,people,game)
for(i in 1:people){
  for(g in 1:game){
    WDL[i,g]<- one.game.score(team.array[,,i],Nout,Nbase,people,Ninning)
  }
}
# 密度に変換する。
Dens<- function(vec){
  a<- rep(0,max(vec)+1)
    for(w in 0:max(vec)){
      a[w+1]<- length(which(vec==w))
    }
  a<- a/length(vec)
  return(a)
}
# 得点分布
Sd<- apply(WDL,1,Dens)
# 2チーム間での試合。行2に対する行1のチームの勝率が出る。
Bcmb<- cbind(battle.cmb,matrix(0,nrow(battle.cmb),3))
  for(k in 1:nrow(battle.cmb)){
    H<- t(t(Sd[[battle.cmb[k,1]]]))%*%Sd[[battle.cmb[k,2]]]
    Bcmb[k,3]<- sum(H[lower.tri(H)])
    Bcmb[k,4]<- sum(diag(H))
    Bcmb[k,5]<- sum(H[upper.tri(H)])
  }
# 9チーム間でのペナント様式に変換する。
res<- Bcmb
c.data<- matrix(0,nr=people,nc=3);dimnames(c.data)<- list(1:people,c("win","draw","lose"))
  for(nr in 3:5){
    for(t in 1:people){
      c.data[t,nr-2]<- sum(res[which(res[,1]==t),nr])
    }
  }
  for(nr in 5:3){
    for(t in 1:people){
      c.data[t,6-nr]<- c.data[t,6-nr]+sum(res[which(res[,2]==t),nr])
    }
  }
# 勝分負の割合
c.data/(people-1)
# 各チームの勝率。ペナント結果。
finalp<- apply(c.data,1,winprob)
# 9*9マスで勝敗表を作成する。
Q<- diag(0,people)
dimnames(Bcmb)<- list(1:nrow(Bcmb),c(1,2,"win","draw","lose"))
QB<- cbind(battle.cmb,apply(Bcmb[,3:5],1,winprob))
  for(x in 1:nrow(QB)){
    Q[QB[x,1],QB[x,2]]<- QB[x,3]
    Q[QB[x,2],QB[x,1]]<- 1-QB[x,3]
  }
persp(Q,zlim=c(0,1),phi=20,theta=50,col=5,
   xlab="Team 1",ylab="Team 2",zlab="Win Prob")

ペナント結果は

finalp
        1         2         3         4         5         6         7         8         9 
0.5579265 0.5443338 0.5317352 0.5231784 0.4937143 0.4812046 0.4791171 0.4593829 0.4293001 

となった。見事に昇順。
ただ、これにχ自乗検定をしてもp値が1…だ…と…?
2チーム間で各々してみたら

chiscore<- rep(0,nrow(battle.cmb))
for(ch in 1:nrow(battle.cmb)){
  chiscore[ch]<- chisq.test(finalp[battle.cmb[ch,]])[3]
}
cbind(battle.cmb,unlist(chiscore))

p値がほとんど1に近くて、しかもエラーでまくり。
(´・ω・`)