Dapatkah saya merekonstruksi distribusi normal dari ukuran sampel, dan nilai min dan maks? Saya bisa menggunakan titik tengah untuk proksi mean

14

Saya tahu ini mungkin sedikit mengikat, secara statistik, tetapi ini adalah masalah saya.

Saya memiliki banyak data rentang, yaitu ukuran minimum, maksimum, dan sampel suatu variabel. Untuk beberapa data ini saya juga memiliki nilai rata-rata, tetapi tidak banyak. Saya ingin membandingkan rentang ini satu sama lain untuk mengukur variabilitas masing-masing rentang, dan juga membandingkan cara. Saya punya alasan kuat untuk mengasumsikan bahwa distribusi simetris di sekitar rata-rata, dan bahwa data akan memiliki distribusi Gaussian. Untuk alasan ini saya berpikir saya dapat membenarkan menggunakan titik tengah distribusi sebagai proksi untuk mean, ketika tidak ada.

Yang ingin saya lakukan adalah merekonstruksi distribusi untuk setiap rentang, dan kemudian menggunakannya untuk memberikan standar deviasi atau kesalahan standar untuk distribusi itu. Satu-satunya informasi yang saya miliki adalah maks dan min yang diamati dari sampel, dan titik tengah sebagai proksi untuk mean.

Dengan cara ini saya berharap dapat menghitung cara tertimbang untuk setiap kelompok, dan juga menghitung koefisien variasi untuk masing-masing kelompok juga, berdasarkan pada rentang data yang saya miliki dan asumsi saya (dari distribusi simetris dan normal).

Saya berencana untuk menggunakan R untuk melakukan ini, jadi bantuan kode apa pun akan dihargai juga.

green_thinlake
sumber
2
Saya bertanya-tanya mengapa Anda mengatakan Anda memiliki data untuk nilai minimum & maksimum & maksimum; kemudian Anda memiliki informasi tentang hanya yang diharapkan minimum & maksimum. Mana itu - diamati atau diharapkan?
Scortchi
Maaf, itu kesalahan saya. Data maksimum dan minimum diamati (diukur dari objek kehidupan nyata). Saya telah mengubah pos.
green_thinlake

Jawaban:

11

Fungsi distribusi kumulatif bersama untuk minimum & maksimum x ( n ) untuk sampel n dari distribusi Gaussian dengan mean μ & standar deviasi σ adalahx(1)x(n)nμσ

F(x(1),x(n);μ,σ)=Pr(X(1)<x(1),X(n)<x(n))=Pr(X(n)<x(n))Pr(X(1)>x(1),X(n)<x(n)=Φ(x(n)μσ)n[Φ(x(n)μσ)Φ(x(1)μσ)]n

di mana adalah CDF Gaussian standar. Diferensiasi sehubungan dengan x ( 1 ) & x ( n ) memberikan fungsi densitas probabilitas gabunganΦ()x(1)x(n)

f(x(1),x(n);μ,σ)=n(n1)[Φ(x(n)μσ)Φ(x(1)μσ)]n2ϕ(x(n)μσ)ϕ(x(1)μσ)1σ2

di mana adalah standar Gaussian PDF. Mengambil log & menjatuhkan istilah yang tidak mengandung parameter memberikan fungsi log-likelihoodϕ()

(μ,σ;x(1),x(n))=(n2)log[Φ(x(n)μσ)Φ(x(1)μσ)]+logϕ(x(n)μσ)+logϕ(x(1)μσ)2logσ

Ini tidak terlihat sangat penurut tapi mudah untuk melihat bahwa itu dimaksimalkan apapun nilai oleh pengaturan μ = μ = x ( n ) + x ( 1 )σ , yaitu titik tengah - istilah pertama dimaksimalkan ketika argumen dari satu CDF adalah negatif dari argumen yang lain; istilah kedua & ketiga mewakili kemungkinan gabungan dari dua varian normal independen.μ=μ^=x(n)+x(1)2

Mengganti μ ke dalam log-likelihood & menulis r = x ( n ) - x ( 1 ) memberikan ( σ ; x ( 1 ) , x ( n ) , μ ) = ( n - 2 ) log [ 1 - 2 Φ ( - rμ^r=x(n)x(1)

(σ;x(1),x(n),μ^)=(n2)log[12Φ(r2σ)]r24σ22logσ

Ungkapan ini harus dimaksimalkan numerik (misalnya dengan optimizedari R statpaket) untuk menemukan σ . (Ternyata σ = k ( n ) r , di mana k adalah konstanta tergantung hanya pada n -mungkin orang yang lebih matematis gesit daripada aku bisa menunjukkan mengapa.)σ^σ^=k(n)rkn

Perkiraan tidak ada gunanya tanpa ukuran presisi yang menyertainya. Informasi Fisher yang diamati dapat dievaluasi secara numerik (misalnya dengan hessiandari numDerivpaket R ) & digunakan untuk menghitung perkiraan kesalahan standar:

I(σ)=-2(σ; μ )

I(μ)=2(μ;σ^)(μ)2|μ=μ^
I(σ)=2(σ;μ^)(σ)2|σ=σ^

Akan menarik untuk membandingkan kemungkinan & metode-saat estimasi untuk dalam hal bias (apakah MLE konsisten?), Varians, & mean-square error. Ada juga masalah estimasi untuk kelompok-kelompok tersebut di mana rata-rata sampel diketahui selain minimum & maksimum.σ

Scortchi - Reinstate Monica
sumber
1
2log(r)σ/rnσ/rnk(n)σ^=k(n)rKisaran studentized .
whuber
@whuber: Terima kasih! Tampak jelas dengan melihat ke belakang. Saya akan memasukkannya ke dalam jawabannya.
Scortchi
1

μσR=x(n)x(1)99.7

μ+3σx(n)

μ3σx(1)

Mengurangkan yang kedua dari yang pertama kita dapatkan

6σx(n)x(1)=R
σ^=16(x¯(n)x¯(1))

Memiliki nilai untuk mean dan untuk standar deviasi sepenuhnya mencirikan distribusi normal.

Alecos Papadopoulos
sumber
3
Itu bukan perkiraan dekat untuk kecil n atau hasil asimptotik untuk yang besar n.
Scortchi
1
@Stortchi Yah, saya tidak mengatakan bahwa ini adalah perkiraan yang bagus -tapi saya percaya bahwa selalu baik untuk memiliki solusi yang mudah diimplementasikan, bahkan sangat kasar, untuk mendapatkan pemahaman kuantitatif dari masalah yang ada, selain itu pendekatan yang canggih dan efisien seperti misalnya yang diuraikan dalam jawaban lain untuk pertanyaan ini.
Alecos Papadopoulos
Saya tidak akan menyerah pada "harapan kisaran sampel ternyata sekitar 6 kali standar deviasi untuk nilai nmulai dari 200 hingga 1000 ". Tetapi apakah saya kehilangan sesuatu yang halus dalam derivasi Anda, atau tidakkah itu bekerja dengan baik untuk membenarkan pembagian rentang dengan angka berapa pun?
Scortchi - Reinstate Monica
@Scortchi Well, the spirit of the approach is "if we expect almost all realizations to fall within 6 sigmas, then it is reasonable to expect that the extreme realizations will be near the border" -that's all there is to it, really. Perhaps I am too used to operate under extremely incomplete information, and obliged to say something quantitative about it... :)
Alecos Papadopoulos
4
I could reply that even more observations would fall within 10σ of the mean, giving a better estimate σ^=R10. I shan't because it's nonsense. Any number over 1.13 will be a rough estimate for some value of n.
Scortchi
1

Sangat mudah untuk mendapatkan fungsi distribusi maksimum dari distribusi normal (lihat "P.max.norm" dalam kode). Dari itu (dengan beberapa kalkulus) Anda bisa mendapatkan fungsi kuantil (lihat "Q.max.norm").

Dengan menggunakan "Q.max.norm" dan "Q.min.norm" Anda bisa mendapatkan median rentang yang terkait dengan N. Menggunakan ide yang disajikan oleh Alecos Papadopoulos (dalam jawaban sebelumnya) Anda dapat menghitung sd.

Coba ini:

N = 100000    # the size of the sample

# Probability function given q and N
P.max.norm <- function(q, N=1, mean=0, sd=1){
    pnorm(q,mean,sd)^N
} 
# Quantile functions given p and N
Q.max.norm <- function(p, N=1, mean=0, sd=1){
    qnorm(p^(1/N),mean,sd)
} 
Q.min.norm <- function(p, N=1, mean=0, sd=1){
    mean-(Q.max.norm(p, N=N, mean=mean, sd=sd)-mean)
} 

### lets test it (takes some time)
Q.max.norm(0.5, N=N)  # The median on the maximum
Q.min.norm(0.5, N=N)  # The median on the minimum

iter = 100
median(replicate(iter, max(rnorm(N))))
median(replicate(iter, min(rnorm(N))))
# it is quite OK

### Lets try to get estimations
true_mean = -3
true_sd = 2
N = 100000

x = rnorm(N, true_mean, true_sd)  # simulation
x.vec = range(x)                  # observations

# estimation
est_mean = mean(x.vec)
est_sd = diff(x.vec)/(Q.max.norm(0.5, N=N)-Q.min.norm(0.5, N=N))

c(true_mean, true_sd)
c(est_mean, est_sd)

# Quite good, but only for large N
# -3  2
# -3.252606  1.981593
Vyga
sumber
2
Melanjutkan pendekatan ini, E(R)=σ1(1Φ(x))nΦ(x)ndx=σd2(n), where R is the range & Φ() the standard normal cumulative distribution function. You can find tabulated values of d2 for small n in the statistical process control literature, numerically evaluate the integral, or simulate for your n.
Scortchi - Reinstate Monica