a. Model Dependensi Spasial

Deskripsi

Analisis regresi spasial digunakan untuk melakukan pemodelan hubungan antara peubah respon dan prediktor dalam bentuk model regresi dengan mempertimbangkan efek spasial. Efek spasial antar lokasi dapat disebabkan oleh dependensi spasial dan heterogenitas spasial (Anselin, 1988). Dependensi spasial adalah ketergantungan antara suatu lokasi dengan lokasi di sekitarnya sedangkan heterogenitas spasial adalah perbedaan pengaruh peubah prediktor terhadap respon di lokasi yang diamati. Eigenvector Spasial Filtering (ESF) adalah metode untuk mengatasi masalah dependensi spasial dengan menambahkan kombinasi linier dari vektor ciri matriks pembobot spasial pada spesifikasi model regresi (Chun & Griffith, 2014). Eigenvector Spasial Filtering (ESF) dapat mengurangi kesalahan spesifikasi spasial, meningkatkan kesesuaian model, serta kenormalan dan homogenitas ragam sisaan. Metode lain yang dapat mengatasi masalah dependensi spasial adalah Random Effect Eigenvector Spatial Filtering (RE-ESF).

Murakami & Griffith (2015) memperluas Eigenvector Spatial Filtering (ESF) yang merupakan pendekatan efek tetap, ke pendekatan efek acak yaitu Random Effect Eigenvector Spatial Filtering (RE-ESF). Random Effect Eigenvector Spatial Filtering (RE-ESF) mempertimbangkan efek spasial sebagai variabel acak yang bervariasi di antara wilayah-wilayah yang diamati. Hughes & Haran (2013) membahas kegunaan Random Effect Eigenvector Spatial Filtering (RE-ESF) untuk analisis dependensi spasial yang meminimumkan nilai ragam yang disebabkan oleh faktor acak antara peubah prediktor dan proses spasial laten.

Eigenvector Spatial Filtering (ESF) yang juga disebut koordinat utama matriks ketetanggaan (Dray et al., 2006) atau pemetaan vektor ciri Moran didasarkan pada koefisien Moran (MC). Koefisien Moran untuk vektor peubah acak \(y\) didefinisikan sebagai berikut:

\[ MC[y] = \frac{n}{(1' W1)} \left( y' MWMy \right) / \left( y' My \right) \tag{1}\]

dengan \(W\) adalah matriks bobot spasial simetris berukuran \(n \times n\) dengan diagonal nol, dan \(M = (I-11'/n)\) adalah matriks pemusatan. Misalkan dekomposisi eigen matriks \(MWM\) sebagai \(E^* \Lambda^* E^{*'}\), dengan \(E^* = [e_1, \ldots, e_n]\) adalah matriks vektor ciri. \(\Lambda^* = \text{diag}(\lambda_1, \ldots, \lambda_n)\) adalah matriks diagonal yang elemennya adalah akar ciri \(\{\lambda_1, \ldots, \lambda_N\}\). Nilai koefisien Moran untuk hasil vektor ciri ke-\(l\) adalah:

\[ MC[e_l] = \frac{n}{(1' W1)} \left( e_l' MWMye_l \right) / \left( e_l' Me_l \right) = \frac{n}{(1' W1)} \left( e_l' E^* \Lambda^* E^{*'} y e_l \right) / \left( e_l' Me_l \right) = \frac{n}{(1' W1)} \lambda_l \tag{2}\]

Fixed Effect Eigenvector Spatial Filtering (ESF)

Kombinasi linier dari akar ciri masih dapat diinterpretasikan dalam istilah MC. Misalkan \(E = [e_1, \ldots, e_L]\) adalah matriks yang terdiri dari vektor ciri \(L (< n)\) di \(E^*\), dan \(\gamma = [\gamma_1, \ldots, \gamma_L]'\) adalah vektor parameter. Model dasar dari ESF adalah:

\[ y = X\beta + E\gamma + \varepsilon, \quad \varepsilon \sim N(0, \sigma^2 I) \tag{3}\]

dengan \(y\) adalah vektor peubah respon berukuran \(n \times 1\), \(E\) adalah matriks berukuran \(n \times L\) yang terdiri dari subset \(L\) (\(L < n\)) vektor ciri dari \(E^*\), \(\varepsilon\) adalah vektor galat berukuran \(n \times 1\), \(\beta\) adalah vektor parameter berukuran \(K \times 1\), \(\gamma = [\gamma_1, \ldots, \gamma_l, \ldots, \gamma_L]'\) adalah vektor parameter berukuran \(L \times 1\), \(\sigma^2\) adalah parameter ragam, dan \(0\) adalah vektor nol berukuran \(n \times 1\).

Pemilihan \(L\) vektor ciri pada \(E\) yang signifikan dilakukan menggunakan metode regresi bertatar. Persamaan 3 identik dengan model regresi linier standar, sehingga linier standar, pemilihan vektor ciri yang signifikan dapat dilakukan dengan metode pemilihan bertatar dengan menggunakan metode kuadrat terkecil. Penduga parameter \(\beta\) dan \(\gamma\) sebagai berikut:

\[ \begin{bmatrix} \hat{\beta} \\ \hat{\gamma} \end{bmatrix} = \left( X^T X \quad X^T E \quad E^T X \quad I \right)^{-1} \tag{4}\]

Random Effect Eigenvector Spatial Filtering (RE-ESF)

Berdasarkan Hughes & Haran (2013), model linier RE-ESF yang didasarkan pada model dengan autoregresif prior bersyarat intrinsik dirumuskan sebagai berikut:

\[ y = X\beta + E\gamma + \varepsilon, \quad \gamma \sim N(0, \sigma_{\gamma}^2 \Lambda), \quad \varepsilon \sim N(0, \sigma^2 I) \tag{5}\]

dengan \(\sigma_{\gamma}^2\) adalah parameter ragam dan \(\Lambda\) adalah matriks diagonal berukuran \(L \times L\) yang elemennya adalah \(L\) pertama akar ciri dari \(MWM\), \(\gamma_1, \ldots, \gamma_L\), sesuai dengan \(L\) vektor ciri di \(E\). Persamaan di atas identik dengan model campuran linier (linear mixed model atau LMM). Model campuran linier merupakan perluasan dari model linear yaitu dengan menambahkan efek acak. Metode ini banyak digunakan dalam rancangan percobaan untuk data yang berkorelasi. Bentuk umum model linear campuran adalah

\[ y = X\beta + Zu + \varepsilon \tag{6}\]

dengan \(X\) adalah matriks desain dari efek tetap yang teramati, \(\beta\) adalah vektor parameter pengaruh efek tetap yang tidak diketahui, \(Z\) adalah matriks desain efek acak yang teramati, \(u\) adalah vektor efek acak yang tidak diketahui, dan \(\varepsilon\) vektor galat acak yang tidak diketahui.

Hughes & Haran (2013) memperkenalkan vektor ciri ke dalam \(E\) yang akar cirinya lebih besar dari nol. Penduga linear tak bias terbaik (best linear unbiased estimator atau BLUE) dari \(\beta\) dan \(\gamma\) diperoleh sebagai berikut (Henderson, 1975),

\[ \begin{bmatrix} \hat{\beta} \\ \hat{\gamma} \end{bmatrix} = \left( X^T X \quad X^T E \quad E^T X \quad I + \frac{\sigma^2}{\sigma_{\gamma}^2} \Lambda^{-1} \right)^{-1} \begin{bmatrix} X^T y \\ E^T y \end{bmatrix} \tag{7}\]

dengan \(\Lambda^{-1}\) adalah matriks diagonal berukuran \(L \times L\) yang elemennya adalah \(1/\lambda_1, \ldots, 1/\lambda_L\). Pada Persamaan 7, diagonal ke-\(l\) (\(l = 1, \ldots, L\)) dari \((\sigma^2/(\sigma_{\gamma}^2)) \Lambda^{-1}\), yang diberikan oleh \((\sigma^2/(\sigma_{\gamma}^2)) \gamma_l^{-1}\), menyusutkan koefisien ke-\(l\) (\(\gamma_l\)) menuju nol. Penyusutan ini signifikan untuk vektor eigen yang berhubungan dengan akar ciri yang lebih kecil, \(\gamma_l\). Penyusutan seperti itu yang berguna untuk menduga model sambil mempertahankan derajat kebebasan yang konstan (Hodges, 2013), memungkinkan menduga model RE-ESF tanpa menerapkan pemilihan vektor ciri bertahap (stepwise eigenvector selection). Ragam dari pendugaan \(\beta\) dan \(\gamma\) diberikan sebagai berikut:

\[ \text{Var}[\hat{\beta}, \hat{\gamma}] = \sigma^2 \left( X^T X \quad X^T E \quad E^T X \quad I + \frac{\sigma^2}{\sigma_{\gamma}^2} \Lambda^{-1} \right)^{-1} \tag{8}\]

Penduga Parameter Model

Persamaan 7 dapat ditulis sebagai berikut:

\[ y = X\beta + EV(\theta)u + \varepsilon, \quad u \sim N(0, \sigma^2 I), \quad \varepsilon \sim N(0, \sigma^2 I) \tag{9}\]

dengan \(V(\theta)u = \gamma\) dengan \(\theta = \left\{ \left( \frac{\sigma_{\gamma}^2}{\sigma^2} \right), \alpha \right\}\), dan \(V(\theta)\) adalah matriks diagonal berukuran \(L \times L\) yang elemen ke-\((l, l)\) adalah \(\left( \frac{\sigma_{\gamma}}{\sigma} \right) \lambda_l \sqrt{\alpha}\). Berdasarkan Bates & M. (2010), Persamaan 10 identik dengan model linear campuran yang pendugaannya menggunakan metode kemungkinan maksimum (maximum likelihood atau ML) dan kemungkinan maksimum terbatas (restricted maximum likelihood atau REML). Loglikelihood pada Persamaan 11 dapat dituliskan sebagai berikut:

\[ \text{log} L(\beta, \theta, \sigma^2) = -\frac{n}{2} \text{log}(2\pi\sigma^2) - \text{log}(|V(\theta)|^2) - \frac{d(\theta)}{\sigma^2} \tag{10}\]

\[ d(\theta) = \lVert y - X\beta + EV(\theta)u \rVert^2 + \lVert u \rVert^2 \tag{11}\]

Penduga linear tak bias terbaik dari \(\beta\) dan \(u\) diperoleh dengan memaksimumkan Persamaan 12 sebagai berikut:

\[ \begin{bmatrix} \hat{\beta} \\ \hat{u} \end{bmatrix} = \left( X^T X \quad X^T EV(\theta) \quad V(\theta)E^T X \quad V(\theta)^2 + I \right)^{-1} \begin{bmatrix} X^T y \\ V(\theta)E^T y \end{bmatrix} \tag{12}\]

dengan \(V(\theta)^2\) adalah matriks diagonal berukuran \(L \times L\) yang elemen ke-\((l, l)\) adalah \(\left( \frac{\sigma_{\gamma}^2}{\sigma^2} \right) \gamma_l \alpha\). Persamaan 11 identik dengan Persamaan 7 dengan penambahan \(\alpha\). \(\theta\) dapat diduga dengan memaksimumkan Persamaan 10. Penduga kemungkinan maksimum dari parameter ragam (\(\theta\) dan \(\sigma^2\)) diketahui berbias (Bates & M., 2010). Oleh karena itu, digunakan pendugaan REML yang memberikan pendugaan parameter ragam yang tidak bias. Residual log-likehood untuk parameter ragam diformulasikan sebagai berikut:

\[ \text{log} L_R(\theta, \sigma^2) = \int \text{log} L(\beta, \theta, \sigma^2) d\beta \tag{13}\]

Dengan memaksimumkan Persamaan 14, penduga tak bias diperoleh sebagai berikut:

\[ \sigma^2 = \frac{\lVert y - X\beta + EV(\theta)u \rVert^2}{(n-K)} \quad (13) \tag{14}\]

Berdasarkan Bates & M. (2010), \(\theta\) diduga dengan memaksimumkan profil residual log-likelihood, \(\text{log} L_R(\theta)\) yang diturunkan dengan langkah (a) mensubstitusikan Persamaan 10 pada Persamaan 12 dan menurunkan fungsi \(\text{log} L_R(\theta, \sigma^2)\) dan (b) mensubstitusikan Persamaan 13 pada \(\text{log} L_R(\theta, \sigma^2)\) dan menurunkan fungsi \(\text{log} L_R(\theta)\). Hasil fungsi profil residual log-likelihood diberikan sebagai berikut:

\[ \text{log} L_R(\theta) = -\frac{1}{2} \text{log} \left| \begin{bmatrix} X^T X & X^T EV(\theta) \\ V(\theta) E^T X & V(\theta)^2 + I \end{bmatrix} \right| - \frac{(n-K)}{2} \left( 1 + \text{log} \left( \frac{2\pi d(\theta)}{(n-K)} \right) \right) \tag{15}\]

Maksimalisasi Persamaan 16 membutuhkan penghitungan determinan pada suku pertama. Beban komputasi dapat diringankan dengan mengurangi jumlah vektor eigen di \(E\). Ragam dari \(\hat{\beta}\) dan \(\hat{\gamma}\) yang diduga adalah sebagai berikut:

\[ \text{Var}[\hat{\beta}\hat{u}] = \text{Var}[\hat{\beta}V(\theta)\hat{u}] = \sigma^2 \left(X'X \quad X'E \atop E'X \quad I + \frac{\sigma^2}{\sigma_{\gamma}^2}\Lambda^{-1} \right)^{-1} \quad \text{(15)} \tag{16}\]

dengan \(V(\theta)^{-2}\) adalah matriks diagonal berukuran \(L \times L\) yang elemen ke-\(l\) adalah \(\left(\frac{\sigma^2}{\sigma_{\gamma}^2}\right) \lambda_l(\alpha)^{-1}\).

Data

Data yang digunakan pada penelitian ini adalah data empiris Stunting per Kabupaten/Kota di Indonesia tahun 2021 dari Badan Pusat Statistika yang terdiri dari 513 Kabupaten/Kota. Data tersebut terdiri dari peubah respon dan peubah penjelas. Peubah-peubah tersebut dapat dilihat pada Tabel 1:

Tabel 1: Keterangan Peubah

Kode	Peubah
\(Y\)	Stunting
\(X_1\)	Tingkat Kemiskinan
\(X_2\)	Akses Layanan Kesehatan
\(X_3\)	Gizi Ibu Hamil
\(X_4\)	Akses Sanitasi
\(X_5\)	Akses Air Bersih

Tahapan Analisis Data

Penelitian ini dilakukan dengan langkah-langkah sebagai berikut:

1. Melakukan analisis eksplorasi data pada data stunting di Indonesia tahun 2021 untuk melihat gambaran umum setiap peubah yang diamati.

2. Melakukan tahap pra-pemrosesan data dengan mengidentifikasi pencilan pada peubah respon dan menanganinya.

3. Memodelkan data stunting menggunakan metode Regresi Linier, memeriksa asumsinya, dan mengevaluasi model.

4. Memodelkan data stunting menggunakan metode Eigenvector Spatial Filtering (ESF), dengan rincian tahapan sebagai berikut:

   1. Menguji data stunting dengan indeks Moran dan Lagrange-Multiplier
   2. Menghitung bobot spasial menggunakan pendekatan matriks invers jarak dan eksponensial.
   3. Memodelkan data stunting menggunakan metode ESF efek tetap (fixed effects) untuk setiap unit spasial.
   4. Memodelkan data stunting menggunakan metode ESF efek acak (random effect) dengan asumsi peubah acak berasal dari sebaran normal.
   5. Memodelkan data stunting menggunakan metode ESF efek acak dengan pendekatan Varying Coefficient (Random Effect-ESF VC).
   6. Melakukan uji hipotesis untuk efek peubah penjelas terhadap peubah respon.
   7. Mengevaluasi model.

5. Membandingkan kinerja model Regresi Linier, Fixed Effect-ESF, Random Effect-ESF, dan Random Effect-ESF VC menggunakan nilai Koefisien Determinasi (R-Sqrd) dan Akaike Information Criterion (AIC).

6. Menguji kelayakan dari model terbaik.

Tahapan Analisis Data dengan R

Package

library(spdep) # Untuk pemodelan dependensi spasial
library(spatialreg) # Untuk pemodelan dependensi spasial
library(stats) # Untuk operasi matematis
library(car) # Untuk menghitung VIF
library(nortest) # Untuk uji kenormalan
library(sf) # Untuk menginput peta
library(ggplot2) # Untuk visualisasi data
library(tidyr) # Untuk mengolah dataframe
library(dplyr) # Untuk mengolah dataframe
library(tidyverse) # Struktur matriks
library(lmtest) # Uji model linier
library(spmoran) # ESF
library(ape)
library(fitdistrplus)
library(readxl)

Data dan Peta

#Input Data
stunting <- read_excel("data/data_stunting_indo.xlsx")

#Import peta SHP
peta <- read_sf("data/shp/BATAS KABUPATEN KOTA DESEMBER 2019 DUKCAPIL.shp")

Modifikasi Data

#Menggabungkan Data ke file SHP
#df = cbind(stunting,KAB_KOTA=shp.jabar$KAB_KOTA)
dfix = left_join(peta, stunting, by = "KAB_KOTA")

dfix = na.omit(dfix)
sf::sf_use_s2(FALSE)
#> Spherical geometry (s2) switched off
longlat = st_coordinates(st_centroid(dfix$geometry))
#> Warning in st_centroid.sfc(dfix$geometry): st_centroid does not
#> give correct centroids for longitude/latitude data

#library("writexl")
#write_xlsx(dfix,"datacek.xlsx")

Plot Spasial

#Pemetaan Kasus Stunting
ggplot(dfix) +
  geom_sf(aes(fill = Y)) +
  geom_text(aes(longlat[, 1], longlat[, 2], label = KAB_KOTA), size = 1) +
  ggtitle("Stunting di Indonesia") +
  xlab("") +
  ylab("") +
  labs(fill = "Stunting")

Bobot Spasial

#Bobot Spasial
D = as.matrix(dist(longlat, method = "euclidean"))
#inverse weight matrix
w = 1 / D
# inverse weight matrix - row-normalized
diag(w) <- 0
rtot <- rowSums(w, na.rm = T)
w_std <- w / rtot
cekrowsum = rowSums(w_std, na.rm = T)
wmoran = mat2listw(w_std, style = "W")

Uji Asumsi

#Uji Moran Indeks (Autokorelasi)
moran.test(dfix$Y, wmoran)
#> 
#>  Moran I test under randomisation
#> 
#> data:  dfix$Y  
#> weights: wmoran    
#> 
#> Moran I statistic standard deviate = 29, p-value <2e-16
#> alternative hypothesis: greater
#> sample estimates:
#> Moran I statistic       Expectation          Variance 
#>          1.84e-01         -1.97e-03          4.06e-05

# Uji Breusch-Pagan (Heteroskesdasitas)
bptest(lm(Y ~ X1 + X2 + X3 + X4 + X5, data = dfix))
#> 
#>  studentized Breusch-Pagan test
#> 
#> data:  lm(Y ~ X1 + X2 + X3 + X4 + X5, data = dfix)
#> BP = 57, df = 63, p-value = 0.7

#LM Test (Model Spsial)
model <- lm(Y ~ X1 + X2 + X3 + X4 + X5, data = dfix)
LM <- lm.LMtests(model, wmoran, test = c('LMerr', 'LMlag', 'RLMerr', 'RLMlag'))
#> 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, data = dfix)
#> test weights: listw
#>  
#>          statistic parameter p.value    
#> RSerr        126.6         1 < 2e-16 ***
#> RSlag        119.4         1 < 2e-16 ***
#> adjRSerr      27.6         1 1.5e-07 ***
#> adjRSlag      20.4         1 6.2e-06 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Fixed ESF

#ESF Fixed
Y = dfix$Y
X = matrix(c(dfix$X1, dfix$X2, dfix$X3, dfix$X4, dfix$X5), nrow(dfix), 5)
meig <- meigen(cmat = w)
#>  149 spatial eigen-pairs
esf <- esf(y = Y, x = X, meig = meig)
#> [1] 50
#>   76/149 eigenvectors are selected
esf
#> Call:
#> esf(y = Y, x = X, meig = meig)
#> 
#> ----Coefficients------------------------------
#>              Estimate       SE t_value  p_value
#> (Intercept)  3.16e+01 2.43e+00  12.976 1.08e-32
#> V1           2.78e-01 4.53e-02   6.132 1.97e-09
#> V2           2.36e-02 2.64e-02   0.896 3.71e-01
#> V3           7.36e-05 4.44e-05   1.658 9.80e-02
#> V4          -1.47e-01 2.06e-02  -7.119 4.65e-12
#> V5           3.50e-03 2.21e-02   0.158 8.74e-01
#> 
#> ----Variance parameter (residuals)---------------
#>                      Estimate
#> SD                      5.581
#> Moran.I/max(Moran.I)    0.114
#> 
#> ----Error statistics--------------------------
#>               stat
#> resid_SE     5.147
#> adjR2        0.648
#> logLik   -1511.474
#> AIC       3188.948
#> BIC       3540.241

Random Effect-ESF

#ESF Random
Y = dfix$Y
X = matrix(c(dfix$X1, dfix$X2, dfix$X3, dfix$X4, dfix$X5), nrow(dfix), 5)
meig <- meigen(cmat = w)
#>  149 spatial eigen-pairs
resf <- resf(y = Y, x = X, meig = meig)
resf
#> Call:
#> resf(y = Y, x = X, meig = meig)
#> 
#> ----Coefficients------------------------------
#>              Estimate       SE t_value  p_value
#> (Intercept)  3.43e+01 2.46e+00  13.963 0.00e+00
#> V1           3.35e-01 5.66e-02   5.914 6.85e-09
#> V2           7.94e-03 2.54e-02   0.312 7.55e-01
#> V3           5.75e-05 4.54e-05   1.267 2.06e-01
#> V4          -1.41e-01 2.07e-02  -6.787 3.82e-11
#> V5          -3.73e-02 2.23e-02  -1.675 9.47e-02
#> 
#> ----Variance parameter------------------------
#> 
#> Spatial effects (residuals):
#>                      (Intercept)
#> random SD: Spatial         2.581
#> Moran.I/max(Moran.I)       0.146
#> 
#> ----Error statistics--------------------------
#>                  stat
#> resid SE        5.086
#> adjR2(cond)     0.653
#> rlogLik     -1675.706
#> AIC          3369.412
#> BIC          3407.504
#> 
#> NULL model: lm( y ~ x )
#>    (r)loglik: -1711.8 ( AIC: 3437.6,  BIC: 3467.227 )
#> 
#> Note: AIC and BIC are based on the restricted/marginal likelihood.
#>       Use method="ml" for comparison of models with different fixed effects (x)

Random Effect ESF - VC

#ESF VC
Y = dfix$Y
X = matrix(c(dfix$X1, dfix$X2, dfix$X3, dfix$X4, dfix$X5), nrow(dfix), 5)
meig <- meigen(coords = longlat)
#>  64 spatial eigen-pairs
resf_vc <- resf_vc(y = Y,
                   x = X,
                   meig = meig,
                   x_sel = FALSE)
#> [1] "-------  Iteration 1  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3591.401"
#> [1] "-------  Iteration 2  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3455.855"
#> [1] "-------  Iteration 3  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3452.081"
#> [1] "-------  Iteration 4  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3451.932"
#> [1] "-------  Iteration 5  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3451.901"
#> [1] "-------  Iteration 6  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3451.895"
#> [1] "-------  Iteration 7  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3451.894"
#> [1] "-------  Iteration 8  -------"
#> [1] "1/6"
#> [1] "2/6"
#> [1] "3/6"
#> [1] "4/6"
#> [1] "5/6"
#> [1] "6/6"
#> [1] "BIC: 3451.894"

resf_vc
#> Call:
#> resf_vc(y = Y, x = X, x_sel = FALSE, meig = meig)
#> 
#> ----Spatially varying coefficients on x (summary)----
#> 
#> Coefficient estimates:
#>   (Intercept)         V1               V2         
#>  Min.   :34.1   Min.   :-0.126   Min.   :-0.0814  
#>  1st Qu.:34.2   1st Qu.: 0.302   1st Qu.:-0.0166  
#>  Median :34.4   Median : 0.394   Median :-0.0023  
#>  Mean   :34.4   Mean   : 0.399   Mean   :-0.0034  
#>  3rd Qu.:34.5   3rd Qu.: 0.489   3rd Qu.: 0.0102  
#>  Max.   :34.8   Max.   : 0.793   Max.   : 0.0634  
#>        V3                  V4               V5         
#>  Min.   :-4.97e-05   Min.   :-0.188   Min.   :-0.0463  
#>  1st Qu.: 5.21e-05   1st Qu.:-0.163   1st Qu.:-0.0410  
#>  Median : 8.63e-05   Median :-0.152   Median :-0.0352  
#>  Mean   : 7.74e-05   Mean   :-0.152   Mean   :-0.0336  
#>  3rd Qu.: 1.02e-04   3rd Qu.:-0.137   3rd Qu.:-0.0259  
#>  Max.   : 2.29e-04   Max.   :-0.125   Max.   :-0.0189  
#> 
#> Statistical significance (adjusted by Benjamini-Yekutieli method):
#>                         Intercept  V1  V2  V3  V4  V5
#> Not significant                 0 378 509 509   0 509
#> Significant (10% level)         0  21   0   0   0   0
#> Significant ( 5% level)         0  75   0   0   0   0
#> Significant ( 1% level)       509  35   0   0 509   0
#> 
#> ----Variance parameters----------------------------------
#> 
#> Spatial effects (coefficients on x):
#>                      (Intercept)    V1     V2       V3     V4
#> random SD: Spatial         0.845 0.281 0.0755 0.000158 0.0294
#> Moran.I/max(Moran.I)       0.783 0.232 0.0909 0.446504 0.6594
#>                          V5
#> random SD: Spatial   0.0195
#> Moran.I/max(Moran.I) 0.7200
#> 
#> ----Estimated probability distribution of y--------------
#>                 Estimates
#> skewness                0
#> excess kurtosis         0
#> 
#> ----Error statistics-------------------------------------
#>                  stat
#> resid SE        5.224
#> adjR2(cond)     0.626
#> rlogLik     -1666.739
#> AIC          3371.478
#> BIC          3451.894
#> 
#> NULL model: lm( y ~ x )
#>    (r)loglik: Inf ( AIC: -Inf,  BIC: -Inf )
#> 
#> Note: AIC and BIC are based on the restricted/marginal likelihood.
#>       Use method="ml" for comparison of models with different fixed effects (x and xconst)

head(resf_vc$b_vc)
#>   (Intercept)    V1       V2       V3     V4      V5
#> 1        34.5 0.553  0.00757 1.14e-04 -0.131 -0.0297
#> 2        34.5 0.705  0.02144 1.18e-04 -0.129 -0.0313
#> 3        34.5 0.365 -0.01972 9.61e-05 -0.136 -0.0297
#> 4        34.5 0.439 -0.00888 1.03e-04 -0.134 -0.0296
#> 5        34.5 0.714  0.01582 1.07e-04 -0.128 -0.0338
#> 6        34.5 0.512  0.01101 9.83e-05 -0.126 -0.0356

Anselin, L. (1988). Spatial econometrics: Methods and models. Kluwer Academic Publisher.

Bates, & M., D. (2010). lme4: Mixed-effects modeling with r. http://lme4.r-forge.r-project.org/book

Chun, Y., & Griffith, D. A. (2014). A quality assessment of eigenvector spatial filtering based parameter estimates for the normal probability model. Spatial Statistics, 10, 1–11.

Dray, S., Legendre, P., & Peres-Neto, P. R. (2006). Spatial modelling: A comprehensive framework for principal coordinate analysis of neighbour matrices (PCNM). Ecological Modelling, 196(3-4), 483–493.

Henderson, C. R. (1975). Best linear unbiased estimation and prediction under a selection model. Biometrics, 423–447.

Hodges, J. S. (2013). Richly parameterized linear models: Additive, time series, and spatial models using random effects. CRC Press.

Hughes, J., & Haran, M. (2013). Dimension reduction and alleviation of confounding for spatial generalized linear mixed models. Journal of the Royal Statistical Society Series B: Statistical Methodology, 75(1), 139–159.

Murakami, D., & Griffith, D. A. (2015). Random effects specifications in eigenvector spatial filtering: A simulation study. Journal of Geographical Systems, 17, 311–331.