Latihan Praktikum 12 - AED

Identifikasi Pencilan Peubah Tunggal dan Penduga Kekar bagi Parameter Pemusatan dan Penyebaran

Library yang Digunakan

library(readxl)

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

library(data.table)

## Warning: package 'data.table' was built under R version 4.1.2

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

Impor Data Peubah y

dt.prak12 <- read_excel("C:/Users/hp/Documents/MY DEPARTEMEN/#4/AED/DATA/Data praktikum 12.xlsx", sheet = 1)
View(dt.prak12)

Uji Formalitas (Uji Normalitas) Peubah y1 dan y2

# uji normalitas y1
qqnorm(dt.prak12$y1)
qqline(dt.prak12$y1, col = "red")

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

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

Berdasarkan uji formal: Kolmogorov Smirnov Test, p-value = 0.4512 > taraf nyata (α) = 0.05 -> tak tolak H0 Jadi, sebaran data peubah tunggal y1 mendekati sebaran normal.

# uji normalitas y2
ks.test(dt.prak12$y2, "pnorm", mean = mean(dt.prak12$y2), sd = sd(dt.prak12$y2))

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

# uji lognormal
ks.test(dt.prak12$y2, "plnorm", meanlog = mean(log(dt.prak12$y2), na.rm = T), sdlog = sd(log(dt.prak12$y2)))

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

# uji chisq
ks.test(dt.prak12$y2, "pchisq", df = mean(dt.prak12$y2))

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

Berdasarkan uji formal: Kolmogorov Smirnov Test, p-value = 0.007139 < taraf nyata (α) = 0.05 -> tolak H0 Jadi, sebaran data peubah tunggal y2 tidak mendekati sebaran normal, melainkan sebaran lognormal dan chi-square.

Pendeteksiaan Pencilan

Pendeteksian pencilan dengan menggunakan dua cara. Pertama, suatu amatan akan disebut pencilan juka berada di luar selang (x bar - 3s, x bar + 3s). Kedua, amatan yang bernilai lebih kecil daripada Batas Bawah atau lebih besar daripada Batas Atas diidentifikasi sebagai pencilan. • Batas Bawah = Q1 – 1.5 * IQR • Batas Atas = Q3 + 1.5 * IQR

###Pendeteksiaan Pencilan pada Peubah y1

#Cara 1
y1 <- dt.prak12$y1
m <- mean(y1) 
s <- sd(y1) 
pencilan <- (y1 > m+3*s) | (y1 < m-3*s) 

#menghitung banyaknya amatan pencilan
sum(pencilan)

## [1] 0

#mengidentifikasi nomor amatan yang menjadi pencilan
which(pencilan)

## integer(0)

#nilai pencilan 
y1[which(pencilan)]

## numeric(0)

#Cara 2
value.outliery1 <- boxplot.stats(y1)$out
which(y1 == value.outliery1)

## integer(0)

value.outliery1

## numeric(0)

#boxplot
boxplot(dt.prak12$y1, ylab = "y1")

Dengan pendeteksiaan pencilan pada data peubah tunggal y1, dapat dilihat bahwa tidak ada pencilan.Selain itu, dapat dilihat juga pada boxplot bahwa data peubah tunggal y1 simetris.

Pendeteksiaan Pencilan pada Peubah y2

#Cara 1
y2 <- dt.prak12$y2
m <- mean(y2) 
s <- sd(y2) 
pencilan <- (y2 > m+3*s) | (y2 < m-3*s) 

#menghitung banyaknya amatan pencilan
sum(pencilan)

## [1] 3

#mengidentifikasi nomor amatan yang menjadi pencilan
which(pencilan)

## [1]  98  99 100

#nilai pencilan 
y2[which(pencilan)]

## [1] 37.56401 43.29616 48.86305

#Cara 2
value.outliery2 <- boxplot.stats(y2)$out
which(y2 == value.outliery2)

## Warning in y2 == value.outliery2: longer object length is not a multiple of
## shorter object length

## integer(0)

value.outliery2

## [1] 37.56401 43.29616 48.86305

#boxplot
boxplot(dt.prak12$y2, ylab = "y2")

Dengan pendeteksiaan pencilan pada data peubah tunggal y2, dapat dilihat bahwa terdapat pencilan pada amatan data ke-98, 99, dan 100. Jadi, terdapat tiga pencilan pada peubah tunggal y2. Selain itu, dapat dilihat juga pada boxplot bahwa data peubah tunggal y2 menjulur ke kanan.

Impor Data Peubah x

dt.prak12.2 <- read_excel("C:/Users/hp/Documents/MY DEPARTEMEN/#4/AED/DATA/Data praktikum 12.xlsx", sheet = 2)
View(dt.prak12.2)

UJi Formalitas (Uji Normalitas)

# uji normalitas x1
qqnorm(dt.prak12.2$x1)
qqline(dt.prak12.2$x1, col = "red")

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

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

Berdasarkan uji formal: Kolmogorov Smirnov Test, p-value = 0.7403 > taraf nyata (α) = 0.05 -> tak tolak H0 Jadi, sebaran data peubah tunggal x1 mendekati sebaran normal.

# uji normalitas x2
qqnorm(dt.prak12.2$x2)
qqline(dt.prak12.2$x2, col = "red")

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

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

Berdasarkan uji formal: Kolmogorov Smirnov Test, p-value = 0.6449 > taraf nyata (α) = 0.05 -> tak tolak H0 Jadi, sebaran data peubah tunggal x2 mendekati sebaran normal.

Pendeteksiaan Pencilan pada Peubah x1

#Cara 1
x1 <- dt.prak12.2$x1
m <- mean(x1) 
s <- sd(x1) 
pencilan <- (x1 > m+3*s) | (x1 < m-3*s) 

#menghitung banyaknya amatan pencilan
sum(pencilan)

## [1] 0

#mengidentifikasi nomor amatan yang menjadi pencilan
which(pencilan)

## integer(0)

#nilai pencilan 
x1[which(pencilan)]

## numeric(0)

#Cara 2
value.outlierx1 <- boxplot.stats(x1)$out
which(x1 == value.outlierx1)

## [1] 10

value.outlierx1

## [1] 24.75342

#boxplot
boxplot(dt.prak12.2$x1, ylab = "x1")

Dengan pendeteksiaan pencilan pada data peubah tunggal x1, dapat dilihat bahwa terdapat pencilan pada amatan data ke-10. Jadi, terdapat satu pencilan pada peubah tunggal x1. Selain itu, dapat dilihat juga pada boxplot bahwa data peubah tunggal x1 menjulur ke kanan.

Pendeteksiaan Pencilan pada Peubah x2

#Cara 1
x2 <- dt.prak12.2$x2
m <- mean(x2) 
s <- sd(x2) 
pencilan <- (x2 > m+3*s) | (x2 < m-3*s) 

#menghitung banyaknya amatan pencilan
sum(pencilan)

## [1] 0

#mengidentifikasi nomor amatan yang menjadi pencilan
which(pencilan)

## integer(0)

#nilai pencilan 
x2[which(pencilan)]

## numeric(0)

#Cara 2
value.outlierx2 <- boxplot.stats(x2)$out
which(x2 == value.outlierx2)

## integer(0)

value.outlierx2

## numeric(0)

#boxplot
boxplot(dt.prak12.2$x2, ylab = "x2")

Dengan pendeteksiaan pencilan pada data peubah x2, dapat dilihat bahwa tidak ada pencilan.Selain itu, dapat dilihat juga pada boxplot bahwa data peubah tunggal x2 simetris.

Penduga Robust bagi Mean

Penduga robust bagi mean dengan menggunakan 2 metode yaitu Trimmed Mean dan Winsored Mean. 1. Trimmed Mean Rataan terpangkas (trimmed mean) merupakan rata-rata dari data yang ada di bagian tengah data, tepatnya data di 1 – 2𝛼 bagian tengah data dengan 0 < 𝛼 < 1. 2. Winsorized Mean Nilai winsorized mean diperoleh dengan menghitung rataan setelah dilakukan penggantian nilai terhadap amatan-amatan terbesar dan terkecil.

Penduga Robust bagi Mean pada Peubah y1

#Trimmed Mean
mean(dt.prak12$y1)

## [1] 0.4906979

mean(dt.prak12$y1, trim=0.05)

## [1] 0.5074825

mean(dt.prak12$y1, trim=0.1)

## [1] 0.5107736

#Winsorized Mean
winsorMEAN <- function(x,probs=c(0.05,0.95))
  {
    xq<-quantile(x,probs=probs)
    x[x < xq[1]]<-xq[1]
    x[x > xq[2]]<-xq[2]
    return(mean(x))
}

#nilai proporsi
wm05 <- winsorMEAN(dt.prak12$y1)                      #nilai peluang 5% dan 95%
wm05

## [1] 0.504838

wm10 <- winsorMEAN(dt.prak12$y1, probs=c(0.10, 0.90)) #nilai peluang 10% dan 90%
wm10

## [1] 0.506197

#menggunakan package psych
winsor.mean(dt.prak12$y1, trim=0.1)

## [1] 0.506197

winsor.mean(dt.prak12$y1, trim=0.05)

## [1] 0.504838

Penduga Robust bagi Mean pada Peubah y2

#Trimmed Mean
mean(dt.prak12$y2)

## [1] 10.82127

mean(dt.prak12$y2, trim=0.05)

## [1] 9.939848

mean(dt.prak12$y2, trim=0.1)

## [1] 9.776771

#Winsorized Mean
winsorMEAN <- function(x,probs=c(0.05,0.95))
  {
    xq<-quantile(x,probs=probs)
    x[x < xq[1]]<-xq[1]
    x[x > xq[2]]<-xq[2]
    return(mean(x))
}

#nilai proporsi
wm05 <- winsorMEAN(dt.prak12$y2)                      #nilai peluang 5% dan 95%
wm05

## [1] 10.11176

wm10 <- winsorMEAN(dt.prak12$y2, probs=c(0.10, 0.90)) #nilai peluang 10% dan 90%
wm10

## [1] 9.989755

#menggunakan package psych
winsor.mean(dt.prak12$y2, trim=0.1)

## [1] 9.989755

winsor.mean(dt.prak12$y2, trim=0.05)

## [1] 10.11176

Penduga robust bagi mean pada data peubah tunggal y1 dan y2 lebih efisien dengan menggunakan trimmed mean dibandingkan dengan winsorized mean karena kedua data tersebut memiliki banyak amatan yang berjumlah 100 observasi. Ketika memangkas beberapa amatan data terbesar dan terkecil, masih banyak amatan lainnya yang dapat mewakili ukuran pemusatan data, sehingga Penduga robust bagi mean pada data peubah tunggal y1 dan y2 masih representatif. Penentuan opsi trim = 0.05 mengindikasikan bahwa rata-rata terpangkas 5% pada masing-masing ujung kiri dan kanan data, sehingga data yang tersisa hanya 90 observasi. Sedangkan, penentuan opsi trim = 0.10 mengindikasikan bahwa rata-rata terpangkas 10% pada masing-masing ujung kiri dan kanan data, sehingga data yang tersisa hanya 80 amatan.

Penduga Robust bagi Mean pada Peubah x1

#Trimmed Mean
mean(dt.prak12.2$x1)

## [1] 21.02891

mean(dt.prak12.2$x1, trim=0.05)

## [1] 21.02891

mean(dt.prak12.2$x1, trim=0.1)

## [1] 20.76255

#Winsorized Mean
winsorMEAN <- function(x,probs=c(0.05,0.95))
  {
    xq<-quantile(x,probs=probs)
    x[x < xq[1]]<-xq[1]
    x[x > xq[2]]<-xq[2]
    return(mean(x))
}

#nilai proporsi
wm05 <- winsorMEAN(dt.prak12.2$x1)                      #nilai peluang 5% dan 95%
wm05

## [1] 20.92647

wm10 <- winsorMEAN(dt.prak12.2$x1, probs=c(0.10, 0.90)) #nilai peluang 10% dan 90%
wm10

## [1] 20.82403

#menggunakan package psych
winsor.mean(dt.prak12.2$x1, trim=0.1)

## [1] 20.82403

winsor.mean(dt.prak12.2$x1, trim=0.05)

## [1] 20.92647

Penduga Robust bagi Mean pada Peubah x2

#Trimmed Mean
mean(dt.prak12.2$x2)

## [1] 2.020437

mean(dt.prak12.2$x2, trim=0.05)

## [1] 2.020437

mean(dt.prak12.2$x2, trim=0.1)

## [1] 1.958662

#Winsorized Mean
winsorMEAN <- function(x,probs=c(0.05,0.95))
  {
    xq<-quantile(x,probs=probs)
    x[x < xq[1]]<-xq[1]
    x[x > xq[2]]<-xq[2]
    return(mean(x))
}

#nilai proporsi
wm05 <- winsorMEAN(dt.prak12.2$x2)                      #nilai peluang 5% dan 95%
wm05

## [1] 2.045956

wm10 <- winsorMEAN(dt.prak12.2$x2, probs=c(0.10, 0.90)) #nilai peluang 10% dan 90%
wm10

## [1] 2.071475

#menggunakan package psych
winsor.mean(dt.prak12.2$x2, trim=0.1)

## [1] 2.071475

winsor.mean(dt.prak12.2$x2, trim=0.05)

## [1] 2.045956

Penduga robust bagi mean pada data peubah tunggal x1 dan x2 lebih efisien dengan menggunakan winsorized mean dibandingkan dengan trimmed mean yang memiliki sedikit amatan data. Hali ini dikarenakan winsorized mean hanya meminimalisasi pengaruh pencilan dalam suatu data dengan menetapkan bobot yang lebih rendah dan menggantikan nilai, sehingga mendekati nilai lain dalam himpunan. Penentuan nilai α = 5% mengindikasikan bahwa nilai proporsi/peluang masing-masing yaitu 0.05 dan 0.95, sehingga 5% data terkecil dan 5% data terbesar diganti dengan nilai pada titik tersebut, yaitu 𝑄(0.05) dan Q(0.95). Sedangkan, α = 10% mengindikasikan bahwa nilai peluang masing-masing yaitu 0.10 dan 0.90, sehingga 10% data terkecil dan 10% data terbesar diganti dengan nilai pada titik tersebut, yaitu 𝑄(0.10) dan Q(0.90).

Penduga Robust bagi Varians/Ragam

Pendugaan robust bagi ragam dapat ditentukan sebagai kuadrat dari pendugaan robust bagi simpangan baku. Pendugaan robust bagi ragam menggunakan 2 metode yaitu MAD (Median Absolute Deviation) dan Gini’s mean difference.

Penduga Robust bagi Varians pada peubah y1

#MAD
value.mad <- mad(dt.prak12$y1, constant=1) 
stdev.y1 <- sd(dt.prak12$y1)
stdev.kekar.y1 <- mad(dt.prak12$y1)
c(value.mad, stdev.y1, stdev.kekar.y1)

## [1] 0.2720740 0.3482005 0.4033769

ragam.y1 <- stdev.y1^2
ragam.kekar.y1 <- (mad(dt.prak12$y1))^2
c(value.mad, ragam.y1, ragam.kekar.y1)

## [1] 0.2720740 0.1212436 0.1627130

#gini's mean difference
#menggunakan package lmomco
gini.mean.diff(dt.prak12$y1)$gini

## [1] 0.3930819

Nilai MAD dari suatu gugus data diperoleh menggunakan fungsi mad() dengan constant = 1. Selain itu, dihitung nilai simpangan baku asli sebesar 0.3482005 dan simpangan baku kekar peubah tunggal y1. Simpangan baku kekar diduga oleh 1.4826 × MAD = 1.4826 × 0.2720740 = 0.4033769. Nilai ragam asli dan ragam kekar dihtung dengan meng-kuadrat-kan simpangan baku asli dan simpangan baku kekar sehingga diperoleh nilainya masing-masing sebesar 0.1212436 dan 0.1627130. Selain itu, metode lain yaitu gini’s mean difference juga dapat digunakan yang menghasilkan output nilai yang cukup berbeda.

Penduga Robust bagi Varians pada peubah y2

#MAD
# 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.5822

value.mad <- mad(dt.prak12$y2, constant=const) 
stdev.y2 <- sd(dt.prak12$y2)
stdev.kekar.y2 <- mad(dt.prak12$y2)
c(value.mad, stdev.y2, stdev.kekar.y2)

## [1] 11.864060  7.081105  4.910295

ragam.y2 <- stdev.y2^2
ragam.kekar.y2 <- (mad(dt.prak12$y2))^2
c(value.mad, ragam.y2, ragam.kekar.y2)

## [1] 11.86406 50.14205 24.11099

#gini's mean difference
#menggunakan package lmomco
gini.mean.diff(dt.prak12$y2)$gini

## [1] 6.403677

Penduga robust bagi ragam dapat dihitung dengan meng-kuadrat-kan penduga robist bagi simpangan baku. Metode MAD (Median Absolute Deviation) pada data peubah tunggal y2 diperlukan menentukan nilai constant dari sebaran lognormal yang menggunakan iterasi karenatidak mendekati sebaran normal. Dengan demikian,constant sebaran lognormal diperoleh sebesar 3.583238 dari rataan terhadap semua nilai constant hasil iterasi. Selanjutnya, nilia MAD diperoleh dari fungsi MAD() dengan menggunakan nilai constant yang sudah ada, didapatkan nilai simpangan baku kekar dari peubah tunggal y2, yaitu 4.910295 dan ragamm kekar sebesar 24.11099. Selain itu, metode lain yaitu gini’s mean difference juga dapat digunakan yang menghasilkan output nilai yang cukup berbeda.

Penduga Robust bagi Varians pada peubah x1

#MAD
value.mad <- mad(dt.prak12.2$x1, constant=1) 
stdev.x1 <- sd(dt.prak12.2$x1)
stdev.kekar.x1 <- mad(dt.prak12.2$x1)
c(value.mad, stdev.x1, stdev.kekar.x1)

## [1] 0.7385415 1.5377346 1.0949616

ragam.x1 <- stdev.x1^2
ragam.kekar.x1 <- (mad(dt.prak12.2$x1))^2
c(value.mad, ragam.x1, ragam.kekar.x1)

## [1] 0.7385415 2.3646277 1.1989409

#gini's mean difference
#menggunakan package lmomco
gini.mean.diff(dt.prak12.2$x1)$gini

## [1] 1.644843

Nilai MAD dari suatu gugus data diperoleh menggunakan fungsi mad() dengan constant = 1. Selain itu, dihitung nilai simpangan baku asli sebesar 1.537734 dan simpangan baku kekar peubah tunggal x1. Simpangan baku kekar diduga oleh 1.4826 × MAD = 1.4826 × 0.7385415 = 1.0949616. Nilai ragam asli dan ragam kekar dihtung dengan meng-kuadrat-kan simpangan baku asli dan simpangan baku kekar sehingga diperoleh nilainya masing-masing sebesar 2.3646277 dan 1.1989409. Selain itu, metode lain yaitu gini’s mean difference juga dapat digunakan yang menghasilkan output nilai yang cukup berbeda.

Penduga Robust bagi Varians pada peubah x2

#MAD
value.mad <- mad(dt.prak12.2$x2, constant=1) 
stdev.x2 <- sd(dt.prak12.2$x2)
stdev.kekar.x2 <- mad(dt.prak12.2$x2)
c(value.mad, stdev.x2, stdev.kekar.x2)

## [1] 0.734091 2.802018 1.088363

ragam.x2 <- stdev.x2^2
ragam.kekar.x2 <- (mad(dt.prak12.2$x2))^2
c(value.mad, ragam.x2, ragam.kekar.x2)

## [1] 0.734091 7.851305 1.184535

#gini's mean difference
#menggunakan package lmomco
gini.mean.diff(dt.prak12.2$x2)$gini

## [1] 3.124344

Nilai MAD dari suatu gugus data diperoleh menggunakan fungsi mad() dengan constant = 1. Selain itu, dihitung nilai simpangan baku asli sebesar 2.802018 dan simpangan baku kekar peubah tunggal x2. Simpangan baku kekar diduga oleh 1.4826 × MAD = 1.4826 × 0.734091 = 1.088363. Nilai ragam asli dan ragam kekar dihtung dengan meng-kuadrat-kan simpangan baku asli dan simpangan baku kekar sehingga diperoleh nilainya masing-masing sebesar 7.851305 dan 1.184535. Selain itu, metode lain yaitu gini’s mean difference juga dapat digunakan yang menghasilkan output nilai yang cukup berbeda.

Kesimpulan: Identifikasi pencilan peubah punggal dan penduga kekar bagi parameter pemusatan dan Penyebaran dat terbagi menjadi tiga hal, yaitu: 1. Pengidentifikasian pencilan peubah tunggal dengan menggunakan pendekatan selang (𝑥̅ – 3s, 𝑥̅ + 3s) dan batas atas serta batas bawah. 2. Penduga robust bagi mean dengan menggunakan trimmed mean dan winsorized mean. 3. penduga robust bagi ragam dapat ditentukan sebagai kuadrat dari penduga robust bagi simpangan baku.