PEMODELAN LOCALLY COMPENSATED RIDGE IMPROVED GEOGRAPHICALLY TEMPORALLY WEIGHTED REGRESSION PADA DATA TINGKAT FERTILITAS (TFR) SULAWESI DAFTAR ISI 0. Pengaturan dan paket 1. Data 2. OLS 3. Bandwidth GWR (pemilihan kernel) 4. Pengujian asumsi (normalitas, heterogenitas spasial, multikolinearitas) 5. GWR 6. Matriks pembobot spasial W, Moran’s I, dan uji LM 7. GWR-SAR: bandwidth, estimator, dan estimasi 8. Perbandingan model

rm(list = ls())
# ==============================================================================
# 0. PENGATURAN DAN PAKET
# ==============================================================================

# ---- 0.1 Lokasi file (sesuaikan) ----
dir_kerja  <- "C:/Users/ASUS/Downloads"
file_func  <- file.path(dir_kerja, "FUNCTION_LCR-IGTWRne.R")
file_data  <- file.path(dir_kerja, "Data_TFR_Sulawesi_2025.xlsx")
sheet_data <- "Data Gabungan"

# ---- 0.2 Pengaturan model ----
alpha  <- 0.05              # taraf signifikansi

# ---- 0.3 Paket ----
# install.packages(c("readxl", "openxlsx", "lmtest", "car", "corrplot",
#                    "spdep", "MASS", "sp", "spgwr", "GWmodel"))   # cukup sekali

paket <- c("readxl",    # membaca data Excel
           "openxlsx",  # menyimpan hasil ke Excel
           "lmtest",    # bptest()
           "car",       # vif()
           "corrplot",  # plot korelasi
           "spdep",     # mat2listw(), moran.test()
           "MASS",      # ginv()
           "sp", "spgwr", "GWmodel")   # gold() untuk bandwidth CV
invisible(lapply(paket, library, character.only = TRUE))
## Warning: package 'readxl' was built under R version 4.5.3
## Warning: package 'openxlsx' was built under R version 4.5.3
## Warning: package 'lmtest' was built under R version 4.5.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.5.3
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Warning: package 'car' was built under R version 4.5.3
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.5.3
## Warning: package 'corrplot' was built under R version 4.5.3
## corrplot 0.95 loaded
## Warning: package 'spdep' was built under R version 4.5.3
## Loading required package: spData
## Warning: package 'spData' was built under R version 4.5.3
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
## Loading required package: sf
## Warning: package 'sf' was built under R version 4.5.3
## Linking to GEOS 3.14.1, GDAL 3.12.1, PROJ 9.7.1; sf_use_s2() is TRUE
## Warning: package 'MASS' was built under R version 4.5.3
## Warning: package 'sp' was built under R version 4.5.3
## Warning: package 'spgwr' was built under R version 4.5.3
## NOTE: This package does not constitute approval of GWR
## as a method of spatial analysis; see example(gwr)
## Warning: package 'GWmodel' was built under R version 4.5.3
## Loading required package: robustbase
## Warning: package 'robustbase' was built under R version 4.5.3
## Loading required package: Rcpp
## Warning: package 'Rcpp' was built under R version 4.5.3
## Welcome to GWmodel version 2.4-2.
# ---- 0.4 Fungsi buatan sendiri ----
source(file_func)

1. DATA

data <- data.frame(read_excel(file_data, sheet = sheet_data))

y     <- data[, 5]              # respon: TFR
x     <- data[, c(6:10)] # prediktor: RLS, TPAKP, TPT, MISKIN
coord <- data[, 3:4]            # koordinat: Longitude, Latitude

n <- nrow(data)
D <- jarak.euclid(coord)        # matriks jarak (dipakai di bagian 6 dan 7)

str(data)
## 'data.frame':    81 obs. of  10 variables:
##  $ prov      : chr  "Sulawesi Utara" "Sulawesi Utara" "Sulawesi Utara" "Sulawesi Utara" ...
##  $ kabkot    : chr  "Bolaang Mongondow" "Minahasa" "Kepulauan Sangihe" "Kepulauan Talaud" ...
##  $ Longitude : num  282361 373638 450241 588654 335866 ...
##  $ Latitude  : num  299607 359465 613607 687412 341272 ...
##  $ TFR       : num  2.16 1.99 2.07 2.16 2.06 2.02 2.19 2.06 2.24 2.27 ...
##  $ RLS       : num  8.49 10.25 8.7 9.97 9.53 ...
##  $ TPAKP     : num  54.2 48.9 54.5 50.8 50 ...
##  $ TPAKL     : num  87.6 80.1 81.2 81 84.4 ...
##  $ TPT       : num  4.97 7.74 2.64 3.4 4.77 6.8 4.07 2.41 3.38 2.91 ...
##  $ MISKIN....: num  6.98 5.88 10.91 8.34 8.26 ...

2. OLS

OLS <- lm(y ~ ., data = data.frame(y = y, x))
summary(OLS)
## 
## Call:
## lm(formula = y ~ ., data = data.frame(y = y, x))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.34952 -0.13151 -0.00616  0.07845  0.50553 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.784685   0.601391   4.630  1.5e-05 ***
## RLS          0.003450   0.024733   0.140  0.88943    
## TPAKP        0.008101   0.002895   2.798  0.00653 ** 
## TPAKL       -0.012853   0.006009  -2.139  0.03570 *  
## TPT         -0.034674   0.018505  -1.874  0.06486 .  
## MISKIN....   0.025759   0.007585   3.396  0.00110 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1793 on 75 degrees of freedom
## Multiple R-squared:  0.4327, Adjusted R-squared:  0.3949 
## F-statistic: 11.44 on 5 and 75 DF,  p-value: 3.182e-08

3. BANDWIDTH GWR (PEMILIHAN KERNEL)

kernel_uji <- c("exponential", "gaussian", "bisquare", "tricube")

bw_kernel <- sapply(kernel_uji, function(k)
  bwd.gwr.twr(y = y, x = x, long.lat = coord, spasial = TRUE, kernel = k))
## Fixed bandwidth: 913788.5 CV score: 2.5205 
## Fixed bandwidth: 564865.3 CV score: 2.340359 
## Fixed bandwidth: 349218.9 CV score: 2.103019 
## Fixed bandwidth: 215942.1 CV score: 1.845772 
## Fixed bandwidth: 133572.5 CV score: 1.64356 
## Fixed bandwidth: 82665.27 CV score: 1.518837 
## Fixed bandwidth: 51202.88 CV score: 1.463511 
## Fixed bandwidth: 31758.06 CV score: 1.688484 
## Fixed bandwidth: 63220.45 CV score: 1.468245 
## Fixed bandwidth: 43775.62 CV score: 1.509792 
## Fixed bandwidth: 55793.18 CV score: 1.460464 
## Fixed bandwidth: 58630.15 CV score: 1.461897 
## Fixed bandwidth: 54039.85 CV score: 1.46072 
## Fixed bandwidth: 56876.81 CV score: 1.460768 
## Fixed bandwidth: 55123.47 CV score: 1.460446 
## Fixed bandwidth: 913788.5 CV score: 2.725857 
## Fixed bandwidth: 564865.3 CV score: 2.526895 
## Fixed bandwidth: 349218.9 CV score: 2.166979 
## Fixed bandwidth: 215942.1 CV score: 1.850646 
## Fixed bandwidth: 133572.5 CV score: 1.637793 
## Fixed bandwidth: 82665.27 CV score: 1.560676 
## Fixed bandwidth: 51202.88 CV score: 2.571048 
## Fixed bandwidth: 102110.1 CV score: 1.552873 
## Fixed bandwidth: 114127.7 CV score: 1.56719 
## Fixed bandwidth: 94682.83 CV score: 1.549737 
## Fixed bandwidth: 90092.53 CV score: 1.55978 
## Fixed bandwidth: 97519.8 CV score: 1.548459 
## Fixed bandwidth: 99273.13 CV score: 1.549094 
## Fixed bandwidth: 96436.17 CV score: 1.548582 
## Fixed bandwidth: 98189.51 CV score: 1.548586 
## Fixed bandwidth: 97105.89 CV score: 1.548456 
## Fixed bandwidth: 913788.5 CV score: 2.229115 
## Fixed bandwidth: 564865.3 CV score: 1.914661 
## Fixed bandwidth: 349218.9 CV score: 1.702733 
## Fixed bandwidth: 215942.1 CV score: 2.246674 
## Fixed bandwidth: 431588.5 CV score: 1.850383 
## Fixed bandwidth: 298311.7 CV score: 1.690326 
## Fixed bandwidth: 266849.3 CV score: 2.568912 
## Fixed bandwidth: 317756.5 CV score: 1.682905 
## Fixed bandwidth: 329774.1 CV score: 1.676433 
## Fixed bandwidth: 337201.3 CV score: 1.686054 
## Fixed bandwidth: 325183.8 CV score: 1.692446 
## Fixed bandwidth: 332611 CV score: 1.680005 
## Fixed bandwidth: 328020.7 CV score: 1.67429 
## Fixed bandwidth: 326937.1 CV score: 1.672986 
## Fixed bandwidth: 326267.4 CV score: 1.694112 
## Fixed bandwidth: 327351 CV score: 1.673481 
## Fixed bandwidth: 326681.3 CV score: 1.672685 
## Fixed bandwidth: 326523.2 CV score: 1.694514 
## Fixed bandwidth: 326779 CV score: 1.6728 
## Fixed bandwidth: 326620.9 CV score: 1.672614 
## Fixed bandwidth: 913788.5 CV score: 2.294018 
## Fixed bandwidth: 564865.3 CV score: 1.961853 
## Fixed bandwidth: 349218.9 CV score: 1.732014 
## Fixed bandwidth: 215942.1 CV score: 2.569846 
## Fixed bandwidth: 431588.5 CV score: 1.938748 
## Fixed bandwidth: 298311.7 CV score: 1.719732 
## Fixed bandwidth: 266849.3 CV score: 3.160788 
## Fixed bandwidth: 317756.5 CV score: 1.718661 
## Fixed bandwidth: 329774.1 CV score: 1.730365 
## Fixed bandwidth: 310329.2 CV score: 1.720281 
## Fixed bandwidth: 322346.8 CV score: 1.721571 
## Fixed bandwidth: 314919.5 CV score: 1.718216 
## Fixed bandwidth: 313166.2 CV score: 1.718566 
## Fixed bandwidth: 316003.2 CV score: 1.718246
bw_kernel
## exponential    gaussian    bisquare     tricube 
##    55123.47    97105.89   326620.91   314919.54
bandwith.GWR <- unname(bw_kernel[1])   # bandwidth kernel terpilih
h_gwr        <- bandwith.GWR

4. PENGUJIAN ASUMSI

# ---- 4.1 Normalitas sisaan ----
sisaan_ols <- residuals(OLS)
print(ks.test(sisaan_ols, "pnorm", mean = mean(sisaan_ols), sd = sd(sisaan_ols)))
## 
##  Exact one-sample Kolmogorov-Smirnov test
## 
## data:  sisaan_ols
## D = 0.10947, p-value = 0.2663
## alternative hypothesis: two-sided
print(shapiro.test(sisaan_ols))   # lebih tepat untuk n kecil dan parameter diduga
## 
##  Shapiro-Wilk normality test
## 
## data:  sisaan_ols
## W = 0.96985, p-value = 0.05296
# ---- 4.2 Heterogenitas spasial (Breusch-Pagan) ----
BPTest_Spasial <- bptest(OLS, studentize = FALSE)
BPTest_Spasial
## 
##  Breusch-Pagan test
## 
## data:  OLS
## BP = 16.376, df = 5, p-value = 0.005848
# ---- 4.3 Multikolinearitas (VIF) ----
Multi <- vif(OLS)
Multi
##        RLS      TPAKP      TPAKL        TPT MISKIN.... 
##   1.879287   1.376575   1.395072   2.403514   1.610396
# ---- 4.4 Korelasi antar variabel ----
Korelasi <- corrplot(cor(data.frame(TFR = y, x)), method = "number")

5. GWR

gwr <- prog.gwr(y = y, x = x, long.lat = coord, bwd = bandwith.GWR, kernel = "exponential")

# ---- 5.1 Estimasi parameter ----
Estimasi_GWR <- data.frame(gwr$beta)
summary(Estimasi_GWR)
##        X1               X2                 X3                  X4           
##  Min.   :-1.360   Min.   :-0.15172   Min.   :-0.009103   Min.   :-0.039108  
##  1st Qu.: 1.218   1st Qu.:-0.02596   1st Qu.:-0.001397   1st Qu.:-0.014609  
##  Median : 1.780   Median : 0.01010   Median : 0.004026   Median :-0.004396  
##  Mean   : 2.074   Mean   : 0.01907   Mean   : 0.003078   Mean   :-0.002695  
##  3rd Qu.: 3.197   3rd Qu.: 0.07061   3rd Qu.: 0.006717   3rd Qu.: 0.010542  
##  Max.   : 5.955   Max.   : 0.15846   Max.   : 0.014776   Max.   : 0.031342  
##        X5                 X6          
##  Min.   :-0.12970   Min.   :-0.03059  
##  1st Qu.:-0.05374   1st Qu.: 0.01044  
##  Median :-0.03359   Median : 0.01877  
##  Mean   :-0.03424   Mean   : 0.02288  
##  3rd Qu.:-0.01170   3rd Qu.: 0.03670  
##  Max.   : 0.10758   Max.   : 0.07435
# ---- 5.2 Matriks jarak dan pembobot ----
Matriks_Jarak      <- data.frame(gwr$jarak)
Matriks_Pembobot_W <- data.frame(gwr$bobot)
# write.csv(Matriks_Pembobot_W, file = "Matriks_Pembobot_W.csv", row.names = FALSE)

# ---- Uji-t lokal GWR biasa ----
inf_gwr <- inferensi.gwr(gwr, y, x)
sig_gwr <- uji.signifikansi(inf_gwr, id = data[, c("prov", "kabkot")], alpha = alpha)
print(sig_gwr$ringkasan, digits = 4)
##    Parameter       Min    Median      Mean     Max  Mean_t N_Signif Prop_Signif
## 1  Intercept -1.359595  1.779968  2.074339 5.95455  1.9803       36      0.4444
## 2        RLS -0.151719  0.010100  0.019073 0.15846  0.2934        9      0.1111
## 3      TPAKP -0.009103  0.004026  0.003078 0.01478  0.8446       13      0.1605
## 4      TPAKL -0.039108 -0.004396 -0.002695 0.03134 -0.2131       17      0.2099
## 5        TPT -0.129701 -0.033587 -0.034245 0.10758 -0.9911       15      0.1852
## 6 MISKIN.... -0.030593  0.018774  0.022879 0.07435  1.4218       33      0.4074
##   N_Positif N_Negatif                        Status
## 1        36         0 Signifikan di sebagian lokasi
## 2         5         4 Signifikan di sebagian lokasi
## 3        13         0 Signifikan di sebagian lokasi
## 4         7        10 Signifikan di sebagian lokasi
## 5         0        15 Signifikan di sebagian lokasi
## 6        30         3 Signifikan di sebagian lokasi
# Buat workbook baru
wb <- createWorkbook()

# Tambahkan sheet untuk Ringkasan
addWorksheet(wb, sheetName = "Ringkasan")
writeData(wb, sheet = "Ringkasan", x = sig_gwr$ringkasan)

# Tambahkan sheet untuk Detail Lokal
addWorksheet(wb, sheetName = "Lokal")
writeData(wb, sheet = "Lokal", x = sig_gwr$lokal)

# Simpan file Excel
saveWorkbook(wb, file = "Uji_Sig_Parameter0.xlsx", overwrite = TRUE)

cat("File 'Uji_Sig_Parameter0.xlsx' berhasil disimpan.\n")
## File 'Uji_Sig_Parameter0.xlsx' berhasil disimpan.

6. MATRIKS PEMBOBOT SPASIAL W, MORAN’S I, DAN UJI LM

kernel<- "exponential"
# ---- 6.1 W lag: kernel, diagonal 0, row-normalized ----
W_lag <- buat.W.lag(D, h_gwr, kernel)

# ---- 6.2 Moran's I pada sisaan GWR ----
sisaan_gwr <- drop((diag(n) - gwr$hat_matrix) %*% as.numeric(y))
lw         <- mat2listw(W_lag, style = "W")

moran_gwr <- moran.test(sisaan_gwr, lw, randomisation = TRUE, alternative = "two.sided")
print(moran_gwr)
## 
##  Moran I test under randomisation
## 
## data:  sisaan_gwr  
## weights: lw    
## 
## Moran I statistic standard deviate = -1.9899, p-value = 0.0466
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      -0.105070194      -0.012500000       0.002164029
# ---- 6.3 Uji LM spasial pada sisaan GWR ----
hasil_lm  <- lm.spatial.tests.gwr(gwr, y, W_lag)
keputusan <- pilih.model.lm(hasil_lm, alpha = alpha)
print(hasil_lm, digits = 5)
##      Uji Statistik df    p.value
## 1  LMerr   4.16291  1 4.1318e-02
## 2  LMlag  43.65755  1 3.9117e-11
## 3 RLMerr   0.28697  1 5.9217e-01
## 4 RLMlag  39.78161  1 2.8401e-10
cat("=> Keputusan LM:", keputusan, "\n")
## => Keputusan LM: SAR

7. GWR-SAR

q <- 2 * ncol(x) + 1      # jumlah kolom instrumen

# ---- 7.1 Bandwidth dan estimator ----
# (a) Awal: WLS
bw.sar    <- bwd.gwr.sar(y, x, coord, W_lag, kernel, "wls")
estimator <- "wls"
ESS_wls   <- hitung.ESS(D, bw.sar$bw, kernel)
cat(sprintf("WLS: h = %.2f, min ESS = %.2f, q = %d\n", bw.sar$bw, min(ESS_wls), q))
## WLS: h = 69670.95, min ESS = 1.36, q = 11
# (b) 2SLS dipakai hanya jika ESS minimum pada bandwidth BARU masih > q
if (min(ESS_wls) > q) {
  bw2  <- bwd.gwr.sar(y, x, coord, W_lag, kernel, "2sls")
  ESS2 <- hitung.ESS(D, bw2$bw, kernel)
  cat(sprintf("2SLS: h = %.2f, min ESS = %.2f\n", bw2$bw, min(ESS2)))
  if (min(ESS2) > q) {
    bw.sar    <- bw2
    estimator <- "2sls"
  } else {
    message("ESS 2SLS <= q pada bandwidth baru -> tetap WLS")
  }
}

# (c) Cek batas grid
if (bw.sar$idx_grid %in% c(1, nrow(bw.sar$grid)))
  warning("Optimum jatuh di batas grid: perlebar rentang bandwidth.")

# ---- 7.2 Estimasi akhir ----
fit <- prog.gwr.sar(y, x, coord, W_lag, bw.sar$bw, kernel, estimator)
cat(sprintf("=> Estimator = %s, h optimum = %.2f (CV = %.5f)\n",
            fit$estimator, bw.sar$bw, bw.sar$cv))
## => Estimator = wls, h optimum = 69670.95 (CV = 1.38291)
# ---- 7.3 Parameter lokal ----
Estimasi_GWRSAR <- data.frame(coord, fit$beta)
print(summary(Estimasi_GWRSAR))
##    Longitude          Latitude            Rho            Intercept       
##  Min.   :-286460   Min.   :-504107   Min.   :-0.6308   Min.   :-0.62847  
##  1st Qu.:-176467   1st Qu.:-270339   1st Qu.: 0.4194   1st Qu.: 0.08699  
##  Median :   8909   Median :-104925   Median : 0.6130   Median : 0.47627  
##  Mean   :  25406   Mean   : -26631   Mean   : 0.5800   Mean   : 0.79280  
##  3rd Qu.: 172975   3rd Qu.: 281036   3rd Qu.: 0.7739   3rd Qu.: 1.49696  
##  Max.   : 588654   Max.   : 687412   Max.   : 1.2689   Max.   : 3.15215  
##       RLS                TPAKP               TPAKL          
##  Min.   :-0.116678   Min.   :-0.009176   Min.   :-0.038878  
##  1st Qu.:-0.040113   1st Qu.:-0.001596   1st Qu.:-0.011760  
##  Median : 0.017348   Median : 0.002685   Median :-0.002065  
##  Mean   : 0.005107   Mean   : 0.002432   Mean   :-0.001311  
##  3rd Qu.: 0.041101   3rd Qu.: 0.006074   3rd Qu.: 0.013714  
##  Max.   : 0.171802   Max.   : 0.011913   Max.   : 0.029306  
##       TPT              MISKIN....      
##  Min.   :-0.133296   Min.   :-0.01798  
##  1st Qu.:-0.051126   1st Qu.: 0.01001  
##  Median :-0.026650   Median : 0.01482  
##  Mean   :-0.030466   Mean   : 0.01890  
##  3rd Qu.:-0.005126   3rd Qu.: 0.02578  
##  Max.   : 0.077838   Max.   : 0.05975
Rho_lokal <- fit$beta[, "Rho"]
cat(sprintf("Proporsi |rho_i| >= 1 (melanggar stasioneritas): %.3f\n",
            mean(abs(Rho_lokal) >= 1)))
## Proporsi |rho_i| >= 1 (melanggar stasioneritas): 0.099
cat(sprintf("Proporsi lokasi signifikan (alpha = %g):\n", alpha))
## Proporsi lokasi signifikan (alpha = 0.05):
print(round(colMeans(fit$pval < alpha), 3))
##        Rho  Intercept        RLS      TPAKP      TPAKL        TPT MISKIN.... 
##      0.556      0.099      0.062      0.210      0.173      0.198      0.259
# Jika Anda tetap ingin melihat Uji t Lokal per wilayah:
uji_sig <- uji.signifikansi(fit, id = data[, c("prov", "kabkot")], alpha = alpha)
print(uji_sig$ringkasan, digits = 4)
##    Parameter       Min    Median      Mean     Max   Mean_t N_Signif
## 1        Rho -0.630806  0.613029  0.579957 1.26887  1.94081       45
## 2  Intercept -0.628469  0.476271  0.792801 3.15215  0.97684        8
## 3        RLS -0.116678  0.017348  0.005107 0.17180 -0.04305        5
## 4      TPAKP -0.009176  0.002685  0.002432 0.01191  0.81033       17
## 5      TPAKL -0.038878 -0.002065 -0.001311 0.02931 -0.17833       14
## 6        TPT -0.133296 -0.026650 -0.030466 0.07784 -0.97919       16
## 7 MISKIN.... -0.017978  0.014823  0.018901 0.05975  1.22429       21
##   Prop_Signif N_Positif N_Negatif                        Status
## 1     0.55556        45         0 Signifikan di sebagian lokasi
## 2     0.09877         8         0 Signifikan di sebagian lokasi
## 3     0.06173         0         5 Signifikan di sebagian lokasi
## 4     0.20988        17         0 Signifikan di sebagian lokasi
## 5     0.17284         3        11 Signifikan di sebagian lokasi
## 6     0.19753         0        16 Signifikan di sebagian lokasi
## 7     0.25926        21         0 Signifikan di sebagian lokasi
# Buat workbook baru
wb <- createWorkbook()

# Tambahkan sheet untuk Ringkasan
addWorksheet(wb, sheetName = "Ringkasan")
writeData(wb, sheet = "Ringkasan", x = uji_sig$ringkasan)

# Tambahkan sheet untuk Detail Lokal
addWorksheet(wb, sheetName = "Lokal")
writeData(wb, sheet = "Lokal", x = uji_sig$lokal)

# Simpan file Excel
saveWorkbook(wb, file = "Uji_Sig_Parameter.xlsx", overwrite = TRUE)

cat("File 'Uji_Sig_Parameter.xlsx' berhasil disimpan.\n")
## File 'Uji_Sig_Parameter.xlsx' berhasil disimpan.

8. PERBANDINGAN MODEL: OLS vs GWR vs GWR-SAR

gwr0  <- prog.gwr.sar(y, x, coord, W_lag, h_gwr, kernel, lag = FALSE)   # GWR tanpa lag (AICc sebanding)
k_ols <- ncol(x) + 1

mse_ols <- mean(residuals(OLS)^2)
mse_gwr <- mean(gwr0$residuals^2)
mse_sar <- mean(fit$residuals^2)

perbandingan <- data.frame(
  Model = c("OLS", "GWR", "GWR-SAR"),
  R2    = c(summary(OLS)$r.squared,     gwr0$R2,    fit$R2),
  R2adj = c(summary(OLS)$adj.r.squared, gwr0$R2adj, fit$R2adj),
  MSE   = c(mse_ols, mse_gwr, mse_sar),
  RMSE  = sqrt(c(mse_ols, mse_gwr, mse_sar)),
  AICc  = c(aicc.fun(sum(residuals(OLS)^2), n - k_ols, k_ols, n), gwr0$AICc, fit$AICc))
print(perbandingan, digits = 5)
##     Model      R2   R2adj       MSE     RMSE    AICc
## 1     OLS 0.43271 0.39489 0.0297525 0.172489 -33.066
## 2     GWR 0.95732 0.74143 0.0022382 0.047310 299.583
## 3 GWR-SAR 0.94098 0.73740 0.0030956 0.055638 189.377