P4 Pembangkitan Bilangan Acak Normal

Indri ramdani (G1401201036)

2024-02-12

Untuk mensimulasikan bilangan acak normal, berbagai teknik dan distribusi dapat digunakan. Berikut beberapa metode yang dapat digunakan untuk membangkitkan peubah acak normal :

  1. Dalil Limit Pusat

  2. Metode Box-Muller

  3. Metode Polar-Marsaglia

1. Dalil Limit Pusat

-Untuk sembarang peubah acak Xi yang saling bebas dan identik, \(Y=\sum_{i=1}^{n}X_i\) akan menyebar normal untuk n cukup besar.

-Untuk \(𝑈_1,𝑈_2, … ,U_n\) ~ \(Uniform (0,1)\) –> \(X=\sum_{i=1}^{n}U_i\) ~ \(Normal (μ,σ^2)\)

Diketahui:

\(E(U_i)=\frac{a+b}{2}=\frac{0+1}{2}=\frac{1}{2}\)

\(Var(U_i)=\frac{(b-a)^2}{12}=\frac{(1-0)^2}{12}=\frac{1}{12}\)

Dengan \(E(X)=μ\) dan \(Var(X)=σ^2\)

Jadi untuk menduga \(X=\sum_{i=1}^{n}U_i-k\)

\(Var(X)=Var(\sum_{i=1}^{n}U_i-k)=\sum_{i=1}^{n}Var(U_i)=\frac{n}{12}=σ^2\)

\(E(N)=E(\sum_{i=1}^{n}U_i-k)=\sum_{i=1}^{n}E(U_i)-k=\frac{n}{2}-k=μ\)

Dari persamaan diatas diperoleh persamaan:

\(n=12σ^2\) dan \(k=\frac{n}{2}-μ\)

#1. Membangkitkan Normal (0,1)
i  = 1000
mu = 0 
sigma2 = 1
n  = 12*sigma2
k  = (n/2) - mu 

set.seed(1234)
U  = runif(n*i) #membangkitkan U ~ uniform (0,1)
Um = matrix (U,i)
X  = apply(Um,1,sum)-k

result = data.frame( parameter = c("miu","sigma2"), 
                    Nilai= c(mean(X),var(X)))
result
hist(X,col = "Light Coral", main = paste ("Histogram of", "X ~ Normal (0,1)"))

#2. Membangkitkan Normal (4,5)
i  = 1000
mu = 4 
sigma2 = 5
n  = 12*sigma2
k  = (n/2) - mu 

set.seed(1234)
U  = runif(n*i) #membangkitkan U ~ uniform (4,5)
Um = matrix (U,i)
X  = apply(Um,1,sum)-k

result2 = data.frame( parameter = c("miu","sigma2"), 
                    Nilai= c(mean(X),var(X)))
result2
hist(X,col = "Light Coral", main = paste ("Histogram of", "X ~ Normal (4,5)"))

2. Metode Box-Muller

Algoritma Box-Muller adalah metode untuk menghasilkan pasangan bilangan acak yang independen, standar, dan terdistribusi normal dari sumber bilangan acak terdistribusi seragam. Algoritma ini melibatkan transformasi dua set bilangan acak dengan distribusi seragam menjadi dua set bilangan acak dengan distribusi normal. Transformasi ini didasarkan pada bentuk polar, yang menggunakan fungsi trigonometri untuk memetakan bilangan acak terdistribusi seragam menjadi bilangan acak terdistribusi normal standar. Algoritma ini banyak digunakan dalam pengambilan sampel statistik dan simulasi Monte Carlo karena kesederhanaan dan efisiensinya.

Untuk dua peubah acak menyebar \(U(0,1)\) yang saling bebas, \(U_1\) dan \(U_2\) , maka : \(N_1=(-2 log_e U_1)^\frac{1}{2}cos(2\pi U_2)\)

\(N_2=(-2 log_e U_1)^\frac{1}{2}sin(2\pi U_2)\)

dengan N1 dan N2 masing-masing menyebar Normal (0,1)

Fungsi di atas dapat dituliskan menjadi :

\(N_1= R cos Θ\)

\(N_2=R sin Θ\)

dengan

\(R=(-2 log_e U_1)^\frac{1}{2}\)~\(Eksponensial(\frac{1}{2})\)

\(Θ=2πU_2\) ~ \(Uniform(0,2π)\)

\((N_1 , N_2)\) ↔︎\((R,π)\)

Koordinat Cartesius ↔︎ Koordinat Polar

i  = 1000
set.seed(1234)
U1 = runif (i)
U2 = runif (i)
R  = sqrt (-2*log(U1))
Theta = 2*pi*U2
N1 = R*cos(Theta)
N2 = R*sin(Theta)

result2 = data.frame( parameter = c("miu N1","sigma2 N1","miu N2","sigma2 N2"), 
                    Nilai= c(mean(N1),var(N1),mean(N2),var(N2)))
result2
par(mfrow = c(1,2))
hist(N1,col = "Light Coral", main = paste ("Histogram of", "N1 ~ Normal (0,1)"))
hist(N2,col = "Light Coral", main = paste ("Histogram of", "N2 ~ Normal (0,1)"))

3. Metode Polar-Marsaglia

Metode pembangkit bilangan acak pseudo untuk menghasilkan sepasang independen peubah acak normal baku.Perbedaan dengan metode Box-Muller yaitu terletak pada penggunaan komputasi sinus dan cosinus. Pada Metode Polar Marsaglia tidak digunakan komputasi sinus dan cosinus.

Dari persamaan metode Box-Muller: \(N_1=(-2 log_e U_1)^\frac{1}{2}cos(2\pi U_2)\)

\(N_2=(-2 log_e U_1)^\frac{1}{2}sin(2\pi U_2)\)

menjadi,

\(N_1=(-2 log_e R^2)^\frac{1}{2}V1(V1^2+V2^2)^\frac{-1}{2}\)

\(N_2=(-2 log_e R^2)^\frac{1}{2}V2(V1^2+V2^2)^\frac{-1}{2}\)

menjadi,

\(N_1=(-2 log_e (V1^2+V2^2))^\frac{1}{2}V1(V1^2+V2^2)^\frac{-1}{2}\)

\(N_2=(-2 log_e (V1^2+V2^2))^\frac{1}{2}V2(V1^2+V2^2)^\frac{-1}{2}\)

sehingga, \(N_1=V1{(-2 log_e W)/W}^\frac{1}{2}\)

\(N_2=V2{(-2 log_e W)/W}^\frac{1}{2}\)

dengan \(W=V1^2+V2^2\)

i  = 1000
set.seed(1234)
U1 = runif (i)
U2 = runif (i)
V1 = 2*U1 - 1
V2 = 2*U2 - 1
W  = V1^2 + V2^2

W1 = W [W<1]
V1 = V1 [W<1]
V2 = V2 [W<1]

i_1 = i-length(W1)
while (i_1>0) {
  U_1 = runif(i_1)
  U_2 = runif(i_1)
  V_1 = 2*U_1-1
  V_2 = 2*U_2-1
  W_1 = V_1^2 + V_2^2
  W1  = c(W1,W_1[W_1<1])
  V1  = c(V1,V_1[W_1<1])
  V2  = c(V2,V_2[W_1<1])
  i_1 = i-length(W1)
}

N1 = V1*sqrt(-2*log(W1)/W1)
N2 = V2*sqrt(-2*log(W1)/W1)

result4 = data.frame( parameter = c("miu N1","sigma2 N1","miu N2","sigma2 N2"), 
                    Nilai= c(mean(N1),var(N1),mean(N2),var(N2)))
result4
par(mfrow = c(1,2))
hist(N1,col = "Light Coral", main = paste ("Histogram of", "N1 ~ Normal (0,1)"))
hist(N2,col = "Light Coral", main = paste ("Histogram of", "N2 ~ Normal (0,1)"))