Principal Component Analysis
Novika SIjabat
3/23/2021
Analisis Komponen Utama untuk jumlah budidaya perikanan menurut jenis budidaya di Indonesia tahun 2000 – 2019 yang terdiri dari 4 variabel jenis budidaya yaitu Tambak (X1), Pembenihan (X2), Air Tawar (X3) dan Laut (X4), dengan sampel sebanyak 20 observasi pada masing – masing variabel. Metode analisis komponen utama dilakukan dengan tujuan mereduksi data yang berasal dari 4 variabel menjadi beberapa komponen utama. Komponen utama yang dipilih lebih sedikit dari jumlah variabel yang tersedia dan sudah cukup merepresentasikan varians dalam data, sehingga analisis data multivariat akan lebih efisien.
Menginput data ke software R Jumlah Perusahaan Budidaya Perikanan Menurut Jenis Budidaya Tahun 2000-2019 yang bersumber dari Badan Pusat Statistik, 2020
#Jumlah Perusahaan Budidaya Perikanan Menurut Jenis Budidaya, 2000-2019
data=read.csv("data pca 1.csv",header=T)
data## x1 x2 x3 x4
## 1 131 67 4 14
## 2 172 85 7 19
## 3 174 85 7 19
## 4 156 89 8 23
## 5 193 104 8 30
## 6 91 30 4 22
## 7 126 54 9 21
## 8 135 59 13 27
## 9 145 54 7 22
## 10 148 51 6 24
## 11 134 54 9 24
## 12 138 55 8 25
## 13 115 64 12 31
## 14 123 63 15 35
## 15 136 69 16 40
## 16 137 75 16 42
## 17 139 75 15 45
## 18 118 82 13 44
## 19 126 73 14 45
## 20 168 70 9 34#MISAL:
tambak=data$x1
pembenihan=data$x2
tawar=data$x3
laut=data$x4Melakukan pengujian asumsi untuk melakukan analisis komponen utama a. Data berdistribusi normal multivariate b. Data memiliki varians yang homogen
#UJI ASUMSI
##UJI NORMALITAS MULTIVARIAT
library(MVN)## Warning: package 'MVN' was built under R version 3.6.3## Registered S3 method overwritten by 'GGally':
## method from
## +.gg ggplot2## sROC 0.1-2 loadedmvn(data=data, mvnTest="mardia")## $multivariateNormality
## Test Statistic p value Result
## 1 Mardia Skewness 19.7327517281494 0.474755076593422 YES
## 2 Mardia Kurtosis -0.939671975764251 0.347385844218097 YES
## 3 MVN <NA> <NA> YES
##
## $univariateNormality
## Test Variable Statistic p value Normality
## 1 Shapiro-Wilk x1 0.9611 0.5669 YES
## 2 Shapiro-Wilk x2 0.9775 0.8984 YES
## 3 Shapiro-Wilk x3 0.9190 0.0950 YES
## 4 Shapiro-Wilk x4 0.9182 0.0914 YES
##
## $Descriptives
## n Mean Std.Dev Median Min Max 25th 75th Skew Kurtosis
## x1 20 140.25 23.485438 136.5 91 193 126.00 150.00 0.325314228 -0.1351834
## x2 20 67.90 16.625599 68.0 30 104 54.75 76.75 -0.002724845 -0.1175834
## x3 20 10.00 3.906809 9.0 4 16 7.00 13.25 0.171054135 -1.4079084
## x4 20 29.30 9.690473 26.0 14 45 22.00 36.25 0.369658272 -1.2981719#UJI VARIANS HOMOGEN
variabel=as.factor(rep(c("tambak","pembenihan","tawar","laut"),each=20))
jenis.budidaya=c(tambak,pembenihan,tawar,laut)
data2=data.frame(variabel,jenis.budidaya)
bartlett.test(jenis.budidaya~variabel,data2)##
## Bartlett test of homogeneity of variances
##
## data: jenis.budidaya by variabel
## Bartlett's K-squared = 48.782, df = 3, p-value = 1.451e-10Memanggil data dalam bentuk matriks
X=cbind(tambak,pembenihan,tawar,laut)
X## tambak pembenihan tawar laut
## [1,] 131 67 4 14
## [2,] 172 85 7 19
## [3,] 174 85 7 19
## [4,] 156 89 8 23
## [5,] 193 104 8 30
## [6,] 91 30 4 22
## [7,] 126 54 9 21
## [8,] 135 59 13 27
## [9,] 145 54 7 22
## [10,] 148 51 6 24
## [11,] 134 54 9 24
## [12,] 138 55 8 25
## [13,] 115 64 12 31
## [14,] 123 63 15 35
## [15,] 136 69 16 40
## [16,] 137 75 16 42
## [17,] 139 75 15 45
## [18,] 118 82 13 44
## [19,] 126 73 14 45
## [20,] 168 70 9 34Menghitung matriks varians kovarians atau matriks korelasi (berdasarkan terpenuhi atau tidaknya asumsi homogenitas varians)
# Function vektor mean
vectorMean <- function(X){
n <- dim(X)[1]
Xbar <- rep(1,n) %*% X / n
return(Xbar)
}
# Run program
print(vectorMean(X))## tambak pembenihan tawar laut
## [1,] 140.25 67.9 10 29.3# Function matriks korelasi
matCor <- function(X){
n <- dim(X)[1] # jumlah baris
p <- dim(X)[2] # jumlah kolom
J <- rep(1,n) # vektor satuan
C <- diag(1,n) - J %*% t(J) / n
D <- c()
for (j in 1:p)
D <- c(D, sd(X[,j])) #diagonal
Xs <- C %*% X %*% solve(diag(D))
matCor <- t(Xs) %*% Xs / (n-1)
# penamaan baris dan kolom
rownames(matCor) <- colnames(matCor) <- colnames(X)
return(matCor)
}
# print out function matriks kovarian X untuk populasi
matrixKor<-print(matCor(X))## tambak pembenihan tawar laut
## tambak 1.0000000 0.7027483 -0.2053566 -0.1865122
## pembenihan 0.7027483 1.0000000 0.2147300 0.2491271
## tawar -0.2053566 0.2147300 1.0000000 0.8619279
## laut -0.1865122 0.2491271 0.8619279 1.0000000Mencari pasangan nilai eigen dan vector eigen dari matriks
##MENCARI EIGEN VALUE DAN EIGEN VECTOR
eigen(matrixKor)## eigen() decomposition
## $values
## [1] 1.9757553 1.6987823 0.1889934 0.1364690
##
## $vectors
## [,1] [,2] [,3] [,4]
## [1,] 0.08345045 -0.72746105 0.67521558 0.08899643
## [2,] -0.26306000 -0.67507158 -0.67977726 -0.11393275
## [3,] -0.67760819 0.09785074 0.27797750 -0.67379588
## [4,] -0.68167631 0.07418905 0.06866811 0.72463652#Eigen Value
eigen(matrixKor)$value## [1] 1.9757553 1.6987823 0.1889934 0.1364690#Eigen Vector
eigen(matrixKor)$vectors## [,1] [,2] [,3] [,4]
## [1,] 0.08345045 -0.72746105 0.67521558 0.08899643
## [2,] -0.26306000 -0.67507158 -0.67977726 -0.11393275
## [3,] -0.67760819 0.09785074 0.27797750 -0.67379588
## [4,] -0.68167631 0.07418905 0.06866811 0.72463652jum.mCor<-sum(diag(matrixKor))
jum.mCor## [1] 4jum.eigen<-sum(c(eigen(matrixKor)$value))
jum.eigen## [1] 4Menghitung proporsi dari total varians komponen utama
cmCor<-c(diag(matrixKor))
cmCor## tambak pembenihan tawar laut
## 1 1 1 1ceigenva<-c(eigen(matrixKor)$value)
ceigenva## [1] 1.9757553 1.6987823 0.1889934 0.1364690ceigenve<-c(eigen(matrixKor)$vectors)
ceigenve## [1] 0.08345045 -0.26306000 -0.67760819 -0.68167631 -0.72746105 -0.67507158
## [7] 0.09785074 0.07418905 0.67521558 -0.67977726 0.27797750 0.06866811
## [13] 0.08899643 -0.11393275 -0.67379588 0.72463652n<-c(1:16)
p<-c(1:4)
q<-rep(1:4,each=4)
q## [1] 1 1 1 1 2 2 2 2 3 3 3 3 4 4 4 4#proporsi
prop<-function(q,w){
prop<-(q/w)
return(prop)
}
proporsi<-c(prop(ceigenva[p],jum.eigen))
proporsi## [1] 0.49393883 0.42469557 0.04724836 0.03411725#Proporsi Kumulatif
propcum=c()
for (i in p){
propcum=c(propcum,sum(proporsi[1:i]))
}
propcum## [1] 0.4939388 0.9186344 0.9658827 1.0000000# Function matriks korelasi untuk populasi
corYX <- function(a,b){
corYnXp<-(a*sqrt(b))
return(corYnXp)
}Menghitung korelasi antara variabel acak asli dan komponen utama
# print out function matriks korelasi antar var acak Z
korelasi<-c(corYX(ceigenve[n],ceigenva[q]))
korelasi## [1] 0.11729926 -0.36976125 -0.95245667 -0.95817488 -0.94815339 -0.87987036
## [7] 0.12753605 0.09669603 0.29353899 -0.29552211 0.12084620 0.02985234
## [13] 0.03287680 -0.04208871 -0.24891173 0.26769313matriks.kor<-matrix(c(korelasi), nrow=3,byrow=T)## Warning in matrix(c(korelasi), nrow = 3, byrow = T): data length [16] is not a
## sub-multiple or multiple of the number of rows [3]matriks.kor## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.1172993 -0.36976125 -0.9524567 -0.9581749 -0.9481534 -0.87987036
## [2,] 0.1275360 0.09669603 0.2935390 -0.2955221 0.1208462 0.02985234
## [3,] 0.0328768 -0.04208871 -0.2489117 0.2676931 0.1172993 -0.36976125Scree Plot
#Scree plot
pev=plot(ceigenva,type="b",main="Scree Plot",xlab="Principal Component",
ylab="Variances")