Tugas Kelompok MSR P10 dan 11

Afris Setiya Intan Amanda

5/4/2023

Simulasi Statistika
Anggota Kelompok 2 Paralel 2 STA1372 :

1. Ervina Dwi Anggrahini       (G1401201017)
2. Afris Setiya Intan Amanda     (G1401201018)
3. Syadiah                 (G1401201042)

Outline

  1. Soal
  2. Jawaban

Soal

Soal 1

Soal 2

Jawaban

Library yang digunakan dalam data yang dianalisis ini adalah sebagai berikut.

library(fitdistrplus)
## Warning: package 'fitdistrplus' was built under R version 4.2.3
## Loading required package: MASS
## Loading required package: survival
## Warning: package 'survival' was built under R version 4.2.2

Soal 1

Soal 1.1, 1.3, & 1.4

Pembangkitan data sebaran simetris dengan menggunakan sebaran normal dengan mengambil contoh dengan \(n\) = 4, 12, 20, 60, 100 dan banyak ulangan (\(k\)) sejumlah 1000 dan ukuran populasi 4 juta.

data.s1<-rnorm(4000000, mean = 0, sd = 1)
hist(data.s1, main="Sebaran Populasi", xlab="Populasi", ylab="Frekuensi", col="coral")

set.seed(17)
k=1000
norm.1 = matrix(sample(data.s1,4*k),k)
norm.1 = apply(norm.1,1,mean)
norm.2 = matrix(sample(data.s1,12*k),k)
norm.2 = apply(norm.2,1,mean)
norm.3 = matrix(sample(data.s1,20*k),k)
norm.3 = apply(norm.3,1,mean)
norm.4 = matrix(sample(data.s1,60*k),k)
norm.4 = apply(norm.4,1,mean)
norm.5 = matrix(sample(data.s1,100*k),k)
norm.5 = apply(norm.5,1,mean)
for(i in 1:5)
  {
cat("Rata-rata sebaran contoh ke-", i, ":")
print(mean(get(paste0("norm.",i))))
cat("Ragam sebaran contoh", i, ":")
print(var(get(paste0("norm.",i))))
cat("\n")
  }
## Rata-rata sebaran contoh ke- 1 :[1] 0.003974997
## Ragam sebaran contoh 1 :[1] 0.2359441
## 
## Rata-rata sebaran contoh ke- 2 :[1] 0.003918884
## Ragam sebaran contoh 2 :[1] 0.07286735
## 
## Rata-rata sebaran contoh ke- 3 :[1] 0.001606943
## Ragam sebaran contoh 3 :[1] 0.04675378
## 
## Rata-rata sebaran contoh ke- 4 :[1] -0.002162009
## Ragam sebaran contoh 4 :[1] 0.0163305
## 
## Rata-rata sebaran contoh ke- 5 :[1] 0.001097989
## Ragam sebaran contoh 5 :[1] 0.01020871

Berdasarkan hasil simulasi di atas, terlihat ketika memperbesar ukuran contoh, maka semakin kecil keragaman data yang mana cenderung semakin mendekati nilai nol.

#Histogram
par(mfrow=c(2,3))
hist(norm.1, xlim=c(-2,2))
hist(norm.2, xlim=c(-2,2))
hist(norm.3, xlim=c(-2,2))
hist(norm.4, xlim=c(-2,2))
hist(norm.5, xlim=c(-2,2))

Visualisasi histogram di atas menunjukan bahwa semakin besar ukuran contoh (\(n\)), bar-bar pada histogram akan semakin sempit yang berarti nilai ragam pada data semakin kecil. Selain itu, QQ-plot digunakan untuk melihat kenormalan data berdasarkan masing-masing \(n\).

#QQ-plot
par(mfrow=c(2,3))
qqcomp(fitdist(norm.1, "norm"), main = "QQ-Plot n = 4")
qqcomp(fitdist(norm.2, "norm"), main = "QQ-Plot n = 12")
qqcomp(fitdist(norm.3, "norm"), main = "QQ-Plot n = 20")
qqcomp(fitdist(norm.4, "norm"), main = "QQ-Plot n = 60")
qqcomp(fitdist(norm.5, "norm"), main = "QQ-Plot n = 100")

Dapat terlihat bahw \(n\) berjumlah kecil berdistribusi sangat mendekati sebaran normal karena amatan-amatan data mengikuti garis tren. Hal ini disebabkan oleh sebaran pada populasi data yang digunakan adalah sebaran simetris, yaitu sebaran normal.

Simpulan

Sebaran rataan mulai mendekati sebaran normal dari mulai \(n\)=4. Oleh karena itu, pada populasi menyebar normal atau simetris, dengan jumlah sampel kecil, sebaran rataan sudah mendekati normal atau simetris.

Soal 1.2.a, 1.3, & 1.4

Data akan dibangkitkan dengan membangkitkan 50% data yang menyebar Normal(0,1) dan 50% data yang menyebar Khi-kuadrat(3) dengan 4 juta amatan data sebagai populasi. Lalu, pengambilan contoh dilakukan dengan \(n\) = 4, 12, 20, 60, 100 dan banyak ulangan (\(k\)) sejumlah 1000.

set.seed(17)
pop.data <- c(rnorm(2000000,mean = 0, sd = 1),rchisq(2000000,3))
hist(pop.data, main="Sebaran Populasi", xlab="Populasi", ylab="Frekuensi", col="coral")

set.seed(17)
k=1000
n1=4
n2=12
n3=20
n4=60
n5=100
nchi1 = matrix(sample(pop.data,n1*k),k)
nchi1 = apply(nchi1,1,mean)
nchi2 = matrix(sample(pop.data,n2*k),k)
nchi2 = apply(nchi2,1,mean)
nchi3 = matrix(sample(pop.data,n3*k),k)
nchi3 = apply(nchi3,1,mean)
nchi4 = matrix(sample(pop.data,n4*k),k)
nchi4 = apply(nchi4,1,mean)
nchi5 = matrix(sample(pop.data,n5*k),k)
nchi5 = apply(nchi5,1,mean)
#Histogram
par(mfrow=c(2,3))
hist(nchi1, xlim=c(0,8))
hist(nchi2, xlim=c(0,8))
hist(nchi3, xlim=c(0,8))
hist(nchi4, xlim=c(0,8))
hist(nchi5, xlim=c(0,8))

visualisasi histogram di atas menunjukan bahwa semakin besar ukuran contoh, semakin menyempit bar-bar histogram yang didapat yang mana seiring semakin kecil ragam pada data. Selain itu, visulaisasi di atas menunjukan bahwa semakin besar ukuran contoh, sebaran data tersebut semakin mendekati sebaran normal.

#QQ-plot
par(mfrow=c(2,3))
qqcomp(fitdist(nchi1, "norm"), main = "QQ-plot n = 4")
qqcomp(fitdist(nchi2, "norm"), main = "QQ-plot n = 12")
qqcomp(fitdist(nchi3, "norm"), main = "QQ-plot n = 20")
qqcomp(fitdist(nchi4, "norm"), main = "QQ-plot n = 60")
qqcomp(fitdist(nchi5, "norm"), main = "QQ-plot n = 100")

Berdasarkan QQ-plot di atas, terlihat bahwa data dengan n = 4 dan n = 12 memiliki sebaran statistik yang tidak terlalu mendekati sebaran normal. Sebaran mulai mendekati sebaran normal ketika n = 20 sampai 100 yang mana ukuran contoh semakin besar.

Simpulan

Sebaran rataan dengan \(n\)=4 dan \(n\)=12 belum mendekati normal, masih cenderung menjulur ke kanan. Sebaran rataan mendekati normal mulai dari \(n\)=20. Oleh karena itu, pada populasi menyebar tidak simetris (gabungan normal dan khi-kuadrat) membutuhkan sampel lebih banyak agar sebaran rataan mendekati normal.

Soal 1.2.b, 1.3, & 1.4

Data akan dibangkitkan dengan membangkitkan 50% data yangyang menyebar Khi-kuadrat(5) dan 50% data yang menyebar Khi-kuadrat(2) dengan 4 juta amatan data sebagai populasi.Lalu, pengambilan contoh dilakukan dengan \(n\) = 4, 12, 20, 60, 100 dan banyak ulangan (\(k\)) sejumlah 1000.

set.seed(17)
popchi <- c(rchisq(2000000,5),rchisq(2000000,2))
hist(popchi, main="Sebaran Populasi", xlab="Populasi", ylab="Frekuensi", col="coral")

set.seed(17)
k=1000
n1=4
n2=12
n3=20
n4=60
n5=100
chi1 = matrix(sample(popchi,n1*k),k)
chi1 = apply(chi1,1,mean)
chi2 = matrix(sample(popchi,n2*k),k)
chi2 = apply(chi2,1,mean)
chi3 = matrix(sample(popchi,n3*k),k)
chi3 = apply(chi3,1,mean)
chi4 = matrix(sample(popchi,n4*k),k)
chi4 = apply(chi4,1,mean)
chi5 = matrix(sample(popchi,n5*k),k)
chi5 = apply(chi5,1,mean)
# Histogram
par(mfrow=c(2,3))
hist(chi1, xlim=c(0,12))
hist(chi2, xlim=c(0,12))
hist(chi3, xlim=c(0,12))
hist(chi4, xlim=c(0,12))
hist(chi5, xlim=c(0,12))

visualisasi histogram di atas menunjukan bahwa semakin besar ukuran contoh, semakin menyempit bar-bar histogram yang didapat yang mana seiring semakin kecil ragam pada data.

par(mfrow=c(2,3))
qqcomp(fitdist(chi1, "norm"), main = "n = 4")
qqcomp(fitdist(chi2, "norm"), main = "n = 12")
qqcomp(fitdist(chi3, "norm"), main = "n = 20")
qqcomp(fitdist(chi4, "norm"), main = "n = 60")
qqcomp(fitdist(chi5, "norm"), main = "n = 100")

Berdasarkan QQ-plot di atas, terlihat bahwa data dengan n = 4, n = 12, dan n =20 memiliki sebaran statistik yang tidak terlalu mendekati sebaran normal.**Sebaran mulai mendekati sebaran normal ketika n = 60 sampai 100 yang mana ukuran contoh semakin *besar**.

Simpulan

Sebaran rataan dengan \(n\)=4, \(n\)=12, dan \(n\)=20 belum mendekati normal, masih cenderung menjulur ke kanan. Sebaran rataan mendekati normal mulai dari \(n\)=60. Oleh karena itu, pada populasi menyebar tidak simetris (gabungan khi-kuadrat dan normal dengan beda parameter) membutuhkan sampel yang cukup banyak agar sebaran rataan mendekati normal.

Soal 1.2.c, 1.3, & 1.4

Data akan dibangkitkan dengan membangkitkan 25% data yangyang menyebar Khi-kuadrat(5), 25% data yang menyebar Khi-kuadrat(2), 25% data yang menyebar Normal(2,5), dan 25% data yang menyebar Normal(5,2) dengan 1 juta amatan data.Lalu, pengambilan contoh dilakukan dengan mengambil Sampel dengan n = 4, 12, 20, 60, 100 dan banyak ulangan (\(k\)) sejumlah 1000 dan ukuran populasi 4 juta.

popnchi <- c(rchisq(1000000, 5), rchisq(1000000, 2),
rnorm(1000000, 2, 5), rnorm(1000000, 5, 2))
hist(popnchi, main="Sebaran Populasi", xlab="Populasi", ylab="Frekuensi", col="coral")

set.seed(17)
k=1000
n1=4
n2=12
n3=20
n4=60
n5=100
nchi1 = matrix(sample(popnchi,n1*k),k)
nchi1 = apply(nchi1,1,mean)
nchi2 = matrix(sample(popnchi,n2*k),k)
nchi2 = apply(nchi2,1,mean)
nchi3 = matrix(sample(popnchi,n3*k),k)
nchi3 = apply(nchi3,1,mean)
nchi4 = matrix(sample(popnchi,n4*k),k)
nchi4 = apply(nchi4,1,mean)
nchi5 = matrix(sample(popnchi,n5*k),k)
nchi5 = apply(nchi5,1,mean)
#Histogram
par(mfrow=c(2,3))
hist(nchi1, xlim=c(0,12))
hist(nchi2, xlim=c(0,12))
hist(nchi3, xlim=c(0,12))
hist(nchi4, xlim=c(0,12))
hist(nchi5, xlim=c(0,12))

visualisasi histogram di atas menunjukan bahwa semakin besar ukuran contoh, semakin menyempit bar-bar histogram yang didapat yang mana seiring semakin kecil ragam pada data.

par(mfrow=c(2,3))
qqcomp(fitdist(nchi1, "norm"), main = "n = 4")
qqcomp(fitdist(nchi2, "norm"), main = "n = 12")
qqcomp(fitdist(nchi3, "norm"), main = "n = 20")
qqcomp(fitdist(nchi4, "norm"), main = "n = 60")
qqcomp(fitdist(nchi5, "norm"), main = "n = 100")

Berdasarkan QQ-plot di atas, terlihat bahwa data dengan n = 4, n = 12, n = 20 memiliki sebaran statistik yang tidak terlalu mendekati sebaran normal.**Sebaran mulai mendekati sebaran normal ketika n = 60 sampai 100 yang mana ukuran contoh semakin *besar**.

Simpulan

Sebaran rataan dengan \(n\)=4, \(n\)=12, dan \(n\)=20 belum mendekati normal, masih cenderung menjulur ke kanan. Sebaran rataan mendekati normal mulai dari \(n\)=60. Oleh karena itu, pada populasi menyebar tidak simetris (gabungan khi-kuadrat dan normal dengan beda parameter) membutuhkan sampel yang cukup banyak agar sebaran rataan mendekati normal.

Soal 2

Untuk mengetahui metode korelasi mana yang lebih baik dalam menangani data yang mengandung pencilan, peneliti perlu membandingkan kinerja korelasi Pearson dan korelasi Spearman untuk setiap ukuran sampel. Untuk melakukan hal ini, peneliti dapat mengikuti langkah-langkah sebagai berikut.

  1. Bangkitkan data berdasarkan distribusi normal ganda dengan vektor rata-rata (10,12) dan korelasi 0.85 serta vektor simpangan baku (2,3).

  2. Bangkitkan penciln pada kedua peubah dengan berdasarkan distribusi 𝑁(25,2) sebanyak 20% dari ukuran contoh pada langkah pertama. Ini bisa dilakukan dengan memilih 20% data dengan nilai paling ekstrim pada kedua peubah.

  3. Ganti nilai amatan dari hasil bangkitan pada langkah pertama dengan hasil bangkitan pada langkah kedua secara acak. Ini akan menciptakan variasi dalam data dan membuatnya lebih realistis.

  4. Bangkitkan data dengan ulangan 1000 dan ukuran sampel yang diteliti adalah 10,20,30,40,50,60,70,80,90,100.

  5. Untuk setiap ulangan, hitung nilai koefisien korelasi Pearson dan Spearman untuk sampel yang diambil dari populasi yang dibangkitkan pada langkah 3.

  6. Hitung rata-rata koefisien korelasi Pearson dan Spearman untuk setiap ukuran sampel.

  7. Bandingkan kinerja korelasi Pearson dan Spearman untuk setiap ukuran sampel dengan membandingkan rata-rata koefisien korelasi masing-masing metode.

Fungsi Perhitungan Korelasi

library(MASS)

# Membuat fungsi untuk menghasilkan data
generate_data <- function(n, mean_vec, sd_vec, corr, outlier_prop) {
  
  # Menghitung matriks kovarians
  cov_mat <- matrix(c(sd_vec[1]^2, sd_vec[1]*sd_vec[2]*corr,
                      sd_vec[1]*sd_vec[2]*corr, sd_vec[2]^2), ncol = 2)
  
  # Menghasilkan data multivariat normal tanpa outlier
  data <- mvrnorm(n = n, mu = mean_vec, Sigma = cov_mat)
  
  # Menghasilkan outlier pada kedua variabel dengan proporsi tertentu
  outlier_count <- round(n * outlier_prop)
  outlier_indices <- sample(1:n, outlier_count)
  outlier_values <- matrix(rnorm(outlier_count*2, mean = 25, sd = 2), ncol = 2)
  data[outlier_indices,] <- outlier_values
  
  return(data)
}

Langkah pertama, adalah menghasilkan data yang sesuai dengan spesifikasi pada soal. Pertama-tama, kita menentukan ukuran sampel n, vektor rata-rata mu, dan matriks kovarians sigma untuk distribusi normal ganda. Kemudian, kita menggunakan fungsi mvrnorm dari paket MASS untuk menghasilkan data berdasarkan distribusi tersebut. Untuk menambahkan pencilan, kita mengambil 20% dari ukuran sampel sebagai pencilan, dan menambahkannya pada kedua peubah dengan menggunakan fungsi rnorm.

# Menjalankan simulasi
n_samp <- c(10, 20, 30, 40, 50, 60, 70, 80, 90, 100) #ukuran sampel yang akan digunakan
n_iter <- 1000 #jumlah ulangan yang akan dilakukan
mean_vec <- c(10, 12) #vektor rata-rata
sd_vec <- c(2, 3) #vektor simpangan baku
corr <- 0.85 #korelasi
outlier_prop <- 0.2 #proporsi pencilan

hasil_pearson <- matrix(NA, nrow = n_iter, ncol = length(n_samp))
hasil_spearman <- matrix(NA, nrow = n_iter, ncol = length(n_samp))

set.seed(17)
for (n in n_samp)
  {
  for (i in 1:n_iter)
    {
    data <- generate_data(n, mean_vec, sd_vec, corr, outlier_prop)
    #menghitung koefisien korelasi Pearson
    hasil_pearson[i, which(n_samp == n)] <- cor(data[,1], data[,2], method = "pearson")
    #menghitung koefisien korelasi Spearman
    hasil_spearman[i, which(n_samp == n)] <- cor(data[,1], data[,2], method = "spearman")
  }
}

Langkah kedua, menghitung koefisien korelasi Pearson dan Spearman dari data yang telah diganti, dan menyimpannya dalam sebuah data frame (hasil).

Tabel

#melihat hasil
hasil_mean_pearson <- apply(hasil_pearson, 2, mean)
hasil_sd_pearson <- apply(hasil_pearson, 2, sd)

hasil_mean_spearman <- apply(hasil_spearman, 2, mean)
hasil_sd_spearman <- apply(hasil_spearman, 2, sd)

hasil <- data.frame(n = n_samp,
                      mean_pearson = hasil_mean_pearson,
                      sd_pearson = hasil_sd_pearson,
                      mean_spearman = hasil_mean_spearman,
                      sd_spearman = hasil_sd_spearman)

hasil
##      n mean_pearson  sd_pearson mean_spearman sd_spearman
## 1   10    0.9553862 0.025083005     0.8752121  0.08971725
## 2   20    0.9483987 0.019626504     0.8927624  0.05318404
## 3   30    0.9464776 0.016632014     0.8980529  0.03975048
## 4   40    0.9454531 0.013932241     0.9003583  0.03437621
## 5   50    0.9444501 0.013612287     0.9011011  0.02948344
## 6   60    0.9444776 0.011589401     0.9040982  0.02694050
## 7   70    0.9443745 0.010866310     0.9045666  0.02436120
## 8   80    0.9441324 0.010071923     0.9059611  0.02202668
## 9   90    0.9432423 0.009586271     0.9053540  0.02170920
## 10 100    0.9438463 0.008854110     0.9054605  0.01980795

Visualisasi

plot(hasil$n, hasil$mean_pearson, type = "b", ylim = c(0.8, 1), xlab = "Sample size", ylab = "Correlation", main = "Perbandingan Korelasi Pearson dan Spearman", col="navy")
lines(hasil$n, hasil$mean_spearman, type = "b", ylim = c(0.8, 1), col="red")
legend("bottomright", c("Pearson", "Spearman"), lty = 1, col = c("navy", "red"))

Simpulan

Berdasarkan hasil perbandingan kinerja korelasi Pearson dan Spearman di atas, metode korelasi yang lebih baik untuk menangani data yang mengandung pencilan adalah metode korelasi Spearman karena nilai korelasinya mendekati nilai korelasi yang telah ditentukan yaitu 0.85. Hal ini sesuai dengan visualisasi plot perbandingan di atas. Vusvitasari et al. (2008) menyatakan bahwa metode korelasi Spearman memang lebih baik digunakan untuk data yang tidak simetris (mengandung pencilan), sedangkan metode korelasi Pearson lebih baik digunakan untuk data normal (simetris).

Daftar Pustaka

Vusvitasari R, Nugroho S, Akbar S. 2008. Kajian hubungan koefisien korelasi Pearson, Spearman-Rho, Kendall-Tau, Gamma, dan Somers. e-Jurnal Statistika. 4(1): 40-54.


TUGAS PRAKTIKUM P11

Soal Tugas 1

Jika data diasumsikan tidak normal bagaimana ukuran sampel berpengaruh kesalahan tipe 1 pada uji hipotesis pada Latihan 1 ? Silahkan buat simulasi berdaasarkan asumsi sebaran berikut:

  1. 50% sebaran normal N(0,1) dan 50% sebaran chi-square χ2(df=1)

  2. 50% sebaran chi-square χ2(df=1) dan 50% dengan chi-square χ2(df=3)

  3. 25% sebaran chi-square χ2(df=1) . 25% dengan chi-square χ2(df=3) , 25% sebaran normal N(0,1) dan 25% sebaran normal N(3,1.5)

Pembahasan Tugas 1

Langkah pertama, yaitu mendefinisikan fungsi untuk mengextract p-value dari uji-t dan uji wilcoxon.

library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.2.2
set.seed(26)
# Fungsi untuk menghitung p-value dari uji t
t_test_p_value <- function(data) {
  t.test(data)$p.value
}
# Fungsi untuk menghitung p-value dari uji wilcoxon
wilcox_test_p_value <- function(data) {
  wilcox.test(data)$p.value
}

Kemudian, kita akan menjalankan simulasi untuk masing-masing kasus, yaitu kasus 1, kasus 2, dan kasus 3 dengan ulangan sebanyak 1000 kali.

T-test

Kasus 1: Gabungan dari 50% sebaran normal N(0,1) dan 50% sebaran chi-square χ2(df=1)

t_test_case1 <- function(sample_size) {
  normal_size <- sample_size %/% 2
  chi_size <- sample_size - normal_size
  p_values <- replicate(1000, {
    normal_data <- rnorm(normal_size)
    chi_squared_data <- rchisq(chi_size, df = 1)
    data <- c(normal_data, chi_squared_data)
    t_test_p_value(data)
  })
  mean(p_values)
}

Kasus 2: Gabungan dari 50% sebaran chi-square χ2(df=1) dan 50% sebaran chi-square χ2(df=3)

t_test_case2 <- function(sample_size) {
  chi_size <- sample_size %/% 2
  p_values <- replicate(1000, {
    chi_squared_1_data <- rchisq(chi_size, df = 1)
    chi_squared_3_data <- rchisq(chi_size, df = 3)
    data <- c(chi_squared_1_data, chi_squared_3_data)
    t_test_p_value(data)
  })
  mean(p_values)
}

Kasus 3: Gabungan dari 25% sebaran chi-square χ2(df=1), 25% sebaran chi-square χ2(df=3), 25% sebaran normal N(0,1), dan 25% sebaran normal N(3,1.5)

t_test_case3 <- function(sample_size) {
  chi_size <- sample_size %/% 4
  normal_size <- sample_size - 3 * chi_size
  p_values <- replicate(1000, {
    chi_squared_1_data <- rchisq(chi_size, df = 1)
    chi_squared_3_data <- rchisq(chi_size, df = 3)
    normal_1_data <- rnorm(normal_size)
    normal_2_data <- rnorm(normal_size, mean = 3, sd = 1.5)
    data <- c(chi_squared_1_data, chi_squared_3_data, normal_1_data, normal_2_data)
    t_test_p_value(data)
  })
  mean(p_values)
}

Mendefinisikan ukuran sampel dan menghitung p-value untuk setiap ukuran sampel.

# Range ukuran sampel yang akan diuji
sample_sizes <- seq(5, 60)

# Menghitung p-value untuk setiap ukuran sampel
t_test_p_values_case1 <- sapply(sample_sizes, t_test_case1)
t_test_p_values_case2 <- sapply(sample_sizes, t_test_case2)
t_test_p_values_case3 <- sapply(sample_sizes, t_test_case3)

Langkah terakhir, menampilkan plot t-test antara nilai-p dan ukuran sampel.

Plot T-test

# Menampilkan plot antara nilai p dan ukuran sampel
plot_data_t_test <- data.frame(
  Sample_Size = rep(sample_sizes, 3),
  P_Value = c(t_test_p_values_case1, t_test_p_values_case2, t_test_p_values_case3),
  Case = rep(c("Kasus Sebaran 1", "Kasus Sebaran 2", "Kasus Sebaran 3"), 
             each = length(sample_sizes))
)

ggplot(plot_data_t_test, aes(x = Sample_Size, y = P_Value, color = Case)) +
  geom_point() +
  geom_hline(yintercept = 0.05, linetype = "dashed", color = "red") +
  labs(title="Uji-t", x = "Ukuran Sampel", y = "Nilai p") +
  scale_color_discrete(name = "Kasus Sebaran") +
  theme_minimal()

Wilcoxon test

Menjalankan simulasi untuk masing-masing kasus, yaitu kasus 1, kasus 2, dan kasus 3 dengan ulangan sebanyak 1000 kali.

Kasus 1: Gabungan dari 50% sebaran normal N(0,1) dan 50% sebaran chi-square χ2(df=1)

wilcox_test_case1 <- function(sample_size) {
  normal_size <- sample_size %/% 2
  chi_size <- sample_size - normal_size
  p_values <- replicate(1000, {
    normal_data <- rnorm(normal_size)
    chi_squared_data <- rchisq(chi_size, df = 1)
    data <- c(normal_data, chi_squared_data)
    wilcox_test_p_value(data)
  })
  mean(p_values)
}

Kasus 2: Gabungan dari 50% sebaran chi-square χ2(df=1) dan 50% sebaran chi-square χ2(df=3)

wilcox_test_case2 <- function(sample_size) {
  chi_size <- sample_size %/% 2
  p_values <- replicate(1000, {
    chi_squared_1_data <- rchisq(chi_size, df = 1)
    chi_squared_3_data <- rchisq(chi_size, df = 3)
    data <- c(chi_squared_1_data, chi_squared_3_data)
    wilcox_test_p_value(data)
  })
  mean(p_values)
}

Kasus 3: Gabungan dari 25% sebaran chi-square χ2(df=1), 25% sebaran chi-square χ2(df=3), 25% sebaran normal N(0,1), dan 25% sebaran normal N(3,1.5)

wilcox_test_case3 <- function(sample_size) {
  chi_size <- sample_size %/% 4
  normal_size <- sample_size - 3 * chi_size
  p_values <- replicate(1000, {
    chi_squared_1_data <- rchisq(chi_size, df = 1)
    chi_squared_3_data <- rchisq(chi_size, df = 3)
    normal_1_data <- rnorm(normal_size)
    normal_2_data <- rnorm(normal_size, mean = 3, sd = 1.5)
    data <- c(chi_squared_1_data, chi_squared_3_data, normal_1_data, normal_2_data)
    wilcox_test_p_value(data)
  })
  mean(p_values)
}

Mendefinisikan ukuran sampel dan menghitung p-value untuk setiap ukuran sampel.

# Range ukuran sampel yang akan diuji
sample_sizes <- seq(5, 60)

# Menghitung p-value untuk setiap ukuran sampel
wilcox_test_p_values_case1 <- sapply(sample_sizes, wilcox_test_case1)
wilcox_test_p_values_case2 <- sapply(sample_sizes, wilcox_test_case2)
wilcox_test_p_values_case3 <- sapply(sample_sizes, wilcox_test_case3)

Menampilkan plot t-test antara nilai-p dan ukuran sampel.

Plot Wilcoxon test

# Menampilkan plot antara nilai p dan ukuran sampel
plot_data_wilcox_test <- data.frame(
  Sample_Size = rep(sample_sizes, 3),
  P_Value = c(wilcox_test_p_values_case1, wilcox_test_p_values_case2, wilcox_test_p_values_case3),
  Case = rep(c("Kasus Sebaran 1", "Kasus Sebaran 2", "Kasus Sebaran 3"), 
             each = length(sample_sizes))
)

plot_data_wilcox_test
##     Sample_Size      P_Value            Case
## 1             5 3.744375e-01 Kasus Sebaran 1
## 2             6 4.153125e-01 Kasus Sebaran 1
## 3             7 3.327969e-01 Kasus Sebaran 1
## 4             8 3.473281e-01 Kasus Sebaran 1
## 5             9 2.741914e-01 Kasus Sebaran 1
## 6            10 2.935605e-01 Kasus Sebaran 1
## 7            11 2.343154e-01 Kasus Sebaran 1
## 8            12 2.550723e-01 Kasus Sebaran 1
## 9            13 2.072598e-01 Kasus Sebaran 1
## 10           14 2.117366e-01 Kasus Sebaran 1
## 11           15 1.688126e-01 Kasus Sebaran 1
## 12           16 1.796122e-01 Kasus Sebaran 1
## 13           17 1.490793e-01 Kasus Sebaran 1
## 14           18 1.616422e-01 Kasus Sebaran 1
## 15           19 1.325496e-01 Kasus Sebaran 1
## 16           20 1.471109e-01 Kasus Sebaran 1
## 17           21 1.209594e-01 Kasus Sebaran 1
## 18           22 1.306836e-01 Kasus Sebaran 1
## 19           23 9.826032e-02 Kasus Sebaran 1
## 20           24 1.085013e-01 Kasus Sebaran 1
## 21           25 9.310902e-02 Kasus Sebaran 1
## 22           26 8.987604e-02 Kasus Sebaran 1
## 23           27 7.254666e-02 Kasus Sebaran 1
## 24           28 8.180967e-02 Kasus Sebaran 1
## 25           29 7.685793e-02 Kasus Sebaran 1
## 26           30 7.963752e-02 Kasus Sebaran 1
## 27           31 6.739056e-02 Kasus Sebaran 1
## 28           32 6.390754e-02 Kasus Sebaran 1
## 29           33 5.889891e-02 Kasus Sebaran 1
## 30           34 5.907191e-02 Kasus Sebaran 1
## 31           35 4.896098e-02 Kasus Sebaran 1
## 32           36 5.062560e-02 Kasus Sebaran 1
## 33           37 4.143760e-02 Kasus Sebaran 1
## 34           38 4.567775e-02 Kasus Sebaran 1
## 35           39 4.045663e-02 Kasus Sebaran 1
## 36           40 3.799917e-02 Kasus Sebaran 1
## 37           41 3.645818e-02 Kasus Sebaran 1
## 38           42 3.604625e-02 Kasus Sebaran 1
## 39           43 2.658114e-02 Kasus Sebaran 1
## 40           44 3.352894e-02 Kasus Sebaran 1
## 41           45 2.692554e-02 Kasus Sebaran 1
## 42           46 2.823673e-02 Kasus Sebaran 1
## 43           47 2.052139e-02 Kasus Sebaran 1
## 44           48 2.154894e-02 Kasus Sebaran 1
## 45           49 2.017880e-02 Kasus Sebaran 1
## 46           50 2.196888e-02 Kasus Sebaran 1
## 47           51 1.672970e-02 Kasus Sebaran 1
## 48           52 2.294338e-02 Kasus Sebaran 1
## 49           53 1.567569e-02 Kasus Sebaran 1
## 50           54 1.610306e-02 Kasus Sebaran 1
## 51           55 1.311487e-02 Kasus Sebaran 1
## 52           56 1.515095e-02 Kasus Sebaran 1
## 53           57 1.169806e-02 Kasus Sebaran 1
## 54           58 1.654558e-02 Kasus Sebaran 1
## 55           59 1.152220e-02 Kasus Sebaran 1
## 56           60 1.161815e-02 Kasus Sebaran 1
## 57            5 1.250000e-01 Kasus Sebaran 2
## 58            6 3.125000e-02 Kasus Sebaran 2
## 59            7 3.125000e-02 Kasus Sebaran 2
## 60            8 7.812500e-03 Kasus Sebaran 2
## 61            9 7.812500e-03 Kasus Sebaran 2
## 62           10 1.953125e-03 Kasus Sebaran 2
## 63           11 1.953125e-03 Kasus Sebaran 2
## 64           12 4.882813e-04 Kasus Sebaran 2
## 65           13 4.882813e-04 Kasus Sebaran 2
## 66           14 1.220703e-04 Kasus Sebaran 2
## 67           15 1.220703e-04 Kasus Sebaran 2
## 68           16 3.051758e-05 Kasus Sebaran 2
## 69           17 3.051758e-05 Kasus Sebaran 2
## 70           18 7.629395e-06 Kasus Sebaran 2
## 71           19 7.629395e-06 Kasus Sebaran 2
## 72           20 1.907349e-06 Kasus Sebaran 2
## 73           21 1.907349e-06 Kasus Sebaran 2
## 74           22 4.768372e-07 Kasus Sebaran 2
## 75           23 4.768372e-07 Kasus Sebaran 2
## 76           24 1.192093e-07 Kasus Sebaran 2
## 77           25 1.192093e-07 Kasus Sebaran 2
## 78           26 2.980232e-08 Kasus Sebaran 2
## 79           27 2.980232e-08 Kasus Sebaran 2
## 80           28 7.450581e-09 Kasus Sebaran 2
## 81           29 7.450581e-09 Kasus Sebaran 2
## 82           30 1.862645e-09 Kasus Sebaran 2
## 83           31 1.862645e-09 Kasus Sebaran 2
## 84           32 4.656613e-10 Kasus Sebaran 2
## 85           33 4.656613e-10 Kasus Sebaran 2
## 86           34 1.164153e-10 Kasus Sebaran 2
## 87           35 1.164153e-10 Kasus Sebaran 2
## 88           36 2.910383e-11 Kasus Sebaran 2
## 89           37 2.910383e-11 Kasus Sebaran 2
## 90           38 7.275958e-12 Kasus Sebaran 2
## 91           39 7.275958e-12 Kasus Sebaran 2
## 92           40 1.818989e-12 Kasus Sebaran 2
## 93           41 1.818989e-12 Kasus Sebaran 2
## 94           42 4.547474e-13 Kasus Sebaran 2
## 95           43 4.547474e-13 Kasus Sebaran 2
## 96           44 1.136868e-13 Kasus Sebaran 2
## 97           45 1.136868e-13 Kasus Sebaran 2
## 98           46 2.842171e-14 Kasus Sebaran 2
## 99           47 2.842171e-14 Kasus Sebaran 2
## 100          48 7.105427e-15 Kasus Sebaran 2
## 101          49 7.105427e-15 Kasus Sebaran 2
## 102          50 7.790492e-10 Kasus Sebaran 2
## 103          51 7.790492e-10 Kasus Sebaran 2
## 104          52 3.607881e-10 Kasus Sebaran 2
## 105          53 3.607881e-10 Kasus Sebaran 2
## 106          54 1.672125e-10 Kasus Sebaran 2
## 107          55 1.672125e-10 Kasus Sebaran 2
## 108          56 7.755146e-11 Kasus Sebaran 2
## 109          57 7.755146e-11 Kasus Sebaran 2
## 110          58 3.599102e-11 Kasus Sebaran 2
## 111          59 3.599102e-11 Kasus Sebaran 2
## 112          60 1.671329e-11 Kasus Sebaran 2
## 113           5 1.479063e-01 Kasus Sebaran 3
## 114           6 1.048672e-01 Kasus Sebaran 3
## 115           7 6.605664e-02 Kasus Sebaran 3
## 116           8 6.029688e-02 Kasus Sebaran 3
## 117           9 4.211914e-02 Kasus Sebaran 3
## 118          10 3.220068e-02 Kasus Sebaran 3
## 119          11 2.166577e-02 Kasus Sebaran 3
## 120          12 1.728271e-02 Kasus Sebaran 3
## 121          13 1.123950e-02 Kasus Sebaran 3
## 122          14 9.391205e-03 Kasus Sebaran 3
## 123          15 6.368767e-03 Kasus Sebaran 3
## 124          16 5.632294e-03 Kasus Sebaran 3
## 125          17 4.500755e-03 Kasus Sebaran 3
## 126          18 3.130062e-03 Kasus Sebaran 3
## 127          19 2.307158e-03 Kasus Sebaran 3
## 128          20 1.768503e-03 Kasus Sebaran 3
## 129          21 1.379830e-03 Kasus Sebaran 3
## 130          22 9.885370e-04 Kasus Sebaran 3
## 131          23 6.731114e-04 Kasus Sebaran 3
## 132          24 5.243133e-04 Kasus Sebaran 3
## 133          25 5.202566e-04 Kasus Sebaran 3
## 134          26 4.085838e-04 Kasus Sebaran 3
## 135          27 1.939399e-04 Kasus Sebaran 3
## 136          28 2.182479e-04 Kasus Sebaran 3
## 137          29 1.798101e-04 Kasus Sebaran 3
## 138          30 9.812944e-05 Kasus Sebaran 3
## 139          31 8.672011e-05 Kasus Sebaran 3
## 140          32 5.214805e-05 Kasus Sebaran 3
## 141          33 5.918846e-05 Kasus Sebaran 3
## 142          34 3.198516e-05 Kasus Sebaran 3
## 143          35 4.066285e-05 Kasus Sebaran 3
## 144          36 2.442193e-05 Kasus Sebaran 3
## 145          37 1.532228e-05 Kasus Sebaran 3
## 146          38 1.321207e-05 Kasus Sebaran 3
## 147          39 6.835600e-06 Kasus Sebaran 3
## 148          40 9.658354e-06 Kasus Sebaran 3
## 149          41 4.974081e-06 Kasus Sebaran 3
## 150          42 3.243838e-06 Kasus Sebaran 3
## 151          43 2.479970e-06 Kasus Sebaran 3
## 152          44 2.840447e-06 Kasus Sebaran 3
## 153          45 1.585269e-06 Kasus Sebaran 3
## 154          46 1.203878e-06 Kasus Sebaran 3
## 155          47 2.668580e-06 Kasus Sebaran 3
## 156          48 9.055343e-07 Kasus Sebaran 3
## 157          49 1.996906e-06 Kasus Sebaran 3
## 158          50 1.613563e-06 Kasus Sebaran 3
## 159          51 1.184803e-06 Kasus Sebaran 3
## 160          52 1.287310e-06 Kasus Sebaran 3
## 161          53 7.563990e-07 Kasus Sebaran 3
## 162          54 6.882095e-07 Kasus Sebaran 3
## 163          55 6.147097e-07 Kasus Sebaran 3
## 164          56 4.108646e-07 Kasus Sebaran 3
## 165          57 3.844906e-07 Kasus Sebaran 3
## 166          58 2.999723e-07 Kasus Sebaran 3
## 167          59 2.230644e-07 Kasus Sebaran 3
## 168          60 1.513368e-07 Kasus Sebaran 3
ggplot(plot_data_wilcox_test, aes(x = Sample_Size, y = P_Value, color = Case)) +
  geom_point() +
  geom_hline(yintercept = 0.05, linetype = "dashed", color = "red") +
  labs(title="Uji Peringkat Bertanda Wilcoxon", x = "Ukuran Sampel", y = "Nilai p") +
  scale_color_discrete(name = "Kasus Sebaran") +
  theme_minimal()

Berdasarkan kedua grafik di atas dapat disimpulkan bahwa semakin besar ukuran sampel, kedua metode uji (t dan wilcoxon) cenderung memiliki alpha empiris yang lebih kecil dibandingkan alpha 5% untuk sebaran ketiga kasus. Hal ini menunjukkan bahwa kedua metode uji cenderung menolak H0 walaupun sebenarnya H0 benar. Kasus 1 nilai alpha < 5% saat ukuran sampel di atas 40, sementara untuk kasus 2 nilai alpha < 5% saat ukuran sampel berada di atas 5.

Soal Tugas 2

Buatlah simulasi statistika untuk menunjukkan pengaruh banyaknya amatan dan tidak terpenuhinya asumsi sebaran normal terhadap pengujian hipotesis di one-way ANOVA! Cobalah dengan sebaran yang ada di Tugas 1.

  1. 50% sebaran normal N(0,1) dan 50% sebaran chi-square χ2(df=1)

  2. 50% sebaran chi-square χ2(df=1) dan 50% dengan chi-square χ2(df=3)

  3. 25% sebaran chi-square χ2(df=1). 25% dengan chi-square χ2(df=3), 25% sebaran normal N(0,1) dan 25% sebaran normal N(3,1.5)

(Hint: model linear padarancangan acak lengkap untuk membangkitkan data dalam simulasi ini)

library(ggplot2)

Pembahasan Tugas 2

50% sebaran normal N(0,1) dan 50% sebaran chi-square χ2(df=1)

set.seed(26)
  1. Menentukan jumlah simulasi yang akan dilakukan, daftar ukuran sampel yang akan digunakan
num_simulations <- 1000
sample_sizes <- seq(10, 300, 10)
  1. Melakukan looping untuk setiap ukuran sampel.Dalam looping ukuran sampel, dilakukan looping untuk jumlah simulasi.

Dalam loop ini akan dibangkitkan data dari sebaran yang telah ditentukan dengan proporsi yang sesuai dan menggunakan model linear pada rancangan acak lengkap. Lalu melakukan pengujian hipotesis dengan one-way ANOVA dan menghasilkan nilai-p untuk setiap ukuran sampel dan nomor simulasi. Setelah didapatkan semua nilai-p, rata-rata nilai-p dihitung untuk setiap ukuran sampel.

p_values <- numeric(length(sample_sizes))

for (i in seq_along(sample_sizes)) {
  p_value_sum <- 0
  for (j in seq_len(num_simulations)) {
    # Bangkitkan data dari sebaran normal dan chi-square
    normal_data <- rnorm(sample_sizes[i]/2, 0, 1)
    chi2_data <- rchisq(sample_sizes[i]/2, df = 1)
    data <- c(normal_data, chi2_data)
    
    # Menggunakan model linear pada rancangan acak lengkap
    group <- rep(1:2, each = sample_sizes[i]/2)
    
    # Melakukan one-way ANOVA
    anova_res <- aov(data ~ group)
    p_value <- summary(anova_res)[[1]][["Pr(>F)"]][1]
    
    p_value_sum <- p_value_sum + p_value
  }
  
  # Menghitung rata-rata nilai-p
  avg_p_value <- p_value_sum / num_simulations
  p_values[i] <- avg_p_value
}
  1. Menampilkan plot rata-rata nilai-p terhadap ukuran sampel
df1 <- data.frame(sample_sizes, p_values)

plot1 <- ggplot(df1, aes(x = sample_sizes, y = p_values)) +
  geom_point(alpha = 0.5, size=2, color="blue") +
  geom_hline(yintercept = 0.05, linetype = "dashed", color = "red") +
  labs(x = "Ukuran Sampel", y = "Nilai-p") +
  ggtitle("Pengaruh Ukuran Sampel terhadap Nilai-p") +
  theme_minimal()
plot1

50% sebaran chi-square χ2(df=1) dan 50% dengan chi-square χ2(df=3)

set.seed(26)
  1. Menentukan jumlah simulasi yang akan dilakukan, daftar ukuran sampel yang akan digunakan
num_simulations <- 1000
sample_sizes <- seq(10, 300, 10)
  1. Melakukan looping untuk setiap ukuran sampel.Dalam looping ukuran sampel, dilakukan looping untuk jumlah simulasi.

Dalam loop ini akan dibangkitkan data dari sebaran yang telah ditentukan dengan proporsi yang sesuai dan menggunakan model linear pada rancangan acak lengkap. Lalu melakukan pengujian hipotesis dengan one-way ANOVA dan menghasilkan nilai-p untuk setiap ukuran sampel dan nomor simulasi. Setelah didapatkan semua nilai-p, rata-rata nilai-p dihitung untuk setiap ukuran sampel.

p_values <- numeric(length(sample_sizes))

for (i in seq_along(sample_sizes)) {
  p_value_sum <- 0
  for (j in seq_len(num_simulations)) {
    # Bangkitkan data dari sebaran chi-square df=1 dan df=3
    chi2_data_df1 <- rchisq(sample_sizes[i]/2, df = 1)
    chi2_data_df3 <- rchisq(sample_sizes[i]/2, df = 3)
    data <- c(chi2_data_df1, chi2_data_df3)
    
    # Menggunakan model linear pada rancangan acak lengkap
    group <- rep(1:2, each = sample_sizes[i]/2)
    
    # Melakukan one-way ANOVA
    anova_res <- aov(data ~ group)
    p_value <- summary(anova_res)[[1]][["Pr(>F)"]][1]
    
    p_value_sum <- p_value_sum + p_value
  }
  
  # Menghitung rata-rata nilai-p
  avg_p_value <- p_value_sum / num_simulations
  p_values[i] <- avg_p_value
}
  1. Menampilkan plot rata-rata nilai-p terhadap ukuran sampel
df2 <- data.frame(sample_sizes, p_values)

plot2 <- ggplot(df2, aes(x = sample_sizes, y = p_values)) +
  geom_point(alpha = 0.5, size=2, color="blue") +
  geom_hline(yintercept = 0.05, linetype = "dashed", color = "red") +
  labs(x = "Ukuran Sampel", y = "Nilai-p") +
  ggtitle("Pengaruh Ukuran Sampel terhadap Nilai-p") +
  theme_minimal()
plot2

25% sebaran chi-square χ2(df=1). 25% dengan chi-square χ2(df=3), 25% sebaran normal N(0,1) dan 25% sebaran normal N(3,1.5)

set.seed(26)
  1. Menentukan jumlah simulasi yang akan dilakukan, daftar ukuran sampel yang akan digunakan
num_simulations <- 1000
sample_sizes <- seq(10, 300, 10)
  1. Melakukan looping untuk setiap ukuran sampel.Dalam looping ukuran sampel, dilakukan looping untuk jumlah simulasi.

Dalam loop ini akan dibangkitkan data dari sebaran yang telah ditentukan dengan proporsi yang sesuai dan menggunakan model linear pada rancangan acak lengkap. Lalu melakukan pengujian hipotesis dengan one-way ANOVA dan menghasilkan nilai-p untuk setiap ukuran sampel dan nomor simulasi. Setelah didapatkan semua nilai-p, rata-rata nilai-p dihitung untuk setiap ukuran sampel.

p_values <- numeric(length(sample_sizes))

for (i in seq_along(sample_sizes)) {
  p_value_sum <- 0
  for (j in seq_len(num_simulations)) {
    # Bangkitkan data dari sebaran chi-square df=1, df=3, normal N(0,1), dan normal N(3,1.5)
    chi2_data_df1 <- rchisq(sample_sizes[i]/4, df = 1)
    chi2_data_df3 <- rchisq(sample_sizes[i]/4, df = 3)
    normal_data1 <- rnorm(sample_sizes[i]/4, 0, 1)
    normal_data2 <- rnorm(sample_sizes[i]/4, 3, 1.5)
    data <- c(chi2_data_df1, chi2_data_df3, normal_data1, normal_data2)
    
    # Menggunakan model linear pada rancangan acak lengkap
    group <- rep(1:4, each = sample_sizes[i]/4)
    
    # Melakukan one-way ANOVA
    anova_res <- aov(data ~ group)
    p_value <- summary(anova_res)[[1]][["Pr(>F)"]][1]
    
    p_value_sum <- p_value_sum + p_value
  }
  
  # Menghitung rata-rata nilai-p
  avg_p_value <- p_value_sum / num_simulations
  p_values[i] <- avg_p_value
}
  1. Menampilkan plot rata-rata nilai-p terhadap ukuran sampel
df3 <- data.frame(sample_sizes, p_values)
plot3 <- ggplot(df3, aes(x = sample_sizes, y = p_values)) +
  geom_point(alpha = 0.5, size=2, color="blue") +
  geom_hline(yintercept = 0.05, linetype = "dashed", color = "red") +
  labs(x = "Ukuran Sampel", y = "Nilai-p") +
  ggtitle("Pengaruh Ukuran Sampel terhadap Nilai-p") +
  theme_minimal()
plot3

Membandingkan Plot

plot1

plot2

plot3

Berdasarkan hasil simulasi statistika yang dilakukan dapat dilihat bahwa dengan meningkatnya ukuran sampel, nilai-p cenderung menurun atau menolak H0. Artinya, semakin besar sample size, semakin baik untuk mendeteksi perbedaan signifikan antara kelompok-kelompok dalam one-way ANOVA.Lalu ketika asumsi sebaran normal tidak terpenuhi namun ukuran sampel besar, hasil pengujian hipotesis cenderung menolak H0. Untuk sebaran dari 25% sebaran chi-square χ2(df=1). 25% dengan chi-square χ2(df=3), 25% sebaran normal N(0,1) dan 25% sebaran normal N(3,1.5), perlu ukuran sampel yang sangat besar agar pengujian hipotesis menghasilkan tolak H0.

Soal Tugas 3

Buatlah langkah-langkah pembangkitan data agar bisa digunakan untuk menyelidiki asumsi ragam homogen dalam model regresi! Realiasikan langkah-langkahmu tersebut dalam R!

Pembahasan Tugas 3

Dengan metode regresi linier sederhana pendekatan MKT menggunakan bantuan software R Studio, agloritma pembangkitan data adalah sebagai berikut.

  1. Menentukan peubah penjelas (\(x\)) dan peubah respon (\(y\)) dengan masing-masing sebaran. Sebaran kedua peubah ini menggunakan Normal(0,1) dengan ukuran contoh sebesar 100.

# Library yang digunakan
library(ggplot2)
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.2.2
## 
## Attaching package: 'dplyr'
## The following object is masked from 'package:MASS':
## 
##     select
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(car)   # Uji Levene 
## Warning: package 'car' was built under R version 4.2.2
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.2.2
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
library(stats) # Uji Bartlett
set.seed(26)  # Untuk memastikan hasil yang konsisten
x <- rnorm(100, mean = 0, sd = 1)
y.duga <- x + rnorm(100, mean = 0, sd = 1)  

  1. Membuat kerangka data dari peubah \(x\) dan \(y\).

data3 <- data.frame(x, y.duga)
head(data3)
##            x     y.duga
## 1 -2.1298417 -2.5703169
## 2  1.1478961  2.0385625
## 3 -0.4895019 -0.5533213
## 4  0.8263438  0.3239139
## 5 -0.4099352 -0.8642651
## 6  0.1487879 -0.4410678

  1. Membuat visualisasi berupa scatter plot untuk melihat hubungan \(x\) dan \(y\).

ggplot(data3, aes(x = x, y = y.duga)) +
  geom_point(shape=18, color="blue") +
  labs(x = "Peubah Penjelas (x)", y = "Peubah Respon (y)") +
  theme_minimal()+
  geom_smooth(method=lm, se=FALSE, linetype="dashed", color="darkred")
## `geom_smooth()` using formula = 'y ~ x'

Berdasarkan scatter plot di atas, dapat disimpulkan bahwa hubungan \(x\) dan \(y\) relatif positif kuat.

  1. Melakukan pembagian data menjadi duah kelompok berdasakan nilai \(x\).

data.new <- data3 %>%
  mutate(Group = ifelse(x < 0, "Kelompok A", "Kelompok B"))

  1. Melakukan visualisasi perbedaan ragam (variability) antara kedua kelompok.

ggplot(data.new, aes(x = Group, y = y.duga, fill = Group)) +
  geom_boxplot() +
  labs(x = " ", y = "Peubah Respon (y)") +
  theme_minimal()

Keragaman data kelompok B lebih tinggi dibandingkan kelompok A.

  1. Melakukan pengujian asumsi ragam homogen menggunakan uji Levene dan uji Bartlett.
    Hipotesis yang diuji adalah sebagai berikut.
    \(H_0\) : Sisaan memiliki ragam homogen
    \(H_1\) : Sisaan tidak memiliki ragam homogen
    Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

levene_test <- leveneTest(y.duga ~ Group, data = data.new)
## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.
bartlett_test <- bartlett.test(y.duga ~ Group, data = data.new)

levene_test
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  1   0.622 0.4322
##       98
bartlett_test
## 
##  Bartlett test of homogeneity of variances
## 
## data:  y.duga by Group
## Bartlett's K-squared = 0.75991, df = 1, p-value = 0.3834

Berdasarkan kedua uji formal di atas, terlihat bahwa sebaran data pada model regresi linier sederhana dapat disimpulkan bahwa sisaan memiliki ragam homogen pada taraf nyata 5%.