library(sp)
library(rgdal) # membaca shapefile
## Please note that rgdal will be retired during October 2023,
## plan transition to sf/stars/terra functions using GDAL and PROJ
## at your earliest convenience.
## See https://r-spatial.org/r/2023/05/15/evolution4.html and https://github.com/r-spatial/evolution
## rgdal: version: 1.6-7, (SVN revision (unknown))
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.6.2, released 2023/01/02
## Path to GDAL shared files: C:/Users/htauf/AppData/Local/R/win-library/4.4/rgdal/gdal
## GDAL does not use iconv for recoding strings.
## GDAL binary built with GEOS: TRUE
## Loaded PROJ runtime: Rel. 9.2.0, March 1st, 2023, [PJ_VERSION: 920]
## Path to PROJ shared files: C:/Users/htauf/AppData/Local/R/win-library/4.4/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.6-0
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading sp or rgdal.
library(spData)
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
library(readxl) # membaca file excel
library(raster) # text pada peta
library(spdep) # pembobot data spasial
## Loading required package: sf
## Linking to GEOS 3.12.2, GDAL 3.9.3, PROJ 9.4.1; sf_use_s2() is TRUE
library(lmtest) # regresi klasik
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(DescTools) # Uji dependensi
library(nortest) # uji normalitas
library(spatialreg) # running model SAR
## Loading required package: Matrix
##
## Attaching package: 'spatialreg'
## The following objects are masked from 'package:spdep':
##
## get.ClusterOption, get.coresOption, get.mcOption,
## get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
## set.coresOption, set.mcOption, set.VerboseOption,
## set.ZeroPolicyOption
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:DescTools':
##
## Recode
library(nortest)
library(GWmodel) # Model GWR
## Loading required package: robustbase
## Loading required package: Rcpp
## Welcome to GWmodel version 2.4-2.
library(ggplot2)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following object is masked from 'package:car':
##
## recode
## The following objects are masked from 'package:raster':
##
## intersect, select, union
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
# untuk melakukan pengaturan working directory
datasulsel <- read_excel("DATA.xlsx")
datasulsel <- datasulsel %>%
mutate(across(c(Y, X1, X2, X3, X4, X5, X6, X7), as.numeric))
data<-datasulsel[3:9]
data
## # A tibble: 24 × 7
## X1 X2 X3 X4 X5 X6 X7
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 72.6 13.6 8.54 11712 8.46 0.371 156
## 2 73.2 13 7.54 9682 10.5 0.389 180
## 3 74.0 13.4 8.26 11392 7.22 0.379 383
## 4 73.9 13.9 8.94 11636 12.7 0.365 130
## 5 74.0 13.7 8.41 10233 7.42 0.376 442
## 6 73.8 12.1 7 9781 13.1 0.34 521
## 7 74.3 13.0 8.93 13451 6.93 0.393 46
## 8 73.6 12.6 8.14 12513 12.7 0.342 45
## 9 73.3 13.4 8.73 10691 12.7 0.341 131
## 10 75.2 15.6 11.6 17889 5.07 0.387 8227
## # ℹ 14 more rows
summary(datasulsel)
## Kabupaten / Kota Y X1 X2
## Length:24 Min. :68.95 Min. :72.57 Min. :12.12
## Class :character 1st Qu.:71.42 1st Qu.:73.27 1st Qu.:12.95
## Mode :character Median :73.27 Median :73.69 Median :13.29
## Mean :73.84 Mean :73.75 Mean :13.39
## 3rd Qu.:74.46 3rd Qu.:74.32 3rd Qu.:13.64
## Max. :84.85 Max. :75.15 Max. :15.61
## X3 X4 X5 X6
## Min. : 7.000 Min. : 9.83 Min. : 5.070 Min. :0.3360
## 1st Qu.: 7.973 1st Qu.: 9756.25 1st Qu.: 7.370 1st Qu.:0.3450
## Median : 8.305 Median :11514.00 Median : 8.725 Median :0.3635
## Mean : 8.519 Mean :10339.58 Mean : 9.332 Mean :0.3618
## 3rd Qu.: 8.633 3rd Qu.:12604.00 3rd Qu.:12.322 3rd Qu.:0.3767
## Max. :11.560 Max. :17889.00 Max. :13.400 Max. :0.3930
## X7
## Min. : 45.0
## 1st Qu.: 146.5
## Median : 216.5
## Mean : 666.8
## 3rd Qu.: 461.5
## Max. :8227.0
dim(datasulsel)
## [1] 24 9
#menampilkan peta
#Peta Dasar Batas-batas Wilayah dari Sulsel
sulsel<-readOGR(dsn=".",layer = "sul-sel")
## Warning in readOGR(dsn = ".", layer = "sul-sel"): OGR support is provided by
## the sf and terra packages among others
## Warning in ogrInfo(dsn = dsn, layer = layer, encoding = encoding, use_iconv =
## use_iconv, : OGR support is provided by the sf and terra packages among others
## Warning in ogrFIDs(dsn = dsn, layer = layer): OGR support is provided by the sf
## and terra packages among others
## Warning in OGRSpatialRef(dsn, layer, morphFromESRI = morphFromESRI, dumpSRS =
## dumpSRS, : OGR support is provided by the sf and terra packages among others
## Warning in ogrListLayers(dsn): OGR support is provided by the sf and terra
## packages among others
## Warning in ogrFIDs(dsn = dsn, layer = layer): OGR support is provided by the sf
## and terra packages among others
## OGR data source with driver: ESRI Shapefile
## Source: "C:\Users\htauf\Documents\Documents\freelancer\C4\revisi", layer: "sul-sel"
## with 24 features
## It has 11 fields
## Integer64 fields read as strings: Y
kode=sulsel$ID
koordinat=coordinates(sulsel)
plot(sulsel)

sulsel$Kabupaten
## [1] "BARRU" "BONE" "BULUKUMBA" "ENREKANG" "GOWA"
## [6] "JENEPONTO" "LUWU TIMUR" "LUWU UTARA" "LUWU" "MAKASSAR"
## [11] "MAROS" "PALOPO" "PANGKEP" "PARE-PARE" "PINRANG"
## [16] "SELAYAR" "SIDRAP" "SINJAI" "SOPPENG" "TAKALAR"
## [21] "TORAJA UTARA" "TORAJA" "WAJO" "BANTAENG"
## eksplor data
colfunc<-colorRampPalette(c("green", "yellow", "orange","red"))
color<-colfunc(16)
sulsel$pertumbuhanUMKM<-datasulsel$Y
spplot(sulsel, "pertumbuhanUMKM", col.regions=color,
main="Peta Sebaran Pertumbuhan UMKM")

sulsel$pertumbuhanUMKM<-datasulsel$X1
spplot(sulsel, "X1", col.regions=color,
main="Peta Sebaran X1")

sulsel$pertumbuhanUMKM<-datasulsel$X2
spplot(sulsel, "X2", col.regions=color,
main="Peta Sebaran X2")

sulsel$pertumbuhanUMKM<-datasulsel$X3
spplot(sulsel, "X3", col.regions=color,
main="Peta Sebaran X3")

sulsel$pertumbuhanUMKM<-datasulsel$X4
spplot(sulsel, "X4", col.regions=color,
main="Peta Sebaran X4")

sulsel$pertumbuhanUMKM<-datasulsel$X5
spplot(sulsel, "X5", col.regions=color,
main="Peta Sebaran X5")

sulsel$pertumbuhanUMKM<-datasulsel$X6
spplot(sulsel, "X6", col.regions=color,
main="Peta Sebaran X6")

sulsel$pertumbuhanUMKM<-datasulsel$X7
spplot(sulsel, "X7", col.regions=color,
main="Peta Sebaran X7")

## plot
ggplot(datasulsel, aes(x = `Kabupaten / Kota` , y = Y, group = 1)) +
geom_line(color = "blue", size = 1) +
geom_point(color = "blue", size = 3) +
geom_text(aes(label = Y), position = 'dodge') +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1)) + # Rotasi teks kabupaten
labs(title = "Plot",
x = "Kabupaten/Kota",
y = "variabel Y")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning: Width not defined
## ℹ Set with `position_dodge(width = ...)`

## Moran's index
w = gw.dist(dp.locat = koordinat)
ww = gw.weight(vdist = w, bw = 2, kernel = 'gaussian', adaptive = T)
#dist = mat2listw(ww)
# Membuat matriks bobot spasial menggunakan Queen contiguity
nb <- poly2nb(sulsel, queen = TRUE) # Tetangga berbasis batas wilayah
dist <- nb2listw(nb, style = "W") # Konversi ke matriks bobot berbobot W
#nb <- poly2nb(sulsel)
#dist <- nb2listw(nb, style = "S")
moran(datasulsel$Y, dist, n=length(dist$neighbours), S0=Szero(dist))
## $I
## [1] 0.007111038
##
## $K
## [1] 5.286209
moran.test(datasulsel$Y, dist, randomisation = F)
##
## Moran I test under normality
##
## data: datasulsel$Y
## weights: dist
##
## Moran I statistic standard deviate = 0.35248, p-value = 0.3622
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.007111038 -0.043478261 0.020598751
moran.plot(datasulsel$Y, dist, labels=datasulsel$`Kabupaten / Kota`)

## regresi klasik
reg.klasik<-lm(Y~X1 + X2 + X3 + X4 + X5 + X6 +X7,datasulsel)
summary(reg.klasik)
##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = datasulsel)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.06155 -0.56936 -0.03539 0.52839 1.57532
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.159e+02 2.701e+01 4.292 0.00056 ***
## X1 -6.733e-01 3.558e-01 -1.892 0.07669 .
## X2 -5.711e-02 5.767e-01 -0.099 0.92235
## X3 2.340e+00 4.174e-01 5.605 3.94e-05 ***
## X4 1.988e-04 5.766e-05 3.448 0.00331 **
## X5 -3.337e-01 8.835e-02 -3.777 0.00165 **
## X6 -2.978e+01 1.373e+01 -2.168 0.04558 *
## X7 3.620e-04 1.599e-04 2.263 0.03788 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9277 on 16 degrees of freedom
## Multiple R-squared: 0.9541, Adjusted R-squared: 0.934
## F-statistic: 47.49 on 7 and 16 DF, p-value: 1.582e-09
## uji asumsi klasik dari model
err.regklasik<-residuals(reg.klasik)
#heterogenitas
bptest(reg.klasik) #digunakan modelnya dalam perhitungan
##
## studentized Breusch-Pagan test
##
## data: reg.klasik
## BP = 2.6302, df = 7, p-value = 0.917
#dependensi
RunsTest(err.regklasik)
##
## Runs Test for Randomness
##
## data: err.regklasik
## runs = 10, m = 12, n = 12, p-value = 0.3009
## alternative hypothesis: true number of runs is not equal the expected number
## sample estimates:
## median(x)
## -0.03539317
#Normalitas
ad.test(err.regklasik)
##
## Anderson-Darling normality test
##
## data: err.regklasik
## A = 0.31394, p-value = 0.523
qqnorm(err.regklasik,datax=T)
qqline(rnorm(length(err.regklasik),mean(err.regklasik),sd(err.regklasik)),datax=T, col="red")

## Spatial Autocorrelation: Uji Moran Index
#Moran Test pada error
moran.test(err.regklasik, dist,randomisation=T, alternative="greater")
##
## Moran I test under randomisation
##
## data: err.regklasik
## weights: dist
##
## Moran I statistic standard deviate = 0.2645, p-value = 0.3957
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.00498414 -0.04347826 0.02118052
## uji efek spasial
#Uji LM
LM<-lm.LMtests(reg.klasik, dist,test=c("LMerr", "LMlag","RLMerr","RLMlag","SARMA"))
## Please update scripts to use lm.RStests in place of lm.LMtests
summary(LM)
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
## data:
## model: lm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data =
## datasulsel)
## test weights: listw
##
## statistic parameter p.value
## RSerr 0.0010114 1 0.9746
## RSlag 0.5969856 1 0.4397
## adjRSerr 0.1402328 1 0.7080
## adjRSlag 0.7362071 1 0.3909
## SARMA 0.7372185 2 0.6917
#Model SEM
sem<-errorsarlm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ,data=datasulsel, dist)
summary(sem)
##
## Call:
## errorsarlm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = datasulsel,
## listw = dist)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.064248 -0.570503 -0.038926 0.526611 1.570754
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.1560e+02 2.2064e+01 5.2393 1.612e-07
## X1 -6.6790e-01 2.9038e-01 -2.3001 0.021440
## X2 -6.6468e-02 4.6961e-01 -0.1415 0.887444
## X3 2.3447e+00 3.4037e-01 6.8886 5.634e-12
## X4 1.9811e-04 4.6987e-05 4.2163 2.484e-05
## X5 -3.3343e-01 7.1964e-02 -4.6334 3.598e-06
## X6 -2.9747e+01 1.1219e+01 -2.6515 0.008013
## X7 3.6161e-04 1.3028e-04 2.7757 0.005508
##
## Lambda: -0.018088, LR test value: 0.002161, p-value: 0.96292
## Asymptotic standard error: 0.26693
## z-value: -0.067765, p-value: 0.94597
## Wald statistic: 0.0045921, p-value: 0.94597
##
## Log likelihood: -27.38695 for error model
## ML residual variance (sigma squared): 0.57366, (sigma: 0.7574)
## Number of observations: 24
## Number of parameters estimated: 10
## AIC: 74.774, (AIC for lm: 72.776)
# Uji Asumsi
err.sem<-residuals(sem)
ad.test(err.sem)
##
## Anderson-Darling normality test
##
## data: err.sem
## A = 0.31268, p-value = 0.5255
moran.test(err.sem,dist,randomisation=T, alternative="greater")
##
## Moran I test under randomisation
##
## data: err.sem
## weights: dist
##
## Moran I statistic standard deviate = 0.27912, p-value = 0.3901
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.002848914 -0.043478261 0.021188896
#Model GSM/SARMA
gsm<-sacsarlm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ,data=datasulsel, dist)
summary(gsm)
##
## Call:
## sacsarlm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = datasulsel,
## listw = dist)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.078060 -0.543789 -0.074447 0.442193 1.525376
##
## Type: sac
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.2076e+02 2.3345e+01 5.1727 2.307e-07
## X1 -5.8819e-01 2.8558e-01 -2.0596 0.039435
## X2 -1.5068e-01 4.5410e-01 -0.3318 0.740018
## X3 2.3966e+00 3.3200e-01 7.2186 5.254e-13
## X4 2.0866e-04 4.8517e-05 4.3007 1.702e-05
## X5 -3.3486e-01 6.9758e-02 -4.8004 1.584e-06
## X6 -3.1383e+01 1.1286e+01 -2.7806 0.005425
## X7 3.3220e-04 1.2992e-04 2.5569 0.010561
##
## Rho: -0.13408
## Asymptotic standard error: 0.13782
## z-value: -0.97286, p-value: 0.33062
## Lambda: -0.10607
## Asymptotic standard error: 0.30354
## z-value: -0.34945, p-value: 0.72675
##
## LR test value: 0.842, p-value: 0.65639
##
## Log likelihood: -26.96704 for sac model
## ML residual variance (sigma squared): 0.54982, (sigma: 0.7415)
## Number of observations: 24
## Number of parameters estimated: 11
## AIC: 75.934, (AIC for lm: 72.776)
#uji asumsi
err.gsm<-residuals(gsm)
ad.test(err.gsm)
##
## Anderson-Darling normality test
##
## data: err.gsm
## A = 0.45677, p-value = 0.2432
moran.test(err.gsm,dist,randomisation=T, alternative="greater")
##
## Moran I test under randomisation
##
## data: err.gsm
## weights: dist
##
## Moran I statistic standard deviate = 0.19239, p-value = 0.4237
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.01550888 -0.04347826 0.02113443
## moran test
moran.test(datasulsel$Y,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$Y
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.37756, p-value = 0.3529
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.007111038 -0.043478261 0.017953794
moran.test(datasulsel$X1,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X1
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 1.1758, p-value = 0.1198
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.12770161 -0.04347826 0.02119626
moran.test(datasulsel$X2,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X2
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.50943, p-value = 0.3052
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.02687401 -0.04347826 0.01907127
moran.test(datasulsel$X3,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X3
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.64342, p-value = 0.26
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.04453481 -0.04347826 0.01871139
moran.test(datasulsel$X4,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X4
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 1.3837, p-value = 0.08323
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.14702262 -0.04347826 0.01895530
moran.test(datasulsel$X5,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X5
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.9175, p-value = 0.1794
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.09174598 -0.04347826 0.02172193
moran.test(datasulsel$X6,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X6
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = -0.098225, p-value = 0.5391
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.05796377 -0.04347826 0.02174814
moran.test(datasulsel$X7,nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X7
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.79781, p-value = 0.2125
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.004067919 -0.043478261 0.002440180
## Autokorelasi spasial pada residual
# Uji autokorelasi spasial dengan Moran's I pada residual model regresi klasik
moran.test(err.regklasik, dist, randomisation = TRUE, alternative = "greater")
##
## Moran I test under randomisation
##
## data: err.regklasik
## weights: dist
##
## Moran I statistic standard deviate = 0.2645, p-value = 0.3957
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.00498414 -0.04347826 0.02118052
moran.plot(err.regklasik, dist, labels=datasulsel$`Kabupaten / Kota`)

# Uji autokorelasi spasial dengan Moran's I pada residual model SEM (Spatial Error Model)
moran.test(err.sem, dist, randomisation = TRUE, alternative = "greater")
##
## Moran I test under randomisation
##
## data: err.sem
## weights: dist
##
## Moran I statistic standard deviate = 0.27912, p-value = 0.3901
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.002848914 -0.043478261 0.021188896
moran.plot(err.sem, dist, labels=datasulsel$`Kabupaten / Kota`)

# Uji autokorelasi spasial dengan Moran's I pada residual model GSM/SARMA (Generalized Spatial Model)
moran.test(err.gsm, dist, randomisation = TRUE, alternative = "greater")
##
## Moran I test under randomisation
##
## data: err.gsm
## weights: dist
##
## Moran I statistic standard deviate = 0.19239, p-value = 0.4237
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.01550888 -0.04347826 0.02113443
moran.plot(err.gsm, dist, labels=datasulsel$`Kabupaten / Kota`)

# Anda bisa menguji juga setiap variabel independen, jika diperlukan
# Moran's Test dan Plot untuk variabel X1
moran.test(datasulsel$X1, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X1
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 1.1758, p-value = 0.1198
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.12770161 -0.04347826 0.02119626
moran.plot(datasulsel$X1, dist, labels=datasulsel$`Kabupaten / Kota`)

# Moran's Test dan Plot untuk variabel X2
moran.test(datasulsel$X2, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X2
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.50943, p-value = 0.3052
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.02687401 -0.04347826 0.01907127
moran.plot(datasulsel$X2, dist, labels=datasulsel$`Kabupaten / Kota`)

# Moran's Test dan Plot untuk variabel X3
moran.test(datasulsel$X3, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X3
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.64342, p-value = 0.26
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.04453481 -0.04347826 0.01871139
moran.plot(datasulsel$X3, dist, labels=datasulsel$`Kabupaten / Kota`)

# Moran's Test dan Plot untuk variabel X4
moran.test(datasulsel$X4, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X4
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 1.3837, p-value = 0.08323
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.14702262 -0.04347826 0.01895530
moran.plot(datasulsel$X4, dist, labels=datasulsel$`Kabupaten / Kota`)

# Moran's Test dan Plot untuk variabel X5
moran.test(datasulsel$X5, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X5
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.9175, p-value = 0.1794
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.09174598 -0.04347826 0.02172193
moran.plot(datasulsel$X5, dist, labels=datasulsel$`Kabupaten / Kota`)

# Moran's Test dan Plot untuk variabel X6
moran.test(datasulsel$X6, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X6
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = -0.098225, p-value = 0.5391
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.05796377 -0.04347826 0.02174814
moran.plot(datasulsel$X6, mat2listw(ww, style="B"), labels = datasulsel$`Kabupaten / Kota`)

# Moran's Test dan Plot untuk variabel X7
moran.test(datasulsel$X7, nb2listw(nb, style = "W"))
##
## Moran I test under randomisation
##
## data: datasulsel$X7
## weights: nb2listw(nb, style = "W")
##
## Moran I statistic standard deviate = 0.79781, p-value = 0.2125
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.004067919 -0.043478261 0.002440180
moran.plot(datasulsel$X7, dist, labels=datasulsel$`Kabupaten / Kota`)

## peta baru (model SARMA)
coef_gsm <- coef(gsm)
# Hitung prediksi Y secara manual dari koefisien SARMA
pred_gsm <- coef_gsm[1] + coef_gsm[2]*datasulsel$X1 + coef_gsm[3]*datasulsel$X2 +
coef_gsm[4]*datasulsel$X3 + coef_gsm[5]*datasulsel$X4 +
coef_gsm[6]*datasulsel$X5 + coef_gsm[7]*datasulsel$X6 +
coef_gsm[8]*datasulsel$X7
# Menyimpan prediksi ke shapefile Sulsel
sulsel$pred_SARMA <- pred_gsm
# Menggunakan 4 warna untuk gradasi
colfunc <- colorRampPalette(c("green", "yellow", "orange", "red"))
color<-colfunc(16)
# Visualisasi peta prediksi SARMA dengan 4 warna
# Menampilkan peta prediksi SARMA
spplot(sulsel, "pred_SARMA", col.regions = color,
main = "Peta Prediksi Pertumbuhan UMKM dengan Model SARMA (4 Warna)",
sp.layout = list("sp.text", coordinates(sulsel), sulsel$Kabupaten, cex = 0.5, col = "black"))
