Latihan Praktikum 12 - AED

Library

library(readxl) 

## Warning: package 'readxl' was built under R version 4.1.3

library(psych)

## Warning: package 'psych' was built under R version 4.1.3

library(lmomco) 

## Warning: package 'lmomco' was built under R version 4.1.3

## 
## Attaching package: 'lmomco'

## The following object is masked from 'package:psych':
## 
##     harmonic.mean

Input Data

datap12 <- read_excel("C:/Users/ACER/Downloads/Data praktikum 12.xlsx") 
datap12.2 <- read_excel("C:/Users/ACER/Downloads/Data praktikum 12.xlsx",sheet=2)
View(datap12)
View(datap12.2)

Uji Sebaran Data (Uji Normalitas)

# uji normalitas y1

qqnorm(datap12$y1)
qqline(datap12$y1, col = "red")

ks.test(datap12$y1, "pnorm", mean = mean(datap12$y1), sd = sd(datap12$y1))

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datap12$y1
## D = 0.085937, p-value = 0.4512
## alternative hypothesis: two-sided

Didapatkan p-value > alpha = 0.05 sehingga tak tolak H0. Artinya, cukup bukti untuk menyatakan bahwa peubah tunggal y1 menyebar normal.

# uji normalitas y2

ks.test(datap12$y2, "pnorm", mean = mean(datap12$y2), sd = sd(datap12$y2))

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datap12$y2
## D = 0.16786, p-value = 0.007139
## alternative hypothesis: two-sided

# uji lognormal

ks.test(datap12$y2, "plnorm", meanlog = mean(log(datap12$y2), na.rm = T), sdlog = sd(log(datap12$y2)))

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datap12$y2
## D = 0.060634, p-value = 0.8558
## alternative hypothesis: two-sided

# uji chisq

ks.test(datap12$y2, "pchisq", df = mean(datap12$y2))

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datap12$y2
## D = 0.095563, p-value = 0.3206
## alternative hypothesis: two-sided

Peubah tunggal y2 tidak menyebar normal. Berdasarkan ketiga uji sebaran teoritis di atas, didapatkan bahwa tak tolak H0 ketika menggunakan sebaran teoritis lognormal dan chi-square.

# uji normalitas x1

qqnorm(datap12.2$x1)
qqline(datap12.2$x1, col = "red")

ks.test(datap12.2$x1, "pnorm", mean = mean(datap12.2$x1), sd = sd(datap12.2$x1))

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datap12.2$x1
## D = 0.20162, p-value = 0.7403
## alternative hypothesis: two-sided

Didapatkan p-value > alpha = 0.05 sehingga tak tolak H0. Artinya, cukup bukti untuk menyatakan bahwa peubah x1 menyebar normal.

# uji normalitas x2

qqnorm(datap12.2$x2)
qqline(datap12.2$x2, col = "red")

ks.test(datap12.2$x2, "pnorm", mean = mean(datap12.2$x2), sd = sd(datap12.2$x2))

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datap12.2$x2
## D = 0.21956, p-value = 0.6449
## alternative hypothesis: two-sided

Didapatkan p-value > alpha = 0.05 sehingga tak tolak H0. Artinya, cukup bukti untuk menyatakan bahwa peubah x2 menyebar normal.

Analisis Pencilan

Peubah y1

y1 <- datap12$y1
m1 <- mean(y1)
s1 <- sd(y1)
pencilan <- (y1 > m1+3*s1) | (y1 < m1-3*s1)

# Menghitung banyaknya amatan pencilan pada y1
sum(pencilan) # tidak ada pencilan 

## [1] 0

# Membuat plot 
boxplot(y1, border = "#2E0249", col = "#3BACB6", horizontal = T) 

plot(y1, col = "#2F8F9D") 

Berdasarkan hasil pendeteksi pencilan di atas, didapatkan bahwa tidak terdeteksi adanya pencilan pada peubah tunggal y1. Selain itu, berdasarkan visualisasi boxplot di atas juga terlihat bahwa peubah tunggal y1 tidak terdapat pencilan.

Peubah y2

y2 <- datap12$y2
m2 <- mean(y2)
s2 <- sd(y2)
pencilan <- (y2 > m2+3*s2) | (y2 < m2-3*s2)

# Menghitung banyaknya amatan pencilan pada y2
sum(pencilan)

## [1] 3

#mengidentifikasi nomor amatan yang menjadi pencilan
which(pencilan) # amatan 98, 99, 100

## [1]  98  99 100

#nilai pencilan
y2[which(pencilan)]

## [1] 37.56401 43.29616 48.86305

# Membuat plot 
boxplot(y2, border = "#2E0249", col = "#3BACB6", horizontal = T) 

plot(y2, col = "#2F8F9D") 

Berdasarkan hasil pendeteksi pencilan di atas, didapatkan bahwa terdapat 3 buah pencilan yang berada pada amatan ke-98, ke-99, dan ke-100. Selain itu, dapat dilihat berdasarkan visualisasi boxplot di atas juga terlihat bahwa peubah tunggal y2 terdapat 3 buah pencilan.

Peubah x1

x1 <- datap12.2$x1
boxplot(datap12.2$x1, xlab = "x1", border = "#2E0249", col = "#3BACB6" , horizontal = T)

plot(x1, col = "#2F8F9D")

pencilan_x1 <- boxplot.stats(datap12.2$x1)$out
pencilan_x1 # amatan ke-10

## [1] 24.75342

Berdasarkan hasil pendeteksi pencilan di atas, didapatkan bahwa terdapat 1 buah pencilan yang berada pada amatan ke-10. Selain itu, dapat dilihat berdasarkan visualisasi boxplot di atas juga terlihat bahwa peubah tunggal x1 terdapat 1 buah pencilan.

Peubah x2

x2 <- datap12.2$x2
boxplot(datap12.2$x2, xlab = "x2", border = "#2E0249", col = "#3BACB6" , horizontal = T)

plot(x2, col = "#2F8F9D")

pencilan_x2 <- boxplot.stats(datap12.2$x2)$out
pencilan_x2 # tidak ada pencilan

## numeric(0)

Berdasarkan hasil pendeteksi pencilan di atas, didapatkan bahwa tidak terdeteksi adanya pencilan pada peubah tunggal x2. Selain itu, berdasarkan visualisasi boxplot di atas juga terlihat bahwa peubah tunggal x2 tidak terdapat pencilan.

Penduga Robust bagi Mean

Peubah y1

mean(datap12$y1)

## [1] 0.4906979

mean(datap12$y1, trim = 0.05) # trimmed mean

## [1] 0.5074825

# winsorized mean
winsor.mean(datap12$y1, trim=0.05) 

## [1] 0.504838

winsor.mean(datap12$y1, trim=0.1)

## [1] 0.506197

Rata-rata dapat dipengaruhi oleh keberadaan pencilan. Pada peubah tunggal y1 tidak terdeteksi adanya pencilan sehingga tidak terdapat permasalahan saat perhitungan rata-rata. Akan tetapi, untuk mendapatkan rata-rata yang kekar atau robust diperlukan untuk melakukan pemangkasan data dan penduga kekar yang digunakan, yaitu trimmed mean. Alasannya, karena jumlah observasi pada data yang tergolong besar (100 observasi). Dengan menggunakan syntax di atas, didapatkan rata-rata asli dari peubah tunggal y1, yaitu 0.4906979 dan rata-rata terpangkas 5%, yaitu 0.5074825. Dugaan rata-rata yang kekar didapatkan dari rataan yang terpangkas tersebut.

Peubah y2

mean(datap12$y2)

## [1] 10.82127

mean(datap12$y2, trim = 0.05) # trimmed mean

## [1] 9.939848

# winsorized mean
winsor.mean(datap12$y2, trim=0.05) 

## [1] 10.11176

winsor.mean(datap12$y2, trim=0.1)

## [1] 9.989755

Rata-rata dapat dipengaruhi oleh keberadaan pencilan. Pada peubah tunggal y2 terdapat 3 buah pencilan atas. Oleh karena itu, dengan menggunakan penduga kekar trimmed mean untuk memangkas (0.05 ×100) = 5 amatan terbesar dan 5 amatan terkecil sehingga pencilan dapat terpangkas. Berdasarkan hasil di atas, melalui pemangkasan 5% menghasilkan dugaan rata-rata yang kekar terhadap pencilan, yaitu 9.939848.

Peubah x1

mean(datap12.2$x1)

## [1] 21.02891

mean(datap12.2$x1, trim = 0.05) # trimmed mean

## [1] 21.02891

# winsorized mean
winsor.mean(datap12.2$x1, trim=0.05) 

## [1] 20.92647

winsor.mean(datap12.2$x1, trim=0.1)

## [1] 20.82403

Rata-rata dapat dipengaruhi oleh keberadaan pencilan. Pada peubah tunggal x1 terdapat 1 buah pencilan atas. Oleh karena itu, dengan menggunakan penduga kekar winsorized mean dan nilai α yang digunakan sebesar 10% untuk memangkas (0.1 ×10) = 1 amatan terbesar dan 1 amatan terkecil sehingga pencilan dapat terpangkas. Alasannya, menggunakan winsorized mean hanya meminimalkan pengaruh pencilan dalam data dengan menetapkan bobot yang lebih rendah pada pencilan dan mengubah nilai sehingga mendekati nilai lain dalam himpunan. Berdasarkan hasil di atas, melalui pemangkasan 10% menghasilkan dugaan rata-rata yang kekar terhadap pencilan, yaitu 20.82403.

Peubah x2

mean(datap12.2$x2)

## [1] 2.020437

mean(datap12.2$x2, trim = 0.05) # trimmed mean

## [1] 2.020437

# winsorized mean
winsor.mean(datap12.2$x2, trim=0.05) 

## [1] 2.045956

winsor.mean(datap12.2$x2, trim=0.1)

## [1] 2.071475

Rata-rata dapat dipengaruhi oleh keberadaan pencilan. Pada peubah tunggal x2 tidak terdeteksi adanya pencilan sehingga tidak terdapat permasalahan saat perhitungan rata-rata. Akan tetapi, untuk mendapatkan rata-rata yang kekar atau robust diperlukan untuk melakukan pemangkasan data dan penduga kekar yang digunakan, yaitu winsorized mean. Alasannya,menggunakan winsorized mean hanya meminimalkan pengaruh pencilan dalam data dengan menetapkan bobot yang lebih rendah pada pencilan dan mengubah nilai sehingga mendekati nilai lain dalam himpunan. Dengan menggunakan syntax di atas, didapatkan rata-rata terpangkas 10%, yaitu 2.071475. Dugaan rata-rata yang kekar didapatkan dari rataan yang terpangkas tersebut.

Penduga Robust bagi Ragam

Peubah y1

## MAD (Median Absolute Deviation)
mad.y1 <- mad(datap12$y1, constant = 1) 
simp.baku.y1 <- sd(datap12$y1) 
simp.baku.kekar.y1 <- mad(datap12$y1) 
c(mad.y1, simp.baku.y1, simp.baku.kekar.y1)

## [1] 0.2720740 0.3482005 0.4033769

ragam.y1 <- var(datap12$y1)
ragam.kekar.y1 <- mad(datap12$y1)^2
c(mad.y1, ragam.y1, ragam.kekar.y1)

## [1] 0.2720740 0.1212436 0.1627130

Penduga kekar bagi ragam dapat dicari dengan mengkuadratkan penduga kekar bagi simpangan baku. Salah satu penduga kekar yang dapat digunakan, yaitu MAD (Median Absolute Deviation). Berdasarkan output di atas, didapatkan nilai ragam dan ragam kekar yang diperoleh dari kuadrat simpangan baku dan simpangan baku kekar. Simpangan baku kekar didapatkan dari 1.4826 × 0.272074 = 0.4033769.

Peubah y2

## MAD (Median Absolute Deviation)

# Menentukan nilai constant
m = 10^5;  n = 1000;  c = numeric(m)
for(i in 1:m) {
  u = rlnorm(n);  s = sd(u);  d = mad(u, const=T)
  c[i] = s/d 
}
const <- mean(c)
const

## [1] 3.581792

# Menentukan nilai sd robust
mad(datap12$y2, constant = const)

## [1] 11.86271

Penduga kekar bagi ragam dapat dicari dengan mengkuadratkan penduga kekar bagi simpangan baku. Salah satu penduga kekar yang dapat digunakan, yaitu MAD (Median Absolute Deviation). Pada peubah y2 diperlukan menentukan nilai constant dari sebaran lognormal menggunakan iterasi karena peubah y2 tidak menyebar secara normal. Selanjutnya, diperoleh constant sebaran lognormal sebesar 3.579393 dari perhitungan rata-rata terhadap seluruh nilai constant hasil iterasi. Kemudian, menggunakan fungsi MAD() dengan menggunakan nilai constant yang sudah diperoleh, didapatkan nilai simpangan baku kekar dari peubah tunggal y2, yaitu 11.85476.

Peubah x1

## MAD (Median Absolute Deviation)
mad.x1 <- mad(datap12.2$x1, constant = 1) 
simp.baku.x1 <- sd(datap12.2$x1) 
simp.baku.kekar.x1 <- mad(datap12.2$x1) 
c(mad.x1, simp.baku.x1, simp.baku.kekar.x1)

## [1] 0.7385415 1.5377346 1.0949616

ragam.x1 <- var(datap12.2$x1)
ragam.kekar.x1 <- mad(datap12.2$x1)^2
c(mad.x1, ragam.x1, ragam.kekar.x1)

## [1] 0.7385415 2.3646277 1.1989409

Penduga kekar bagi ragam dapat dicari dengan mengkuadratkan penduga kekar bagi simpangan baku. Salah satu penduga kekar yang dapat digunakan, yaitu MAD (Median Absolute Deviation). Berdasarkan output di atas, didapatkan nilai ragam dan ragam kekar yang diperoleh dari kuadrat simpangan baku dan simpangan baku kekar sebesar 2.3646277 dan 1.1989409 . Didapatkan juga nilai MAD, yaitu 0.7385415. Simpangan baku kekar didapatkan dari 1.4826 × 0.7385415 = 1.0949616 (1.4826 × MAD) pada data yang menyebar normal.

Peubah x2

## MAD (Median Absolute Deviation)
mad.x2 <- mad(datap12.2$x2, constant = 1) 
simp.baku.x2 <- sd(datap12.2$x2) 
simp.baku.kekar.x2 <- mad(datap12.2$x2) 
c(mad.x2, simp.baku.x2, simp.baku.kekar.x2)

## [1] 0.734091 2.802018 1.088363

ragam.x2 <- var(datap12.2$x2)
ragam.kekar.x2 <- mad(datap12.2$x2)^2
c(mad.x2, ragam.x2, ragam.kekar.x2)

## [1] 0.734091 7.851305 1.184535

Penduga kekar bagi ragam dapat dicari dengan mengkuadratkan penduga kekar bagi simpangan baku. Salah satu penduga kekar yang dapat digunakan, yaitu MAD (Median Absolute Deviation). Berdasarkan output di atas, didapatkan nilai ragam dan ragam kekar yang diperoleh dari kuadrat simpangan baku dan simpangan baku kekar sebesar 7.851305 dan 1.184535. Didapatkan juga nilai MAD, yaitu 0.734091. Simpangan baku kekar didapatkan dari 1.4826 × 0.734091 = 1.088363 (1.4826 × MAD) pada data yang menyebar normal.