library(dplyr)
library(ggplot2)
library(vcd)
library(stringr)
library(lubridate)
library(forcats)
library(readr)
library(tidyr)
library(purrr)
library(sf)
library(RColorBrewer)
library(grid)
library(ggalluvial)
library(ggridges)
library(tigris)
library(patchwork)
df <- read.csv("Parks_Special_Events_20251027.csv")
parks <- read.csv("Parks_Properties_20251209.csv")
parks <- parks |>
mutate(
geometry = st_as_sfc(multipolygon)
) |>
st_as_sf(crs = 4326)theme_edav <- function(base_size = 14, base_family = "Helvetica") {
theme_minimal(base_size = base_size, base_family = base_family) +
theme(
panel.grid = element_blank(),
axis.line = element_line(color = "black", linewidth = 0.6),
axis.title = element_text(face = "bold"),
plot.title = element_text(face = "bold", size = base_size + 2),
legend.key.height = unit(0.6, "cm"),
legend.key.width = unit(0.6, "cm"),
legend.title = element_text(face = "bold"),
panel.border = element_blank()
)
}
theme_set(theme_edav())
options(tigris_use_cache = TRUE)
ny_counties <- counties(state = "NY", cb = TRUE, year = 2020) |>
st_transform(4326)
nyc_outline <- ny_counties |>
filter(NAME %in% c("New York", "Kings", "Queens", "Bronx", "Richmond"))
df <- df |>
mutate(
datetime = mdy_hms(Date.and.Time),
year = year(datetime),
month = month(datetime, label = TRUE, abbr = TRUE),
dow = wday(datetime, label = TRUE, abbr = TRUE),
hour = hour(datetime),
doy = yday(datetime),
Attendance_num = parse_number(as.character(Attendance)),
Attendance_num = ifelse(Attendance_num < 0, NA, Attendance_num)
) |>
filter(year != 2002 & year != 2017 & LocationType == "Park") |>
mutate(
park_clean = Location |>
tolower() |>
str_replace_all("[^a-z0-9 ]", " ") |>
str_squish()
)
# head(df)events_by_month <- df |>
count(month, name = "n") |>
arrange(month)
ggplot(events_by_month, aes(x = month, y = n)) +
geom_col() +
labs(
title = "Number of Events by Month",
x = "Month",
y = "Number of Events"
) +
theme_edav()
From Figure 3.1 we can clearly see the seasonal pattern of the number of events held in NYC parks. The events are heavily concentrated in late spring and summer, peaking between June and August, Droping sharply during winter. It hightlights the strong dependence of oourdoor public activities on weather and daylight.
top8_cat <- df |>
count(Category, sort = TRUE) |>
slice_head(n = 8) |>
pull(Category)
df_top <- df |>
filter(Category %in% top8_cat)
ggplot(df_top, aes(x = doy, y = Category, fill = Category)) +
geom_density_ridges(alpha = 0.7, scale = 1.5) +
labs(
title = "Seasonal Distribution of Events by Category (Top 8)",
x = "Day of Year",
y = "Category"
) +
theme_edav()
Then we use ridge plot to decompose the overall seasonal pattern into categories. We can see from Figure 3.2 that fitness, sport and KIM related programs show sharply peaked distributions centered around the mid of the year, showing the preference on summer conditions. In contrast, arts and culture as well as family festival events exhibit broader spreads, suggesting greater flexibility. These patterns suggest that the seasonal rhythm of park usage is not uniform, but is structured by the functional nature of events.
df_hex <- df |>
filter(!is.na(Attendance_num), !is.na(doy)) |>
mutate(
Attendance_capped = pmin(
Attendance_num,
quantile(Attendance_num, 0.99, na.rm = TRUE)
)
)
ggplot(df_hex, aes(x = doy, y = Attendance_capped)) +
geom_hex(bins = 20, color = NA) +
scale_y_log10() +
scale_fill_gradientn(
colours = c("#e0f3f8", "#abd9e9", "#74add1", "#4575b4", "#313695"),
# trans = "log10",
name = "Count"
) +
labs(
title = "Attendance vs Day-of-Year",
x = "Day of Year",
y = "Attendance",
fill = "Count"
)
The hexbin heatmap of attendance versus day-of-year(Figure 3.3) shows that the seasonal increase in park activity is driven primarily by a large number of small- to mid-sized events. The highest concentration of events occurs during the summer months with attendance levels between approximately 5 and 20 participants. In contrast, high-attendance events are relatively rare and distributed sparsely across the year.
This pattern suggests that peak park usage reflects frequent, lower-density gatherings rather than occasional large-scale events, highlighting the role of parks as flexible spaces for decentralized public activity.
df_year <- df |>
filter(!is.na(year), year >= 2019, !is.na(Category))
top_cat <- df_year |>
count(Category, sort = TRUE) |>
slice_head(n = 7) |>
pull(Category)
df_year_filtered <- df_year |>
filter(Category %in% top_cat)
summary_df_year <- df_year_filtered |>
group_by(year, Category) |>
summarise(
total_events = n(),
.groups = "drop"
) |>
mutate(
year = factor(year),
Category = fct_infreq(Category)
)ggplot(
summary_df_year,
aes(
x = year,
stratum = Category,
alluvium = Category,
y = total_events,
fill = Category
)
) +
geom_flow(
stat = "alluvium",
lode.guidance = "frontback",
color = "grey60",
alpha = 0.8
) +
geom_stratum(width = 0.3, color = "grey30") +
scale_fill_brewer(palette = "Set2", name = "Category") +
labs(
title = "NYC Event Categories Over Time",
x = "Year",
y = "Number of Events"
) +
theme_edav()
Base on our common knowledge, activities should be effected greatly during the pandemic. But surprisingly, we can see from Figure 3.4 that with the immediate impact of COVID-19 restrictions, the number of events in 2020 is higher compared to the other years. Looking at the poportion over different categories, there’s a huge decrease on the number of Mobile Unit-related events.
df_year <- df |>
filter(!is.na(year), year >= 2019, !is.na(Category))
top_cat <- df_year |>
count(Category, sort = TRUE) |>
slice_head(n = 7) |>
pull(Category)
df_year_filtered <- df_year |>
filter(Category %in% top_cat)
summary_df_year_att <- df_year_filtered |>
group_by(year, Category) |>
summarise(
total_attendance = sum(Attendance_num, na.rm = TRUE),
.groups = "drop"
) |>
mutate(
year = factor(year),
Category = fct_infreq(Category)
)ggplot(
summary_df_year_att,
aes(
x = year,
stratum = Category,
alluvium = Category,
y = total_attendance,
fill = Category
)
) +
geom_flow(
stat = "alluvium",
lode.guidance = "frontback",
color = "grey60",
alpha = 0.8
) +
geom_stratum(width = 0.3, color = "grey30") +
scale_fill_brewer(palette = "Set2", name = "Category") +
labs(
title = "NYC Event Categories by Total Attendance Over Time",
x = "Year",
y = "Total Attendance"
) +
theme_edav()
To study the strange pattern in the first alluvial graph, we plot the total attendance instead of number of events(shown in Figure 3.5). When examining total attendance, a different recovery pattern emerges. While event counts is high during pandemics, total attendance is obviously lower than in any other years, particularly for mobile unit and family festival programs. This indicates that those events happened during pandemics is mostly made up of small scale events.
mobile_df <- df |>
filter(Category == "Mobile Unit",
!is.na(Classification),
str_trim(Classification) != "" )
top_mobile_class <- mobile_df |>
count(Classification, sort = TRUE) |>
slice_head(n = 10) |>
pull(Classification)
mobile_year_class <- mobile_df |>
filter(Classification %in% top_mobile_class) |>
count(year, Classification, name = "n") |>
filter(year >= 2017) |>
mutate(
year = factor(year),
Classification = fct_infreq(Classification)
)
ggplot(mobile_year_class,
aes(x = year, y = n, fill = Classification)) +
geom_col(position = "fill") +
scale_fill_brewer(palette = "Set2", name = "Classification") +
labs(
title = "Mobile Unit Classification Trends Over Time",
x = "Year",
y = "Proportion",
fill = "Classification"
) +
theme_edav()
Since mobile unit is a very important category of the events, and the category name itself doesn’t represent the function, we dive deeper into the detailed calssification to see the patterns.
From Figure 3.6, we can see the proportional composition of mobile unit classifications shows a clear structural shift over time. Prior to the pandemic, mobile units were dominated by movie-related programming, reflecting their role as mobile entertainment platforms. In 2020, this structure collapses into a single, pandemic-specific classification, indicating an emergency-driven reorganization rather than gradual change.
In the post-pandemic period, mobile units increasingly shift toward M.U.T.S, play, and fitness-related functions. This new composition stabilizes after 2022 and does not return to pre-pandemic patterns, suggesting that mobile units evolved from cultural entertainment vehicles into flexible public service infrastructure.
parks <- parks |>
mutate(
park_clean = NAME311 |>
tolower() |>
str_replace_all("[^a-z0-9 ]", " ") |>
str_squish()
)
theme_hex_small <- function() {
theme_minimal(base_size = 8) +
theme(
axis.text = element_blank(),
axis.title = element_blank(),
axis.ticks = element_blank(),
panel.grid = element_blank(),
plot.title = element_text(size = 10, face = "bold", hjust = 0.5),
legend.title = element_text(size = 8),
legend.text = element_text(size = 7),
legend.key.height = unit(0.3, "in"),
legend.key.width = unit(0.15, "in")
)
}
plot_hex_for_year <- function(year_sel, pal) {
park_attendance_year <- df |>
filter(year == year_sel, !is.na(park_clean)) |>
group_by(park_clean) |>
summarise(total_attendance = sum(Attendance_num, na.rm = TRUE)) |>
ungroup()
parks_joined_year <- parks |>
left_join(park_attendance_year, by = "park_clean") |>
mutate(total_attendance = replace_na(total_attendance, 0))
parks_valid_year <- parks_joined_year |>
st_make_valid()
parks_pts_year <- parks_valid_year |>
st_centroid(of_largest_polygon = TRUE) |>
mutate(
lon = st_coordinates(geometry)[, 1],
lat = st_coordinates(geometry)[, 2]
) |>
st_drop_geometry()
ggplot() +
stat_summary_hex(
data = parks_pts_year,
aes(x = lon, y = lat, z = total_attendance),
bins = 20,
fun = "sum",
color = NA
) +
scale_fill_gradientn(
colours = pal,
name = "Total Attendance"
) +
geom_sf(
data = nyc_outline,
fill = NA,
color = "grey70",
linewidth = 0.4
) +
coord_sf() +
labs(
title = paste("NYC Events Attendance in", year_sel),
x = "Longitude",
y = "Latitude"
) +
theme_hex_small()
}pal_blue <- c("#e0f3f8", "#abd9e9", "#74add1", "#4575b4", "#313695")
pal_red <- c("#fee5d9", "#fcae91", "#fb6a4a", "#de2d26", "#a50f15")
p2019 <- plot_hex_for_year(2019, pal_blue)
p2020 <- plot_hex_for_year(2020, pal_red)
p2021 <- plot_hex_for_year(2021, pal_blue)
p2022 <- plot_hex_for_year(2022, pal_blue)
(p2019 + p2020) / (p2021 + p2022)
Figure 3.7 shows a clear transformation over time. In 2019, park event participation was highly concentrated in a small number of central locations, reflecting a centralized pattern of public gathering. In 2020, attendance levels dropped sharply and spatial hotspots largely disappeared, resulting in a more uniform distribution across the city.
In 2021 and 2022, attendance gradually recovered, with high-attendance areas re-emerging. However, these hotspots were more spatially distributed than in 2019, suggesting a shift toward a more decentralized pattern of park usage.