Email         : vanessasupit0910@gmail.com
RPubs       : https://rpubs.com/vanessasupit/
Department : Business Statistics
Address : ARA Center, Matana University Tower
   Jl. CBD Barat Kav, RT.1, Curug Sangereng, Kelapa Dua, Tangerang, Banten 15810.

Berikut merupakan library-library yang akan digunakan dalam pemodelan spasial

library(rgdal)           
library(raster)    
library(spdep)
library(sp)
library(raster)
library(rgeos)
library(rspatial)
library(latticeExtra)
library(RColorBrewer)
library(spatialreg)
library(DT)

1 Input Data

Link untuk file bisa diambil di sini

data = readOGR(dsn = "data", layer = "NCVACO")

## OGR data source with driver: ESRI Shapefile 
## Source: "C:\Users\User_Pro\Documents\TSSpatialModel\data", layer: "NCVACO"
## with 234 features
## It has 49 fields
## Integer64 fields read as strings:  FIPS qtystores PCI COUNTMXBV DC GA KY MD SC TN WV VA COUNTBKGR TOTALPOP POP18OV LABFORCE HHOLDS POP25OV POP16OV

-Daftar variabel dataset

names(data)

##  [1] "GEO_ID"     "STATE"      "COUNTY"     "NAME"       "LSAD"      
##  [6] "CENSUSAREA" "FIPS2"      "Lon"        "Lat"        "FIPS"      
## [11] "qtystores"  "SALESPC"    "PCI"        "COMM15OVP"  "COLLENRP"  
## [16] "SOMECOLLP"  "ARMEDP"     "NONWHITEP"  "UNEMPP"     "ENTRECP"   
## [21] "PUBASSTP"   "POVPOPP"    "URBANP"     "FOREIGNBP"  "BAPTISTSP" 
## [26] "ADHERENTSP" "BKGRTOMIX"  "COUNTMXBV"  "MXBVSQM"    "BKGRTOABC" 
## [31] "MXBVPPOP18" "DUI1802"    "FVPTHH02"   "DC"         "GA"        
## [36] "KY"         "MD"         "SC"         "TN"         "WV"        
## [41] "VA"         "AREALANDSQ" "COUNTBKGR"  "TOTALPOP"   "POP18OV"   
## [46] "LABFORCE"   "HHOLDS"     "POP25OV"    "POP16OV"

Dataset ini adalah beberapa data dari beberapa penelitian yang dilakukan Mark Burkey tentang permintaan minuman keras menggunakan data dari sekitar tahun 2003. Secara khusus, dia melihat negara bagian Virginia dan Carolina Utara. Dataset ini terkait dengan, tetapi tidak sama dengan data yang digunakan untuk hibah NIH dan diterbitkan dalam makalah:

Burkey, Mark L. Geographic Access and Demand in the Market for Alcohol. The Review of Regional Studies, 40(2), Fall 2010, 159-179

Unit analisis: kabupaten di Virginia dan Carolina Utara

Deskripsi Variabel:

Lon Lat Bujur dan Lintang dari County Centroid
FIPS Kode FIPS untuk Wilayah (Standar Pemrosesan Informasi Federal)
qtystores #Toko Minuman Keras di County
SALESPC Penjualan Minuman Keras per kapita per tahun, $
PCI Pendapatan per kapita
COMM15OVP % bepergian selama 15 menit ke kantor
COLLENRP % orang yang saat ini terdaftar di perguruan tinggi
SOMECOLLP % orang yang kuliah atau berpendidikan lebih tinggi
ARMEDP % di angkatan bersenjata
NONWHITEP % bukan kulit putih
UNEMPP % pengangguran
ENTRECP % dipekerjakan di bidang hiburan atau rekreasi (perwakilan untuk bidang pariwisata)
PUBASSTP % tentang bantuan publik
POVPOPP % dalam kemiskinan
URBANP % tinggal di daerah perkotaan
FOREIGNBP % lahir di luar negeri
BAPTISTSP % baptis selatan (secara historis anti-alkohol)
ADHERENTSP % penganut agama apapun
BKGRTOMIX wtd. jarak rata-rata dari kelompok blok ke bar terdekat yang menjual minuman keras
COUNTMXBV hitungan bar yang menjual minuman keras
MXBVSQM bar per mil persegi
BKGRTOABC jarak dari grup blok ke outlet ritel minuman keras terdekat (“toko ABC”)
MXBVPPOP18OV Bar per 1.000? orang berusia 18 tahun ke atas
DUI1802 Penangkapan DUI per 1.000 orang 18+
FVPTHH02 Pelanggaran terhadap keluarga dan anak (kekerasan dalam rumah tangga) per 1.000 rumah tangga
DC GA KY MD SC TN WV VA Variabel tiruan untuk negara bagian yang berbatasan dengan negara bagian lain
AREALANDSQMI Luas wilayah
COUNTBKGR jumlah kelompok blok di daerah
TOTALPOP Penduduk daerah
POP18OV orang berusia 18+ di daerah
LABFORCE jumlah angkatan kerja di daerah
HHOLDS # rumah tangga di daerah
POP25OV Pop 25+ di daerah
POP16OV Pop 16+ di daerah

2 Pemetaan Distribusi Spasial

Sebelumnya kita akan mentransformasi struktur variabel PCI (Pendapatan per kapita) dari faktor ke numerik

data$PCI <- as.numeric(levels(data$PCI))[data$PCI]

Sekarang kita akan memvisualisasikan penjualan minuman keras per kapita

spplot(data, "SALESPC")

Dari hasil visualisasi di atas kita bisa melihat bahwa rata-rata penjualan minuman keras di tiap daerah berkisar di angka 50-150. Namun ada satu daerah dengan penjualan kapita yang sangat tinggi yang mencapai300 per kapita, yaitu area berwarna kuning.

3 Pemodelan Regresi Spasial

Di bagian sini kita akan membangun model linier untuk data kita. Model spasial secara global terdiri dari:

SEM: $y=XB+u$
SAR: $y=\rho Wy+XB+\varepsilon$
SARMA: $y=\rho Wy+XB+u$

Sebelum memulai pemodelan kita perlu membuat pemberat dengan fungsi berikut

queen.nb <- poly2nb(data)             # Menciptakan relasi
listw <- nb2listw(queen.nb)           # Mengkonversi nb ke listw

Selanjutnya kita akan membuat persamaan model linier kita

lm <-DUI1802 ~ SALESPC + COLLENRP + BKGRTOABC + BAPTISTSP + BKGRTOMIX + ENTRECP

3.1 Model SEM

reg1 <- errorsarlm(lm, data=data, listw)
summary(reg1)

## 
## Call:errorsarlm(formula = lm, data = data, listw = listw)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -7.21199 -2.66820 -0.85841  1.57794 25.85910 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##               Estimate Std. Error z value  Pr(>|z|)
## (Intercept)  5.1218025  1.4443139  3.5462 0.0003909
## SALESPC      0.0111587  0.0076284  1.4628 0.1435286
## COLLENRP     0.0533827  0.0617209  0.8649 0.3870915
## BKGRTOABC    0.1533586  0.1454500  1.0544 0.2917121
## BAPTISTSP    0.0326804  0.0394447  0.8285 0.4073797
## BKGRTOMIX   -0.0232354  0.0832068 -0.2792 0.7800541
## ENTRECP      0.1067644  0.1449288  0.7367 0.4613242
## 
## Lambda: 0.40342, LR test value: 20.075, p-value: 7.4467e-06
## Asymptotic standard error: 0.075449
##     z-value: 5.347, p-value: 8.9438e-08
## Wald statistic: 28.59, p-value: 8.9438e-08
## 
## Log likelihood: -694.1478 for error model
## ML residual variance (sigma squared): 21.247, (sigma: 4.6095)
## Number of observations: 234 
## Number of parameters estimated: 9 
## AIC: 1406.3, (AIC for lm: 1424.4)

3.2 Model SAR

reg2 <- lagsarlm(lm, data = data, listw)
summary(reg2)

## 
## Call:lagsarlm(formula = lm, data = data, listw = listw)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -7.18539 -2.59660 -0.77605  1.61369 25.93228 
## 
## Type: lag 
## Coefficients: (asymptotic standard errors) 
##              Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5421139  1.3012975  1.1851  0.23599
## SALESPC     0.0138185  0.0077171  1.7906  0.07335
## COLLENRP    0.0327744  0.0603469  0.5431  0.58706
## BKGRTOABC   0.1164671  0.1358375  0.8574  0.39122
## BAPTISTSP   0.0392792  0.0332872  1.1800  0.23800
## BKGRTOMIX   0.0092150  0.0700728  0.1315  0.89538
## ENTRECP     0.1532086  0.1359122  1.1273  0.25963
## 
## Rho: 0.39828, LR test value: 22.924, p-value: 1.6858e-06
## Asymptotic standard error: 0.07512
##     z-value: 5.3019, p-value: 1.1459e-07
## Wald statistic: 28.11, p-value: 1.1459e-07
## 
## Log likelihood: -692.7235 for lag model
## ML residual variance (sigma squared): 21.012, (sigma: 4.5839)
## Number of observations: 234 
## Number of parameters estimated: 9 
## AIC: 1403.4, (AIC for lm: 1424.4)
## LM test for residual autocorrelation
## test value: 13.047, p-value: 0.0003037

3.3 Model SARMA

reg3 <- sacsarlm(lm,data = data, listw)
summary(reg3)

## 
## Call:sacsarlm(formula = lm, data = data, listw = listw)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -7.01945 -2.34197 -0.78281  1.41007 23.48056 
## 
## Type: sac 
## Coefficients: (asymptotic standard errors) 
##               Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.8353870  0.9351383 -0.8933  0.37168
## SALESPC      0.0171848  0.0067612  2.5417  0.01103
## COLLENRP     0.0138969  0.0505268  0.2750  0.78328
## BKGRTOABC    0.1404632  0.1045255  1.3438  0.17901
## BAPTISTSP    0.0211294  0.0225887  0.9354  0.34958
## BKGRTOMIX   -0.0214055  0.0491538 -0.4355  0.66321
## ENTRECP      0.1409664  0.1061888  1.3275  0.18434
## 
## Rho: 0.74152
## Asymptotic standard error: 0.069065
##     z-value: 10.737, p-value: < 2.22e-16
## Lambda: -0.60359
## Asymptotic standard error: 0.13376
##     z-value: -4.5125, p-value: 6.4082e-06
## 
## LR test value: 32.473, p-value: 8.8832e-08
## 
## Log likelihood: -687.9487 for sac model
## ML residual variance (sigma squared): 16.534, (sigma: 4.0662)
## Number of observations: 234 
## Number of parameters estimated: 10 
## AIC: 1395.9, (AIC for lm: 1424.4)

Teknik Sampling

UAS Nomor 4

Vanessa Supit

May 27, 2021