Estimasi Distribusi dan Parameter Model

Teorema Limit Pusat

• Rata-rata dari contoh acak yang berasal dari sebaran apapun memiliki sebaran normal jika ukuran contohnya sangat besar

• Jika contoh acak diambil dari populasi dengan mean \(μ\) dan ragam \(σ^2\), maka semakin besar ukuran contoh, sebaran dari \(\bar{x}\) akan semakin mendekati sebaran normal dengan mean \(μ\) dan ragam \(σ^2/n\)

Algoritma

1. Tentukan ukuran contoh (\(n\))
2. Tentukan sebaran data
3. Ulang \(k\) kali

• Ambil \(n\) contoh acak dari sebaran data yang sudah ditentukan

• Hitung rataannya lalu simpan

4. Periksa sebaran \(k\) rataan

Aplikasi di R

Hampiran Normal Terhadap Geometrik

library(probs)

## 
## Attaching package: 'probs'

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, union

par(mfrow = c(3,1))
set.seed(123)
populasi = rgeom(20, 0.1)

n1 = 2
contoh_geo1 = urnsamples(populasi, size = 2, replace = F, ordered = F)
mean_geo1 = matrix(apply(contoh_geo1, 1, mean))

n2 = 5
contoh_geo2 = urnsamples(populasi, size = 5, replace = F, ordered = F)
mean_geo2 = matrix(apply(contoh_geo2, 1, mean))

n3 = 10
contoh_geo3 = urnsamples(populasi, size = 10, replace = F, ordered = F)
mean_geo3 = matrix(apply(contoh_geo3, 1, mean))

hist(mean_geo1, main = paste("Hampiran Normal terhadap Geometrik (n = 2)"), xlab = "xbar")
hist(mean_geo2, main = paste("Hampiran Normal Terhadap Geometrik (n = 5)"), xlab = "xbar")
hist(mean_geo3, main = paste("Hampiran Normal Terhadap Geometrik (n = 10)"), xlab = "xbar")

menunjukkan bagaimana distribusi rata-rata sampel \(\bar{x}\) dari distribusi geometrik mendekati distribusi normal ketika ukuran sampel (\(n\)) meningkat.

Hampiran Normal Terhadap Eksponensial

par(mfrow=c(3,1))
library(probs)
set.seed(123)
populasi = rexp(20)

n = 2
contoh_exp1 = urnsamples(populasi, size = 2, replace = F, ordered = F)
mean_exp1 = matrix(apply(contoh_exp1, 1, mean))

n2 = 5
contoh_exp2 = urnsamples(populasi, size = 5, replace = F, ordered = F)
mean_exp2 = matrix(apply(contoh_exp2, 1, mean))

n3 = 10
contoh_exp3 = urnsamples(populasi, size = 10, replace = F, ordered = F)
mean_exp3 = matrix(apply(contoh_exp3, 1, mean))
 
hist(mean_exp1, main = paste("Hampiran Normal Terhadap Eksponensial (n = 2)"), xlab = "xbar")
hist(mean_exp2, main = paste("Hampiran Normal Terhadap Eksponensial (n = 5)"), xlab = "xbar")
hist(mean_exp3, main = paste("Hampiran Normal Terhadap Eksponensial (n = 10)"), xlab = "xbar")

Hampiran Normal Terhadap Seragam

par(mfrow = c(3,1))
library(probs)
set.seed(123)
populasi = runif(20)

n1 = 2
contoh_unif1 = urnsamples(populasi, size = 2, replace = F, ordered = F)
mean_unif1 = matrix(apply(contoh_unif1, 1, mean))

n2 = 5
contoh_unif2 = urnsamples(populasi, size = 5, replace = F, ordered = F)
mean_unif2 = matrix(apply(contoh_unif2, 1, mean))

n3 = 10
contoh_unif3 = urnsamples(populasi, size = 10, replace = F, ordered = F)
mean_unif3 = matrix(apply(contoh_unif3, 1, mean))

hist(mean_unif1,main = paste("Hampiran Normal Terhadap Seragam (n = 2)"),xlab = "xbar")
hist(mean_unif2,main = paste("Hampiran Normal Terhadap Seragam (n = 5)"),xlab = "xbar")
hist(mean_unif3,main = paste("Hampiran Normal Terhadap Seragam (n = 10)"),xlab = "xbar")

kesimpulan : Semakin besar ukuran contoh, maka sebaran rata-rata dari contoh acak yang berasal dari sebaran geometrik, eksponensial, maupun uniform akan mendekati sebaran normal. Hal ini ditunjukkan dari histogram yang mana ketika n seakin besar akan semakin cenderung memberntuk kurva normal

Sebaran Percontohan Sebaran Normal

par(mfrow = c(3,1))
library(probs)
set.seed(123)
populasi = rnorm(20,5, sqrt(12)) # Membangkitkan bil, acak ~ Normal (miu = 5, sigma2 = 12)

n1 = 3
contoh_norm1 = urnsamples(populasi, size = 3, replace = F, ordered = F)
mean_norm1 = matrix(apply(contoh_norm1, 1, mean))
mean_xbar1 = mean(mean_norm1)
var_xbar1 = var(mean_norm1)

n2 = 4
contoh_norm2 = urnsamples(populasi, size = 4, replace = F, ordered = F)
mean_norm2 = matrix(apply(contoh_norm2, 1, mean))
mean_xbar2 = mean(mean_norm2)
var_xbar2 = var(mean_norm2)

n1 = 15
contoh_norm3 = urnsamples(populasi, size = 15, replace = F, ordered = F)
mean_norm3 = matrix(apply(contoh_norm3, 1, mean))
mean_xbar3 = mean(mean_norm3)
var_xbar3 = var(mean_norm3)

hist(mean_norm1, main = paste("(n = 3)"),xlab = "xbar")
hist(mean_norm2,main = paste("(n = 4)"),xlab = "xbar")
hist(mean_norm3,main = paste("(n = 15)"),xlab = "xbar")

hasil = data.frame("."=c("mean", "varian"), "Populasi" =c(5,12), "n=3"=c(mean_xbar1, var_xbar1), "n=4"=c(mean_xbar2,var_xbar2),"n=15"=c(mean_xbar3,var_xbar3))
hasil

##        . Populasi      n.3      n.4      n.15
## 1   mean        5 5.490599 5.490599 5.4905992
## 2 varian       12 3.219489 2.271056 0.1892278

Kesimpulan:

• Berdasarkan output di atas, contoh acak yang diambil dari populasi dengan mean = \(μ\) dan ragam = \(σ^2\), maka semakin besar ukuran contoh, mean dari \(\bar{x}\) akan semakin mendekati \(μ\) dan ragamnya semakin mendekati \(σ^2/n\)

• Dilihat dari histogramnya pun, semakin besar ukuran contoh, kurva mendekati sebaran normal.

Ketakbiasan Penduga Parameter

\(\bar{x}\) adalah penduga tak bias bagi \(μ\) , jika E(\(\bar{x}\)) = \(μ\) \(s^2\) adalah penduga tak bias bagi \(σ^2\) , jika E(\(s^2\)) = \(σ^2\)

Maka dalam hal ini kita akan membuktikan apakah benar nilai harapan penduga parameter sama dengan nilai parameternya

Algoritme:

1. Tentukan sebaran yang akan digunakan
2. Ulangi sebanyak k kali (bangkitkan n buah data dari sebaran yang sudah ditentukan, hitung nilai \(\bar{x}\) dan \(s^2\))
3. Hitung rata-rata dari \(\bar{x}\) dan \(σ^2\)

Aplikasi di R

Penerapan ketakbiasan penduga (mean)

library(probs)
set.seed(123)

# Populasi terhingga

# 1. Sebaran Normal
n = 10
populasi1 = rnorm(20)
mean_pop1 = mean(populasi)
sampel_normal1 = urnsamples(populasi, size = 10, replace = F, ordered = F)
mean_normal1 = matrix(apply(sampel_normal1, 1, mean))
median_normal1 = matrix(apply(sampel_normal1, 1, mean))
harapan_mean_norm1 = mean(mean_normal1)
harapan_median_norm1 = mean((median_normal1))

# 2. Sebaran Eksponensial
n = 10
populasi2 = rexp(20)
mean_pop2 = mean(populasi2)
sampel_exp1 = urnsamples(populasi2, size = 10, replace = F, ordered = F)
mean_exp1 = matrix(apply(sampel_exp1, 1, mean))
median_exp1 = matrix (apply(sampel_exp1, 1, median))
harapan_mean_exp1 = mean(mean_exp1)
harapan_median_exp1 = mean(median_exp1)

# 3. Uniform
n = 10
populasi3 = runif(20)
mean_pop3 = mean(populasi3)
sampel_unif1 = urnsamples(populasi3, size = 10, replace = F, ordered = F)
mean_unif1 = matrix(apply(sampel_unif1, 1, mean))
median_unif1 = matrix(apply(sampel_unif1, 1, mean))
harapan_mean_unif1 = mean(mean_unif1)
harapan_median_unif1 = mean(median_unif1)

hasil = data.frame("Hasil" = c("mean_populasi", "harapan_mean_contoh", "harapan_median_contoh"), "Sebaran Normal"=c(mean_pop1, harapan_mean_norm1, harapan_median_norm1), "Sebaran Eksponensial" =c(mean_pop2, harapan_mean_exp1, harapan_median_exp1), "Sebaran Seragam" = c(mean_pop3, harapan_mean_unif1, harapan_median_unif1))

hasil

##                   Hasil Sebaran.Normal Sebaran.Eksponensial Sebaran.Seragam
## 1         mean_populasi       5.490599             1.113617        0.443685
## 2   harapan_mean_contoh       5.490599             1.113617        0.443685
## 3 harapan_median_contoh       5.490599             1.132277        0.443685

Kesimpulan:

• Berdasarkan output di atas, dengan populasi terhingga maupun tak hingga serta tiga sebaran yang berbeda, nilai harapan median contoh tetap berbeda dengan \(μ\) dan nilai harapan rataan populasi \(μ\) (pada populasi terhingga) sehingga penduga tak bias bagi \(μ\) adalah (\(\bar{x}\))

• Pada populasi terhingga, percontohan bersifat unik artinya tidak ada percontohan yang berulang sehingga dapat dipastikan kombinasi contoh hanya muncul satu kali sehingga nilai parameter dan nilai harapan penduga parameter yang tak bias sama persis.

• Pada populasi tak hingga, percontohan yang terambil secara acak merupakan sebagian dari keseluruhan kemungkinan percontohan yang ada sehingga nilai parameter dan nilai harapan penduga parameter yang tak bias tidak sama persis, namun sangat mendekati.

Penerpaan ketakbiasan penduga (Ragam)

library(probs)
set.seed(999)
# Populasi Terhingga

# Sebaran Normal
n        = 10
populasi = rnorm(20) 
sigma2   = var(populasi)*(20-1)/20 #fungsi var pada R adalah varian contoh (penyebut n-1) sehingga perlu dikali (n-1)/n

sampel = urnsamples(populasi, size = 10, replace =F, ordered = F)

## Pembagi (n-1)
s2.n1 = matrix(apply(sampel, 1, var))
E.s2.n1 = mean(s2.n1)

## Pembagi (n)
s2.n = s2.n1*(10-1)/10
E.s2.n = mean(s2.n)

# Sebaran Eksponensial
n = 10
populasi2 = rexp(20)
sigma2.exp = var(populasi2)*(20-1)/20
  
sampel_exp = urnsamples(populasi2, size = 10, replace = F, ordered = F)

## Pembagi (n-1)
s2.n1.exp = matrix(apply(sampel_exp, 1, var))
E.s2.n1.exp = mean(s2.n1.exp)

## Pembagi (n)
s2.n.exp    = s2.n1.exp*(10-1)/10
E.s2.n.exp  = mean(s2.n.exp)

hasil = data.frame( "."  = c("ragam populasi","nilai harapan ragam contoh (n-1)","nilai harapan ragam contoh (n)"), 
                    "Sebaran Normal" = c(sigma2, E.s2.n1, E.s2.n),"Sebaran Eksponensial" = c(sigma2.exp, E.s2.n1.exp, E.s2.n.exp))
hasil

##                                  . Sebaran.Normal Sebaran.Eksponensial
## 1                   ragam populasi      0.9669853            0.8550756
## 2 nilai harapan ragam contoh (n-1)      1.0178793            0.9000796
## 3   nilai harapan ragam contoh (n)      0.9160914            0.8100716

Selang Kepercayaann

Apa arti dari SK 95%? - SK 95% bagi \(θ\) : Kita percaya 95% bahwa selang a sampai b memuat nilai parameter \(θ\) yang sebenarnya - SK 95%: Jika kita melakukan 100 kali percontohan acak dan setiap percontohan acak dibuat selang kepercayaannya, maka dari 100 SK yang terbentuk, ada 95 SK yang mencakup parameter sedangkan sisanya sebanyak 5 SK tidak mencakup parameter.

Algoritme

1. Tentukan sebaran yang akan digunakan
2. Ulangi sebanyak k kali

• Bangkitkan n buah data dari sebaran yang sudah ditentukan
• Hitung nilai \(\bar{x}\) dan \(s^@\)
• Hitung \(\bar{σ}_{x}^2\) dan buat selang kepercayaan (1 - \(α\))%

3. Hitung proporsi benyaknya selang kepercayaann yang memuat \(μ\) dibandingkan dengan (1 - α)

Aplikasi di R

n1 = 10
k = 100 #ulangan
alpha = 0.05
mu = 50
std = 10
set.seed(123)
sampel_norm1 = matrix(rnorm(n1*k, mu, std),k)
xbar.norm1 = apply(sampel_norm1, 1, mean)
s.norm1 = apply(sampel_norm1, 1, mean)
SE.norm1 = s.norm1/sqrt(n1)
z.norm1 = qnorm(1-alpha/2)
SK.norm1 = (xbar.norm1-z.norm1 < mu & mu < xbar.norm1+z.norm1*SE.norm1)
x.norm1 = sum(SK.norm1)/k # Proporsi banyaknya SK yang memuat mu

n2 = 30
k = 100 #ulangan
alpha = 0.05
mu = 50
std = 10
set.seed(123)
sampel.norm2 = matrix(rnorm(n2*k, mu, std),k)
xbar.norm2   = apply(sampel.norm2,1,mean)
s.norm2      = apply(sampel.norm2,1,sd)
SE.norm2     = s.norm2/sqrt(n2)
z.norm2      = qnorm(1-alpha/2)
SK.norm2     = (xbar.norm2-z.norm2*SE.norm2 < mu & mu < xbar.norm2+z.norm2*SE.norm2)
x.norm2      = sum(SK.norm2)/k #proporsi banyaknya SK yang memuat mu

n3     = 100
k      = 100 #ulangan
alpha  = 0.05
mu     = 50
std    = 10
set.seed(123)
sampel.norm3 = matrix(rnorm(n3*k,mu,std),k)
xbar.norm3   = apply(sampel.norm3,1,mean)
s.norm3      = apply(sampel.norm3,1,sd)
SE.norm3     = s.norm3/sqrt(n3)
z.norm3      = qnorm(1-alpha/2)
SK.norm3     = (xbar.norm3-z.norm3*SE.norm3 < mu & mu < xbar.norm3+z.norm3*SE.norm3)
x.norm3      = sum(SK.norm3)/k #proporsi banyaknya SK yang memuat mu

hasil = data.frame("n" =c(10,30,100),"Ketepatan SK Sebaran Normal"=c(x.norm1, x.norm2, x.norm3))
hasil

##     n Ketepatan.SK.Sebaran.Normal
## 1  10                        0.74
## 2  30                        0.93
## 3 100                        0.96

matplot(rbind (xbar.norm2-z.norm2*SE.norm2, xbar.norm2+z.norm2*SE.norm2), rbind(1:k,1:k), col=ifelse(SK.norm2,"blue","red"), type = "l", lty = 1,main='Selang Kepercayaan 95% (n=100)', xlab='SK', ylab='banyak ulangan')
abline(v=mu)

• Gambar ini adalah hasil dari simulasi selang kepercayaan 95% untuk rata-rata populasi (\(μ\)=50) dengan ukuran sampel \(n\)=100

• Simulasi dilakukan sebanyak \(k\)=100 kali, sehingga ada 100 selang kepercayaan yang dihasilkan.

• Tujuan gambar ini adalah untuk memvisualisasikan seberapa sering selang kepercayaan berhasil menangkap nilai rata-rata populasi (\(μ\))

• Garis vertikal di \(x\)=50 mengartikan nilai rata-rata populasi (\(μ\)=50). Selang kepercayaan yang berhasil menangkap μ akan melintasi garis ini.

• Garis Horizontal mengartikan bahwa setiap garis mewakili selang kepercayaan dari satu sampel. Jika garis tersebut melintasi garis vertikal di x=50, artinya selang kepercayaan tersebut berhasil menangkap \(μ\)

• Semakin besar ukuran contoh (\(n\)), maka proporsi SK yang memuat nilai parameter semakin mendekati kebenaran (1-alpha)

# Interal Kepercayaan
library(car)

## Loading required package: carData

data("Prestige")

# Menghitung rata-rata
m <- mean(Prestige$income)
m

## [1] 6797.902

p <- dim(Prestige) [1]
se <- sd(Prestige$income)/sqrt(p)
se

## [1] 420.4089

# Menghitung standar error
p <- dim(Prestige)[1]
se <- sd(Prestige$income)/sqrt(p)
se

## [1] 420.4089

# Menghitung nilai kritis t
tval <- qt(0.975, df=p-1)

# Menghitung interval kepercayaan
cat(paste("KI: [", round(m-tval*se, 2),",",round(m+tval*se,2),"]"))

## KI: [ 5963.92 , 7631.88 ]

Artinya, dengan tingkat kepercayaan 95%, rata-rata pendapatan populasi berada dalam rentang 5963.92 hingga 7631.88.