Pembangkitan Bilangan Acak

Pendahuluan

Bilangan acak (random) adalah bilangan yang kemunculannya terjadi secara acak dan diharapkan memenuhi sebaran statistik tertentu (pdf/pmf, cdf). Pembangkitan bilangan acak merupakan alat yang diperlukan dalam komputasi statistik, umumnya untuk simulasi. Semua metode pembangkitan bilangan acak tergantung dari pembangkitan bilangan acak uniform.

Pseudo-Random Number Generator merupakan algoritma untuk mendapatkan bilangan secara acak dan tidak dapat diprediksi kemunculan selanjutnya. Bilangan acak yang dibangkitkan oleh komputer disebut bilangan acak semu atau pseudo random number. Disebut pseudo karena kemunculan bilangan tersebut sangat bergantung padaseed atau nilai awal. Jika nilaiseed yang digunakan adalah bilangan yang benar-benar acak, maka bilangan yang dihasilkan juga akan acak sempurna (tidak memiliki pola). Nilai seed yang sama akan menghasilkan bilangan acak yang sama pula.

Fungsi sebaran statistik di R :

  • dnamasebaran : density or probability mass function (pmf)

  • pnamasebaran : cumulative distribution function (cdf)

  • qnamasebaran : quantiles

  • rnamasebaran : random number generation

Berikut adalah beberapa sebaran statistik dan namanya dalam R :

Distribution Name
Beta beta
Binomial binom
Chisquare chisq
Exponential exp
F f
Gamma gamma
Gemotric geom
Lognormal lnorm
Logistic logis
Negative Binomial nbinom
Normal (Gaussian norm
Poisson pois
Student's t t
Uniform unif
Weibull weibull

Latihan

#1. Membangkitkan 15 bilangan acak dari sebaran normal baku (0,1)
rnorm(15)
 [1]  1.75118872 -0.05664911  0.99862382  0.30116267 -1.75655153 -0.45812103
 [7] -0.03765437 -0.11450039 -1.42671114 -0.72084220 -0.91716536  1.14438401
[13] -0.07784343 -0.40088263  1.50461925
#2. Membangkitkan 15 bilangan acak dari sebaran normal (10,2)
rnorm(15, mean = 10, sd = 2)
 [1]  9.863149 10.386806 11.104595 10.013181 10.207742 12.971959  9.470645
 [8]  9.573681  7.291151  7.071614 10.452450  8.815495  8.922989  8.902545
[15] 11.878963
#3. Membangkitkan 15 bilangan acak dari sebaran seragam (0,1)
runif(15) #runif(n, min = 0, max = 1)
 [1] 0.54847951 0.66400776 0.36554399 0.41851871 0.47909992 0.88658965
 [7] 0.95513099 0.29239513 0.11537893 0.03455096 0.86056947 0.50079681
[13] 0.41783628 0.04045636 0.08033844
#4. Membangkitkan 15 bilangan acak dari sebaran seragam (5,20)
runif(15,5,20) 
 [1] 14.453216  9.835738 19.855227 16.758759  8.173334 18.642049  5.259067
 [8] 19.291619  7.538851 19.582394 19.192940 13.110540  7.752578  6.133072
[15]  7.083736
#5. Plot Bilangan Acak Sebaran Normal
x = 1:15
y = rnorm (15,10,2) # Membangkitkan 15 bilangan acak dari sebaran normal (10,2)
plot(x,y,type="p",xlab="Sumbu x",ylab="Sumbu y",
     main="Bilangan Acak Normal",col=2,pch=16)

Teknik Pembangkitan Bil. Acak

Apabila fungsi pembangkit bilangan acak belum ada karena sebarannya belum ada, teknik yang umum digunakan untuk membangkitkan bilangan acak adalah sebagai berikut :

  • Inverse-transform method -> dapat digunakan untuk peubah acak yang berdistribusi kontinu maupun diskrit
  • Acceptance-rejection method -> dapat digunakan jika sebuah distribusi memiliki fungsi massa dan metode lain tidak efisien digunakan
  • Other Special techniques :
    • Direct Transformation -> menggunakan teori dari sebaran (sebaran dari suatu peubah acak merupakan sebaran dari peubah acak lain)
    • Convolution -> untuk beberapa distribusi yang mungkin dapat dinyatakan dalam jumlah \((X\) = \(Y_1\) + \(Y_2\) + ... + \(Y_m)\)
    • dan metode yang lain

1. Inverse-Transform Method

#Sebaran Kontinu

Langkah :

  • Untuk fungsi cdf : f = F(x)
  • Bangkitkan u dari U(0,1)
  • Maka, X = F^(-1)(u)

Latihan

Diketahui : \[ pdf : f(x) = 3x^2 , 0 < x < 1\] \[F(x) = x^3 , 0 < x < 1 \] \[F^{-1}(u) = u^{1/3}\]

n = 1000 # misalkan membangkit 1000 bilangan acak
u = runif(n)
x = u^(1/3)
hist(x, prob = TRUE, col = "Light Coral")
y = seq(0,1,.01)
lines(y, 3*y^2)

Diketahui : \[ X \sim Exp (\lambda) , x > 0 \] \[F(x) = 1 - e^{- \lambda x}\] \[F^{-1}(u) = -\frac{1}{\lambda} {log(1-u)}\]

set.seed(1234)
exp = function (n,lambda){
  u = runif(n)
  x = -log(u)/lambda
  return (x)
}

#metode invers
exp1 = exp(1000,5)

#fungsi bawaan R
exp2 = rexp(1000, rate = 3)

par(mfrow = c(1,2))

hist(exp1, prob = TRUE, main = " Exp : Metode Invers-Transform" , col = "Light Coral")
hist(exp2, prob = TRUE, main = " Exp : Fungsi Rexp" , col = "Light Coral")

#Sebaran Diskrit

Jika X ~ p.a diskrit dan \(... < x_{i-1} < x_i < x_{i+1}< ...\) adalah titik tidak kontinu dari F(x), maka transformasi inversnya adalah : \(F^{-1}(u) = x_i\) dimana \(F(x_{i-1}) < u < F(x_i)\)

Langkah :

  • Bangkitkan U(0,1)
  • Tentukan \(x_i\) dimana \(F(x_{i-1}) < u < F(x_i)\)

Latihan

#X ~ Bernoulli (0.4) 
#F(0) = f(0) = 1-p dan F(1) = 1
#F^(-1)(u) = 1 jika u > 0.6
#F^(-1)(u) = 0 jika u <= 0.6

n = 1000
p = 0.4
u = runif(n)
x = as.integer (u > 0.6)
barplot (table(x) , col = "Thistle")

2. Acceptance-Rejection Method

Pada metode sebelumnya, bilangan acak dapat dibangkitkan dari suatu sebaran jika mana sebaran tersebut memiliki invers dari sebaran kumulatifnya. Jika sebarang kumulatif tidak bisa diperoleh atau sulit diperoleh, maka dapat menggunakan metode Acceptance-Rejection Method.

Misalkan X dan Y adalah peubah acak dengan pdf/pmf f dan g dan terdapat konstanta c sehingga : \(\frac{f(t)}{g(t)} \leq c\) untuk semua \(t:f(t)>0\)

Langkah :

  1. Tetapkan peubah acak Y dengan density g yang memenuhi \(\frac{f(t)}{g(t)} \leq c\) untuk semua \(t:f(t)>0\)

  2. Untuk setiap satu bilangan acak :

    • Bangkitkan y acak dari sebaran dengan density g

    • Bangkitkan u acak dari sebaran uniform (0,1)

    • Jika u < \(\frac{f(y)}{cg(y)}\) terima y dan x = y ; selainnya tolak y dan ulangi langkah 2 awal

Latihan

Diketahui : X ~ Beta (= 2 , = 2) \[pdf : f(x) = 6x(1-x) , 0 < x < 1\]

#1. Ambil g(x) dari sebaran uniform (0,1)
#2. Maka f(x)/g(X) <= 6 untuk 0 < x < 1
#3. Sebuah x acak dari g(x) diterima jika : f(x)/[c g(x)] = 6x(1-x) / [6(1)] = x(1-x) > u

n = 1000
k = 0
j = 0
y = numeric(n)
while (k < n) {
u = runif(1)
j = j+1
x = runif(1)
if (x*(1-x) > u) {
k = k+1 # terima x
y[k] = x
}
}

hist(y, prob = TRUE , col = "Misty Rose")

3. Direct Transformation Method

Beberapa transformasi dari transformasi invers sebaran dapat digunakan untuk membangkitkan bilangan acak.

Berikut beberapa contoh peubah acak transformasi :

https://www.ynufe.edu.cn/pub/tsxy/jhzhao/teach/CS/notes/lec3.pdf

https://www.ynufe.edu.cn/pub/tsxy/jhzhao/teach/CS/notes/lec3.pdf

Latihan

Diketahui :

Jika \[U \sim Gamma (r,\lambda) \space \space dan \space \space V \sim Gamma (s,\lambda)\] kemudian \[ X = \frac{U}{U + V} \sim Beta (r,s)\]

#1. Bangkitkan u dari Gamma (a,1)
#2. Bangkitkan v dari Gamma (b,1)
#3. Tentukan x = u/(u+v)

n = 1000
a = 3
b = 2
u = rgamma (n , shape =a , rate =1);
v = rgamma (n , shape =b , rate =1);
x = u/( u + v )

q = qbeta ( ppoints ( n ) , a , b )
qqplot (q, x , cex =0.25 , xlab =" Beta (3 , 2)", ylab =" Sample ")

Diketahui :

\(U, V\) peubah acak saling bebas menyebar \(U (0,1)\), kemudian : \[ X = [1 + \frac{log(V)}{log(1-(1-\theta)^U)}]\]

# Kombinasikan  Sebaran logarithmic (0.5) dengan
# sampel untuk transformasi
n = 1000 
theta = 0.5
u = runif ( n ) # membangkitkan angka logarithmic
v = runif ( n )
x = floor (1 + log( v ) / log (1 - (1 - theta )^ u ))
k = 1: max( x ) # hitung peluang logarithmic
p = -1 / log (1 - theta ) * theta ^ k / k
se = sqrt ( p*(1 - p )/n )
p.hat = tabulate ( x )/n
print ( round ( rbind (p.hat , p , se) , 3))
       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]
p.hat 0.712 0.177 0.062 0.030 0.013 0.004 0.001 0.001
p     0.721 0.180 0.060 0.023 0.009 0.004 0.002 0.001
se    0.014 0.012 0.008 0.005 0.003 0.002 0.001 0.001

4. Convolution

Jika \(X_1, ..., X_n\) i.i.d. dengan \(X_j\sim X\) dan \(S = X_1 + ... + X_n\). Fungsi sebaran dari penjumlahan S disebut n-fold convolution. Sangat mudah untuk mensimulasikan konvolusi dengan langsung membangkitkan \(X_1, ..., X_n\) dan menghitung penjumlahan.

Contoh :

  • Jika v > 0 adalah bilangan bulat menyebar chi-square dengan db v adalah konvolusi dari v iid kuadrat peubah normal baku

  • Sebaran negative binomial , NegBin ( r , p ) adalah konvolusi dari r iid peubah acak Geom ( p )

  • Konvolusi dari r (independen) peubah acak Exp ( \(\lambda\) ) mempunyai sebaran Gamma ( r, \(\lambda\) )

Latihan

Bangkitkan peubah acak \(X^2(v)\). Jika \(Z_1 , ... , Z_v\) iid peubah acak \(N(0,1)\) kemudian \(V = Z^2_1 + ... + Z^2_v\) menyebar \(X^2(v)\)

#1. Isi matriks n x v dengan p.a v dari N(0,1)
#2. Hitung jumlah baris dari kuadrat normal
#3. Tentukan vektor jumlah baris

n  = 1000
nu = 2 # v
X  = matrix ( rnorm (n*nu ) , n , nu )^2 # matriks dari sq. normals
# Menjumlahkan kuadrat normal di setiap baris : method 1
y  = rowSums (X )
y  = apply (X , MARGIN =1 , FUN = sum ) # vektor length n, method 2
mean(y)
[1] 1.957253
mean(y^2)
[1] 7.549724

Referensi

  1. Firdaus F. 2015. Pseudo-random number generator menggunakan waktu lokal, konsep efek kupu-kupu silang.Makalah. Sekolah Teknik Elektro dan Informatika ITB.

  2. Jianhua Z. Chapter 3 : Methods for Generating Random Variables. Department of Statistics Yunnan University of Finance and Economics.

  3. Rahmi A. 2021. Pembangkitan Bilangan Acak. [internet]. [diakses 2021 Sept 12]. Tersedia pada : https://rpubs.com/annisarahmi/749913

  4. Soleh AM.2019. STA571 : Pemrograman Statistika: Pembangkitan Bilangan Acak dan Resampling Departemen Statistika Institut Pertanian Bogor.

"Semua kehidupan adalah percobaan. Semakin banyak percobaan yang kamu lakukan, semakin baik." -- Ralph Waldo Emerson