Jumlah variabel acak terpotong normal

8

Misalkan saya punya n variabel acak normal independen

X1N(μ1,σ12)X2N(μ2,σ22)XnN(μn,σn2)

dan . Bagaimana saya menandai kerapatan jika distribusi setiap masing -masing terpotong ke dalam ? Dengan kata lain, saya mengambil sampel dari distribusi normal independen, membuang sampel yang tidak berada dalam dari setiap rata-rata, dan menjumlahkannya.Y=X1+X2++XnYXsaya(μsaya-2σsaya,μsaya+2σsaya)n2σsaya

Saat ini, saya melakukan ini dengan kode R di bawah ini:

x_mu <- c(12, 18, 7)
x_sd <- c(1.5, 2, 0.8)
a <- x_mu - 2 * x_sd
b <- x_mu + 2 * x_sd

samples <- sapply(1:3, function(i) {
  return(rtruncnorm(100000, a[i], b[i], x_mu[i], x_sd[i]))
})

y <- rowSums(samples)

Apakah ada metode untuk menghasilkan kerapatan secara langsung?Y

Devin
sumber
2
Pertanyaan Anda menyiratkan Anda tahu semuaσsaya. Apakah itu benar-benar terjadi atau apakah Anda memperkirakannya ? Ada perbedaan besar! Karena penasaran, mengapa Anda membuang data seperti itu? Tergantung pada tujuan Anda, saya curiga ada (banyak) prosedur yang lebih baik.
whuber
Saya tahu semua sarana dan SD untuk data saya, ya.
Devin
7
Saya percaya bahwa Anda dapat menggambarkannya sebagai "berantakan". Makalah ini, jstor.org/stable/2236545 , meneliti masalah ini, dengan kekakuan yang lebih ilmiah.
Alecos Papadopoulos
2
Di luar perkiraan melalui CLT, ini relatif rumit. Saya kira jikancukup kecil Anda dapat mencoba konvolusi numerik.
Glen_b -Reinstate Monica
2
@Silverfish Bergantung pada implementasi, platform, dan seberapa halus sebuah grid dapat ditoleransi, ratusan harus baik-baik saja (mungkin lebih); selain kecepatan, dengan persyaratan yang cukup, Anda harus lebih berhati-hati tentang detail implementasi atau sejumlah masalah numerik dapat mulai muncul.
Glen_b -Reinstate Monica

Jawaban:

2

Anda bisa menggunakan pendekatan dengan metode saddlepoint, untuk jumlah normals terpotong. Saya tidak akan memberikan detailnya sekarang, Anda dapat melihat jawaban saya untuk jumlah umum distribusi Gamma untuk petunjuk. Yang kita butuhkan adalah menemukan fungsi penghasil momen untuk normal terpotong, yang mudah. Saya akan melakukannya di sini untuk standar normal terpotong±2, yang memiliki kepadatan

f(x)={1Cϕ(x),|x|20,|x|>2
dimana C=Φ(2)-Φ(-2) sini ϕ(x),Φ(x) adalah kepadatan dan cdf untuk standar normal, masing-masing.

Fungsi pembangkit momen dapat dihitung sebagai

M.(t)=EetX=1C-22etxϕ(x)dx=1Ce12t2[Φ(2-t)-Φ(-2-t)]
dan kemudian kita bisa menggunakan pendekatan saddlepoint.
kjetil b halvorsen
sumber
-3

Saya ingin tahu mengapa, tapi ya, ada cara sederhana untuk menghasilkan pdf dari jumlah distribusi ini:

## install.packages("truncnorm")
## install.packages("caTools")
library(truncnorm)

x.mu <- c(12, 18, 7)
x.sd <- c(1.5, 2, 0.8)
x.a <- x.mu - 2*x.sd
x.b <- x.mu + 2*x.sd

dmulti <- function(x, a, b, mu, sd)
  rowSums(
    sapply(1:length(mu),
           function(idx)
             dtruncnorm(x, a=a[idx], b=b[idx], mean=mu[idx], sd=sd[idx])))/length(mu)
pmulti <- function(q, a, b, mu, sd)
  rowSums(
    sapply(1:length(mu),
           function(idx)
             ptruncnorm(q, a=a[idx], b=b[idx], mean=mu[idx], sd=sd[idx])))/length(mu)

pointrange <- range(c(x.a, x.b))
pointseq <- seq(pointrange[1], pointrange[2], length.out=100)
## Plot the probability density function
plot(pointseq, dmulti(pointseq, x.a, x.b, x.mu, x.sd),
     type="l")

## Plot the cumulative distribution function
plot(pointseq, pmulti(pointseq, x.a, x.b, x.mu, x.sd),
     type="l")
Bill Denney
sumber
Jika saya membaca kode ini dengan benar, Anda tampaknya menerapkan sesuatu seperti campuran daripada penjumlahan. Plot yang dihasilkan kode ini sangat salah. Ini bahkan bukan fungsi kepadatan probabilitas yang valid!
whuber
@whuber, terima kasih atas tangkapannya. Saya menormalkan pdf dan menambahkan cdf.
Bill Denney
3
Terima kasih. Namun, kesalahan dasar tetap ada: Anda menghitung distribusi campuran daripada jumlahnya.
whuber