library(fitdistrplus)
## Warning: package 'fitdistrplus' was built under R version 4.4.2
## Loading required package: MASS
## Loading required package: survival
library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Input Data
data <- read.csv(file.choose(),header=TRUE,sep=",",dec=".")
str(data)
## 'data.frame':    175 obs. of  15 variables:
##  $ Perairan                  : chr  "Laut Bali" "Laut Flores bagian barat" "Laut Flores bagian timur" "Laut Flores bagian utara" ...
##  $ Tinggi.gelombang.rata.rata: num  0.875 0.875 1.875 0.875 0.875 ...
##  $ Jenis.Gelombang           : chr  "rendah" "rendah" "sedang " "rendah" ...
##  $ J1                        : int  1 1 0 1 1 0 0 1 0 0 ...
##  $ J2                        : int  0 0 1 0 0 1 1 0 1 1 ...
##  $ Rata.Rata.Kecepatan.Angin : num  12 11.5 20 15 11.5 12.5 10.5 12 12.5 12 ...
##  $ Cuaca                     : chr  "Cerah Berawan" "Berawan" "Hujan Sedang" "Hujan Ringan" ...
##  $ CC                        : int  1 1 0 0 0 1 1 1 1 1 ...
##  $ C1                        : int  0 1 0 0 0 1 0 1 0 1 ...
##  $ C2                        : int  0 0 0 0 0 0 1 0 1 0 ...
##  $ C3                        : int  1 0 0 0 0 0 0 0 0 0 ...
##  $ C4                        : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ C5                        : int  0 0 0 1 1 0 0 0 0 0 ...
##  $ C6                        : int  0 0 1 0 0 0 0 0 0 0 ...
##  $ Potensi.Daya..Y.          : num  103077 103077 150889 103077 103077 ...
Y <- data$Potensi.Daya..Y.
X1 <- data$Rata.Rata.Kecepatan.Angin
X2 <- data$Tinggi.gelombang.rata.rata
J1 <- data$J1 ## Kategori Jenis Gelombang
J2 <- data$J2 ## Kategori Jenis Gelombang
C1 <- data$C1 ## Kategori Cuaca
C2 <- data$C2 ## Kategori Cuaca
C3 <- data$C3 ## Kategori Cuaca
C4 <- data$C4 ## Kategori Cuaca
C5 <- data$C5 ## Kategori Cuaca
C6 <- data$C6 ## Kategori Cuaca
## Pemodelan Regresi Linier
model1 <- lm(Y ~ X1 + J1 + J2 + C1 + C2 + C3 + C4 + C5 + C6 , data=data)
summary(model1)
## 
## Call:
## lm(formula = Y ~ X1 + J1 + J2 + C1 + C2 + C3 + C4 + C5 + C6, 
##     data = data)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -1.302e-10 -1.799e-11 -4.880e-12  5.110e-12  2.187e-09 
## 
## Coefficients:
##               Estimate Std. Error    t value Pr(>|t|)    
## (Intercept)  5.564e+04  9.384e-11  5.930e+14   <2e-16 ***
## X1           4.286e-12  4.795e-12  8.940e-01    0.373    
## J1           4.743e+04  3.916e-11  1.211e+15   <2e-16 ***
## J2           9.524e+04  4.806e-11  1.982e+15   <2e-16 ***
## C1          -1.043e-11  9.211e-11 -1.130e-01    0.910    
## C2          -9.825e-13  1.013e-10 -1.000e-02    0.992    
## C3           6.615e-11  9.419e-11  7.020e-01    0.483    
## C4          -2.178e-11  1.527e-10 -1.430e-01    0.887    
## C5          -1.523e-11  9.171e-11 -1.660e-01    0.868    
## C6          -1.654e-11  9.356e-11 -1.770e-01    0.860    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.744e-10 on 165 degrees of freedom
## Multiple R-squared:      1,  Adjusted R-squared:      1 
## F-statistic: 6.147e+29 on 9 and 165 DF,  p-value: < 2.2e-16
model2 <- lm(Y ~ J1 + J2, data=data)
summary(model2)
## 
## Call:
## lm(formula = Y ~ J1 + J2, data = data)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -9.238e-10 -1.470e-11 -1.470e-11  1.148e-11  2.350e-09 
## 
## Coefficients:
##              Estimate Std. Error   t value Pr(>|t|)    
## (Intercept) 5.564e+04  3.660e-11 1.520e+15   <2e-16 ***
## J1          4.743e+04  4.161e-11 1.140e+15   <2e-16 ***
## J2          9.524e+04  4.654e-11 2.047e+15   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.971e-10 on 172 degrees of freedom
## Multiple R-squared:      1,  Adjusted R-squared:      1 
## F-statistic: 2.167e+30 on 2 and 172 DF,  p-value: < 2.2e-16
## Cullen Frey
descdist((model2$residuals))

## summary statistics
## ------
## min:  -9.238048e-10   max:  2.349669e-09 
## median:  -1.469963e-11 
## mean:  5.198997e-26 
## estimated sd:  1.959429e-10 
## estimated skewness:  9.301391 
## estimated kurtosis:  126.2099
## Pemodelan Regresi Gamma
model3 <- glm(Y ~ X1 + J1 + J2 + C1 + C2 + C3 + C4 + C5 + C6 ,family=Gamma("inverse"), data=data)
summary(model3)
## 
## Call:
## glm(formula = Y ~ X1 + J1 + J2 + C1 + C2 + C3 + C4 + C5 + C6, 
##     family = Gamma("inverse"), data = data)
## 
## Coefficients:
##               Estimate Std. Error    t value Pr(>|t|)    
## (Intercept)  1.797e-05  3.467e-21  5.183e+15   <2e-16 ***
## X1          -1.323e-22  1.336e-22 -9.900e-01    0.323    
## J1          -8.270e-06  2.075e-21 -3.986e+15   <2e-16 ***
## J2          -1.134e-05  2.149e-21 -5.278e+15   <2e-16 ***
## C1          -9.141e-24  3.087e-21 -3.000e-03    0.998    
## C2          -1.535e-22  3.222e-21 -4.800e-02    0.962    
## C3          -1.045e-21  3.112e-21 -3.360e-01    0.738    
## C4           4.572e-22  5.028e-21  9.100e-02    0.928    
## C5          -2.426e-22  3.060e-21 -7.900e-02    0.937    
## C6           8.516e-22  3.146e-21  2.710e-01    0.787    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 3.43536e-31)
## 
##     Null deviance: 1.6452e+01  on 174  degrees of freedom
## Residual deviance: 5.5109e-15  on 165  degrees of freedom
## AIC: -2265.6
## 
## Number of Fisher Scoring iterations: 1
## Karena variabel cuaca tidak signifikan maka coba hilangkan variabel cuaca
model4 <- glm(Y~ J1 + J2, family=Gamma("inverse"), data=data)
summary(model4)
## 
## Call:
## glm(formula = Y ~ J1 + J2, family = Gamma("inverse"), data = data)
## 
## Coefficients:
##               Estimate Std. Error    t value Pr(>|t|)    
## (Intercept)  1.797e-05  3.603e-21  4.987e+15   <2e-16 ***
## J1          -8.270e-06  3.754e-21 -2.203e+15   <2e-16 ***
## J2          -1.134e-05  3.752e-21 -3.024e+15   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 1.165984e-30)
## 
##     Null deviance: 1.6452e+01  on 174  degrees of freedom
## Residual deviance: 2.0358e-15  on 172  degrees of freedom
## AIC: -2453.9
## 
## Number of Fisher Scoring iterations: 1
## Pembentukan model null
modelNull <- glm(Y~1, family=Gamma("inverse"), data=data)
summary(modelNull)
## 
## Call:
## glm(formula = Y ~ 1, family = Gamma("inverse"), data = data)
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 9.254e-06  2.014e-07   45.96   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.08285843)
## 
##     Null deviance: 16.452  on 174  degrees of freedom
## Residual deviance: 16.452  on 174  degrees of freedom
## AIC: 4129.8
## 
## Number of Fisher Scoring iterations: 5
## Uji Simultan
# Melakukan uji likelihood ratio antara model penuh dan model null
lr_test <- lrtest(model4, modelNull)
print(lr_test)
## Likelihood ratio test
## 
## Model 1: Y ~ J1 + J2
## Model 2: Y ~ 1
##   #Df  LogLik Df  Chisq Pr(>Chisq)    
## 1   4  1231.0                         
## 2   2 -2062.9 -2 6587.7  < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Uji RMSE
dataa <- data.frame(Y, J1, J2)
predicted_values <- predict(model4, newdata = dataa, type="response")# Menghitung RMSE
rmse <- sqrt(mean((dataa$Y - predicted_values)^2))
print(paste("RMSE:", rmse))
## [1] "RMSE: 1.54814089365581e-10"