先日ついにtidyrパッケージで縦横変換をやってみた。今までずっと覚えなきゃなと思ってて結局億劫でやってなかったんだけれど手伝いで仕方なく。

データハンドリングで，横に広がる（列数の多い）wide-formatのdataを分析に便利なlong-formatに変換するというのは，分析をやるまえのデータ整形の過程で必ず通らなければならない道。Excel上で地道にやるのもあるにはあるだろうけど，Rのtidyrを使えばすぐできるというお話。これ系の話はググればネット上にゴロゴロ記事が転がっている。なので，今回は具体的に外国語教育研究系のデータ分析ではどんなときに役に立つかなーという話です。特に無駄に変数が多い時なんかは結構役に立つと思います。

sample data

適当にこんな感じでデータが入っているとします（数値は乱数）。事前・事後・遅延事後のテストがあって，それぞれで流暢さ（fluency），正確さ（accuracy），複雑さ（complexity）の指標があるみたいな。1列目には被験者番号があって，それぞれの列名は”test.measure”という感じでドット区切りにしてあります。で，これを，「事前・事後・遅延事後」の3つのカテゴリカル変数からなる”test”列と，3つの指標（CAF）の”measure”列の2列に分解しましょうという話し。

R上に上のようなデータを読み込みます。

read.table関数で読み込んだ場合（※MacなのでちょっとWindowsと違います）

まずはパッケージの準備

install.packages(“tidyr”) #パッケージのインストール

library(tidyr) #tidyrパッケージの呼び出し

datに入っているデータフレームをまずは，変数名の列（variable）とその値（value）の列に集約します。key=変数名の列，value=値の列という引数指定です。

#%>%はパイプ演算子。以下の関数で扱うデータフレームを指定するということです

dat %>%

#-subjectでsubjectの列を除外。変換したdataをdat2に入れる

gather(key = variable,value = value,-subject) -> dat2

long formatのデータに変換されました

これで変数名を1列に集約できました。ただし，この列にはpre, post, delayedというテスト実施時期の情報と，fluency, accuracy, complexityという指標の情報が混ざってしまっています。これを，2列に分けます。

dat2 %>%

#colで分けたい列名，intoで分けた後の列名を指定。sepでセパレータを指定しますが，デフォルトはあらゆる非アルファベットになっているので，ドットなら指定しなくても大丈夫

separate(col = variable, into = c(‘test’, ‘measure’)) -> dat 3

さて，dat3の中身を確認してみると…

long formatのデータに変換されました

おおおおー！！！！！

というわけで，横に変数がたくさんある場合もこのtidyrのgatherとseparateを使えば簡単にlong型に変換できそうです。

参考サイト：http://meme.biology.tohoku.ac.jp/students/iwasaki/rstats/tidyr.html

なにをゆう　たむらゆう

おしまい。

ストレスが溜まっているのかなんなのか帰り道にセブン-イレブンでチーズとワインを買って飲むとかいう珍しいことをしているという言い訳を最初に書いておこうと思います。それから，自戒もたっぷり込めていることも申し添えておきたいと思います。

データ分析の方法の相談についてです。色んな所で起こっていることだと思うし色んな人がきっと同じようなことを考えているだろうと思います。一言でいうと，データを取る前に相談してみましょうということです。データを取ったあとに，「このデータってどうやって分析すればいいですか？」っていう質問は結構つらいところがあったりします（恐れ多くもなぜだかこういう質問をされる立場に私もなってしまいました）。なぜなら，本当にその質問者の方が調べたいこと，明らかにしたいことが，そのデータ収集方法で明らかにすることができるのかちょっと怪しかったりすることがあるからです。

そういう状況にならないために重要なことは，研究課題を見つけて，それを具体的な実験計画なりに落としこむ際に，その先のことまでちょっと考えてからデータ収集をしてみることだと思います。つまり，

こういう方法を使うと，こういうデータを収集することができる。それをこういう分析にかければ，こういう結果が得られるはず。もしこういう結果が得られた場合，それはこういうことを意味していると解釈できるだろう。

みたいなことを考えてから実験を組んで，データを取っていただきたいなということです。そうすれば，データを取ってからどうしたらいいかという質問が来ることはなくなるからです。上述の事でイマイチ自分がうまく説明できないとなったら，その時点で相談をすればよいと思います。そうすれば，実はその方法で得られたデータで，質問者の方の主張をするのはちょっと違うので，こういうデータ収集をしないといけませんよ，なんてアドバイスができたりするかもしれないからです。これって仮説検証型の研究でも探索的な研究でも同じだと思います。たとえ探索的な研究だったとしても，得られた結果に対してなんらかの解釈を加える必要はあるはずですから。研究課題がwh疑問文で始まる形であったとしても，結果に対しては，きっとこういう理由でこうだったのだろうということを言えなければいけないわけですし。

データをとって，分析をしてみて，これであってるのか不安だから相談してみようということもあるのかもしれませんが，それって別にデータ取る前にもできることなんですよね。

こういう研究課題をもっていて，これを明らかにするためにこういう分析をしてこういう結果が得られたらこういう解釈ができると思うのですが，これで合ってますか？

ってデータ取る前に相談すればすむと思うのです。なのにどうして，「データを取ったあとにこれってどうやって分析すればいいんですか？」ってなってしまうのでしょう。

これはあくまで持論なのですが，研究って，研究のデザインとか実験計画ですべて決まると私は思っています。面白いなと思う研究はやっぱりデザインがしっかりしていて，分析の結果から得られたことの解釈がとってもシンプル（注1)なんです（そういう研究が良い研究だと思う私からすると，discussionを長々と書いている論文やそれが良しとされる風潮，また「discussionが短い」みたいなコメントってなんだかなと思ったりします）。それって，デザインを練る段階で色んなことを考えているからなんですよね。研究の計画を立てる段階で，データの分析方法まで明確にイメージできない場合は，データ収集に入るべきではないと私は考えています。私の場合は基本的に，ある研究で行われることを考えることは，R上で分析を行う際のデータセットの形をイメージすることとほぼ同義です。変数のデータフレームがイメージできていて，それをどうやって分析すれば何をどうしていることになるのか，こういう結果になったらそれはどう解釈するのか，また別の結果がでたらそれをどう解釈するのかまで考えられて初めて研究のデータ収集ができるのではないでしょうか。

こういった考えが当たり前になれば「統計屋さん」（私は違いますけど）とか言われる人たちもずいぶん楽になるのではないかなと思っています。

なにをゆう　たむらゆう

おしまい。

*注1: だからといって，t検定一発で明らかになることってそんなにはないとは思います。少なくとも今なら昔はt検定一発だったこともGLMMするでしょう。例えば正答率の平均値をt検定というのはもはや時代遅れになりつつあるわけですし。

もう2ヶ月くらい前の話ですが，Nagoya.R #12で発表しました。

今更って感じなのですが，ちょっと今お手伝いでLMEをやっていて，自分でも理解があやふやな点がぼろぼろと出てきたのでメモ的に。

僕は発表中にlmer関数と，モデル比較に使えるstep関数の使い方を説明したのですが，どうやらstep関数はlmer関数で出力したもののみに対応している模様。というか，正規分布を仮定したモデルにしか適用できないみたいです。なので，lmerでfamily指定を使ってポアソン分布(poisson) や二項分布(binomial)を指定した場合には出力がされません。また，lmer関数でfamily指定すると警告メッセージでglmer関数を使うようにと言われます。なので，正規分布以外でやるときはglmer関数を使うほうがいいかもしれません。回帰の式の入力はほぼ一緒です。

ただし，glmer関数は結果の出力の解釈が実はちょっと難しくて，一発だけじゃ全水準の多重比較までみれないんですよね。なので，ダミー変数にいれたものの順番を入れ替えて計算を回していくようなのですが，それに使う引数がstartってやつっぽいのですよね。そんないちいちダミー変数入れ替えるなんてめちゃめんどくさいわけで，これで引数指定して一番最初にいれる水準を指定できるようになっているみたいなのです。ただしちょっと使い方がまだよく理解できていなくてですね…

実際に自分の研究でちゃんと使えるようになるにはここは避けて通れないわけなのですが，そっちばっかりに手を回しているわけにもいかず，RのヘルプやマニュアルとにらめっこしてはGoogle検索して…とかやってたらなんか1日終わっているみたいな幸せなのか不幸せなのかわからない昨日今日です。

いろいろ終わってません

なにをゆう　たむらゆう

おしまい。

論文に出てる記述統計からシュミレーションをしてみようと思った。いろいろググってみたのだけれど，なかなかうまい方法がなかったのでメモ。

Rでは，様々な分布に従う乱数を発生させることができる。一番有名なのはrnormというやつ。正規分布に従う乱数を発生させる。平均と標準偏差を指定して，こんなかんじで。

>dat <-rnorm(100, mean=50, sd=5) #平均値が50で標準偏差が5の乱数を100個

> dat
[1] 47.49772 48.21094 56.25741 49.96995 44.78259 48.71722 50.17432 55.17608
[9] 44.92424 46.83734 48.06748 53.36903 46.53421 50.15734 40.72731 48.05790
[17] 55.66313 53.06009 50.12512 46.79404 54.13993 45.70641 45.05168 51.28307
[25] 50.83705 48.71708 48.01274 49.00765 56.62082 49.58216 52.19713 63.85808
[33] 48.90386 56.36748 48.44979 53.96921 43.93527 50.40558 46.51421 48.81157
[41] 39.57205 51.43287 50.73080 52.36210 44.21936 59.07660 52.05478 54.53502
[49] 52.26897 54.84681 45.80431 53.86563 53.97093 55.28294 55.68036 57.21981
[57] 48.21196 42.46615 44.65935 46.11213 52.64849 47.62270 46.36596 52.57540
[65] 47.27981 63.84164 52.44592 43.90385 50.75210 51.00874 57.80445 50.17152
[73] 47.95725 50.36983 46.72289 47.57410 44.48939 53.98110 59.12215 47.21243
[81] 42.14203 40.78201 50.64352 48.90662 53.63576 55.65404 46.23628 48.17049
[89] 51.76097 52.69611 49.41447 53.25395 52.08408 41.06224 53.22397 41.33806
[97] 58.23008 44.21629 53.06959 53.20374

ただしこのやり方では，SDが大きくなると負の数値が出てしまったりする。また，例えばテストの点数とかの乱数を作りたい場合，100より上もいらないのに100を超える数値も出る。

> dat1 <- rnorm(100, mean=60, sd=30)
> dat1
[1] 44.494760 89.170093 27.292283 66.450161 62.407602 66.882294
[7] 61.546427 67.186071 49.042632 63.244362 84.543721 60.576172
[13] 51.848288 83.655269 41.251925 116.468853 64.299882 88.146247
[19] 18.110042 47.436176 56.263604 96.849951 89.739825 23.088249
[25] 34.150477 2.516926 4.948293 56.340174 25.665198 95.580383
[31] 83.448675 113.616310 -3.453262 108.500719 56.467913 100.086726
[37] 66.637288 43.201111 59.247994 45.427133 39.968471 37.975361
[43] 16.056423 22.321241 29.502743 39.415555 68.695397 94.841580
[49] 49.784221 82.842063 53.112870 7.554488 34.891446 38.734127
[55] 105.762073 39.075296 89.386547 56.279559 28.455354 84.684350
[61] 84.719927 86.950250 72.399872 48.875458 73.908410 41.047765
[67] 38.112460 38.006715 25.781026 64.453294 63.488881 62.025456
[73] 88.706781 53.570347 71.571211 37.240862 32.510384 20.557917
[79] 48.065483 94.325533 94.798381 64.708053 108.530045 26.543721
[85] 58.603449 52.491463 85.727697 106.653768 47.484039 102.625380
[91] 49.518054 36.609159 27.157635 59.725756 39.526326 82.433037
[97] 90.007551 97.304846 29.513355 80.113184

これと同じようなことは，Excelを使うと次のような式でできる。

=NORMINV(RAND(),60, 30)

ただしこれでも先ほどのように負の値と100位上の値ができてしまう。うーむ。数学的なことがわからないので困った。Excelを使う場合，ある程度の数値ならばここに書いてあるような方法を使ってできるみたい。

http://www016.upp.so-net.ne.jp/sige-lab/before2002/1995make_rand.pdf

ただしこのやり方もSDが広すぎると対応できない。と，ぐぐっていたらこんなのをみつけた。

R言語について。正規分布に従う疑似乱数を発生させる「rnorm()」関数で、正の値だ… http://detail.chiebukuro.yahoo.co.jp/qa/question_detail/q1177484151?fr=pc_tw_share_q #

ここで回答者さんが平均値を指定して，負の値を出さないようにする乱数発生のRスクリプトを書いてくださっている。引用します。

>zz<-rnorm(10000,mean=3.34)

>while(any(zf<-zz<0)){#ひたすら粘る(お勧めできない)#
>zz<-rnorm(10000,mean=3.34)}

>abs(zz)#原点で折り返す#

>zz[zz>=0]#負値は省く#

>(zz>0)*zz#0に圧縮#

>ifelse(zz<0,3.34,zz)*zz#平均を強制#

>while(any(zf<-zz<0)){#正値で上書き#
zz[zf]<-rnorm(sum(zf),mean=3.34)}

これのrnormの引数でSDも指定してやればいいのではないか，と。

>zz<-rnorm(100,mean=60. sd=30)

>while(any(zf<-zz<0)){#ひたすら粘る(お勧めできない)#
zz<-rnorm(100,mean=60, sd=30)}

>abs(zz)#原点で折り返す#

>zz[zz>=0]#負値は省く#

>(zz>0)*zz#0に圧縮#

>ifelse(zz<0,60,zz)*zz#平均を強制#

>while(any(zf<-zz<0)){#正値で上書き#
zz[zf]<-rnorm(sum(zf),mean=60, sd=30)}

ただしこれだけでは100以上の値が含まれるので，

>ifelse(zz>100, 60, zz)

として0以下の値の処理と同じように平均に置き換えてあげる。シュミレーションした平均値と標準偏差を出すと，こんな感じ。

> mean(zz)
[1] 65.13294
> sd(zz)
[1] 29.33622

平均が5点ほど上にずれたけどsdはまあまあ。このやり方は，範囲の外を取る値を平均値で置き換えるやり方だけれど，それらの値を境界の0と100で置き換えるということもありなんじゃないだろうか。やってみる。

> zz<-rnorm(100,mean=60,sd=30)

while(any(zf<-zz<0)){#ひたすら粘る(お勧めできない)#
zz<-rnorm(100,mean=60, sd=30)}

abs(zz)#原点で折り返す#

zz[zz>=0]#負値は省く#

(zz>0)*zz#0に圧縮#

ifelse(zz<0,0,zz)*zz#負の値を0に置き換え#

while(any(zf<-zz<0)){#正値で上書き#
zz[zf]<-rnorm(sum(zf),mean=60,sd=30)}

ifelse(zz>100, 100, zz) #100以上の値を100に置き換え#

> mean(zz)
[1] 60.1755
> sd(zz)
[1] 27.46259

標準偏差は少し小さくなったが平均はさきほどよりも近くなった。もともとこういうやり方では最初に指定したものとぴったり一致することはないにせよ，このやり方だとまあまあの精度でシュミレーションができそう。ただまあサンプルサイズがこれより少なくなると推定の精度は下がる。逆に増やせば上がる。

最後に。ここまでやってきたやり方は勝手に標本分布を正規分布と仮定してやっているのでその点は注意が必要。t検定なんかやってるんだったらまあ正規性を仮定してのことなので良いけれど，ポアソン分布とかカイ二乗分布とかだったらこのやりかたはできない。そういう分布では平均と標準偏差の意味するものが違ってくる。そういう場合も分布のパラメータがわかればシュミレーションはできるのかな。またなにか機会があったらそっちでも遊んでみよう。

そんな感じで昨日の基礎研ではシュミレーションしたデータを使って，論文のRQに答える検定を自分で考えてやってみるということをやろうとしたのでした。バタバタと準備したのでなんかつまらなくなってしまったけれど。もうちょっとちゃんと準備してやればもっと勉強にもなったよなと思うので，また基礎研で僕の番が回ってきたら同じようなことにチャレンジしてみたい。今度は，テキストファイルのデータをExcelにまとめてそれをRでもlangtestでもSPSSでなんでもいいから使って結果を解釈するところまで。

なんかあの準備不足で色々あたふたする感じが僕の普段の授業がどうなっているかを表している気がしないでもなかったり…（遠い目

というわけで，そうちゃんと戯れてこようと思います（じゅるり

なにをゆう　たむらゆう

おしまい