関数解読histprior

Bayesian Computation with R: Second Edition (Use R!)

（ヒストグラムのように）均等な区間ごとに確率が定義された確率分布の密度を計算する

使用法

histprior(p,midpts,prob)

引数
p:連続型確率変数において，対応する密度を計算する確率変数の値をまとめたベクトル
　以下の例のように，できるだけ細かく設定する。
midpts:各区間の中点をまとめたベクトル
prob:各区間の確率をまとめたベクトル

値
各pに対応する確率密度をまとめたベクトル

histprior = function (p, midpts, prob)
#関数histpriorの中身
{
    binwidth = midpts[2] - midpts[1]
     #ベクトルmidptsの第2要素と第1要素の差から，
     #1つの区間の長さを求めて，binwidthに代入。
browser()
    lo = round(10000 * (midpts - binwidth/2))/10000
     # 10000 * (midpts - binwidth/2)の小数点以下を四捨五入
     # その後，10000で割り，loに代入。
     #
     # 結局，各区間の下限を求めている。

    val = 0 * p # ベクトルpと同じ長さの0ベクトルを作成（空箱作り）。

    for (i in 1:length(p)) {
        val[i] = prob[sum(p[i] >= lo)]
         # p[i] >= lo
         #  p[i]が各区間の下限以上かどうかを判定
         #  p[i]が第1区間に含まれるなら
         #  TRUE FALSE FALSE.....となる。
         #  p[i]が第2区間に含まれるなら
         #  TRUE TRUE FALSE.....となる。
         # sum(p[i] >= lo)
         #  p[i] >= loが成立する区間の数を集計
         #  p[i]が第1区間に含まれるなら，TRUEは1個だから1
         #  p[i]が第2区間に含まれるなら，TRUEは2個だから2
         # prob[sum(p[i] >= lo)]
         #  引数のprobの第"sum(p[i] >= lo)"要素を出力。
         #  p[i]が第1区間に含まれるなら，prob[1]
         #  p[i]が第2区間に含まれるなら，prob[2]
         #
         # 結局，pの各要素（密度を計算する各点）について，
         # その点を含む区間の確率を密度と定めるということ。
    }

    return(val)#valを出力する。
}

例

#各区間の中点を設定。
midpts=c(.1,.3,.5,.7,.9)
#各区間の確率を設定。
prob=c(.2,.2,.4,.1,.1)
#密度を計算する各点を設定。
p=seq(.01,.99,by=.0001)

#密度を保存。
histden = histprior(p,midpts,prob)

#(p,histden)を図示。
plot(p,histden,type="l")