ダメ出し：提示する統計グラフに注意 その2

前の記事（元記事も）での religion = unknown の記号の色が gray で，バックグラウンドの色と近くて，目に入りにくいということもあるし，プログラムがごちゃごちゃしているので religion 別の散布図を描いた方がよさそう。

ということで，描いてみた。シンプルに plot を使って（余分な情報はいらない）。

df.split <- split(df.merged, df.merged$religion)
layout(matrix(1:4, 2, byrow=TRUE))
par(mgp=c(1.6, 0.6, 0), mar=c(3, 3, 0.5, 0.8))
lapply(df.split, function(d) {
    plot(babies~income, data=d, xlim=c(0, 80000), ylim=c(0, 7),]
         pch=19, col="#00330060", cex=0.7)
    text(35000, 6, d[1, "religion"])
    })
layout(1)

• Eastern religions は国の数も少ないが，babies が 3 ぐらいまでで，他とはちょっと違う。
• income < 15000 では babies > 3 の国がかなりある。income > 15000 では babies はほぼ一定。

ダメ出し：提示する統計グラフに注意

"Getting data from GapMinder.org" のグラフについて

当初は，このグラフ（図1）を見て，「記号の大きさが人口に比例していないなあ。困ったものだ」と思っていた。おまけに，記号の大きさは人口そのままではなく，log(人口+1e7) になっている。

図1

ggplot は記号の大きさを，最大値と最小値から決めるようで，筆者が 1e7 を採用したのは，記号の大きさのバランスがちょうど好みにあったのだろう。population をそのまま使った図2と比較すればよい。

図2

そして，このふたつの図を比べると，受ける印象が全く違うことにびっくりする。どちらかの図が不適切なわけだ。だって，どちらも適切なわけがない。
筆者がなぜ babies とincome, population, religion を取り上げたのか真意は不明だが，このデータを分析してみる。
まず，religion 別に population, babies, income 相互間の相関を見る。
人口が極端に大きいふたつの国(インドと中国）の影響を除くために，スピアマンの順位相関係数を計算した。

> df.split <- split(df.merged, df.merged$religion)
> lapply(df.split, function(d) round(cor(d[,3:5], use="pair", method="spearman"), 3))
$Christian
           population babies income
population      1.000 -0.029 -0.028
babies         -0.029  1.000 -0.721
income         -0.028 -0.721  1.000

$`Eastern religions`
           population babies income
population      1.000  0.104 -0.275
babies          0.104  1.000 -0.871
income         -0.275 -0.871  1.000

$Muslim
           population babies income
population      1.000  0.023 -0.098
babies          0.023  1.000 -0.757
income         -0.098 -0.757  1.000

$unknown
           population babies income
population      1.000  0.043  0.028
babies          0.043  1.000 -0.689
income          0.028 -0.689  1.000

この結果を見ると，religion に係わらず，3変数は同じような相関関係にあることがわかる。Eastern religions で，population と babies, income の間の相関が他と比べると若干高いが。babies も income も population で調整されているので，こういうことをしなくても population の要因は除外してよいだろうということがわかる。とすれば，当初描かれた図は，特に人口サイズが強調された図は不適切ということであろう。ということで，図3を得る。

図3

（つづく）

ダメ出し：そんなことやっちゃダメだってば!

「haploid から diploidへの座標変換」 の，「座標で表す」の項だけど，

むちゃくちゃしますね。びっくりします。ほんとに。

p1 <- c(0.3,0.7)
p2 <- c(0.4,0.6)
f <- 0.6 # r^2の値
# 連鎖平衡の下でのハプロタイプ頻度
h.LE <- p1 %*% t(p2)
# 連鎖不平衡の下でのハプロタイプ頻度
delta <- f^(1/2) * prod(h.LE)^(1/4)
h.LD <- h.LE + delta * matrix(c(1,-1,-1,1),2,2)
# r^2は以下で求まる統計量
chisq.test(h.LE,correct=FALSE)
chisq.test(h.LD,correct=FALSE)

> h.LE
     [,1] [,2]
[1,] 0.12 0.18
[2,] 0.28 0.42

> h.LD
          [,1]        [,2]
[1,] 0.2938965 0.006103479
[2,] 0.1061035 0.593896521

こんなものを chisq.test に渡しちゃダメダメ!!
引数は，集計表!
こんなこと，初心者だってやらないよ。