2010年10月18日月曜日

イラレの破線とMS Word

イラレの破線は「線分」と「間隔」に数値を入れて線と間隔の長さを調整するんだが、線分に数値を入れただけで破線が作れてしまう。
というか線の設定パネルにある破線のチェックボタンを入れると、「線分」のみに数値の入った状態で破線ができる。
これがデフォなんだろうね。
そのままイラレで印刷するんなら問題ないのだが、MS WordにPDFで貼り付けて印刷すると実線として出てきてしまう。
試行錯誤したあげく、イラレで間隔にも数値を入力してやったものをWordにもってくるとちゃんと破線として印刷されることが判明。
ちなみにプレビューに貼りつけてもWordと同じ出方をした。
イラレCS4 (Mac) がアフォってことでしょうかね。

2010年9月17日金曜日

統計よりも一人のストーリー

この記事がかなり面白かった。
いわく、ある災害などで被害にあった人が◯◯人、と言うよりも、そのうちの1人にスポットを当ててどういういきさつでどういう被害にあったか、というストーリーを伝えたほうが、寄付金が多く集まるのだという。
これを「身元のわかる被害者効果(identifiable victim effect)」と言うらしい。
その通りだなあと同意していたら、話にはまだ続きがあった。

この、統計データと一人のストーリーとでの効果の違いには、個人差があるらしい。
ものごとについて、分析的な傾向を持つ人は、この差が小さいそうだ。
これって、要は、理系的な人と文系的な人ってことじゃん。
理系的、文系的って表現はいかにもミーハーな感じがしていたけれども、自分の中で急に説得力が増してしまった。

とかく研究に携わる人は、「分析的」もとい「理系的」でなければならない。
でも完全に理系的な人っていないよなあ。
論文のイントロだって、こういうケースでこういう条件で調べられた研究はない、とか、「一人のストーリー」的な論旨に持ち込んで研究の動機とする論文も山ほどあるわけで。
そもそも、文系的にアピール力がない研究はつまらないような気もする。

なんだかとりとめもない文章だけど、ここで思考が終了したので、記事もここまで。

2010年8月21日土曜日

自分用反復測定データにおけるGLMMまとめ

現時点ではlme4パッケージのlmer関数を用いるのが一般的と思われる。
測定日が一定間隔でなかったり、測定日によって繰り返し数が異なる(missing dataがある)ようなデータセットの場合、repeated measures ANOVAはあまり有用ではない。
そこでGLMMの登場。
繰り返し(rep)をランダム効果にすることで、missing dataがある場合でも手元にある全データを活用して統計できる。
固定効果として、測定日(day)、生物種(sp)の2つあったとする。
どんな固定構造のモデルが最も当てはまりが良いか探索することにする。

> attach(dataset)
> sp=factor(sp)
> rep=factor(rep)
> library(lme4)

(1)最適モデルの探索
> model0 <- lmer(Y ~ 1 + (1|rep), data=dataset, REML=FALSE)
#固定効果なしモデル
> model1 <- lmer(Y ~ sp*day+(1|rep), data=dataset, REML=FALSE)
#フルモデル
> model2 <- lmer(Y ~ sp+day+(1|rep), data=dataset, REML=FALSE)
#交互作用なしモデル
> model3 <- lmer(Y ~ sp+(1|rep), data=dataset, REML=FALSE)
#種のみモデル
> model4 <- lmer(Y ~ day+(1|rep), data=dataset, REML=FALSE)
#調査日のみモデル

REML=FALSEにしたのは、固定構造(fixed structure)の異なるモデルを比較する場合、REMLでなくMLで推定したモデルでなければならないためである。
lmerではFALSEにするとML、TRUEにするとREMLを使った推定となる。
(参考文献) Zuur et al. 2009 "Mixed Effects Models and Extensions in Ecology with R"
Pinheiro & Bates 2000 "Mixed Effects Models in S and S-Plus"

(1-1)パターン1 手動でAIC比較
> AIC(model0, model1, model2, model3, model4)

(1-2)パターン2 自動でAIC比較
フルモデル(model1)のみあればよい。
> library(MuMIn)
> dredge(model1)
これで固定効果の全組み合わせを自動的に計算して、AICやらAICcを一覧にしてくれる。

(1-3)パターン3 尤度比検定
> anova(model0, model1, model2, model3, model4)
ただし3つ以上のモデルを比較する場合、多重比較の問題が出てくるので、この後bonferroni補正などが必要。
具体的にどうやるかは、勉強中。
AICによるモデル選択に惹かれつつあるのでこのまま放置かも。

(2) 固定効果の多重比較: どれとどれが異なるのか?
種(sp)のみが固定効果として残り(model3)、どの種とどの種が異なるのか検定したい場合について。

(2-1)多重比較法
> library(multcomp)
> summary(glht(model3, linfct=mcp(sp="Tukey")))
これはTukey検定を行う場合。他の方法も使える。多分。

(2-2)多重比較的モデル選択(勝手に命名)
A, B, Cという3種の比較を行うとする。
データセットに列を追加し、ダミー変数を入れて、モデル選択を行う。
例えばAだけがB, Cと異なり、BとCは差がない、というモデルなら、種Aの行には「A」というダミー変数を、種BとCの行には「BC」というダミー変数を入れる。
同様に、A=B≠C、A≠B≠Cなど、あらゆるパターンの分だけ列を追加してダミー変数を入れていく。
rep sp day Y ABC A_BC AB_C AC_B A_B_C
1 A 1 Y ABC A AB AC A
2 A 1 Y ABC A AB AC A
~
1 A X Y ABC A AB AC A
2 A X Y ABC A AB AC A
~
1 B X Y ABC BC AB B B
2 B X Y ABC BC AB B B
~
1 C X Y ABC BC C AC C
2 C X Y ABC BC C AC C
~

このあと(1)と同様にモデル選択する。
> modelABC <- lmer(Y ~ ABC+(1|rep), data=dataset, REML=FALSE)
#3種で差がないモデル
> modelA_BC <- lmer(Y ~ A_BC+(1|rep), data=dataset, REML=FALSE)
#A≠B=Cのモデル
> modelAB_C <- lmer(Y ~ AB_C+(1|rep), data=dataset, REML=FALSE)
#A=B≠Cのモデル
> modelAC_B <- lmer(Y ~ AC_B+(1|rep), data=dataset, REML=FALSE)
#A=C≠Bのモデル
> modelA_B_C <- lmer(Y ~ A_B_C+(1|rep), data=dataset, REML=FALSE)
#A≠C≠Bのモデル

> AIC(modelABC, modelA_BC, modelAB_C, modelAC_B, modelA_B_C)

2010年8月20日金曜日

ML vs. REML

Mixed Model勉強中。
今トレンド(?)のlmerですが、デフォルトはREMLになってます。
しかしいろいろ調べてみると、固定効果部分が異なるモデルを比較する場合、REMLを使っちゃいかんようです。
REMLは固定構造を考慮せずランダム構造のみの尤度を最大化するのがその理由らしい。
というか固定構造が変わるとREMLの評価基準も変わってしまうようだ。
固定構造を考慮しない推定法で、固定構造の異なるモデルの比較を行うこと自体ナンセンス、ということだろうか。
ここの解釈がまだアヤシイ。

いろいろ見てたらTREEにこんな論文を見つけた。
Bolker et al. 2008 "Generalized linear mixed models: a practical guide for ecology and evolution"

この本もいい。ほすぃ。
Zuur et al. 2009 "Mixed Effects Models and Extensions in Ecology with R"