[過去ログ] 高校数学の質問スレ(医者・東大卒専用) Part438 (1002レス)
上下前次1-新
このスレッドは過去ログ倉庫に格納されています。
次スレ検索 歴削→次スレ 栞削→次スレ 過去ログメニュー
842: 05/24(土)02:35 ID:VetM3rz7(1/5) AAS
LearnBayes::beta.selectをoptimを使って算出
beta.optim <- function(x1, p1, x2, p2, verbose = TRUE) {
# -------------------------
# モーメント近似による初期値推定
# -------------------------
mu0 <- (x1 + x2) / 2 # 仮の平均
sigma2 <- ((x2 - x1) / 4)^2 # 仮の分散(中間50%幅から)
v <- mu0 * (1 - mu0) / sigma2 - 1
a0 <- mu0 * v
b0 <- (1 - mu0) * v
start <- c(a0, b0)
# -------------------------
# 最適化対象の誤差関数定義
# -------------------------
objective <- function(params) {
a <- params[1]
b <- params[2]
# pbeta による累積確率との差を2乗誤差で評価
err1 <- pbeta(x1, a, b) - p1
err2 <- pbeta(x2, a, b) - p2
return(err1^2 + err2^2)
}
# -------------------------
# 最適化(境界付き)
# -------------------------
result <- optim(
par = start,
fn = objective,
method = "L-BFGS-B",
lower = c(1e-4, 1e-4) # ベータ分布は a, b > 0
)
# -------------------------
# 結果の取り出し
# -------------------------
a_hat <- result$par[1]
b_hat <- result$par[2]
if (verbose) {
cat("推定されたパラメータ: a =", round(a_hat, 4),
", b =", round(b_hat, 4), "\n")
cat("pbeta(x1) =", round(pbeta(x1, a_hat, b_hat), 4),
"(目標:", p1, ")\n")
cat("pbeta(x2) =", round(pbeta(x2, a_hat, b_hat), 4),
"(目標:", p2, ")\n")
}
# -------------------------
# 結果を返す
# -------------------------
return(list(a = a_hat, b = b_hat, ss_value = result$value))
}
上下前次1-新書関写板覧索設栞歴
あと 160 レスあります
スレ情報 赤レス抽出 画像レス抽出 歴の未読スレ AAサムネイル
ぬこの手 ぬこTOP 0.012s