Dalam statistik, kita sering ingin melakukan resampling untuk menguji
keakuratan estimasi sampel kita. Ini disebut bootstrapping—sebuah tes
yang didasarkan pada pengambilan sampel acak berulang dengan
pengembalian (replacement). Dalam pembelajaran kali ini, kita
akan:
(1) menghasilkan 1000 sampel acak baru dari sebuah vektor
(dengan pengembalian),
(2) menghitung estimasi tertentu (misalnya,
mean, median, dll.) untuk setiap sampel tersebut, dan
(3)
menghitung standar deviasi dari distribusi estimasi tersebut.
Kita akan membuat sebuah vektor berisi angka yang terdistribusi secara normal dengan nama myData.
Fungsi rnorm() memungkinkan kita untuk menghasilkan urutan angka acak dari distribusi normal hipotetis dengan rata-rata (mean) 20 dan standar deviasi 4.5. Di sini, kita memilih untuk menghasilkan 2000 observasi.
set.seed(300) # Setting the seed for replication purposes
myData <- rnorm(2000,20,4.5) # Creating a random normal distribution (n=300, mean=20, sd=4.5)
Penjelasan:
- set.seed(300) digunakan untuk
menetapkan seed (angka acak) agar hasil yang dihasilkan oleh fungsi acak
(seperti rnorm) dapat direplikasi. Ini memastikan bahwa setiap kali kode
dijalankan, hasilnya akan sama.
- rnorm(2000, 20, 4.5)
menghasilkan 2000 observasi dari distribusi normal dengan mean 20 dan
standar deviasi 4.5. Hasilnya disimpan dalam vektor myData.
Untuk memastikan bahwa vektor myData telah dibuat dengan benar, kita
memverifikasi bahwa:
1. Terdapat 2000 observasi (menggunakan
length()),
2. Nilai rata-rata (mean) mendekati 20 (menggunakan
mean()), dan
3. Standar deviasi mendekati 4.5 (menggunakan
sd()).
Catatan: Kita tidak seharusnya mengharapkan mean dan standar deviasi memiliki nilai yang persis sama dengan yang kita masukkan ke dalam rnorm(), karena ada unsur acak yang terlibat saat menggunakan fungsi ini (rnorm menghasilkan distribusi normal acak di sekitar mean dan standar deviasi yang diberikan).
length(myData) # How many observations?
## [1] 2000
length(myData) menghitung jumlah observasi dalam vektor myData. Harusnya menghasilkan 2000 karena kita membuat 2000 observasi.
mean(myData) # What is the mean?
## [1] 20.25773
mean(myData) menghitung rata-rata dari vektor myData. Harusnya mendekati 20 karena kita menetapkan mean = 20 saat membuat data.
sd(myData) # What is the standard deviation?
## [1] 4.590852
sd(myData) menghitung standar deviasi dari vektor myData. Harusnya mendekati 4.5 karena kita menetapkan standar deviasi = 4.5 saat membuat data.
Untuk tujuan visual, berikut adalah grafik dari vektor myData. Garis putus-putus merah menunjukkan nilai mean (20.25773) dari distribusi ini.
# Membuat histogram dari myData
hist(myData, breaks = 30, col = "lightblue", main = "Distribusi myData", xlab = "Nilai", ylab = "Frekuensi")
# Menambahkan garis vertikal untuk mean
abline(v = mean(myData), col = "red", lwd = 2, lty = 2)
# Menambahkan legenda
legend("topright", legend = paste("Mean =", round(mean(myData), 5)), col = "red", lty = 2, lwd = 2)
Penjelasan:
Apa itu bootstrapping? Bootstrapping adalah metode untuk memperhitungkan ketidakpastian saat kita mengukur estimasi sampel. Dengan menggunakan pengambilan sampel acak dengan pengembalian (random sampling with replacement), kita menghitung estimasi kita untuk berbagai sampel acak yang diambil dari distribusi asli kita.
Analogi dengan “Mean of Means” pada Dadu Perhatikan kejadian saat menghitung “rata-rata dari rata-rata” (mean of means) untuk lemparan dadu, kita memperkirakan seberapa besar kemungkinan kita akan mengamati estimasi sampel tertentu (misalnya, mean, standar deviasi). Jika kita mengambil sampel ulang secara acak dari distribusi dasar ini dengan penggantian sejumlah besar kali. Alasan kita mengambil sampel ulang dengan penggantian adalah untuk memungkinkan angka-angka tertentu dipilih beberapa kali, yang menciptakan varians yang lebih besar dalam distribusi yang dihasilkan.
Mengapa Menggunakan Pengembalian (Replacement)? Alasan kita menggunakan pengembalian adalah untuk memungkinkan angka tertentu dipilih lebih dari satu kali, yang menciptakan variasi yang lebih besar dalam distribusi yang dihasilkan.
Langkah-langkah Bootstrapping: 1. Ambil sampel acak
dengan pengembalian dari distribusi asli.
Secara sederhana, yang kita lakukan di sini adalah menghitung mean dari
1000 sampel yang masing-masing terdiri dari 2000 observasi dari myData
menggunakan perulangan for(i in x). Mean dari setiap sampel ini disimpan
dalam sebuah vektor (bootstrap.results).
Tujuan: Kita tertarik untuk mengukur ketidakpastian dalam distribusi estimasi kita (dalam hal ini, mean). Untuk melakukannya, kita menghitung standar deviasi dari setiap mean yang dihitung menggunakan pengambilan sampel acak dengan pengembalian (random sampling with replacement).
set.seed(200) # Setting the seed for replication purposes
sample.size <- 2000 # Sample size
n.samples <- 1000 # Number of bootstrap samples
bootstrap.results <- c() # Creating an empty vector to hold the results
for (i in 1:n.samples) { obs <- sample(1:sample.size, replace=TRUE) bootstrap.results[i] <- mean(myData[obs]) # Mean of the bootstrap sample }
Penjelasan:
length(bootstrap.results) # Sanity check: this should contain the mean of 1000 different samples
## [1] 1000
length(bootstrap.results) memeriksa panjang vektor bootstrap.results. Harusnya 1000 karena kita mengambil 1000 sampel.
summary(bootstrap.results) # Sanity check
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 19.92 20.19 20.26 20.26 20.33 20.57
summary(bootstrap.results) memberikan ringkasan statistik dari vektor bootstrap.results, termasuk min, max, mean, dan kuartil.
sd(bootstrap.results) # Checking the standard deviation of the distribution of means (this is what we are interested in!)
## [1] 0.1021229
sd(bootstrap.results) menghitung standar deviasi dari distribusi mean yang dihasilkan oleh bootstrapping.
par(mfrow=c(2,1), pin=c(3,0.98)) # Combining plots (2 rows, 1 column) and setting the plots size
hist(bootstrap.results, # Creating an histogram
col="#d83737", # Changing the color
xlab="Mean", # Giving a label to the x axis
main=paste("Means of 1000 bootstrap samples from myData")) # Giving a title to the graph
hist(myData, # Creating an histogram
col="#37aad8", # Changing the color
xlab="Value", # Giving a label to the x axis
main=paste("Distribution of myData")) # Giving a title to the graph
Penjelasan:
Demikian pula, kita juga dapat mengambil 1000 sampel acak dari proses pembuatan data asli.
set.seed(200) # Setting the seed for replication purposes
sample.size <- 2000 # Sample size
n.samples <- 1000 # Number of bootstrap samples
bootstrap.results <- c() # Creating an empty vector to hold the results
for (i in 1:n.samples)
{
bootstrap.results[i] <- mean(rnorm(2000,20,4.5)) # Mean of the bootstrap sample
}
Penjelasan:
for (i in 1:n.samples) melakukan perulangan sebanyak 1000
kali.
length(bootstrap.results) # Sanity check: this should contain the mean of 1000 different samples
## [1] 1000
length(bootstrap.results) memeriksa panjang vektor bootstrap.results. Harusnya menghasilkan 1000 karena kita mengambil 1000 sampel bootstrap.
summary(bootstrap.results) # Sanity check
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 19.64 19.93 20.00 20.00 20.07 20.32
summary(bootstrap.results) memberikan ringkasan statistik dari vektor bootstrap.results, termasuk nilai minimum, kuartil pertama (Q1), median (Q2), mean, kuartil ketiga (Q3), dan nilai maksimum. Ini membantu kita memahami distribusi dari mean yang dihasilkan oleh bootstrapping.
sd(bootstrap.results) # Checking the standard deviation of the distribution of means (this is what we are interested in!)
## [1] 0.1041927
sd(bootstrap.results) menghitung standar deviasi dari distribusi mean yang dihasilkan oleh bootstrapping. Standar deviasi ini mengukur variabilitas atau ketidakpastian dari estimasi mean yang kita peroleh melalui bootstrapping.
par(mfrow=c(2,1), pin=c(3,1)) # Combining plots (2 rows, 1 column) and setting the plots size
hist(bootstrap.results, # Creating an histogram
col="#d83737", # Changing the color
xlab="Mean", # Giving a label to the x axis
main=paste("Means of 1000 bootstrap samples from the DGP")) # Giving a title to the graph
hist(myData, # Creating an histogram
col="#37aad8", # Changing the color
xlab="Value", # Giving a label to the x axis
main=paste("Distribution of myData")) # Giving a title to the graph
Penjelasan:
hist(bootstrap.results, col=“#d83737”, xlab=“Mean”, main=paste(“Means
of 1000 bootstrap samples from the DGP”)) membuat histogram dari
distribusi mean hasil bootstrapping.
hist(myData, col=“#37aad8”, xlab=“Value”, main=paste(“Distribution of
myData”)) membuat histogram dari data asli myData.