Pertama, install dan load terlebih dahulu package yang digunakan, yaitu sf, tigris, ggplot2, dan dplyr.

library(sf)
## Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1
library(ggplot2)
library(tigris)
## To enable 
## caching of data, set `options(tigris_use_cache = TRUE)` in your R script or .Rprofile.
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

Data yang akan digunakan adalah data yang berisi SHP Indonesia per kecamatan. Masukkan data SHP tersebut dengan melakukan akses pada file SHP.

Admin3Kecamatan <- "C:/Users/Asus/Downloads/Admin3Kecamatan"
glimpse(Admin3Kecamatan)
##  chr "C:/Users/Asus/Downloads/Admin3Kecamatan"

Untuk membaca file tersebut gunakan fungsi st_read, lalu dapat dilihat bahwa data terdiri dari 17 kolom dan 7.069 baris.

Admin3<-st_read(Admin3Kecamatan)
## Reading layer `idn_admbnda_adm3_bps_20200401' from data source 
##   `C:\Users\Asus\Downloads\Admin3Kecamatan' using driver `ESRI Shapefile'
## Simple feature collection with 7069 features and 16 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 95.01079 ymin: -11.00762 xmax: 141.0194 ymax: 6.07693
## Geodetic CRS:  WGS 84
glimpse(Admin3)
## Rows: 7,069
## Columns: 17
## $ Shape_Leng <dbl> 0.2798656, 0.7514001, 0.6900061, 0.6483629, 0.2437073, 1.35~
## $ Shape_Area <dbl> 0.003107633, 0.016925540, 0.024636382, 0.010761277, 0.00116~
## $ ADM3_EN    <chr> "2 X 11 Enam Lingkung", "2 X 11 Kayu Tanam", "Abab", "Abang~
## $ ADM3_PCODE <chr> "ID1306050", "ID1306052", "ID1612030", "ID5107050", "ID7471~
## $ ADM3_REF   <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,~
## $ ADM3ALT1EN <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,~
## $ ADM3ALT2EN <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,~
## $ ADM2_EN    <chr> "Padang Pariaman", "Padang Pariaman", "Penukal Abab Lematan~
## $ ADM2_PCODE <chr> "ID1306", "ID1306", "ID1612", "ID5107", "ID7471", "ID9432",~
## $ ADM1_EN    <chr> "Sumatera Barat", "Sumatera Barat", "Sumatera Selatan", "Ba~
## $ ADM1_PCODE <chr> "ID13", "ID13", "ID16", "ID51", "ID74", "ID94", "ID94", "ID~
## $ ADM0_EN    <chr> "Indonesia", "Indonesia", "Indonesia", "Indonesia", "Indone~
## $ ADM0_PCODE <chr> "ID", "ID", "ID", "ID", "ID", "ID", "ID", "ID", "ID", "ID",~
## $ date       <date> 2019-12-20, 2019-12-20, 2019-12-20, 2019-12-20, 2019-12-20~
## $ validOn    <date> 2020-04-01, 2020-04-01, 2020-04-01, 2020-04-01, 2020-04-01~
## $ validTo    <date> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA~
## $ geometry   <MULTIPOLYGON [°]> MULTIPOLYGON (((100.2811 -0..., MULTIPOLYGON (~

Data kuantitatif yang akan digunakan adalah jumlah penduduk kabupaten Karawang per kecamatan pada tahun 2020. Masukkan data dengan kolom yang sama seperti data dbf pada Admin3Kecamatan ditambah dengan kolom data penduduk dalam bentuk csv.

data_penduduk_karawang <- read.csv("/Users/Asus/data_penduduk_karawang.csv",header=TRUE, sep=";")

Berikut adalah dimensi dari data, yaitu memiliki 30 baris dan 5 kolom.

dim(data_penduduk_karawang)
## [1] 30  5

Melihat sekilas isi dari data_penduduk_karawang

head(data_penduduk_karawang)
##    Shape_Leng  Shape_Area        ADM3_EN ADM3_PCODE   DATA
## 1 0,393308277  0,00455583      Banyusari  ID3215071  56833
## 2 0,456589633  0,00616456       Batujaya  ID3215190  83944
## 3 0,630162331 0,009592045        Ciampel  ID3215020  43840
## 4 0,627465558    0,009628        Cibuaya  ID3215170  54211
## 5 0,282410175 0,003000352       Cikampek  ID3215050 119230
## 6 0,419235428 0,005525726 Cilamaya Kulon  ID3215082  66597
glimpse(data_penduduk_karawang)
## Rows: 30
## Columns: 5
## $ Shape_Leng <chr> "0,393308277", "0,456589633", "0,630162331", "0,627465558",~
## $ Shape_Area <chr> "0,00455583", "0,00616456", "0,009592045", "0,009628", "0,0~
## $ ADM3_EN    <chr> "Banyusari", "Batujaya", "Ciampel", "Cibuaya", "Cikampek", ~
## $ ADM3_PCODE <chr> "ID3215071", "ID3215190", "ID3215020", "ID3215170", "ID3215~
## $ DATA       <int> 56833, 83944, 43840, 54211, 119230, 66597, 83904, 43914, 80~

Merge data dilakukan dengan membandingkan primary key antara spatial data(Admin3) dan data frame(data_penduduk_karawang)

merged_Karawang <- geo_join(spatial_data=Admin3, 
                             data_frame=data_penduduk_karawang, by_sp="ADM3_PCODE", 
                             by_df="ADM3_PCODE", how = "inner")
## Warning: We recommend using the dplyr::*_join() family of functions instead.

Untuk menunjukkan perbandingan banyaknya penduduk tiap kecamatan ditandai dengan gradasi warna, atur warna yang ingin digunakan dengan syntax berikut:

mycol <- c("green", "yellow", "brown1", "brown3","red3", "red4")

Warna urutan pertama akan menjadi penanda kecamatan dengan jumlah penduduk yang relatif kecil, kemudian dilanjut dengan warna lain di urutan setelahnya bersamaan dengan meningkatnya jumlah penduduk dalam suatu kecamatan. Sehingga kecamatan yang berwarna “red4” merupakan daerah yang memiliki jumlah penduduk sangat tinggi di Karawang.

Berikut ini adalah syntax untuk menampilkan Plot Peta Choropleth:

pDATA<-ggplot()+
  geom_sf(data=merged_Karawang,aes(fill=DATA))+
  scale_fill_gradientn(colours=mycol)+
  labs(title="Jumlah Penduduk Kabupaten Karawang")
pDATA

Kemudian akan dibuat plot yang menunjukkan daerah kecamatan saya yaitu Kotabaru.

data_penduduk_karawang_kotabaru <- read.csv("/Users/Asus/data_penduduk_karawang_kotabaru.csv",header=TRUE, sep=";")
dim(data_penduduk_karawang_kotabaru)
## [1] 30  5
head(data_penduduk_karawang_kotabaru)
##    Shape_Leng  Shape_Area        ADM3_EN ADM3_PCODE DATA
## 1 0,393308277  0,00455583      Banyusari  ID3215071   NA
## 2 0,456589633  0,00616456       Batujaya  ID3215190   NA
## 3 0,630162331 0,009592045        Ciampel  ID3215020   NA
## 4 0,627465558    0,009628        Cibuaya  ID3215170   NA
## 5 0,282410175 0,003000352       Cikampek  ID3215050   NA
## 6 0,419235428 0,005525726 Cilamaya Kulon  ID3215082   NA
glimpse(data_penduduk_karawang_kotabaru)
## Rows: 30
## Columns: 5
## $ Shape_Leng <chr> "0,393308277", "0,456589633", "0,630162331", "0,627465558",~
## $ Shape_Area <chr> "0,00455583", "0,00616456", "0,009592045", "0,009628", "0,0~
## $ ADM3_EN    <chr> "Banyusari", "Batujaya", "Ciampel", "Cibuaya", "Cikampek", ~
## $ ADM3_PCODE <chr> "ID3215071", "ID3215190", "ID3215020", "ID3215170", "ID3215~
## $ DATA       <int> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 133367,~
merged_Karawang_kotabaru <- geo_join(spatial_data=Admin3, 
                             data_frame=data_penduduk_karawang_kotabaru, by_sp="ADM3_PCODE", 
                             by_df="ADM3_PCODE", how = "inner")
## Warning: We recommend using the dplyr::*_join() family of functions instead.
pDATA<-ggplot()+
  geom_sf(data=merged_Karawang_kotabaru,aes(fill=DATA))+
  scale_fill_gradientn(colours=mycol)+
  labs(title="Jumlah Penduduk Kecamatan Kotabaru")
pDATA

Kemudian akan dibuatkan plot jumlah penduduk Kabupaten Karawang dan di plot pada daerah Jawa Barat

data_penduduk_karawang_jabar <- read.csv("/Users/Asus/data_penduduk_karawang_jabar.csv",header=TRUE, sep=";")

Terdapat 630 kecamatan di Jawa Barat.

dim(data_penduduk_karawang_jabar)
## [1] 630   5
head(data_penduduk_karawang_jabar)
##      Shape_Leng    Shape_Area     ADM3_EN ADM3_PCODE DATA
## 1 0,34759297543 0,00521506222    Arjasari  ID3204150   NA
## 2 0,43579048121 0,00327566089   Baleendah  ID3204140   NA
## 3 0,35030473799 0,00318069469    Banjaran  ID3204160   NA
## 4 0,27441385318 0,00229243718 Bojongsoang  ID3204280   NA
## 5 0,27409098899 0,00189792016   Cangkuang  ID3204161   NA
## 6 0,43809655336 0,00355622343  Cicalengka  ID3204100   NA
glimpse(data_penduduk_karawang_jabar)
## Rows: 630
## Columns: 5
## $ Shape_Leng <chr> "0,34759297543", "0,43579048121", "0,35030473799", "0,27441~
## $ Shape_Area <chr> "0,00521506222", "0,00327566089", "0,00318069469", "0,00229~
## $ ADM3_EN    <chr> "Arjasari", "Baleendah", "Banjaran", "Bojongsoang", "Cangku~
## $ ADM3_PCODE <chr> "ID3204150", "ID3204140", "ID3204160", "ID3204280", "ID3204~
## $ DATA       <int> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,~
merged_Karawang_jabar <- geo_join(spatial_data=Admin3, 
                             data_frame=data_penduduk_karawang_jabar, by_sp="ADM3_PCODE", 
                             by_df="ADM3_PCODE", how = "inner")
## Warning: We recommend using the dplyr::*_join() family of functions instead.
pDATA<-ggplot()+
  geom_sf(data=merged_Karawang_jabar,aes(fill=DATA))+
  scale_fill_gradientn(colours=mycol)+
  labs(title="Jumlah Penduduk Kabupaten Karawang")
pDATA