Spatial Data
Pendahuluan
Jenis Data Spasial
Data titik (point patterns): Lokasi suatu kejadian
Data kontinu (Continus data): kejadian-kejadian yang bersifat kontinu dalam suatu ruang
Data Area : kejadian dipisahkkan dalam zona
Tipe Data Spasial Data spasial terdiri dari data vector dan data raster. Data vektor terdiri dari point, line, polygon.
PACKAGES
Packages yang digunakan adalah sp, rgdal, spData, readxl, raster, dan spdep
library(sp)
library(rgdal) # membaca shape file## rgdal: version: 1.5-23, (SVN revision 1121)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.2.1, released 2020/12/29
## Path to GDAL shared files: C:/Users/hendr/OneDrive/Documents/R/win-library/4.0/rgdal/gdal
## GDAL binary built with GEOS: TRUE
## Loaded PROJ runtime: Rel. 7.2.1, January 1st, 2021, [PJ_VERSION: 721]
## Path to PROJ shared files: C:/Users/hendr/OneDrive/Documents/R/win-library/4.0/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.4-5
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading rgdal.
## Overwritten PROJ_LIB was C:/Users/hendr/OneDrive/Documents/R/win-library/4.0/rgdal/projlibrary(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.8.0, GDAL 3.0.4, PROJ 6.3.1DATA
Dataset package R
Untuk melakukan uji coba dapat memanfaatkan dataset yang tersedia di dalam package R
ls() # melihat dataset ## character(0)Menggunakan data columbus di package spData
data(columbus)
ls() # memanggil list untuk melihat berhasil atau tiDak## [1] "bbs" "col.gal.nb" "columbus" "coords" "polys"Untuk melihat struktur data dapat menggunakan fungsi str
str(columbus)## 'data.frame': 49 obs. of 22 variables:
## $ AREA : num 0.3094 0.2593 0.1925 0.0838 0.4889 ...
## $ PERIMETER : num 2.44 2.24 2.19 1.43 3 ...
## $ COLUMBUS. : int 2 3 4 5 6 7 8 9 10 11 ...
## $ COLUMBUS.I: int 5 1 6 2 7 8 4 3 18 10 ...
## $ POLYID : int 1 2 3 4 5 6 7 8 9 10 ...
## $ NEIG : int 5 1 6 2 7 8 4 3 18 10 ...
## $ HOVAL : num 80.5 44.6 26.4 33.2 23.2 ...
## $ INC : num 19.53 21.23 15.96 4.48 11.25 ...
## $ CRIME : num 15.7 18.8 30.6 32.4 50.7 ...
## $ OPEN : num 2.851 5.297 4.535 0.394 0.406 ...
## $ PLUMB : num 0.217 0.321 0.374 1.187 0.625 ...
## $ DISCBD : num 5.03 4.27 3.89 3.7 2.83 3.78 2.74 2.89 3.17 4.33 ...
## $ X : num 38.8 35.6 39.8 36.5 40 ...
## $ Y : num 44.1 42.4 41.2 40.5 38 ...
## $ AREA : num 10.39 8.62 6.98 2.91 16.83 ...
## $ NSA : num 1 1 1 1 1 1 1 1 1 1 ...
## $ NSB : num 1 1 1 1 1 1 1 1 1 1 ...
## $ EW : num 1 0 1 0 1 1 0 0 1 1 ...
## $ CP : num 0 0 0 0 0 0 0 0 0 0 ...
## $ THOUS : num 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 ...
## $ NEIGNO : num 1005 1001 1006 1002 1007 ...
## $ PERIM : num 2.44 2.24 2.19 1.43 3 ...Data tersebut merupakan class data.frame terdiri dari 49 observasi dengan 22 peubah.
Diberikan penjelasan untuk masing-masing variabel berupa tipe datanya. Untuk angka angka yang muncul adalah cuplikan data dari masing-masing peubah.
Untuk melihat deskripsi dari dataset dapat menggunakan fungsi help
help(columbus)## starting httpd help server ... donePenjelasan dataset
Columbus OH spatial analysis data set. The columbus data frame has 49 rows and 22 columns. Unit of analysis: 49 neighbourhoods in Columbus, OH, 1980 data. In addition the data set includes a polylist object polys with the boundaries of the neighbourhoods, a matrix of polygon centroids coords, and col.gal.nb, the neighbours list from an original GAL-format file. The matrix bbs is DEPRECATED, but retained for other packages using this data set.
Data Eksternal
Selain menggunakan data yang berada di dalam package R, kita juga dapat memasukkan data dari luar package. Pada contoh ini akan digunakan data gizi buruk berupa file excel.
setwd("F:/RProject/Spatial with R")
data.giziburuk<-read_excel("data giziburuk jabar.xlsx", sheet="Sheet1")Untuk melihat dimensi dari data dapat menggunakan fungsi dim
dim(data.giziburuk)## [1] 27 5Terlihat bahwa dimensi dari data tersebut adalah 27x5
Untuk melihat struktur data dapat menggunakan fungsi str. Terlihat bahwa data gizi buruk ini terdiri dari 27 observasi dengan 5 peubah.
str(data.giziburuk)## tibble [27 x 5] (S3: tbl_df/tbl/data.frame)
## $ kode kabko: num [1:27] 3201 3202 3203 3204 3205 ...
## $ KABKOT : chr [1:27] "BOGOR" "SUKABUMI" "CIANJUR" "BANDUNG" ...
## $ Giziburuk : num [1:27] 211 117 133 186 132 110 104 111 159 129 ...
## $ Long : num [1:27] 107 107 107 108 108 ...
## $ Lat : num [1:27] -6.54 -7.07 -7.05 -7.07 -7.34 ...Untuk melihat data gizi buruk tersebut dalam bentuk tabel dapat menggunakan fungsi view. Bentuk tabel tersebut menyerupai tabel pada excel
View(data.giziburuk)PLOT
Terdapat beberapa alternatif dalam R untuk membuat peta, diantaranya menggunakan base plot atau menggunakan levelplot yang dapat diakses dengan fungsi spplot. Untuk membuat peta diperlukan shape file.
Shape file merupakan digital vector dari peta tersebut, berisi beberapa file yaitu .shp, .shx, dan .dbf. Untuk membaca peta pada R menggunakan fungsi readOGR yang berada di library rgdal.
Menampilkan peta dalam bentuk dasar berupa garis-garis wilayah
#Plot peta columbus
columbusmap <- readOGR(system.file("shapes/columbus.shp", package="spData"))## OGR data source with driver: ESRI Shapefile
## Source: "C:\Users\hendr\OneDrive\Documents\R\win-library\4.0\spData\shapes\columbus.shp", layer: "columbus"
## with 49 features
## It has 20 fields
## Integer64 fields read as strings: COLUMBUS_ COLUMBUS_I POLYIDplot(columbusmap)
Untuk mewarnai peta dengan data yang ada dalam data columbus tersebut dapat menggunakan fungsi spplot. Peta yang hanya terdiri dari batas wilayah saja bisa memiliki arti ketika diberikan perbedaan warna sesuai dengan data yang berada di dalam peta tersebut.
Sebagai contoh, dibuat peta yang menampilkan informasi harga perumahan. Gradasi warna akan menunjukkan harga harga perumahan pada data tersebut.
# HOVAL adalah salah satu peubah di dalam peta
spplot(columbusmap,"HOVAL",sub="Map of HOVAL") 
Impor Data SHP dari External
Selanjutnya kita juga dapat menggunakan file shp yang berada di luar package. Langkah pertama, kita lakukan setwd terlebih dahulu (jika belum dilakukan), lalu menuliskan dsn-nya. dsn adalah nama folder dari peta tersebut. layer adalah nama shp dari peta yang akan ditampilkan
#Plot data shp peta Jabar
jabar<-readOGR(dsn="petajabar27", layer="Peta Jabar 27")## OGR data source with driver: ESRI Shapefile
## Source: "F:\RProject\Spatial with R\petajabar27", layer: "Peta Jabar 27"
## with 27 features
## It has 2 fieldsMenambahkan text pada Peta
Untuk menambahkan text di dalam peta, dapat menggunakan fungsi text. Untuk mengatur ukuran dari text yang ditampilkan dapat menggunakan fungsi cex.
plot(jabar)
text(jabar,'KABKOT',cex=0.4) #menambahkan nama wilayah pada peta
Memberi warna pada peta
Untuk memberikan warna dapat menggunakan fungsi pallete. Fungsi ini dapat membatasi warna apa yang akan dipakai.
palette(rainbow(6))
#memberi warna pada peta
plot(jabar,col=jabar$IDKAB2) 
Mengisi Peta dengan Data
Untuk menampilan peta yang berisi data mengenai peta tersebut maka dapat menggunakan fungsi sebagai berikut:
#Memberi warna peta berdasarkan data
colfunc<-colorRampPalette(c("green", "yellow","red"))
jabar$gizi<-data.giziburuk$Giziburuk
spplot(jabar, "gizi", col.regions=colfunc(6), cuts = 5, main="Peta Sebaran")
SPATIAL WEIGHTS
Pembobot spasial terdiri dari pembobot area/poligon dan titik.
Pembobot untuk area/poligon adalah ketetanggaan (contiguity) terdiri dari beberapa cara bishop, rook, dan queen.
Pembobot untuk titik adalah jarak (distance), terdiri dari beberapa cara seperti K-nn, Inverse, dan Eksponential.
Package yang digunakan untuk membuat penimbang spasial adalah spdep.
Penimbang berdasarkan contiguity
Queen’s Case
columbus.map <- st_read(system.file("shapes/columbus.shp", package="spData"), quiet=TRUE)
queen.w <- poly2nb(as(columbus.map, "Spatial"),queen = T)
summary(queen.w)## Neighbour list object:
## Number of regions: 49
## Number of nonzero links: 236
## Percentage nonzero weights: 9.829238
## Average number of links: 4.816327
## Link number distribution:
##
## 2 3 4 5 6 7 8 9 10
## 5 9 12 5 9 3 4 1 1
## 5 least connected regions:
## 1 6 42 46 47 with 2 links
## 1 most connected region:
## 20 with 10 linksQueen weight - matriks
#Untuk melihat bentuk matriks bobot
queen.w1<-nb2mat(queen.w,style = "B")
View(queen.w1)Plot the links between the polygons
plot(st_geometry(columbus.map), border="white",col='gray')
coords <- st_coordinates(st_centroid(st_geometry(columbus.map)))
plot(queen.w, coords, add = TRUE, col = "red")
Rook’s Case
plot(st_geometry(columbus.map), border="white",col='gray')
rook.w <- poly2nb(as(columbus.map, "Spatial"), queen=FALSE)
coords <- st_coordinates(st_centroid(st_geometry(columbus.map)))
plot(rook.w, coords, add = TRUE, col = "red")
Plot the links between the polygons
Plot Queen dan Rook WMatrices
plot(st_geometry(columbus.map), border="white",col='grey')
plot(queen.w, coords, add = TRUE, col = "red")
plot(rook.w, coords, add = TRUE, col = "blue")
Penimbang berdasarkan Distance
Tahap pertama, akan mengekstrak koordinat, kemudian dibuat matriks berisi penimbang spasial.
#Ekstrak koordinat
#coordinates(columbus) <-~X+Y #x=long; y=lat
D <- as.matrix(dist(coordinates(columbus), method="euclidean"))
head(D)## 1005 1001 1006 1002 1007 1008 1004 1003
## 1005 0.00000 37.07994 56.56666 53.56437 68.59988 54.10403 24.17551 51.30504
## 1001 37.07994 0.00000 24.49558 25.86193 42.83320 24.66782 40.33898 25.55160
## 1006 56.56666 24.49558 0.00000 16.72850 23.95452 9.76021 59.63632 17.32145
## 1002 53.56437 25.86193 16.72850 0.00000 27.90103 20.07869 55.01842 12.78083
## 1007 68.59988 42.83320 23.95452 27.90103 0.00000 27.32557 72.99820 23.37633
## 1008 54.10403 24.66782 9.76021 20.07869 27.32557 0.00000 55.55063 20.01434
## 1018 1010 1038 1037 1039 1040 1009 1036
## 1005 42.02125 31.75656 98.19590 94.44389 80.82717 86.14240 75.32251 94.86906
## 1001 36.97166 59.52493 84.27394 79.81817 75.25821 81.17488 48.88561 79.79545
## 1006 36.60188 71.85081 66.40066 62.29815 64.33849 69.55758 28.71374 62.15548
## 1002 40.27883 67.13648 71.39838 68.10726 67.99644 72.87520 30.96508 67.38333
## 1007 41.51472 76.64141 58.10954 55.83231 62.21652 64.26446 21.17025 55.11953
## 1008 31.22484 68.67832 65.80145 61.47626 61.59284 67.49313 30.04163 60.61093
## 1011 1042 1041 1017 1043 1019 1012 1035
## 1005 53.91820 82.48801 93.80881 42.01996 99.11878 67.19808 49.74838 74.70231
## 1001 38.73196 84.73298 84.23977 58.95921 86.78280 49.34614 39.69784 73.70668
## 1006 30.58184 79.68763 70.13787 72.82125 72.82245 36.75166 41.12922 68.19724
## 1002 29.78647 83.27702 74.87221 74.73122 79.84725 41.36249 43.81283 70.72295
## 1007 33.87101 80.22679 64.89317 84.71421 68.22090 38.80258 50.89159 68.11019
## 1008 25.96574 77.31596 67.78798 66.37720 69.13314 31.85043 33.67985 63.67520
## 1032 1020 1021 1031 1033 1034 1045 1013
## 1005 98.33540 80.66349 78.80169 94.21592 93.01180 104.54333 96.68331 63.82410
## 1001 82.55910 60.50248 65.00193 79.84132 82.70211 91.28023 88.96477 51.32032
## 1006 64.84498 44.18285 52.36817 63.54707 68.63469 75.39718 81.30780 49.07368
## 1002 69.14075 47.73862 55.74638 67.29627 71.26517 79.31005 87.06008 52.02196
## 1007 54.95650 41.38094 47.17642 55.12256 60.48024 65.28459 85.33352 56.20647
## 1008 63.15269 40.91610 50.88765 60.95906 66.02982 73.84909 76.26345 41.84404
## 1022 1044 1023 1046 1030 1024 1047 1016
## 1005 85.06426 99.35983 83.73819 99.38709 82.67482 94.95609 101.99862 65.00571
## 1001 68.16005 90.32309 68.39659 93.64101 75.89641 79.28130 97.62375 65.75159
## 1006 54.10367 81.03969 56.24367 87.68411 68.09134 64.87389 91.80842 70.45546
## 1002 57.12011 86.70628 60.10677 93.17669 71.29517 67.39760 96.43631 72.72182
## 1007 51.21922 82.64499 54.21017 92.09241 65.90887 58.19684 94.49518 77.50555
## 1008 49.96878 75.99260 51.66738 82.37208 63.45836 61.84989 86.39894 63.88464
## 1014 1049 1029 1025 1028 1048 1015 1027
## 1005 73.45247 104.52995 95.40871 86.84862 93.81061 101.1269 81.46998 96.66114
## 1001 63.70802 102.20867 82.36111 75.88744 81.49586 109.2405 72.72558 84.11485
## 1006 61.98752 97.74379 71.21834 68.82206 72.44941 110.2554 69.25111 75.29042
## 1002 64.29584 103.14758 74.71309 72.26271 76.10283 112.2191 70.68618 78.33910
## 1007 67.97006 100.92889 69.88803 70.75961 72.41371 114.2364 71.03826 76.61968
## 1008 55.11663 91.92015 66.26125 62.71664 66.82763 104.7452 62.96282 69.89376
## 1026
## 1005 91.53183
## 1001 82.31526
## 1006 76.72837
## 1002 79.92787
## 1007 79.00506
## 1008 70.37637Invers Distance Weights
#inverse weight matrix
w=1/D
# inverse weight matrix - row-normalized
diag(w)<-0
rtot<-rowSums(w, na.rm =T)
w_std<-w/rtot
rowSums(w_std, na.rm=T)## 1005 1001 1006 1002 1007 1008 1004 1003 1018 1010 1038 1037 1039 1040 1009 1036
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 1011 1042 1041 1017 1043 1019 1012 1035 1032 1020 1021 1031 1033 1034 1045 1013
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 1022 1044 1023 1046 1030 1024 1047 1016 1014 1049 1029 1025 1028 1048 1015 1027
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 1026
## 1invers.w=mat2listw(w_std, style="W") #untuk melihat matriks W
summary(invers.w)## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 49
## Number of nonzero links: 2352
## Percentage nonzero weights: 97.95918
## Average number of links: 48
## Link number distribution:
##
## 48
## 49
## 49 least connected regions:
## 1005 1001 1006 1002 1007 1008 1004 1003 1018 1010 1038 1037 1039 1040 1009 1036 1011 1042 1041 1017 1043 1019 1012 1035 1032 1020 1021 1031 1033 1034 1045 1013 1022 1044 1023 1046 1030 1024 1047 1016 1014 1049 1029 1025 1028 1048 1015 1027 1026 with 48 links
## 49 most connected regions:
## 1005 1001 1006 1002 1007 1008 1004 1003 1018 1010 1038 1037 1039 1040 1009 1036 1011 1042 1041 1017 1043 1019 1012 1035 1032 1020 1021 1031 1033 1034 1045 1013 1022 1044 1023 1046 1030 1024 1047 1016 1014 1049 1029 1025 1028 1048 1015 1027 1026 with 48 links
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 49 2401 49 2.716383 197.1685#inverse weight matrix - plot
plot(st_geometry(columbus.map), border="white", col='gray' , sub="Inverse Distance")
plot(invers.w,coords,add=TRUE,col="red")
Eksponential Distance Weights
#Eksponential weight matrix
alpha=2
w2=exp(-alpha*D)
#Eksponential weight matrix - row-normalized
diag(w2)<-0
rtot<-rowSums(w2, na.rm =T)
rtot## 1005 1001 1006 1002 1007 1008
## 1.003269e-21 1.002146e-21 3.329584e-09 7.922652e-12 4.155522e-19 3.329581e-09
## 1004 1003 1018 1010 1038 1037
## 1.003269e-21 7.920604e-12 1.429511e-23 2.609779e-28 3.809716e-07 4.627262e-07
## 1039 1040 1009 1036 1011 1042
## 3.934315e-10 3.935495e-10 4.090544e-19 8.186012e-08 4.713461e-19 2.736346e-27
## 1041 1017 1043 1019 1012 1035
## 1.329262e-13 3.024577e-34 5.482211e-25 5.901717e-14 3.564140e-14 4.797072e-24
## 1032 1020 1021 1031 1033 1034
## 7.133463e-09 2.160092e-11 3.356299e-22 7.216211e-09 8.276463e-11 4.810802e-15
## 1045 1013 1022 1044 1023 1046
## 2.951847e-09 3.837424e-13 2.197172e-08 1.485294e-09 2.195000e-08 1.907799e-09
## 1030 1024 1047 1016 1014 1049
## 1.341360e-19 1.778688e-13 4.681446e-10 7.055729e-21 1.998307e-12 2.713481e-11
## 1029 1025 1028 1048 1015 1027
## 1.308403e-10 1.478249e-08 5.026807e-09 7.520383e-31 1.650206e-12 4.081598e-09
## 1026
## 1.392096e-08w_std<-w2/rtot
rowSums(w_std, na.rm=T)## 1005 1001 1006 1002 1007 1008 1004 1003 1018 1010 1038 1037 1039 1040 1009 1036
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 1011 1042 1041 1017 1043 1019 1012 1035 1032 1020 1021 1031 1033 1034 1045 1013
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 1022 1044 1023 1046 1030 1024 1047 1016 1014 1049 1029 1025 1028 1048 1015 1027
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 1026
## 1eksp.w=mat2listw(w_std, style="W") #untuk melihat matriks W
summary(eksp.w)## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 49
## Number of nonzero links: 2352
## Percentage nonzero weights: 97.95918
## Average number of links: 48
## Link number distribution:
##
## 48
## 49
## 49 least connected regions:
## 1005 1001 1006 1002 1007 1008 1004 1003 1018 1010 1038 1037 1039 1040 1009 1036 1011 1042 1041 1017 1043 1019 1012 1035 1032 1020 1021 1031 1033 1034 1045 1013 1022 1044 1023 1046 1030 1024 1047 1016 1014 1049 1029 1025 1028 1048 1015 1027 1026 with 48 links
## 49 most connected regions:
## 1005 1001 1006 1002 1007 1008 1004 1003 1018 1010 1038 1037 1039 1040 1009 1036 1011 1042 1041 1017 1043 1019 1012 1035 1032 1020 1021 1031 1033 1034 1045 1013 1022 1044 1023 1046 1030 1024 1047 1016 1014 1049 1029 1025 1028 1048 1015 1027 1026 with 48 links
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 49 2401 49 68.55016 227.8029#Eksponential weight matrix - plot
plot(st_geometry(columbus.map), border="white", col='gray' , sub="Eksponential Distance Weights")
plot(eksp.w, coords, add = TRUE, col = "red")
k-Nearest Neighbor Weights
#K-nn weight matrix
columbus <- st_read(system.file("shapes/columbus.shp", package="spData")[1], quiet=TRUE)
coords <- st_centroid(st_geometry(columbus), of_largest_polygon=TRUE)
w3=knearneigh(coords, k=1)
knn.w=nb2listw(knn2nb(w3))
summary(knn.w)## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 49
## Number of nonzero links: 49
## Percentage nonzero weights: 2.040816
## Average number of links: 1
## Non-symmetric neighbours list
## Link number distribution:
##
## 1
## 49
## 49 least connected regions:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 with 1 link
## 49 most connected regions:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 with 1 link
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 49 2401 49 75 230#K-nn weight matrix - plot
plot(st_geometry(columbus.map), border="white", col='gray' , sub="1-nn")
plot(knn.w, coords, add = TRUE, col = "red")