Musisz dopasować te skumulowane dane do jakiegoś modelu dystrybucyjnego, ponieważ jest to jedyny sposób ekstrapolacji na górny kwartyl.
Wzór
Z definicji taki model daje funkcja cadlag rosnąca od 0 do 1 . Prawdopodobieństwo, które przypisuje do dowolnego przedziału ( a , b ], to F ( b ) - F ( a ) . Aby dopasować, musisz założyć rodzinę możliwych funkcji indeksowanych parametrem (wektorowym) θ , { F θ } Zakładając, że próba podsumowuje zbiór osób wybranych losowo i niezależnie od populacji opisanej przez niektóre konkretne (ale nieznane) F θF01(a,b]F(b)−F(a)θ{Fθ}FθPrawdopodobieństwo próbki (lub prawdopodobieństwa , ) jest iloczynem z poszczególnych prawdopodobieństwa. W tym przykładzie byłby równyL
L(θ)=(Fθ(8)−Fθ(6))51(Fθ(10)−Fθ(8))65⋯(Fθ(∞)−Fθ(16))182
51Fθ(8)−Fθ(6)65Fθ(10)−Fθ(8)
Dopasowanie modelu do danych
θLL
θ=(μ,σ)
F(μ,σ)(x)=12π−−√∫(log(x)−μ)/σ−∞exp(−t2/2)dt.
LR
log(L(θ))log(L)Llog(L)
logL <- function(thresh, pop, mu, sigma) {
l <- function(x1, x2) ifelse(is.na(x2), 1, pnorm(log(x2), mean=mu, sd=sigma))
- pnorm(log(x1), mean=mu, sd=sigma)
logl <- function(n, x1, x2) n * log(l(x1, x2))
sum(mapply(logl, pop, thresh, c(thresh[-1], NA)))
}
thresh <- c(6,8,10,12,14,16)
pop <- c(51,65,68,82,78,182)
fit <- optim(c(0,1), function(theta) -logL(thresh, pop, theta[1], theta[2]))
θ=(μ,σ)=(2.620945,0.379682)fit$par
Sprawdzanie założeń modelu
F
predict <- function(a, b, mu, sigma, n) {
n * ( ifelse(is.na(b), 1, pnorm(log(b), mean=mu, sd=sigma))
- pnorm(log(a), mean=mu, sd=sigma) )
Jest on stosowany do danych w celu uzyskania dopasowanych lub „przewidywanych” populacji bin:
pred <- mapply(function(a,b) predict(a,b,fit$par[1], fit$par[2], sum(pop)),
thresh, c(thresh[-1], NA))
Możemy narysować histogramy danych i prognozę, aby porównać je wizualnie, pokazane w pierwszym rzędzie tych wykresów:
Aby je porównać, możemy obliczyć statystyki chi-kwadrat. Jest to zwykle określane jako rozkład chi-kwadrat w celu oceny istotności :
chisq <- sum((pred-pop)^2 / pred)
df <- length(pop) - 2
pchisq(chisq, df, lower.tail=FALSE)
0.00876−8630.40
Wykorzystanie dopasowania do oszacowania kwantyli
63(μ,σ)(2.620334,0.405454)F75th
exp(qnorm(.75, mean=fit$par[1], sd=fit$par[2]))
18.066317.76
Te procedury i ten kod można stosować ogólnie. Teorię maksymalnego prawdopodobieństwa można dalej wykorzystać do obliczenia przedziału ufności wokół trzeciego kwartylu, jeśli jest to interesujące.