Общая процедура для нормального распределения

Предположим, у нас есть Нормальное распределение X ~ N(mu, sigma) с неизвестным средним mu и неизвестным стандартным отклонением sigma. И мы стремимся решить для mu и sigma, учитывая два квантильных уравнения:

Pr(X < q1) = alpha1
Pr(X < q2) = alpha2

Мы рассматриваем стандартизацию: Z = (X - mu) / sigma, так что

Pr(Z < (q1 - mu) / sigma) = alpha1
Pr(Z < (q2 - mu) / sigma) = alpha2

Другими словами,

(q1 - mu) / sigma = qnorm(alpha1)
(q2 - mu) / sigma = qnorm(alpha2)

RHS явно известен, и мы определяем beta1 = qnorm(alpha1), beta2 = qnorm(alpha2). Теперь это упрощается до системы двух линейных уравнений:

mu + beta1 * sigma = q1
mu + beta2 * sigma = q2

Эта система имеет матрицу коэффициентов:

1  beta1
1  beta2

с определителем beta2 - beta1. Единственная ситуация для сингулярности - beta2 = beta1. Пока система не является единственной, мы можем использовать solve, чтобы найти mu и sigma.

Подумайте, что означает ситуация сингулярности. qnorm строго монотонный для нормального распределения. Итак, beta1 = beta2 то же самое, что и alpha1 = alpha2. Но этого можно легко избежать, так как это предусмотрено вашей спецификацией, поэтому в дальнейшем я не буду проверять особенность.

Оберните вышеупомянутую функцию оценки:

est <- function(q, alpha) {
  beta <- qnorm(alpha)
  setNames(solve(cbind(1, beta), q), c("mu", "sigma"))
  }

Давайте проведем тест:

x <- est(c(158, 168), c(0.025, 0.975))
#        mu      sigma 
#163.000000   2.551067 

## verification
qnorm(c(0.025, 0.975), x[1], x[2])
# [1] 158 168

Мы также можем сделать что-нибудь произвольное:

x <- est(c(1, 5), c(0.1, 0.4))
#      mu    sigma 
#5.985590 3.890277 

## verification
qnorm(c(0.1, 0.4), x[1], x[2])
# [1] 1 5