Libraries

library(ggplot2)
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
library(tidyr)
library(forcats)
library(hrbrthemes)
library(viridis)
## Loading required package: viridisLite
library(ggridges)
library(igraph)
## 
## Attaching package: 'igraph'
## The following object is masked from 'package:tidyr':
## 
##     crossing
## The following objects are masked from 'package:dplyr':
## 
##     as_data_frame, groups, union
## The following objects are masked from 'package:stats':
## 
##     decompose, spectrum
## The following object is masked from 'package:base':
## 
##     union
library(ggraph)
library(sf)
## Linking to GEOS 3.14.1, GDAL 3.12.1, PROJ 9.7.1; sf_use_s2() is TRUE
library(rnaturalearth)
library(geodata)
## Loading required package: terra
## terra 1.9.50
## 
## Attaching package: 'terra'
## The following objects are masked from 'package:igraph':
## 
##     blocks, compare
## The following object is masked from 'package:tidyr':
## 
##     extract

Data

Q1, Q2, and Q3 use the same dataset from Github, while Q4 and Q5 will use separate data that will be detailed in their respective sections.

Github link: https://raw.githubusercontent.com/zonination/perceptions/master/probly.csv

The data contains series of probability values linked to several adverbs of probability.

data <- read.table("https://raw.githubusercontent.com/zonination/perceptions/master/probly.csv",
                   header = TRUE, sep = ",")
head(data)
##   Almost.Certainly Highly.Likely Very.Good.Chance Probable Likely Probably
## 1               95            80               85       75     66       75
## 2               95            75               75       51     75       51
## 3               95            85               85       70     75       70
## 4               95            85               85       70     75       70
## 5               98            95               80       70     70       75
## 6               95            99               85       90     75       75
##   We.Believe Better.Than.Even About.Even We.Doubt Improbable Unlikely
## 1         66               55         50       40         20       30
## 2         51               51         50       20         49       25
## 3         80               60         50       30         10       25
## 4         80               60         50       30         10       25
## 5         65               60         50       10         50        5
## 6         80               65         50        7         15        8
##   Probably.Not Little.Chance Almost.No.Chance Highly.Unlikely
## 1           15            20                5              25
## 2           49             5                5              10
## 3           25            20                1               5
## 4           25            20                1               5
## 5           20             5                1               2
## 6           15             5                1               3
##   Chances.Are.Slight
## 1                 25
## 2                  5
## 3                 15
## 4                 15
## 5                 10
## 6                 20

Q1. Horizontal Violin Plot with Boxplot

Since the first visualisation requires a single series, we have to first combine the values across different adverbs of probability.

long <- data %>%
  select(where(is.numeric)) %>%
  pivot_longer(everything(), names_to = "phrase", values_to = "probability")

The result below shows the violin plot of the combined data. The image shows the density of each values on the Y axis and the probability percentage on the X axis. Based on the image, the density of low probability in the data peaks at approximately 10%, then declines up untill its lowest density point at 40%. Past the 40% point, the density increased again with some fluctuations untill in finally rescinded past the 80% point. This means that most of the adverbs of probability in the data used, covers low probability, with medium probability having the lowest coverage.

ggplot(long, aes(x = "all", y = probability)) +
  geom_violin(fill = "#B4D9D2", color = "grey25", trim = TRUE) +
  geom_boxplot(width = 0.12, fill = "#FAFAFA", color = "black",
               outlier.size = 1) +
  coord_flip() +
  scale_y_continuous(limits = c(0, 100), breaks = seq(0, 100, 20)) +
  labs(x = "y", y = "probability (%)") +
  theme_light(base_size = 13)

Q2. Ridgeline Plot

On the image below, we can see the density of probability values given to each adverbs of probability. With phrases like “Almost Certainly” and “Highly Likely” showing an overwhelmingly high density in high probability values, while phrases like ““Almost No Chance” and “Highly Unlikely” showing an overwhelmingly high density in low probability values. While this aligns with our day-to-day use of the phrases, the image gives an insight into how each phrases might be understood differently from one person to another.

long %>%
  mutate(phrase = fct_reorder(phrase, probability, .fun = median, na.rm = TRUE)) %>%
  ggplot(aes(x = probability, y = phrase, fill = phrase)) +
  geom_density_ridges(scale = 1.8, color = "black", linewidth = 0.3,
                      alpha = 1, na.rm = TRUE) +
  scale_x_continuous(limits = c(0, 100), breaks = seq(0, 100, 20)) +
  labs(x = "Probability (%)", y = "Phrase") +
  theme_minimal(base_size = 13) +
  theme(legend.position = "none",
        panel.grid.major.y = element_blank())
## Picking joint bandwidth of 3.43

# Q3. Histogram with Error Bar

In this section, we selected a few phrases. Based on the image, it can be inferred that phrases with a high certainty have a higher average probability values attributed to them, instead of thos with lower certainty. As the previous section also suggests, this aligns with our day-to-day understanding of these phrases. However, on the image below we can see the error rate more clearly, which supports the hipothesis that each phrases may register to different probability values. The phrases “Highly Likely”, “Probable”, “Unlikely”, and “Highly Unlikely” shows a lower error bar, which means that these phrase may be registered at different values on more observations. Meanshile, the phrase “About Even” has a minimal error bar, which shows that the prabability values assigned to the phrase is more homogenous.

chosen <- c("Highly.Likely", "Probable", "About.Even", "Unlikely", "Highly.Unlikely")

summ <- long %>%
  filter(phrase %in% chosen) %>%
  group_by(phrase) %>%
  summarise(mean = mean(probability, na.rm = TRUE),
            se   = sd(probability, na.rm = TRUE) / sqrt(sum(!is.na(probability))),
            .groups = "drop") %>%
  mutate(phrase = factor(phrase, levels = chosen),
         lower = mean - 1.96 * se,
         upper = mean + 1.96 * se)

ggplot(summ, aes(x = phrase, y = mean)) +
  geom_col(fill = "#B4D9D2", alpha = 0.6, width = 0.6) +
  geom_errorbar(aes(ymin = lower, ymax = upper),
                color = "#5FB09A", width = 0.1, linewidth = 1) +
  geom_point(color = "#5FB09A", size = 1.5) +
  scale_y_continuous(limits = c(0, 100), breaks = seq(0, 100, 20)) +
  labs(x = "phrase", y = "mean probability (%)") +
  theme_light(base_size = 13)

Q4. Trade Network (Dummy Data)

Dummy Data

countries <- c("Indonesia", "Singapore", "Malaysia", "Thailand",
               "Philippines", "Viet Nam", "Cambodia")
dummy_matrix <- matrix(c(
  0.0, 18.5, 12.3,  9.8,  7.4,  6.1,  0.9,   # Indonesia
  22.7,  0.0, 28.4, 14.2,  9.6, 11.3,  1.8,   # Singapore
  10.4, 26.9,  0.0, 15.7,  5.2,  8.8,  1.1,   # Malaysia
  8.1, 12.6, 19.3,  0.0,  4.7, 10.5,  6.4,   # Thailand
  3.2,  5.9,  2.8,  3.6,  0.0,  2.4,  0.3,   # Philippines
  5.5,  9.1,  6.7, 11.8,  3.1,  0.0,  4.2,   # Viet Nam
  0.2,  0.7,  0.4,  3.9,  0.2,  2.6,  0.0    # Cambodia
), nrow = 7, byrow = TRUE, dimnames = list(countries, countries))
trade <- as.data.frame(as.table(dummy_matrix), stringsAsFactors = FALSE) %>%
  setNames(c("from", "to", "value")) %>%
  filter(from != to)

Visualisation

Based on the image below, Indonesia has a high trade volume towards Singapore and vice-versa, while only having a low trade volume with Vietnam. Cambodia is shown to only do imports, while not actively trading in exports. Singapore seems to be the main trade hub, having major imports and exports with several countries, such as Malaysia, Indonesia, and Thailand.

trade_top <- trade %>%
  slice_max(value, n = 25) %>%
  mutate(involves_idn = if_else(from == "Indonesia" | to == "Indonesia",
                                "Involves Indonesia", "Other flows"))

g <- graph_from_data_frame(trade_top, directed = TRUE)
V(g)$focus <- V(g)$name == "Indonesia"

ggraph(g, layout = "circle") +
  geom_edge_arc(aes(edge_width = value, edge_colour = involves_idn),
                strength = 0.2, alpha = 0.85,
                arrow = arrow(length = unit(3, "mm"), type = "closed"),
                end_cap = circle(6, "mm")) +
  geom_node_point(aes(fill = focus), shape = 21, size = 10,
                  colour = "black", stroke = 0.8, show.legend = FALSE) +
  geom_node_text(aes(label = name), repel = TRUE, fontface = "bold", size = 3.8) +
  scale_fill_manual(values = c(`TRUE` = "#FF7F00", `FALSE` = "#0B3C5D")) +
  scale_edge_colour_manual(values = c("Involves Indonesia" = "#D00000",
                                      "Other flows"        = "#003F7F"),
                           name = NULL) +
  scale_edge_width(range = c(0.3, 3.5), name = "Trade volume (USD bn)") +
  theme_void(base_size = 13) +
  theme(legend.position = "bottom")

Q5. Indonesian

Data (BPS 2025)

csv_file <- "C:/Users/LENOVO/Downloads/Jumlah dan Persentase Penduduk Miskin Menurut Provinsi 2025.csv"
poor <- read.csv(csv_file, check.names = FALSE, fileEncoding = "UTF-8")
 
poverty <- poor %>%
  transmute(province = Provinsi,
            n_poor   = `Jumlah Penduduk Miskin - September (ribu) (Ribu)`,
            pct      = `Persentase Penduduk Miskin - September`) %>%
  # change both columns above to the "- Maret" versions to map March instead
  filter(province != "Indonesia")

Visualisation

The image shows that based on 2025 data from Badan Pusat Statistik (BPS), the eastern part of Indonesia has a higher poverty rate than other regions. Papua provinces have a significantly higher rate of poverty, between 15% to 20% compared to other regions, while the western part of Indonesia has significantly lower poverty rate. Provinces in Kalimantan is observed to have the lowest poverty rate compared to other major Indonesian islands. Thus, the data shows that even in 2025, poverty across indonesia is still a significant problem the government needs to resolve, with higher magnitude on several parts of the nation.

provinces <- gadm("IDN", level = 1, path = tempdir()) %>% st_as_sf()
## Cached as: C:\Users\LENOVO\AppData\Local\Temp\RtmpWwYPb2/gadm/gadm41_IDN_1_pk.rds
parent <- c("Papua Selatan"    = "Papua",
            "Papua Tengah"     = "Papua",
            "Papua Pegunungan" = "Papua",
            "Papua Barat Daya" = "Papua Barat")
 
poverty_map <- poverty %>%
  mutate(province = recode(province, !!!parent),
         pop = n_poor / (pct / 100)) %>%
  group_by(province) %>%
  summarise(pct = sum(n_poor) / sum(pop) * 100, .groups = "drop") %>%
  # GADM spellings differ from BPS for three provinces
  mutate(province = recode(province,
                           "DKI Jakarta"                = "Jakarta Raya",
                           "DI Yogyakarta"              = "Yogyakarta",
                           "Kepulauan Bangka Belitung"  = "Bangka Belitung"))
 
map_data <- provinces %>% left_join(poverty_map, by = c("NAME_1" = "province"))
 
# Anything printed here did not match: fix the spelling in the recode() calls above
setdiff(poverty_map$province, provinces$NAME_1)
## character(0)
map_data$NAME_1[is.na(map_data$pct)]
## character(0)
ggplot(map_data) +
  geom_sf(aes(fill = pct), colour = "grey20", linewidth = 0.2) +
  scale_fill_viridis_c(option = "magma", direction = -1,
                       name = "Poverty Rate",
                       guide = guide_colourbar(title.position = "left",
                                               barwidth = unit(10, "cm"),
                                               barheight = unit(0.5, "cm"))) +
  theme_void(base_size = 13) +
  theme(legend.position = "bottom")