Package

library(readxl)
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.2
## 
## 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
library(corrplot)
## corrplot 0.92 loaded
library(plotly)
## Warning: package 'plotly' was built under R version 4.3.2
## Loading required package: ggplot2
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(ggplot2)
library(ggcorrplot)
## Warning: package 'ggcorrplot' was built under R version 4.3.3
library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.3.2
## Warning: package 'lubridate' was built under R version 4.3.2
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ lubridate 1.9.3     ✔ tibble    3.2.1
## ✔ purrr     1.0.2     ✔ tidyr     1.3.0
## ✔ readr     2.1.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ plotly::filter() masks dplyr::filter(), stats::filter()
## ✖ dplyr::lag()     masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(reshape2)
## 
## Attaching package: 'reshape2'
## 
## The following object is masked from 'package:tidyr':
## 
##     smiths
library(ggforce)
## Warning: package 'ggforce' was built under R version 4.3.3
library(sf)
## Warning: package 'sf' was built under R version 4.3.3
## Linking to GEOS 3.11.2, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUE
library(ggplot2)
library(tigris)
## Warning: package 'tigris' was built under R version 4.3.3
## To enable caching of data, set `options(tigris_use_cache = TRUE)`
## in your R script or .Rprofile.
library(dplyr)
library(maps)
## Warning: package 'maps' was built under R version 4.3.3
## 
## Attaching package: 'maps'
## 
## The following object is masked from 'package:purrr':
## 
##     map
library(viridis)
## Warning: package 'viridis' was built under R version 4.3.3
## Loading required package: viridisLite
## 
## Attaching package: 'viridis'
## 
## The following object is masked from 'package:maps':
## 
##     unemp
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## 
## The following object is masked from 'package:dplyr':
## 
##     combine
library(ggspatial)
## Warning: package 'ggspatial' was built under R version 4.3.3
library("rnaturalearthdata")
## Warning: package 'rnaturalearthdata' was built under R version 4.3.3
library("rnaturalearth")
## Warning: package 'rnaturalearth' was built under R version 4.3.3
## 
## Attaching package: 'rnaturalearth'
## 
## The following object is masked from 'package:rnaturalearthdata':
## 
##     countries110
library(indonesia)
devtools::install_github("rasyidstat/indonesia")
## Skipping install of 'indonesia' from a github remote, the SHA1 (5005a6dc) has not changed since last install.
##   Use `force = TRUE` to force installation

Hubungan Antar Peubah Numerik

Data (Materi 8)

datajabar <- read_excel("C:/Users/user/OneDrive/Dokumen/semester 4/visdat/p10/JabarData(gabung).xlsx")
view(datajabar)
datajabar.baru <- datajabar %>% select(j.miskin16,AHH2016,MYS2016,EXP2016)
view(datajabar.baru)

Keterangan Data

  1. Jumlah Penduduk Miskin (ribu jiwa) 2016 ditunjukkan dengan kode j.miskin16
  2. Angka Harapan Hidup tahun 2016 ditunjukkan dengan kode AHH2016
  3. Rata-rata lama sekolah tahun 2016 ditunjukkan dengan kode MYS2016
  4. Pengeluaran per kapita riil yang disesuaikan dengan formula Atkinson (Tahun 2016) ditunjukkan dengan kode EXP2016

Visualisasi Korelasi

ggplot(datajabar.baru, aes(x = MYS2016, y = j.miskin16)) +
  geom_point() +   
  geom_smooth(method = "lm", se = FALSE, col = "red", lwd = 1.5) +  
  labs(title = "Jumlah Penduduk Miskin vs Rata-Rata Lama Sekolah Tahun 2016", x = "rata-rata lama sekolah (tahun)", y = "jumlah penduduk miskin (ribu jiwa)") +
  theme_bw()
## `geom_smooth()` using formula = 'y ~ x'

Visualisasi antar peubah numerik yang ditampilkan dengan scatter plot di atas menggambarkan hubungan antara jumlah penduduk miskin dalam ribu jiwa dengan rata-rata lama sekolah dalam tahun. Berdasarkan analisis plot tersebut, hubungan antara jumlah penduduk miskin dengan rata-rata lama sekolah memiliki hubungan linear negatif karena plot cenderung memiliki arah dari kiri atas ke kanan bawah. Interpretasinya, ketika rata-rata lama sekolah meningkat, maka jumlah penduduk miskin menurun. Begitu pula sebaliknya, ketika rata-rata lama sekolah menurun, maka jumlah penduduk miskin meningkat.

Visualisasi Matrix Plot : Heatmap/Correlogram

dt_numerik <- select_if(datajabar.baru, is.numeric)
str(dt_numerik)
## tibble [26 × 4] (S3: tbl_df/tbl/data.frame)
##  $ j.miskin16: num [1:26] 491 199 261 273 299 ...
##  $ AHH2016   : num [1:26] 70.7 70.1 69.4 73.1 70.8 ...
##  $ MYS2016   : num [1:26] 7.83 6.74 6.61 8.5 6.88 6.94 7.55 7.34 6.41 6.89 ...
##  $ EXP2016   : num [1:26] 9537 8077 7074 9580 7079 ...
dt_melt <- cor(dt_numerik[sapply(dt_numerik,is.numeric)])

dt_melt <- melt(dt_melt) 

ggplot(dt_melt, aes(x = Var1, y = Var2, fill = value)) +
  geom_tile() +
  geom_text(aes(label = round(value, 2)), color = "white", size = 3) +
  labs(title = "Correlation Heatmap",
       x = "variabel 1",
       y = "variabel 2")

ggplot(dt_melt, aes(Var1, Var2, fill = value)) +
  geom_tile(color = "white") +
  geom_text(aes(label = round(value, 2)), color = "white", size = 3) +
  scale_fill_gradient2(low = "brown", mid = "green", high = "blue", midpoint = 0, limits = c(-1,1), name="Korelasi") +
  labs(title = "Correlogram",  
       x = "variabel 1",
       y = "variabel 2") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 50, vjust = 1, hjust=1))

Visualisasi antar peubah numerik yang ditampilkan di atas menggunakan correlation heatmap dan correlogram. Hubungan antara satu peubah bebas dengan peubah bebas lain yaitu dari angka harapan hidup tahun 2016, rata-rata lama sekolah tahun 2016 dan pengeluaran per kapita riil tahun 2016 yang disesuaikan dengan formula Atkinson memiliki hubungan linear positif yaitu sebesar 0.76 (tanpa hubungan satu peubah bebas dengan dirinya sendiri). Untuk peubah bebas yang berhubungan dengan dirinya sendiri memiliki hubungan linear positif sempurna yaitu sebesar 1. Sementara itu, masing-masing peubah bebas tersebut memiliki hubungan linear negatif dengan peubah terikat (jumlah penduduk miskin tahun 2016). Hubungan linear antara jumlah penduduk miskin dengan rata-rata lama sekolah memiliki hubungan linear negatif tertinggi yaitu sebesar -0.54 sedangkan hubungan linear antara jumlah penduduk miskin dengan angka harapan hidup memiliki hubungan linear negatif terendah yaitu sebesar -0.3.

Visualisasi Locally Estimated Scatter Plot Smoothing (LOESS)

Non linear Outlier

ggplot(datajabar.baru, aes(x = AHH2016, y = MYS2016)) +
  geom_point(color = "blue", size = 3, alpha = 0.6) +
  geom_smooth(method = "loess", formula = y ~ x, color = "red", linetype = "dashed", size = 1.5) +
  labs(
    x = "angka harapan hidup (tahun)", 
    y = "rata-rata lama sekolah (tahun)", 
    title = "Visualisasi LOESS dari Rata-Rata Lama Sekolah vs Angka Harapan Hidup",
    subtitle = "Scatterplot yang dihaluskan dengan kurva LOESS"
  ) +
  theme_minimal() 
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Visualisasi antar peubah numerik yang ditampilkan dengan scatter plot di atas menampilkan hubungan antara angka harapan hidup dalam tahun dengan rata-rata lama sekolah dalam tahun yang menunjukkan hubungan cenderung positif. Scatter plot tersebut telah disesuaikan dengan kurva LOESS. Penyesuaian ini bertujuan untuk menghaluskan pola hubungan antara kedua variabel, sehingga memungkinkan untuk melihat tren atau pola secara lebih jelas. Penyesuaian dengan kurva LOESS membantu menyajikan pola hubungan ini dengan lebih halus daripada menggunakan garis regresi linear. Kurva LOESS memungkinkan model untuk menyesuaikan diri dengan pola lokal dalam data, sehingga lebih mampu menangkap hubungan yang mungkin tidak terlihat dengan jelas dalam scatter plot biasa.

Data Deret Waktu (Time Series)

Data (Materi 9)

datasaham <- read_excel("C:/Users/user/OneDrive/Dokumen/semester 4/visdat/p10/saham wika dan pp 2022.xlsx")
str(datasaham)
## tibble [12 × 3] (S3: tbl_df/tbl/data.frame)
##  $ Tanggal    : POSIXct[1:12], format: "2022-01-01" "2022-02-01" ...
##  $ ClosePPROJK: num [1:12] 5.3 5.4 5.2 5 5.1 5 5 5.3 5.1 5 ...
##  $ CloseWIKAJK: num [1:12] 103.5 101 99.5 95 96.5 ...
view(datasaham)

Keterangan Data

  1. Data menunjukkan harga penutupan saham perbulan setiap tanggal 1 pada tahun 2022
  2. Harga penutupan saham PT Wijaya Karsa (Persero) Tbk (WIKA.JK) ditunjukkan dengan kode CloseWIKAJK
  3. Harga penutupan saham PT Pembangunan Perumahan Properti Tbk (PPRO.JK) ditunjukkan dengan kode ClosePPROJK

Visualisasi Scatter Plot & Time Series

ggplot(datasaham, aes(x = Tanggal, y = CloseWIKAJK)) +
  geom_point() +
  labs(title = "Scatter Plot dari Data Time Series Saham WIKA.JK Tahun 2022",
       x = "tanggal",
       y = "harga (puluh juta)")

ggplot(datasaham, aes(x = Tanggal, y = CloseWIKAJK)) +
  geom_point() + 
  geom_line() + 
  labs(title = "Scatter Plot yang Terhubung Garis",
  subtitle = "Data Time Series Saham WIKA.JK Tahun 2022",
       x = "tanggal",
       y = "harga (puluh juta)")

ggplot(datasaham, aes(x = Tanggal, y = CloseWIKAJK)) +
  geom_line() +
  labs(title = "Random Walk Time Series Plot",
  subtitle = "Data Time Series Saham WIKA.JK Tahun 2022",
       x = "tanggal",
       y = "harga (puluh juta)")

calculate_moving_average <- function(datasaham, window_size) {
  ma_values <- zoo::rollmean(datasaham$CloseWIKAJK, k = window_size, align = "center", fill = NA)
  datasaham$ma <- ma_values
  return(datasaham)
}

window_size <- 3
datasaham <- calculate_moving_average(datasaham, window_size)

ggplot(datasaham, aes(x = Tanggal)) +
  geom_line(aes(y = CloseWIKAJK), color = "white", size = 1, na.rm = TRUE) +  
  geom_line(aes(y = ma), color = "red", linetype = "dashed", size = 1, na.rm = TRUE) +
  geom_ribbon(aes(ymin = -Inf, ymax = ma), fill = "red", alpha = 0.2) + 
    labs(title = "Data Time Series dengan Moving Average (Ukuran Window : 10)",
       x = "tanggal",
       y = "harga (puluh juta)") +
  theme_minimal()

data <- datasaham$CloseWIKAJK

data.ts <- ts(data)

plot(data.ts, xlab ="tanggal", ylab = "harga (puluh juta)", col="red", main = "Plot Data Saham WIKA.JK Tahun 2022")
points(data.ts)

Plot diatas menunjukkan bahwa data tidak stasioner. Nilai harga saham PT Wijaya Karsa (Persero) Tbk (WIKA.JK) berubah secara fluktuatif setiap bulannya pada tahun 2022. Namun, secara keseluruhan mengalami tren penurunan yang signifikan. Titik-titik puncak dari plot menunjukkan periode dalam bulan dimana harga saham berada dalam level tertinggi. Saham WIKA.JK terlihat menurun dari awal Januari hingga April, Juni hingga Juli, September hingga Oktober dan November hingga Desember. Sementara itu, peningkatan salam terjadi pada bulan Juli hingga Agustus dan September hingga Oktober. Hal ini mengindikasikan bahwa adanya pola tren penurunan pada periode tahun 2022.

Visualisasi Dua Peubah Data Time Series

ggplot(datasaham, aes(x = Tanggal, y = CloseWIKAJK)) +
  geom_line() +
  labs(title = "Random Walk Time Series Plot",
  subtitle = "Data Time Series Saham WIKA.JK Tahun 2022",
       x = "tanggal",
       y = "harga (puluh juta)")

ggplot(datasaham, aes(x = Tanggal, y = ClosePPROJK)) +
  geom_line() +
  labs(title = "Random Walk Time Series Plot",
  subtitle = "Data Time Series Saham PPRO.JK Tahun 2022",
       x = "tanggal",
       y = "harga (puluh juta)")

ggplot() +
  geom_line(data = datasaham, aes(x = Tanggal, y = CloseWIKAJK, color = "WIKAJK")) +
  geom_line(data = datasaham, aes(x = Tanggal, y = ClosePPROJK, color = "PPROJK")) +
  scale_color_manual(values = c("blue", "red")) +
  labs(title = "Perbandingan Data Time Series (WIKAJK vs PPROJK)",
       x = "tanggal",
       y = "harga (puluh juta)",
       color = "keterangan") +
  theme_minimal() 

Plot diatas menunjukkan bahwa antara dua variabel data memiliki nilai yang terlampau jauh. Harga penutupan saham WIKAJK berada di atas harga penutupan saham PPROJK pada tahun 2022. Tren penurunan dari WIKAJK dapat terlihat jelas. Namun, untuk tren dari PPROJK perlu dilihat secara terpisah untuk mengetahui tren yang ada. Setelah dilakukan pemisahan, terlihat bahwa tren penurunan dan peningkatan PPROJK terlihat lebih signifikan daripada WIKAJK. Namun, untuk harga penutupan saham menunjukkan bahwa dominasi dari WIKAJK lebih besar daripada PPROJK di tahun 2022.

Data Geospasial

Data (Materi 10)

produksikopi <- read.csv("C:/Users/user/OneDrive/Dokumen/semester 4/visdat/p10/produksi kopi januari 2021.csv", sep=";")
View(produksikopi)
indonesia_provinsi <- id_map("indonesia", "provinsi")
nama_provinsi <- unique(indonesia_provinsi$nama_provinsi)
print(nama_provinsi)
##  [1] "Aceh"                 "Bali"                 "Kep. Bangka Belitung"
##  [4] "Banten"               "Bengkulu"             "Gorontalo"           
##  [7] "Papua Barat"          "DKI Jakarta"          "Jambi"               
## [10] "Jawa Barat"           "Jawa Tengah"          "Jawa Timur"          
## [13] "Kalimantan Barat"     "Kalimantan Selatan"   "Kalimantan Tengah"   
## [16] "Kalimantan Timur"     "Kalimantan Utara"     "Kep. Riau"           
## [19] "Lampung"              "Maluku Utara"         "Maluku"              
## [22] "Nusa Tenggara Barat"  "Nusa Tenggara Timur"  "Papua"               
## [25] "Riau"                 "Sulawesi Barat"       "Sulawesi Selatan"    
## [28] "Sulawesi Tengah"      "Sulawesi Tenggara"    "Sulawesi Utara"      
## [31] "Sumatera Barat"       "Sumatera Selatan"     "Sumatera Utara"      
## [34] "DI Yogyakarta"
merged_data <- merge(indonesia_provinsi, produksikopi, by.x = "nama_provinsi", by.y = "Provinsi", all.x = TRUE)
## old-style crs object detected; please recreate object with a recent sf::st_crs()

Keterangan Data

Jumlah produksi kopi dalam ton tahun 2021 ditunjukkan dengan kode Produksi

Visualisasi Geospasial

ggplot() +
  geom_sf(data = merged_data, aes(fill = Produksi), color = "black") +
  scale_fill_gradient(low = "brown", high = "green")  +
  theme_void() +
  theme(legend.position = "right") +
  ggtitle("Sebaran Jumlah Produktivitas Kopi di Indonesia Tahun 2021") +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Peta geospasial di atas menggambarkan sebaran jumlah produktivitas kopi di Indonesia tahun 2021. Rentang spektrum warna cokelat mengindikasikan bahwa jumlah produksi kopi di suatu provinsi memiliki jumlah sedikit. Rentang spektrum warna hijau mengindikasikan bahwa jumlah produksi kopi di suatu provinsi memiliki jumlah banyak. Semakin gelap warna cokelat, semakin sedikit jumlah produksi kopi yang dihasilkan di suatu provinsi. Semakin terang warna hijau, semakin banyak jumlah produksi yang dihasilkan di suatu provinsi. Wilayah Sumatra memiliki jumlah produksi kopi terbanyak, khususnya di Provinsi Sumatra Selatan dengan jumlah mencapai 1936 ton di tahun 2021. Wilayah dengan rata-rata produksi kopi terendah berada di Pulau Kalimantan, Maluku dan Papua. Sementara itu, wilayah yang tidak memproduksi kopi ada di Provinsi Bangka Belitung, Kepulauan Riau, DKI Jakarta dan Maluku Utara.