#install package
library(ggplot2)
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.2 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ lubridate 1.9.2 ✔ tibble 3.2.1
## ✔ purrr 1.0.1 ✔ tidyr 1.3.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(corrplot)
## corrplot 0.92 loaded
library(caret)
## Loading required package: lattice
##
## Attaching package: 'caret'
##
## The following object is masked from 'package:purrr':
##
## lift
library(tidyr)
library(dplyr)
library(knitr)
library(readxl)
library(MASS)
##
## Attaching package: 'MASS'
##
## The following object is masked from 'package:dplyr':
##
## select
library(ResourceSelection)
## ResourceSelection 0.3-6 2023-06-27
library(lattice)
library(ggplot2)
library(caret)
#input data
dataDBDRSI <- read_excel("D:/ADK/Ordinal - DBD/Data DBD RSI Jemursari.xlsx")
head(dataDBDRSI)
#keterangan
str(dataDBDRSI)
## tibble [28 × 2] (S3: tbl_df/tbl/data.frame)
## $ Y : num [1:28] 3 2 1 2 1 2 1 2 1 1 ...
## $ X1: num [1:28] 7.16 2.58 10.11 4.42 6.27 ...
#mengubah variabel ke factor
dataDBDRSI$Y<- as.factor(dataDBDRSI$Y)
str(dataDBDRSI)
## tibble [28 × 2] (S3: tbl_df/tbl/data.frame)
## $ Y : Factor w/ 3 levels "1","2","3": 3 2 1 2 1 2 1 2 1 1 ...
## $ X1: num [1:28] 7.16 2.58 10.11 4.42 6.27 ...
library(plotrix)
mytable <- table(dataDBDRSI$Y)
lbls <- c("Derajat I","Derajat II", "Derajat III", "Derajat IV")
pct<-round(mytable/sum(mytable)*100)
lbls<-paste(lbls,pct)
lbls<-paste(lbls,"%",sep = "")
pie(mytable, labels = lbls,
main="DBD",
col=c("brown","coral", "grey"))
hist(dataDBDRSI$X1,
main = "Histogram Usia",
xlab = "Tingkat Usia",
ylab = "jumlah",
col = "lightblue")
library(RColorBrewer)
par(bty="o")
boxplot(X1~Y, data=dataDBDRSI, col=brewer.pal(n = 4, name = "RdBu") , xlab="Kategori DBD dengan Usia", horizontal=TRUE)
# Uji asumsi ## Uji multikolineritas
m1 <- glm(as.factor(dataDBDRSI$Y)~X1 , family = binomial(link = "logit"), data = dataDBDRSI)
summary(m1)
##
## Call:
## glm(formula = as.factor(dataDBDRSI$Y) ~ X1, family = binomial(link = "logit"),
## data = dataDBDRSI)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.68505 0.63117 1.085 0.2778
## X1 -0.07198 0.04335 -1.660 0.0968 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 38.243 on 27 degrees of freedom
## Residual deviance: 33.471 on 26 degrees of freedom
## AIC: 37.471
##
## Number of Fisher Scoring iterations: 5
VIF<10
#model
library(MASS)
modelDBDRSI <- polr(Y~X1, data=dataDBDRSI)
summary(modelDBDRSI)
##
## Re-fitting to get Hessian
## Call:
## polr(formula = Y ~ X1, data = dataDBDRSI)
##
## Coefficients:
## Value Std. Error t value
## X1 -0.07189 0.0419 -1.716
##
## Intercepts:
## Value Std. Error t value
## 1|2 -0.6914 0.6190 -1.1170
## 2|3 0.6726 0.6533 1.0295
##
## Residual Deviance: 49.59032
## AIC: 55.59032
##uji serentak
library(pscl)
## Classes and Methods for R developed in the
## Political Science Computational Laboratory
## Department of Political Science
## Stanford University
## Simon Jackman
## hurdle and zeroinfl functions by Achim Zeileis
pR2(modelDBDRSI)
## fitting null model for pseudo-r2
## llh llhNull G2 McFadden r2ML r2CU
## -24.79516084 -27.27174613 4.95317056 0.09081139 0.16213555 0.18909292
##chisquare tabel
qchisq(0.90, 1)
## [1] 2.705543
Berdasarkan output tersebut diperoleh informasi bahwa pada taraf signifikansi 10%, H0 ditolak karena G2 = 4,95317056 > 2,705543. Artinya, ada pengaruh variabel independen terhadap variabel dependen secara bersama-sama. ##partial test
coefmodel <- c(-modelDBDRSI$coefficients)
koefisien<-coef(summary(modelDBDRSI))
##
## Re-fitting to get Hessian
#menghitung p-value
p<-pnorm(abs(koefisien[,"t value"]),lower.tail = FALSE)*2
(ctabel<-cbind <-cbind(round(koefisien,4),"p-value"=round(p,4)))
## Value Std. Error t value p-value
## X1 -0.0719 0.0419 -1.7157 0.0862
## 1|2 -0.6914 0.6190 -1.1170 0.2640
## 2|3 0.6726 0.6533 1.0295 0.3033
##confussion matrix
prediksi.test <- predict(modelDBDRSI, dataDBDRSI, type="class")
dataDBDRSI$Y <- as.factor(dataDBDRSI$Y)
confusionMatrix(as.factor(prediksi.test), dataDBDRSI$Y, positive="1")
## Confusion Matrix and Statistics
##
## Reference
## Prediction 1 2 3
## 1 16 7 5
## 2 0 0 0
## 3 0 0 0
##
## Overall Statistics
##
## Accuracy : 0.5714
## 95% CI : (0.3718, 0.7554)
## No Information Rate : 0.5714
## P-Value [Acc > NIR] : 0.579
##
## Kappa : 0
##
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: 1 Class: 2 Class: 3
## Sensitivity 1.0000 0.00 0.0000
## Specificity 0.0000 1.00 1.0000
## Pos Pred Value 0.5714 NaN NaN
## Neg Pred Value NaN 0.75 0.8214
## Prevalence 0.5714 0.25 0.1786
## Detection Rate 0.5714 0.00 0.0000
## Detection Prevalence 1.0000 0.00 0.0000
## Balanced Accuracy 0.5000 0.50 0.5000
#odds ratio
data.frame(coef(modelDBDRSI), exp(coef(modelDBDRSI)))