Praktikum 3 : Eksplorasi Data Areal

Author

Rosita Ria Rusesta - Statistika dan Sains Data

Review Data Areal/Lattice

Data areal (lattice data) merupakan bentuk data spasial yang bersifat diskrit karena merepresentasikan informasi dalam batas wilayah tertentu. Data ini muncul ketika suatu fenomena atau jumlah kejadian digabungkan ke dalam area administratif, misalnya jumlah kasus penyakit dalam epidemiologi spasial atau jumlah kelahiran pada suatu kabupaten. Dengan pendekatan ini, peneliti dapat mengidentifikasi pola geografis, memahami distribusi fenomena, dan menelusuri faktor risiko yang mungkin berhubungan dengan kondisi lingkungan maupun sosial.

Dalam praktik analisis, data areal umumnya tersedia dalam bentuk objek vektor dengan format shapefile. Format ini menjadi standar paling umum untuk menyimpan data spasial karena mampu memuat informasi geometri (seperti poligon batas wilayah) sekaligus atribut nonspasial yang menyertainya. Di R, data spasial tersebut dapat dikelola menggunakan paket sf yang dirancang khusus untuk memudahkan manipulasi dan analisis data berbasis spasial. Sementara itu, data nonspasial biasanya disimpan dalam format umum seperti file CSV atau .xlsx. Dengan kombinasi kedua jenis data ini, analisis spasial dapat dilakukan secara lebih komprehensif, baik dari sisi geografis maupun dari sisi atribut yang melekat pada setiap area.

knitr::opts_chunk$set(warning = FALSE, message = FALSE)
library(sf)
library(sp)
library(mapview)
library(spData)
library(ggplot2)
library(spatstat)
library(tidyverse)
library(spData)
library(readxl)
list.files(system.file("shapes", package = "spData"))
 [1] "auckland.gpkg"          "boston_tracts.gpkg"     "columbus.gpkg"         
 [4] "cycle_hire_osm.geojson" "cycle_hire.geojson"     "eire.gpkg"             
 [7] "NY8_bna_utm18.gpkg"     "NY8_utm18.gpkg"         "sids.gpkg"             
[10] "wheat.gpkg"             "world.gpkg"            
d2 <- st_read(system.file("shapes/columbus.gpkg",
                         package = "spData"), quiet = TRUE)
st_crs(d2) <- NA
names(d2)
 [1] "AREA"       "PERIMETER"  "COLUMBUS_"  "COLUMBUS_I" "POLYID"    
 [6] "NEIG"       "HOVAL"      "INC"        "CRIME"      "OPEN"      
[11] "PLUMB"      "DISCBD"     "X"          "Y"          "NSA"       
[16] "NSB"        "EW"         "CP"         "THOUS"      "NEIGNO"    
[21] "geom"      
print(d2)
Simple feature collection with 49 features and 20 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 5.874907 ymin: 10.78863 xmax: 11.28742 ymax: 14.74245
CRS:           NA
First 10 features:
       AREA PERIMETER COLUMBUS_ COLUMBUS_I POLYID NEIG  HOVAL    INC     CRIME
1  0.309441  2.440629         2          5      1    5 80.467 19.531 15.725980
2  0.259329  2.236939         3          1      2    1 44.567 21.232 18.801754
3  0.192468  2.187547         4          6      3    6 26.350 15.956 30.626781
4  0.083841  1.427635         5          2      4    2 33.200  4.477 32.387760
5  0.488888  2.997133         6          7      5    7 23.225 11.252 50.731510
6  0.283079  2.335634         7          8      6    8 28.750 16.029 26.066658
7  0.257084  2.554577         8          4      7    4 75.000  8.438  0.178269
8  0.204954  2.139524         9          3      8    3 37.125 11.337 38.425858
9  0.500755  3.169707        10         18      9   18 52.600 17.586 30.515917
10 0.246689  2.087235        11         10     10   10 96.400 13.598 34.000835
       OPEN    PLUMB DISCBD     X     Y NSA NSB EW CP THOUS NEIGNO
1  2.850747 0.217155   5.03 38.80 44.07   1   1  1  0  1000   1005
2  5.296720 0.320581   4.27 35.62 42.38   1   1  0  0  1000   1001
3  4.534649 0.374404   3.89 39.82 41.18   1   1  1  0  1000   1006
4  0.394427 1.186944   3.70 36.50 40.52   1   1  0  0  1000   1002
5  0.405664 0.624596   2.83 40.01 38.00   1   1  1  0  1000   1007
6  0.563075 0.254130   3.78 43.75 39.28   1   1  1  0  1000   1008
7  0.000000 2.402402   2.74 33.36 38.41   1   1  0  0  1000   1004
8  3.483478 2.739726   2.89 36.71 38.71   1   1  0  0  1000   1003
9  0.527488 0.890736   3.17 43.44 35.92   1   1  1  0  1000   1018
10 1.548348 0.557724   4.33 47.61 36.42   1   1  1  0  1000   1010
                             geom
1  POLYGON ((8.624129 14.23698...
2  POLYGON ((8.25279 14.23694,...
3  POLYGON ((8.653305 14.00809...
4  POLYGON ((8.459499 13.82035...
5  POLYGON ((8.685274 13.63952...
6  POLYGON ((9.401384 13.5504,...
7  POLYGON ((8.037741 13.60752...
8  POLYGON ((8.247527 13.58651...
9  POLYGON ((9.333297 13.27242...
10 POLYGON ((10.08251 13.03377...
ggplot(d2) + geom_sf(aes(fill = CRIME)) + theme_minimal()

Membuat objek sf

Objek sf (simple feature) di R adalah struktur data yang digunakan untuk menyimpan data spasial vektor. Objek ini menggabungkan data atribut (seperti data frame) dengan geometri spasial.

  1. Data Atribut:

    Mirip dengan data frame biasa, bagian ini berisi kolom-kolom yang menyimpan informasi non-spasial (atribut) tentang setiap fitur. Misalnya, nama, ID, populasi, dll.

  2. Kolom Geometri:

    Kolom khusus yang menyimpan informasi geometri dari setiap fitur. Kolom ini biasanya dinamai geometry atau geom, tetapi bisa diberi nama lain. Geometri ini bisa berupa titik, garis, poligon, multipoint, atau multipoligon.

  3. Sistem Referensi Koordinat (CRS):

    Informasi tentang sistem referensi koordinat yang digunakan untuk geometri. CRS menentukan bagaimana koordinat dipetakan ke permukaan bumi.

  4. Atribut Metadata:

    Informasi tambahan tentang objek sf, seperti nama kolom geometri, CRS, dan tipe geometri.

Kita dapat menggunakan st_sf() fungsi tersebut untuk membuat objek sf dengan menyediakan dua elemen, yaitu data frame dengan atribut setiap fitur dan kolom daftar geometri fitur sederhana (simple feature geometri list-column/sfc) yang berisi geometri fitur sederhana (simple feature geomteries/sfg).

# Membuat data frame dengan atribut
df <- data.frame(id = 1:3,
                 name = c("A", "B", "C"))

# Membuat geometri (titik)
p1 <- st_point(c(7.35, 52.42))
p2 <- st_point(c(7.22, 52.18))
p3 <- st_point(c(7.44, 52.19))
geom <- st_sfc(list(p1, p2, p3), crs = 'EPSG:4326')

# Menggabungkan data frame dan geometri menjadi objek sf
sf_obj <- st_sf(df, geometry = geom)

# Menampilkan objek sf
print(sf_obj)
Simple feature collection with 3 features and 2 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 7.22 ymin: 52.18 xmax: 7.44 ymax: 52.42
Geodetic CRS:  WGS 84
  id name           geometry
1  1    A POINT (7.35 52.42)
2  2    B POINT (7.22 52.18)
3  3    C POINT (7.44 52.19)

  • POINT Satu titik. Contoh: lokasi sensor atau titik GPS.

  • MULTIPOINT Sekumpulan titik dalam satu fitur. Contoh: sebaran pohon dalam satu plot.

  • LINESTRING Garis tunggal yang terdiri atas segmen-segmen. Contoh: satu ruas jalan atau sungai.

  • MULTILINESTRING Beberapa garis dalam satu fitur. Contoh: jaringan jalan yang dimiliki satu entitas.

  • POLYGON Satu poligon dengan cangkang luar (boleh ada lubang/holes). Contoh: batas satu kelurahan.

  • MULTIPOLYGON Sekumpulan poligon dalam satu fitur. Contoh: kabupaten yang terdiri atas beberapa pulau terpisah.

  • GEOMETRYCOLLECTION Kumpulan berbagai tipe geometri (titik, garis, poligon) dalam satu fitur.

Implementasi Data Area : Data Jawa Barat

Data yang digunakan dalam praktik ini merupakan data sekunder berupa data kependudukan yang diperoleh dari Badan Pusat Statistik (BPS) dan peta shapefile berisi peta Provinsi Jawa Barat. Berdasarkan data tersebut, tahap awal yaitu proses menggabungkan data bentuk geometri wilayah (peta shp Jawa Barat) dengan data statistik atribut dalam format xlsx. Dataset mencakup informasi spasial dan tabular tingkat kabupaten/kota di Provinsi Jawa Barat untuk tahun 2015 dan 2016. Praktik ini melibatkan variabel utama, yaitu %Kemiskinan tahun 2015. Data dapat didownload pada link berikut https://bit.ly/data_praktikum3

setwd("/Users/rositariarusesta/Documents/Regresi Spasial")
shp_jabar <- st_read("petaJabar2/Jabar2.shp")
Reading layer `Jabar2' from data source 
  `/Users/rositariarusesta/Documents/Regresi Spasial/PetaJabar2/Jabar2.shp' 
  using driver `ESRI Shapefile'
Simple feature collection with 26 features and 7 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 106.3705 ymin: -7.823398 xmax: 108.8338 ymax: -5.91377
CRS:           NA
data_jabar <- read_excel("data_jabar.xlsx")
data_jabar$KABKOTNO <- as.character(data_jabar$KABKOTNO)

# join
data_spasial <- shp_jabar %>%
  left_join(data_jabar, by = "KABKOTNO") %>%
  na.omit()
#Cek semua nama kolom
names(data_spasial)
 [1] "PROVNO.x"   "KABKOTNO"   "KODE2010.x" "PROVINSI.x" "KABKOT.x"  
 [6] "SUMBER"     "IDSP2010.x" "PROVNO.y"   "KODE2010.y" "PROVINSI.y"
[11] "KABKOT.y"   "IDSP2010.y" "Long"       "Lat"        "p.miskin15"
[16] "p.miskin16" "j.miskin15" "j.miskin16" "AHH2015"    "AHH2016"   
[21] "EYS2015"    "EYS2016"    "MYS2015"    "MYS2016"    "EXP2015"   
[26] "EXP2016"    "APM.SD15"   "APM.SMP15"  "APM.SMA15"  "APM.PT15"  
[31] "APK.SD15"   "APK.SMP15"  "APK.SMA15"  "APK.PT15"   "APS.USIA15"
[36] "APS.USIA2"  "APS.USIA3"  "APS.USIA4"  "geometry"  
#Memilih data yang mengandung variabel saja
data_filter<-data_spasial[,15:37]
#Visualisasi data Jabar
plot(data_filter)

library(ggplot2)

plot_kemiskinan <- ggplot(data = data_spasial) +  
  geom_sf(aes(fill = p.miskin15), color = "white", linewidth = 0.2) +  
  # Menambahkan label nama kabupaten/kota:
  geom_sf_text(aes(label = KABKOT.x), size = 2, color = "black", check_overlap = TRUE) +
  scale_fill_viridis_c(option = "magma", name = "Persentase (%)") + 
  theme_minimal() +  
  labs(
    title = "%Kemiskinan Jawa Barat",       
    subtitle = "Tahun 2015",       
    x = "Longitude", 
    y = "Latitude"
  ) +  
  theme(legend.position = "bottom")

plot_kemiskinan

Explorasi Data Area

Import Dara Area

Kita dapat membaca objek sf dengan fungsi st_read(). Misalnya, di sini kita membaca shapefile ncyang berisi kabupaten-kabupaten di North Carolina, AS, beserta nama kabupaten, jumlah kelahiran, dan jumlah kematian bayi mendadak pada tahun 1974 dan 1979.

imp <- system.file("shape/nc.shp", package = "sf")
nc <- st_read(imp, quiet = TRUE)
class(nc)
[1] "sf"         "data.frame"

Selanjutnya kita dapat mengidentifikasi CRS yang digunakan pada data tersebut dengan fungsi st_crs().

st_crs(nc)
Coordinate Reference System:
  User input: NAD27 
  wkt:
GEOGCRS["NAD27",
    DATUM["North American Datum 1927",
        ELLIPSOID["Clarke 1866",6378206.4,294.978698213898,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["latitude",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["longitude",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4267]]
plot(nc)

Kita dapat membuat subset set fitur dengan menggunakan notasi tanda kurung siku dan menggunakan argumen drop untuk menghilangkan geometri.

plot(nc[9])

Menghapus Beberapa Poligon

map <- nc[-which(nc$FIPS %in% c("37125", "37051")), ]
ggplot(nc) + geom_sf(aes(fill = SID79))

Menghitung Jumlah Titik dalam Poligon

Kita dapat menggunakan st_intersects() untuk menghitung jumlah titik dalam poligon suatu objek sfdan st_sample() menghasilkan titik sampel acak. Objek yang dikembalikan adalah daftar dengan id fitur yang berpotongan di setiap poligon.

# Peta (objek sf)
map <- read_sf(system.file("shape/nc.shp", package = "sf"))

# Titik di atas peta (kolom daftar geometri fitur sederhana sfc)
points <- st_sample(map, size = 100)

# Peta titik di dalam poligon
ggplot() + geom_sf(data = map) + geom_sf(data = points)

# Interseksi (argumen pertama peta, kemudian titik)
inter <- st_intersects(map, points)

# Menambahkan jumlah titik ke setiap poligon
map$count <- lengths(inter)

# Peta jumlah titik di dalam poligon
ggplot(map) + geom_sf(aes(fill = count))