※現在、ブログ記事を移行中のため一部表示が崩れる場合がございます。
順次修正対応にあたっておりますので何卒ご了承いただけますよう、お願い致します。
2017年 05月 26日
linkFunction <- function(b1, b2, x) { function(r) { 1 / (1 + exp(-b1 - b2 * x - r)) } } l <- linkFunction(-0.5, 0.3, 2) cat("r\ty=0\t\ty=1\t\ty=2\t\ty=3\n") sum <- vector(length = 4) for (r in seq(-5, 5, by = 0.5)) { q <- l(r) cat(r, "\t") for (y in 0:3) { b <- dbinom(y, 3, q) n <- dnorm(r, mean = 0, sd = 2) mixture <- b * n sum[y + 1] <- sum[y + 1] + mixture * 0.5 cat(mixture, "\t") } cat("\n") } cat("sum\t", sum[1], "\t", sum[2], "\t", sum[3], "\t", sum[4], "\n")リンク関数は `x` と `r` を変数に含んでいます。 `x` は `2` に固定した場合について計算しました。 なので、あとは `r` を与えれば、二項分布している事象の発生確率 `q` が求まります。 `r` は正規分布する変数です。 今回は`r` は平均 `0`、標準偏差 `2` の正規分布から得られるとしました。 事象の試行回数は `3` にして、`0` 回から `3` 回事象が発生する場合の確率を計算して出力しています。 `r` は `-5` から `5` まで、`0.5` 刻みで計算しました。 コードを実行すると、次の出力が得られます。
r y=0 y=1 y=2 y=3 -5 0.008571241 0.0001914794 1.425867e-06 3.539279e-09 -4.5 0.01529937 0.0005635068 6.918364e-06 2.831304e-08 -4 0.02542036 0.00154367 3.124683e-05 2.108318e-07 -3.5 0.03909264 0.003913947 0.0001306212 1.453086e-06 -3 0.05514584 0.009102904 0.0005008711 9.186514e-06 -2.5 0.07038014 0.01915423 0.001737632 5.254481e-05 -2 0.07963943 0.03573468 0.005344786 0.0002664708 -1.5 0.07772425 0.05749969 0.01417925 0.00116552 -1 0.06325714 0.0771553 0.031369 0.004251228 -0.5 0.04148674 0.08342818 0.05592358 0.01249557 0 0.02138051 0.07088734 0.07834263 0.02886066 0.5 0.008601664 0.04701976 0.08567559 0.05203704 1 0.002741934 0.02471168 0.07423798 0.07434107 1.5 0.0007137071 0.01060504 0.05252712 0.08672284 2 0.0001570974 0.003848653 0.03142876 0.08555086 2.5 3.018185e-05 0.001219082 0.0164134 0.07366188 3 5.187402e-06 0.0003454491 0.007668262 0.0567399 3.5 8.116416e-07 8.911395e-05 0.003261413 0.03978732 4 1.169644e-07 2.1173e-05 0.001277585 0.02569661 4.5 1.564146e-08 4.66824e-06 0.0004644167 0.01540073 5 1.950226e-09 9.596395e-07 0.0001574019 0.008605787 sum 0.2548242 0.2235203 0.2303399 0.2828235`r` の増加に従って `q` も増加するので、`r` が小さいうちは事象の発生回数は少ない方が確率が高く、`r` が大きくなると多い方が確率が高くなるのが分かります。 混合分布では、`r` についての積分がとられるので、最後の行にざっくりとした積分の近似値を出力しておきました。 二項分布と正規分布の混合分布は、だいたいこういう数字になります。