平均への回帰をわかりやすく例える
なんか平均への回帰ってことばはよく知られているものの，外国語教育研究界隈では，ときに難しいとされていて，そして研究実践においてもものすごく大事なことなのにあんまり重要視されていない気がする。まあ確かに，頑張って教科書読んでも，「相関係数が1ではないということに他ならない」とか書いてあるし。
さて，平均への回帰というのは，たとえば2回テストするとして，1回目のテストで成績がよくないひとは，2回目のテストではよい成績をだすっていうだけの話なんだけど，これってピンとこないかな。
別にこれって正規分布に限ったことじゃないので，サイコロの話をする。
サイコロを30人が2回投げる。
- 一回目投げて目が1だったひとにとって，次の目が今よりも大きい確率は5/6。
- 一回目投げて目が2だったひとにとって，次の目が今よりも大きい確率は4/6。
- 一回目投げて目が3だったひとにとって，次の目が今よりも大きい確率は3/6。
- 一回目投げて目が4だったひとにとって，次の目が今よりも大きい確率は2/6。
- 一回目投げて目が5だったひとにとって，次の目が今よりも大きい確率は1/6。
- 一回目投げて目が6だったひとにとって，次の目が今よりも大きい確率は0。
このとき，2回目の目-1回目の目の期待値は，サイコロ投げの期待値が常に3.5なので，3.5 - 一回目の目。なので，グラフにすると，
これで，1回目1の人だけを集めて，2回目と1回目の目の差の平均値をもとめたら2.5なわけね。これが0ではないかと検定したりしたら，標本サイズによって当然有意になるでしょう。でもこれは，新しいイカサマサイコロに持ち替えたわけではない。
一回目の期待値，というかサイコロの期待値はもちろん，
うん。これだけの話。これ，別の分布でも同じ。
基本的にこれは完全に独立なサイコロを2回投げるからシンプルなんだけど，基本的には相関係数が1じゃない限りこうなる。
外国語教育研究が例外ではなくて，むしろ同じようなことが色々と起きる。
事前―事後でTOEICのテストの点を取ったとする。このとき実は差得点の平均値が0とする。
事前の成績を並べて下から33%のひとだけに検定したら，そりゃ有意になるだろう。このサイコロと全く同じことが起こる。
事前の成績で並べなくても，事前または事後の成績と相関が高い別の変数で並べても似たようなことが起きるのは当たりまえ。
ある読解テストを事前事後に実施して，読解テストと相関が高いTOEICスコアをカテゴリカルにして層別に検定してもそう。
これで成績下位層の子に特に効果がある指導法！というハッピーな指導法が世の中にたくさんになったらハッピーなのかな。
その度毎に，研究仮説に「熟達度の効果がどうなるかわからないから調べる」とか後付けで足して（HARKingして）そういう研究がいっぱいになればハッピーなのかな。
Rで混同行列から感度・特異度などの指標
例えば，何かの指標なりモデルを用いて，
二値的に何かの状態がある（陽性）かない（陰性）かを検査する方法を作るとする。
そのパフォーマンスは以下のような真の状態と検査結果の混同行列を見るとわかる。
この場合の混同行列は，
| 真の状態が陽性 | 真の状態が陰性 | |
| 検査結果が陽性 | 真陽性 | 偽陽性 |
| 検査結果が陰性 | 偽陰性 | 真陰性 |
っていうかんじで，ここから，
- 適中率/Accuracy（全体の内，真陽性と真陰性が占める割合）
- エラー率/Error Rate（1-的中率）
- 感度/Sensitivity （真の状態が陽性であるもののうち，陽性と判断できた割合）
- 特異度/Specificity （真の状態が陰性であるもののうち，陰性と判断できた割合）
- 陽性適中率 （陽性と判断したもののうち，真の状態が陽性である割合）
- 陰性適中率 （陰性と判断したもののうち，真の状態が陰性である割合）
- 陽性尤度比 （真の状態が陽性であるひとが陰性であるひとよりも何倍検査結果が陽性になるか）
- 陰性尤度比 （真の状態が陰性であるひとが陽性であるひとよりも何倍検査結果が陰性になるか）
といった指標が与えられる。
当たり前だけど，完全によい検査ってのは実現しにくくい。
感度をあげようとすると特異度が下がるような場合がほとんど。
これらの指標のバランスを実務的状況から判断して適切に保たなければならない。
私には実務的にこれをやらなければならない案件がある。もちろんRにはこれらの指標を与えたり，さらにROC曲線を描いてくれたりするパッケージがいっぱいあるんだけど，これらを調べるよりも自分でコード書くほうが絶対楽なんで，自分で書いてみた。混同行列を与えたら以上のものを全部返す。
my.indices<-function(CM){ Rate<-CM/sum(CM) TP<-CM[1,1] TPr<-Rate[1,1] FP<-CM[1,2] FPr<-Rate[1,2] FN<-CM[2,1] FNr<-Rate[1,1] TN<-CM[2,2] TNr<-Rate[2,2] Acc<-(TP+TN)/sum(CM) Err<-1-Acc Sens<-TP/(TP+FN) Spec<-TN/(FP+TN) PPV<-TP/(TP+FP) NPV<-TN/(FN+TN) LRP<-Sens/(1-Spec) LRN<-(1-Sens)/Spec Result<-list( "Confusion Matrix"=CM, "Confusion Rates"=Rate, "Accuracy"=Acc, "Error Rate"=Err, "Sensitivity"=Sens, "Sepcificity"=Spec, "PPV"=PPV, "NPV"=NPV, "LR+"=LRP, "LR-"=LRN ) print(Result,digits=3) }
こんな感じ。
my.indices(matrix(c(50,8,321,650),2,2)) $`Confusion Matrix` [,1] [,2] [1,] 50 321 [2,] 8 650 $`Confusion Rates` [,1] [,2] [1,] 0.04859 0.312 [2,] 0.00777 0.632 $Accuracy [1] 0.68 $`Error Rate` [1] 0.32 $Sensitivity [1] 0.862 $Sepcificity [1] 0.669 $PPV [1] 0.135 $NPV [1] 0.988 $`LR+` [1] 2.61 $`LR-` [1] 0.206
うむ。
ゼロ過剰ポアソン分布をデータにフィットさせる
ゼロ過剰ポアソン分布では2つのプロセスが考えられている。まず，最初に0の確率がσである二項分布で，値が0でないときのカウントデータがμに従うという。統計モデリングにおいて「ゼロが多いときに使うといい」とはよく聞くものの，厳密にいえば，このような二段階の過程が十分に考えられる事象を対象にするといいと思う。そして，殆どの場合は，先んじている0かそうでないかということにはあまり関心がないときに。
#使うパッケージ library(fitdistrplus) library(gamlss) #ある実測のデータ d<-c(0,5,6,0,1,8,0,2,0,1,4,2,0,3,8,9,1,1,0,10,1,2,0,2,11,3,0,8,5,11,3,8,8,8,14,7,8,6,5,1,8,7,2,0,6,8,3,3,0,6,2,11,2,0,0,0,0,0,0) #可視化 hist(d) #0の割合 length(d[d==0])/length(d) #まずはポアソンを普通にデータにフィット fit.pois<-fitdist(d,"pois") summary(fit.pois) #次にゼロ過剰ポアソン分布をデータにフィット fit.ZIP<-fitdist(d,"ZIP",start=1) summary(fit.ZIP) #２つのモデルを比較してみる gofstat(list(fit.pois,fit.ZIP),fitnames=c("Poisson","ZIP")) #これではσを出してくれない。gamlssを使う fit.ZIP2<-gamlss(d~1,family=ZIP) fit.ZIP2
うむ。
わーい！どんなときでも有意差を見つけられるフレンズなんだね！
「思ったように有意差が出なかったんですけどなにかこのデータから言えることはありませんか？」
「どんなデザインですか？」
「処置群・統制群，事前・事後，成果変数は1です」
「そうですね…まずは処置群を恣意的に何パターンかに分割してみましょう！そうするとそれぞれのグループの標本サイズが減るし検定が多重になるので，事前事後の得点差について第一種の過誤の可能性が高まります！まあ4-5グループに分割したら，有意水準を補正しない限りどこか有意になるでしょう！これが王道です」
「はい！」
「その後に，さもそのパターン自体に最初から興味があったように研究仮説自体を改定するのです！これはHARKingと呼ばれる有名な技法です！」
「先生，大変です！性別，学年，事前の得点，学校，恋人の有無，そういった変数で分けたグループでも，どこにも有意はありません！」
「安心してください！まだ手はあります！今度は事前と事後で伸びたひと，伸びなかったひとの数を，手持ちの2値のデモグラフィック変数すべての組み合わせについてカイ二乗検定を別々にするのです！別々にやるのが大事ですよ！」
「先生ぇ…どこも有意になりません」
「大丈夫！まだ手はあります！外れ値だとして，事後得点の下側のひとりずつを外していき，それぞれのステップで何度も検定をして，有意になったときにやめてそれを報告するのです！」
「先生ぇ…ダメですぅ」
「そうですか，ではまず10個のデモグラフィック変数をすべて投入してクラスタリングをしましょう。まずは階層的クラスタリングを何回も繰り返して，有意になる階層でやめましょう。それでダメならクラスタリング方法を変えましょう。それでダメなら，k-means法，それでダメなら混合分布モデル，それでダメならクラスタリングに入れるデモグラフィック変数を減らしていきましょう。方法や，入れる変数の組み合わせで考えると無数にパターンができますから，どこかで必ず第一種の過誤が起きるはずです」
「先生ぇ…どうしましょう」
「（まじかそろそろ有意水準いじろっかな…）いえ！まだまだ手はあります。そもそも成果変数はなんですか？」
「英語のテストの点です」
「そうですか。項目数は？」
「30です」
「よかった，では項目分析をして，弁別力がなかったり信頼性を下げたりする項目をひとつずつ減らしてその度毎に検定をするのです！そして有意になったらやめるのです。なあに，その項目はなかったことにするか，materialのところにテストの質の観点から予め外したと書けばいいのですよ！もうなんなら全項目別々に検定したっていいくらいだ！」
「先生，外していったら項目がなくなりました…」
「（おいまじかよ）こんなことは数学的に絶対無いはずですが…ううむ」
「わたしはもうダメでしょうか…」
「大丈夫です！その被験者にまた来てもらって新しくテストをやりましょう！遅延テストといえばいいのです」
「…それはできなさそうです」
「では，どんなデータでもいいから被験者について知っていることはありませんでしたか」
「ええと，実験ノートに実験中あくびをした人をメモっておきました。靴下の色もメモってますよ！」
「素晴らしいですね！そのあくびをした数人を外して分析すると？」
「あ！あ！先生！有意です！やっと有意になりました！先生ありがとうございます！私の研究に有意差がありました！」
「わーい！有意差だ！先生はどんなときでも有意差を見つけられる先生なんですね！」
「たーのし―！」
＊この記事はQRP（quetionable research practice）およびHARKing（研究仮説の改定），そして（本来するべきではない）事後的分析において無理に有意差を見出すことの悪質さについての理解を（私自身が）深めるために書いたものです。
モデルの中で何が捨象できるかを語らない科学
数理モデルというものは，その記述の仕方の形式性の割には，数理モデルということばに親しみを感じないほとんどのひとが思うより，本来結果主義的で効用主義的なものだ。
モデルは，もちろん現象それ自体ではないし，その現象を大幅に捨象していて，しかしそこから得られる予測や知見が有益だと見込まれているものだ。この「有益だ」という考え方は一部の学術分野にはないこともあるが。
人文社会系のほとんどの数理モデルは，世界が数字に支配されていて，その世界の斉一的な決定則を表すものだなんてことを意味しない。せいぜいが「観測がうまく当てはまる」，「うまい数理的な近似になっている」という程度の含意である。しかし，そのモデル（世界の決定則それ自体ではない）について考えることで，人が適切に意思決定をできたり，そして個人間の合意が得られ，判断の公共性が発生する（たとえばある種のエビデンスになる）というに考えられている。
たとえば，TOEICのスコアは大学が自前でやっている単語テストの成績から予測できるとする。
y = ax + b
という簡単なモデルを考えて，
TOEIC = 8.5×単語テストの得点 - 120
とか，そんなふうに。
このモデルは，それがうまくTOEICのスコアを予測できるとか，つまりこの現象のいい数理的な近似になっていれば，まあいいモデルだといえる。
ただし，このとき，TOEICは単語の力だけではないのだから，このモデルは世界の斉一的な決定則の記述として不完全だ！というのはちょっと勇み足だと思う。
TOEICは単語のテストだけではない！とか。
もっというなら，TOEICのときの部屋の温度をモデルに取り入れるべきかだとか，世界の斉一的な決定則を考えれば無限にそんな要因はありえる。
数理モデルを扱う人がもつ，結果主義というか，もっと大きくいってプラグマティズム的な考え方の下では，単にこのモデルはそれらを捨象しているのである。線形モデルなら，それら全部をひっくるめて誤差（残差）として考える。
ここで問題は，何が捨象できるかそして何を積極的に捨象するべきかということになる。
…外国語教育研究やその関連分野である第二言語習得研究が扱ってきた変数や現象は莫大であり，確かに複雑極まりない。もはや還元主義の傾向が強いテーマに関しては指数的に用語が増えていっていて，収拾付かない例もある。
だからこそ，なにが捨象できるか，そしてなにを捨象するべきか，という考えが今後重要になってくると思う。
そして個人的には，大局的に人間の意思決定や判断の合意形成や公共性の創出という点から見れば，捨象できるものはずっとずっと多いと思っている。
…世界の斉一的な決定則それ自体を考えるならば，むしろなにも捨象できないとも思っている。
Rで日付データの処理
自分用のメモ。
#日付クラスへ変換 d<-"2016-1-1" d2<-as.Date(d) class(d2) #日付データの足し引き（日付クラスだとこれができるようになるのが最高） d2-1 d2+1 d2-1000 d2+1000 #基準日から1日毎にログイン回数を累積計算 #datは時間とログイン回数のデータフレーム #基準日を設定 start<-as.Date("2016-04-08") r<-numeric(130) #計算 for(i in 1:130){ r[i]<-sum(dat[dat[,1]<start+i,2]) }
