2山のヒストグラムの間のしきい値は・・・？

２山のヒストグラムにおいて、
その２山が完全に離れているとします。
その時、２山間にしきい値を引くとしたらどの位置が妥当なところなのでしょうか？
妥当とは、
山(1)がしきい値を超える確率と山(2)がしきい値を超える確率の和が、
最も小さい状態のことです。確率はエクセルのNORMDIST関数で求めました。

ちなみに、
山(1)の平均値μ1：236
山(1)の標準偏差σ1：17.35
山(2)の平均値μ2：467
山(2)の標準偏差σ2：26.05
です。

最初は（μ1+a*σ1）=（μ2-a*σ2）という式をたて、
aについて解き、
μ1+a*σ1をしきい値としたのですが、
山(1)がしきい値を超える確率と、
山(2)がしきい値を超える確率を足してみると、
μ1+a*σ1より少しずれたところに確率が最も小さくなるところがありそうです。
ちなみにa=5.34248…で、
実際はa'=5.16056…のようです。
この5.16056…を計算で出せる方法はないでしょうか？

よろしくお願いいたします。

No.3 ベストアンサー 回答者： quaestio 回答日時： 2009/06/26 21:58

> お手数ですが、"二つの正規分布の確率密度が同じになる点、つまり赤い線と青い線が交差する点を求める式"を教えていただけると助かります。

これについては、添付画像をご覧ください。

No.2への返答で理解しました。
それなら一致しませんね。
こちらでも計算してみたのですがまだ一致しませんが、おそらく精度の違いのためでしょうか？

どの方法が最も妥当であるかは、考え方により変わると思いますが、得られたデータがどちらの山に属するのかということを考えた場合、私なら確率密度が同じになる点337.4185を閾値とするでしょう。
例えば337.42という値が得られた場合、その確率密度から山(2)の方がより出やすいので山(2)から得られた可能性が高いと考えます。

ところで、それぞれの値の求めるにはExcelでもできるでしょうが、Rというソフトウェアを紹介しておきます。
フリーの統計ソフトでCUIのためとっつきにくいところはありますが、慣れると使いやすいです。
ダウンロードはこちらから　http://www.r-project.org/
日本語による使い方　　　　http://cse.naro.affrc.go.jp/takezawa/r-tips/r2.h …

下手な関数ですが、3通りの閾値を求める関数のソースも載せておきます。


###### ここから ######

# 確率密度が一致する点を求める関数
threshold1 <- function(m1, s1, p1, m2, s2, p2)
{
m <- m2 - m1
s <- s2 / s1
p <- p2 / p1
if (s1 == s2) # 本当はこういう比較はよくない
m1 + m2 + 2 * s1^2 * log(p / s) / m
else
(-m / s1 + s * sqrt((m / s1)^2 + 2 * (s^2 - 1) * log(s / p))) / (s^2 - 1) * s1 + m1
}

# 確率密度の和が最も小さくなる点を求める関数
threshold2 <- function(m1, s1, p1, m2, s2, p2)
{
optimize(function(x) p1 * dnorm(x, m1, s1) + p2 * dnorm(x, m2, s2), interval = c(m1, m2))
}

# 山(1)で閾値を超える確率と山(2)で閾値より低くなる確率の和がもっとも小さくなる点を求める関数
threshold3 <- function(m1, s1, p1, m2, s2, p2)
{
optimize(function(x) p1 * pnorm(x, m1, s1, lower.tail = FALSE) + p2 * pnorm(x, m2, s2), interval = c(m1, m2))
}

###### ここまで ######

###### 実行例 ######
> threshold1(236, 17.35, 0.9863, 467, 26.05, 0.0133)
[1] 337.4185
> threshold2(236, 17.35, 0.9863, 467, 26.05, 0.0133)
$minimum
[1] 338.529

$objective
[1] 1.657723e-09

> threshold3(236, 17.35, 0.9863, 467, 26.05, 0.0133)
$minimum
[1] 337.4186

$objective
[1] 6.845272e-09

No.4 回答者： quaestio 回答日時： 2009/06/26 22:37

訂正

threshold1の中で
> if (s1 == s2) # 本当はこういう比較はよくない
> m1 + m2 + 2 * s1^2 * log(p / s) / m
に間違いがありました。
以下が正しいソースです。

# 確率密度が一致する点を求める関数
threshold1 <- function(m1, s1, p1, m2, s2, p2)
{
m <- m2 - m1
s <- s2 / s1
p <- p2 / p1
if (s1 == s2) # 本当はこういう比較はよくない
(m1 + m2) / 2 - s1^2 * log(p / s) / m
else
(-m / s1 + s * sqrt((m / s1)^2 + 2 * (s^2 - 1) * log(s / p))) / (s^2 - 1) * s1 + m1
}

この回答へのお礼

本当にご丁寧にありがとうございます。
計算式は恥ずかしながら、ほとんど理解できなさそうです。
logとか平方根とかをちゃんと勉強しないとダメですね。
高校の数学で習った気がしますが、必要性を感じなかったため、
真剣にやらなかったツケが今になってくるとは・・・orz
確率密度の概念も勉強します。

ありがとうございました。

お礼日時：2009/06/27 18:33

No.2 回答者： quaestio 回答日時： 2009/06/25 23:03

> 山(1)のデータは252,973件　(98.63%)

> 山(2)のデータは 　3,419件 ( 1.33%)

1：1の分布ではなかったのですね。
上の割合で再計算してみると、
谷間の最も確率密度が低くなる点は338.5289で、二つの正規分布の確率密度が同じになる点は、谷間のところでは337.4185となるようです。
割合を考えてみてもa'=5.16056…(μ1+a'*σ1=325.5357)と一致しませんでした。
あなたの詳しい計算方法を教えていただければ、一致しない理由がわかるかもしれません。

今回、図を添付してみました。
赤い線が山(1)の正規分布、青い線が山(2)の正規分布、黒い線が山(1)及び(2)の混合分布となります。
確率密度は図でいうと縦軸の高さ(density)を意味します。
確率は横軸の線とグラフで囲まれた面積となります。
二つの正規分布の確率密度が同じになる点、つまり赤い線と青い線が交差する点は、求める式があります。必要であれば記載しますがどうしましょう？
谷間の最も確率密度が低くなる点、つまり黒い線の谷間のところで一番低くなる点は、コンピュータ等で数値計算する必要があります。

この回答へのお礼

丁寧な回答ありがとうございます。
もう一度計算したら、しきい値337.7999999で確率の総和が最も低そうです。しかし、確率密度が最小になる点：338.5289、確率密度が同じ値になる点：337.4185ともずれているようです。
確率はエクセルのNORMDIST関数で求めました。
山(1)のしきい値以上になる確率は、
=(1-NORMDIST(X,μ1,σ1,TRUE))×0.9863・・・・I
山(2)のしきい値以下になる確率は、
=NORMDIST(X,μ2,σ2,TRUE)×0.0133・・・・II
で求めました。Xはしきい値です。
しっかりした計算方法があればよいのですが、
自分はエクセル上でこつこつしきい値を変えて、
IとIIの総和が小さくなる点を探しています。
そこが、不良の発生確率が最も小さいところだと思うので・・・

お手数ですが、"二つの正規分布の確率密度が同じになる点、つまり赤い線と青い線が交差する点を求める式"を教えていただけると助かります。
何卒、よろしくお願いいたします。

お礼日時：2009/06/26 09:59

No.1 回答者： quaestio 回答日時： 2009/06/24 22:38

いくつか疑問点があります。

まず一つ目、

> 山(1)がしきい値を超える確率と山(2)がしきい値を超える確率の和が、最も小さい状態のことです。

とのことですが、確率ではなく確率密度ではないですか？
確率ですと+∞で閾値を超える確率が最も小さくなるのですが…

二つ目、その2つの山は1：1の割合と考えていいのでしょうか？

三つ目、

> 山(1)の平均値μ1：236
> 山(1)の標準偏差σ1：17.35
> 山(2)の平均値μ2：467
> 山(2)の標準偏差σ2：26.05

これらは既知の値なのでしょうか？
ヒストグラムから推定したのであれば、どう推定したのでしょうか？

とりあえず、N(236, 17.35^2)とN(467, 26.05^2)の1:1の混合分布とすると、谷間の最も確率密度が低くなる点は329.9901で、二つの正規分布の確率密度が同じになる点は、谷間のところでは329.1408となるようです。
あなたが計算されたa'=5.16056…(μ1+a'*σ1=325.5357)と一致していないのは、計算方法が異なるからでしょうか、それとも数値の丸めが原因なのでしょうか？

この回答へのお礼

回答ありがとうございます。
すみません。情報が不足していました。
全データ数は256,478件有り、内訳は以下の通りです。
山(1)のデータは252,973件　(98.63%)
山(2)のデータは 　3,419件 ( 1.33%)
山(3)のデータは 　　　77件 ( 0.03%)
山(4)のデータは 　　　　9件 ( 0.004%)
まず、山(1)と山(2)の間にしきい値を引きたくて質問しました。
平均と標準偏差はデータから計算した値です。

確率密度について調べようとしましたが、
積分の概念が難しくて理解できませんでした。
上記の山の割合と確率密度を使えば、下の値は得られますか？
a'=5.16056…(μ1+a'*σ1=325.5357)

>> 確率ですと+∞で閾値を超える確率が最も小さくなるのですが…
確かに1つの山で考えると＋∞か-∞のところが小さくなりますが、
2山の間に引くしきい値となると、2山の間のどこかになると思うのです。この考えが間違っているということですか？

よろしくお願いいたします。

お礼日時：2009/06/25 09:50