ベータ分布に従うデータが与えられたとき、2つのパラメータ \(\alpha\), \(\beta\) を簡単に求めたい。
理論的には、ベータ分布に従う確率変数 \(X\) の期待値 \({\rm E}[X]\) および分散 \({\rm var}[X]\) は次のようになる。
\[ \begin{eqnarray} {\rm E}[X] &=& \frac{\alpha}{\alpha + \beta} \\ {\rm var}[X] &=& \frac{\alpha\beta}{(\alpha+\beta)^2(\alpha+\beta+1)} \end{eqnarray} \]
ここで、
\[ \begin{eqnarray} 1 - {\rm E}[X] &=& \frac{\beta}{\alpha + \beta} \\ \alpha + \beta &=& \frac{\alpha}{{\rm E}[X]} \end{eqnarray} \]
を用いると、
\[ \begin{eqnarray} {\rm var}[X] &=& \frac{\alpha}{\alpha + \beta}\frac{\beta}{\alpha + \beta}\frac{1}{\alpha+\beta+1} \\ &=& {\rm E}[X](1-{\rm E}[X])\frac{1}{\frac{\alpha}{{\rm E}[X]}+1} \\ \frac{\alpha}{{\rm E}[X]}+1 &=& {\rm E}[X](1-{\rm E}[X])\frac{1}{{\rm var}[X]} \\ \alpha+{\rm E}[X] &=& {\rm E}[X]^2(1-{\rm E}[X])\frac{1}{{\rm var}[X]} \\ \alpha &=& {\rm E}[X]^2(1-{\rm E}[X])\frac{1}{{\rm var}[X]} - {\rm E}[X]\\ \end{eqnarray} \]
となり、期待値と分散で \(\alpha\) を表すことができる。
\(\beta\) については、
\[ \begin{eqnarray} \alpha + \beta &=& \frac{\alpha}{{\rm E}[X]} \\ \beta &=& \frac{\alpha}{{\rm E}[X]} - \alpha \end{eqnarray} \]
となる。
したがって、次のような関数によりベータ分布のパラメータを推定することができる。
estimate_beta_params <- function(values) {
mean_ <- mean(values)
var_ <- var(values)
alpha_ <- mean_^2 * (1-mean_) / var_ - mean_
beta_ <- alpha_ / mean_ - alpha_
c(alpha=alpha_, beta=beta_)
}実際にやってみる。
test_estimation <- function(alpha, beta, n = 1000) {
data <- rbeta(n, alpha, beta)
estimate_beta_params(data)
}
test_estimation(alpha = 1, beta = 1)## alpha beta
## 1.009365 1.050909test_estimation(alpha = 1, beta = 2)## alpha beta
## 1.032552 1.962597test_estimation(alpha = 10, beta = 2)## alpha beta
## 9.877167 1.969278Enjoy!