Catchment Households and Deprivation, Brighton and Hove

Author

Adam Dennett

In the map below, every circle represents a postcode in Brighton and Hove, sized according to the estimated number of households with dependent children living in that postcode, and coloured according to the proportion of all children in that postcode aged 0-15, living in income deprived families.

Hover the mouse over each dot to see the postcode and estimated counts of households with dependent children.

The boundaries on the map represent the current (2024/25) school catchment boundaries (blue) and the proposed 2026/27 boundaries (red) taken from the Catchment Area Postcode List published by Brighton and Hove District Council

The graphs below the map

These estimates use postcode level population and household estimates derived from the 2021 Census, obtained from NOMIS - https://www.nomisweb.co.uk/sources/census_2021_pc

Details of households with dependent children taken from the 2021 Census - Table TS003 Household Composition - https://www.nomisweb.co.uk/sources/census_2021_bulk

Deprivation data is taken from the 2019 Index of Multiple Deprivation - https://imd-by-geo.opendatacommunities.org - In this case I have used the Income Deprivation Affecting Children Index - IDACI.

Counts of households with dependent children have been estimated using household proportions at postcode level to distribute output area level numbers of households with dependent children proportionally to each postcode in each output area.

Code

library(sf)          # For spatial data manipulation
library(dplyr)       # For data manipulation
library(ggplot2)     # For plotting
library(RColorBrewer) # For color palettes

# 1. Perform spatial join (use your existing data)
joined_data_city <- st_join(bn_postcodes_pop1, optionZ, join = st_within)

# 2. Aggregate the data for the entire city
aggregated_data_city <- joined_data_city %>%
  group_by(idaci_decile) %>%
  summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>%
  ungroup()

# 3. Create a bar graph for the entire city
ggplot(aggregated_data_city, aes(x = factor(idaci_decile), y = total_dep_ch, fill = idaci_decile)) +
  geom_bar(stat = "identity") +
  scale_fill_gradientn(
    colors = brewer.pal(11, "RdYlBu"),
    limits = c(0, 10),  # Ensure 0 maps to red and 10 maps to blue
    breaks = 0:10,      # Optional: discrete legend ticks for each decile
    name = "IDACI Decile"
  ) +
  labs(
    title = "Count of Households with Dependent Children by IDACI Decile - Brighton and Hove",
    x = "IDACI Decile",
    y = "Count of Households with Dependent Children"
  ) +
  theme_minimal() +
  theme(legend.position = "bottom")  # Adjust the legend position

Rayshader Plot

Code

# library(basemaps)
# 
# # Define bounding box for your map area
# # Calculate bounding box explicitly
# bbox <- st_bbox(bn_postcodes_pop1, crs = 4326)
# 
# # set defaults for the basemap
# set_defaults(ext = bbox, map_service = "esri", map_type = "world_light_gray_base")
# 
# # Overlay with your height plot or main plot
# osm_basemap <- ggplot() +
#   basemap_gglayer(ext = bbox) +  # Pass bbox explicitly here
#   scale_fill_identity() +
#   coord_sf()
# 
# 
# osm_basemap
# 
# # Render the rayshader plot with the basemap
# osm_basemap %>%
#   plot_gg(
#     ggobj_height = height_plot,
#     multicore = TRUE,
#     width = 5,
#     height = 5,
#     scale = 50,  # Reduce height scaling
#     windowsize = c(800, 800),
#     shadow_intensity = 0.5  # Optional: adjust shadow for realism
#   )

Code

height_plot <- ggplot(rayshader_data, aes(x = x, y = y, colour = z)) +
  geom_point(size = 0.05) +  # Raster fills for continuous height mapping
  scale_colour_gradientn(
    colors = brewer.pal(9, "YlGnBu"),  # Adjust as needed
    limits = c(min(rayshader_data$z, na.rm = TRUE), max(rayshader_data$z, na.rm = TRUE)),
    name = "Height"
  ) +
  theme_minimal(base_size = 10) +
  labs(
    title = "Households with Dependent Children",
    colour = "IDACI Decile"
  ) +
  theme(
    panel.background = element_rect(fill = "white", color = NA), # White background
    plot.background = element_rect(fill = "white", color = NA),  # White plot background
    legend.background = element_rect(fill = "white", color = NA),
    legend.position = "none"  # Disable all legends# White legend background
  )# Minimal theme for height mapping


#height_plot

main_plot <- ggplot(rayshader_data, aes(x = x, y = y, colour = w)) +
  geom_point(size = 0.05) +
  scale_colour_gradientn(
    colors = brewer.pal(11, "RdYlBu"),
    limits = c(0, 10),
    name = "IDACI Decile"
  ) +
  theme_minimal(base_size = 10) +
  labs(
    title = "IDACI Decile",
    colour = "IDACI Decile"
  ) +
  theme(
    panel.background = element_rect(fill = "white", color = NA), # White background
    plot.background = element_rect(fill = "white", color = NA),  # White plot background
    legend.background = element_rect(fill = "white", color = NA),
    legend.position = "none"  # Disable all legends# White legend background
  )

#main_plot

main_plot %>%
  plot_gg(
    ggobj_height = height_plot,  # Specify height mapping plot
    multicore = TRUE,
    width = 5,
    height = 5,
    scale = 50,  # Adjust scale for height exaggeration
    windowsize = c(1200, 1200), 
    sunangle=225
  )

# Add 3D rendering
render_camera(theta = 10, phi = 30, zoom = 0.4)
rgl::rglwidget()

Code

# 
# bn_bbox = st_bbox(bn_postcodes_pop1)
# 
# bn_roads <- osmdata::opq(bn_bbox) %>% 
#   osmdata::add_osm_feature("highway") %>% 
#   osmdata::osmdata_sf() 
# 
# bn_roads <- bn_roads$osm_lines
# 
# ggplot(bn_roads, aes(color = highway)) + 
#   geom_sf() +
#   theme(legend.position = "none") +
#   labs(title = "Open Street Map `highway` attribute in Brighton")
# 
# rayshader_data_sf <- rayshader_data %>%
#   st_as_sf(coords = c("x", "y"), crs = 4326)
# 
# ggplot(rayshader_data_sf,aes(color=z)) + 
#   geom_sf() +
#   theme(legend.position = "none") 
# 
# map_plot <- ggplot() +
#   geom_sf(data = rayshader_data_sf, aes(color = z)) +
#   geom_sf(data = bn_roads, aes(fill = highway), size = 0.5) +
#   theme_minimal() +
#   theme(
#     panel.background = element_rect(fill = "white", color = NA), # White background
#     plot.background = element_rect(fill = "white", color = NA),  # White plot background
#     legend.background = element_rect(fill = "white", color = NA),
#     legend.position = "none") 
# 
# map_plot %>%
#   plot_gg(
#     ggobj_height = height_plot,  # Specify height mapping plot
#     multicore = TRUE,
#     width = 5,
#     height = 5,
#     scale = 50,  # Adjust scale for height exaggeration
#     windowsize = c(1200, 1200), 
#     sunangle=225
#   )# Disable all legends# White legend

Code

library(sf)
library(dplyr)
library(ggplot2)
library(RColorBrewer)

# Perform spatial join
joined_data_city <- st_join(bn_postcodes_pop1, option0, join = st_within)

# Drop spatial attributes before further operations
joined_data_city <- st_drop_geometry(joined_data_city)

# Aggregate for city-level proportions
city_proportions <- joined_data_city %>%
  group_by(idaci_decile) %>%
  summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>%
  mutate(city_prop = total_dep_ch / sum(total_dep_ch)) %>%
  ungroup()

# Aggregate for catchment-level proportions
catchment_proportions <- joined_data_city %>%
  group_by(catchment, idaci_decile) %>%
  summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>%
  mutate(catchment_prop = total_dep_ch / sum(total_dep_ch)) %>%
  ungroup()

`summarise()` has grouped output by 'catchment'. You can override using the
`.groups` argument.

Code

# Merge city and catchment proportions
# comparison_data <- catchment_proportions %>%
#   left_join(city_proportions %>% select(idaci_decile, city_prop), by = "idaci_decile")

# # Plot proportions comparison
# ggplot(comparison_data, aes(x = factor(idaci_decile))) +
#   geom_bar(aes(y = catchment_prop, fill = "Catchment"), stat = "identity", position = "dodge") +
#   geom_line(aes(y = city_prop, group = 1, color = "City"), size = 1.2) +
#   scale_fill_manual(values = c("Catchment" = "steelblue")) +
#   scale_color_manual(values = c("City" = "darkred")) +
#   labs(
#     title = "Proportions of Households by IDACI Decile",
#     subtitle = "Comparison of Catchment Proportions with City-Level Proportions",
#     x = "IDACI Decile",
#     y = "Proportion of Households",
#     fill = "Legend",
#     color = "Legend"
#   ) +
#   facet_wrap(~ catchment, scales = "free_y") +
#   theme_minimal() +
#   theme(
#     legend.position = "bottom",
#     axis.text.x = element_text(angle = 45, hjust = 1)
#   )

Code

# Initialize an empty data frame to hold the results
comparison_data <- data.frame()

# Get the unique IDACI deciles from city_proportions
unique_deciles <- unique(city_proportions$idaci_decile)

# Loop through each catchment in catchment_proportions
for (catchment in unique(catchment_proportions$catchment)) {
  # Subset data for the current catchment
  catchment_data <- catchment_proportions %>%
    filter(catchment == !!catchment)
  
  # Loop through each IDACI decile
  for (decile in unique_deciles) {
    # Check if the current decile exists in the catchment data
    if (decile %in% catchment_data$idaci_decile) {
      # If it exists, get the proportion for that decile
      catchment_prop <- catchment_data$catchment_prop[catchment_data$idaci_decile == decile]
    } else {
      # If it does not exist, assign a proportion of 0
      catchment_prop <- 0
    }
    
    # Add the current row to the comparison_data
    comparison_data <- rbind(
      comparison_data,
      data.frame(
        catchment = catchment,
        idaci_decile = decile,
        catchment_prop = catchment_prop
      )
    )
  }
}


comparison_data <- comparison_data %>%
  left_join(city_proportions, by = "idaci_decile")

# Add city-level proportions to the data
# comparison_data <- comparison_data %>%
#   left_join(city_proportions %>% select(idaci_decile, city_prop), by = "idaci_decile")

comparison_data <- comparison_data %>%
  mutate(
    diff = catchment_prop - city_prop
  )

Code

ggplot(comparison_data, aes(x = factor(idaci_decile), y = diff, fill = idaci_decile)) +
  geom_bar(stat = "identity") +
  scale_fill_gradientn(
    colors = brewer.pal(11, "RdYlBu"),
    limits = c(0, 10),  # Ensure 0 maps to red and 10 maps to blue
    breaks = 0:10,      # Optional: discrete legend ticks for each decile
    name = "IDACI Decile"
  ) +
  facet_wrap(~ catchment, scales = "fixed") +
  labs(
    title = "Deviation from City Proportion of Deprivation - Current Catchments",
    x = "IDACI Decile",
    y = "Proportion of Households in Catchment"
  ) +
  theme_minimal() +
  theme(legend.position = "bottom")

 [1] "#808080" "#808080" "#808080" "#808080" "#808080" "#808080" "#808080"
 [8] "#808080" "#808080" "#808080"

--- title: "Catchment Households and Deprivation, Brighton and Hove" author: "Adam Dennett" --- In the map below, every circle represents a postcode in Brighton and Hove, sized according to the estimated number of households with dependent children living in that postcode, and coloured according to the proportion of all children in that postcode aged 0-15, living in income deprived families. Hover the mouse over each dot to see the postcode and estimated counts of households with dependent children. The boundaries on the map represent the current (2024/25) school catchment boundaries (blue) and the proposed 2026/27 boundaries (red) taken from the [Catchment Area Postcode List](https://www.brighton-hove.gov.uk/schools-and-learning/apply-school/proposed-admission-arrangements-brighton-hove-schools-2026-2027#tab--catchment-area-postcode-list) published by Brighton and Hove District Council The graphs below the map These estimates use postcode level population and household estimates derived from the 2021 Census, obtained from NOMIS - <https://www.nomisweb.co.uk/sources/census_2021_pc> Details of households with dependent children taken from the 2021 Census - Table TS003 Household Composition - <https://www.nomisweb.co.uk/sources/census_2021_bulk> Deprivation data is taken from the 2019 Index of Multiple Deprivation - <https://imd-by-geo.opendatacommunities.org> - In this case I have used the [Income Deprivation Affecting Children Index](https://opendatacommunities.org/def/concept/general-concepts/imd/idaci) - IDACI. Counts of households with dependent children have been estimated using household proportions at postcode level to distribute output area level numbers of households with dependent children proportionally to each postcode in each output area. ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE, include=FALSE, results='hide'} library(tidyr) library(dplyr) library(writexl) library(tidyverse) library(here) library(janitor) library(sf) library(usethis) library(tmap) library(readxl) library(r5r) library(RColorBrewer) library(accessibility) library(data.table) library(ggplot2) library(interp) library(h3jsr) library(h3r) library(osmextract) library(stplanr) library(od) library(tidytransit) library(ggrepel) library(DBI) library(RPostgreSQL) library(RPostgres) library(rpostgis) ``` ```{r Connect-to-Database, warning=FALSE, message=FALSE, echo=FALSE} #connect to a local database to store some data later on con <- dbConnect(RPostgres::Postgres(), dbname = 'gisdb2', host = 'localhost', port = 5432, user = 'postgres', password = 'postgres', sslmode = 'disable') #check connection #dbListTables(con) ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE, include=FALSE} bn_postcodes <- st_read(con, "bn_postcode_centroids") ew_oa <- st_read(con, "oa_2021_ew_bfc_v8") %>% st_set_crs(27700) # Filter the data frame bn_oa <- ew_oa %>% filter(grepl("Brighton and Hove", lsoa21nm)) bh_hh_comp <- read_csv("data/bh_hh_comp.csv") bh_hh_oa_map <- left_join(bh_hh_comp, bn_oa, by = join_by(oa21 == oa21cd)) bh_hh_oa_map <- left_join(bn_oa, bh_hh_comp, by = join_by(oa21cd == oa21)) optionZ <- st_read("data/optionZ_Mar25.geojson") option0 <- st_read("data/current_catchments_fix.geojson") bn_all_deprivation <- read_csv("data/bn_all_deprivation.csv") %>% clean_names() bn_dep_ch_ethnicity_hrp <- read_csv("data/bh_dep_ch_ethnicity_hrp.csv") %>% clean_names() pcds_p003 <- read_csv("data/pcds_p003.csv") %>% clean_names() pcd_hh_count <- read_csv("data/pcd_p002.csv") %>% clean_names() pcd_person_count <- read_csv("data/pcd_p001.csv") %>% clean_names() bh_pcd_hh_count <- pcd_hh_count %>% filter(grepl("^BN", postcode)) bh_pcd_person_count <- pcd_person_count %>% filter(grepl("^BN", postcode)) # Summing count for each unique postcode bh_agg_pcd_person_count <- bh_pcd_person_count %>% group_by(postcode) %>% summarise(total_count = sum(count, na.rm = TRUE)) bn_postcodes_pop <- bn_postcodes %>% dplyr::select(c("pcd", "pcd2", "pcds", "oa21", "lsoa21", "msoa21", "oseast1m", "osnrth1m", "geometry")) %>% merge(bh_pcd_hh_count, by.x = "pcds", by.y = "postcode") %>% merge(bh_agg_pcd_person_count, by.x = "pcds", by.y = "postcode") bn_postcodes_pop <- bn_postcodes_pop %>% rename( household_count = count, person_count = total_count ) # Calculate the proportion of households in each postcode bn_postcodes_pop1 <- bn_postcodes_pop %>% group_by(oa21) %>% mutate(pcd_oa_hh_prop = household_count / sum(household_count)) %>% mutate(pcd_oa_pers_prop = person_count / sum(person_count)) %>% ungroup() #merge households with dependent children bn_postcodes_pop1 <- bn_postcodes_pop1 %>% merge(bh_hh_comp, by.x = "oa21", by.y = "oa21") #merge dependent children by ethnicity bn_postcodes_pop1 <- bn_postcodes_pop1 %>% merge(bn_dep_ch_ethnicity_hrp, by.x = "oa21", by.y = "oa21") bn_postcodes_pop1 <- bn_postcodes_pop1 %>% mutate(pcd_dep_ch_hh_count = pcd_oa_hh_prop * all_dependent_ch_hh) bn_postcodes_pop1 <- bn_postcodes_pop1 %>% mutate(pcd_dep_ch_hh_count_round = round(pcd_oa_hh_prop * all_dependent_ch_hh)) bn_postcodes_pop1 <- bn_postcodes_pop1 %>% mutate(pcd_dep_ch_0_18_total_count = round(pcd_oa_pers_prop * total_dep_ch_0_18)) %>% mutate(pcd_dep_ch_0_18_white_count = round(pcd_oa_pers_prop * white_dep_ch_0_18)) %>% mutate(pcd_dep_ch_0_18_non_white_count = round(pcd_oa_pers_prop * non_white_dep_ch_0_18)) bn_postcodes_pop1 <- merge( bn_postcodes_pop1, bn_all_deprivation, by.x = "pcds", by.y = "postcode") library(sf) library(dplyr) # Ensure the CRS (Coordinate Reference System) is the same for spatial operations bn_postcodes_pop1 <- st_transform(bn_postcodes_pop1, st_crs(option0)) option0 <- st_transform(option0, st_crs(bn_postcodes_pop1)) optionZ <- st_transform(optionZ, st_crs(bn_postcodes_pop1)) # Add catchment_2024 based on intersections with option0 bn_postcodes_pop1 <- bn_postcodes_pop1 %>% mutate(catchment_2024 = sapply( st_intersects(geometry, option0), function(idx) { if (length(idx) > 0) { option0$catchment[idx] # Extract the "catchment" value from the intersecting polygon } else { NA # No intersection, assign NA } } )) # Add catchment_2026 based on intersections with optionZ bn_postcodes_pop1 <- bn_postcodes_pop1 %>% mutate(catchment_2026 = sapply( st_intersects(geometry, optionZ), function(idx) { if (length(idx) > 0) { optionZ$catchment[idx] # Extract the "catchment" value from the intersecting polygon } else { NA # No intersection, assign NA } } )) # Check the structure of the updated dataframe str(bn_postcodes_pop1) write_csv(bn_postcodes_pop1, "bn_postcodes_pop1.csv") ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} tmap_mode("view") tmap_options(check_and_fix = T) selection <- bn_postcodes_pop1 %>% dplyr::select(c("pcds", "pcd_dep_ch_hh_count_round", "pcd_dep_ch_0_18_non_white_count","total_hh", "idaci_decile")) %>% rename( Postcode = pcds, `Households with Dependent Children` = pcd_dep_ch_hh_count_round, `Total non-white children 0-18` = pcd_dep_ch_0_18_non_white_count, `Total Households` = total_hh, `IDACI Decile` = idaci_decile ) # # Filter the data # filtered_data <- bn_postcodes_pop1 %>% # filter(pcd_dep_ch_hh_count_round > 5) # # # Adjust the data to emphasize small values # tm_shape(optionZ) + # tm_polygons(alpha = 0, # border.col = "red") + # tm_shape(option0) + # tm_polygons(alpha = 0) + # tm_shape(selection) + # tm_dots( # size = "Households with Dependent Children", # alpha = 0.5, # fill = "IDACI Decile", # col = "IDACI Decile", # col.alpha = 0.2, # lwd = 0, # palette = "brewer.rd_yl_bu" # ) + # tm_layout( # frame = FALSE, # legend.outside = TRUE # ) ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} library(leaflet) optionZ <- st_transform(optionZ, 4326) option0 <- st_transform(option0, 4326) selection <- st_transform(selection, 4326) bn_postcodes_pop1 <- st_transform(bn_postcodes_pop1, 4326) leaflet() %>% # Add a less intrusive basemap addProviderTiles(providers$CartoDB.Positron) %>% # Light, minimalistic basemap # Add the `optionZ` polygon layer with a red border addPolygons( data = optionZ, color = "red", opacity = 0.5, fillOpacity = 0 ) %>% # Add the `option0` polygon layer addPolygons( data = option0, color = NULL, opacity = 0.5, fillOpacity = 0 ) %>% # Add the `selection` point layer with labels addCircleMarkers( data = selection, radius = ~`Households with Dependent Children` / 2.5, color = ~colorFactor(palette = "RdYlBu", domain = selection$`IDACI Decile`)(`IDACI Decile`), fillColor = ~colorFactor(palette = "RdYlBu", domain = selection$`IDACI Decile`)(`IDACI Decile`), fillOpacity = 0.5, opacity = 0.2, stroke = FALSE, label = ~paste( Postcode, "Households Dependent Children:", `Households with Dependent Children`, "Total Households:", `Total Households`, "IDACI Decile:", `IDACI Decile` ) # Tooltips for mouseover interactivity ) %>% # Add a legend for IDACI Decile addLegend( pal = colorFactor(palette = "RdYlBu", domain = selection$`IDACI Decile`), values = selection$`IDACI Decile`, title = "IDACI Decile", opacity = 1, position = "bottomright" ) ``` ```{r} library(sf) # For spatial data manipulation library(dplyr) # For data manipulation library(ggplot2) # For plotting library(RColorBrewer) # For color palettes # 1. Perform spatial join (use your existing data) joined_data_city <- st_join(bn_postcodes_pop1, optionZ, join = st_within) # 2. Aggregate the data for the entire city aggregated_data_city <- joined_data_city %>% group_by(idaci_decile) %>% summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>% ungroup() # 3. Create a bar graph for the entire city ggplot(aggregated_data_city, aes(x = factor(idaci_decile), y = total_dep_ch, fill = idaci_decile)) + geom_bar(stat = "identity") + scale_fill_gradientn( colors = brewer.pal(11, "RdYlBu"), limits = c(0, 10), # Ensure 0 maps to red and 10 maps to blue breaks = 0:10, # Optional: discrete legend ticks for each decile name = "IDACI Decile" ) + labs( title = "Count of Households with Dependent Children by IDACI Decile - Brighton and Hove", x = "IDACI Decile", y = "Count of Households with Dependent Children" ) + theme_minimal() + theme(legend.position = "bottom") # Adjust the legend position ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} library(sf) # For spatial data manipulation library(dplyr) # For data manipulation library(ggplot2) # For plotting # 1. Perform spatial join joined_data1 <- st_join(bn_postcodes_pop1, optionZ, join = st_within) # 2. Aggregate the data aggregated_data1 <- joined_data1 %>% group_by(id, catchment, idaci_decile) %>% summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>% ungroup() # Add the required library for palette library(RColorBrewer) # 3. Create a faceted bar graph with custom colors ggplot(aggregated_data1, aes(x = factor(idaci_decile), y = total_dep_ch, fill = idaci_decile)) + geom_bar(stat = "identity") + scale_fill_gradientn( colors = brewer.pal(11, "RdYlBu"), limits = c(0, 10), # Ensure 0 maps to red and 10 maps to blue breaks = 0:10, # Optional: discrete legend ticks for each decile name = "IDACI Decile" ) + facet_wrap(~ catchment, scales = "fixed") + labs( title = "Count of Households with Dependent Children by IDACI Decile - New Catchments", x = "IDACI Decile", y = "Count Households with Dependent Children" ) + theme_minimal() + theme(legend.position = "bottom") # Adjust the legend position ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} library(sf) # For spatial data manipulation library(dplyr) # For data manipulation library(ggplot2) # For plotting # 1. Perform spatial join joined_data2 <- st_join(bn_postcodes_pop1, option0, join = st_within) # 2. Aggregate the data aggregated_data2 <- joined_data2 %>% group_by(id, catchment, idaci_decile) %>% summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>% ungroup() # Add the required library for palette library(RColorBrewer) # 3. Create a faceted bar graph with custom colors ggplot(aggregated_data2, aes(x = factor(idaci_decile), y = total_dep_ch, fill = idaci_decile)) + geom_bar(stat = "identity") + scale_fill_gradientn( colors = brewer.pal(11, "RdYlBu"), limits = c(0, 10), # Ensure 0 maps to red and 10 maps to blue breaks = 0:10, # Optional: discrete legend ticks for each decile name = "IDACI Decile" ) + facet_wrap(~ catchment, scales = "fixed") + labs( title = "Count of Households with Dependent Children by IDACI Decile - Current Catchments", x = "IDACI Decile", y = "Count Households with Dependent Children" ) + theme_minimal() + theme(legend.position = "bottom") # Adjust the legend position ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} library(dplyr) # Step 1: Drop geometry and create the unique ID column aggregated_data1 <- aggregated_data1 %>% st_set_geometry(NULL) %>% # Drop geometry mutate(unique_id = paste(catchment, idaci_decile, sep = "_")) aggregated_data2 <- aggregated_data2 %>% st_set_geometry(NULL) %>% # Drop geometry mutate(unique_id = paste(catchment, idaci_decile, sep = "_")) # Step 2: Perform a full join on the unique ID column aggregated_data_diff <- aggregated_data1 %>% full_join(aggregated_data2, by = "unique_id", suffix = c("_1", "_2")) %>% mutate( total_dep_ch_1 = replace_na(total_dep_ch_1, 0), # Set missing values to zero total_dep_ch_2 = replace_na(total_dep_ch_2, 0), total_dep_ch_diff = total_dep_ch_1 - total_dep_ch_2 # Calculate the difference ) %>% mutate( catchment = if_else(!is.na(catchment_1), catchment_1, catchment_2), # Resolve `catchment` idaci_decile = if_else(!is.na(idaci_decile_1), idaci_decile_1, idaci_decile_2) # Resolve `idaci_decile` ) %>% select( catchment, idaci_decile, total_dep_ch_diff ) ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} # # 3. Create a faceted bar graph with custom colors # ggplot(aggregated_data_diff, aes(x = factor(idaci_decile), y = total_dep_ch_diff, fill = idaci_decile)) + # geom_bar(stat = "identity") + # scale_fill_gradientn( # colors = brewer.pal(11, "RdYlBu"), # limits = c(0, 10), # Ensure 0 maps to red and 10 maps to blue # breaks = 0:10, # Optional: discrete legend ticks for each decile # name = "IDACI Decile" # ) + # facet_wrap(~ catchment, scales = "fixed") + # labs( # title = "Change in Households with Dependent Children by IDACI Decile - Current Catchments", # x = "IDACI Decile", # y = "Difference in Count Households with Dependent Children" # ) + # theme_minimal() + # theme(legend.position = "bottom") # Adjust the legend position ``` ```{r osm-download, warning=FALSE, message=FALSE, echo=FALSE} # #get some output areas # brighton_oa <- st_read(here("data","oa_brighton.geojson")) # #get some roads # roads_w_sussex = oe_get(c("West Sussex"), stringsAsFactors = FALSE, quiet = TRUE) # roads_e_sussex = oe_get(c("East Sussex"), stringsAsFactors = FALSE, quiet = TRUE) # #names(roads_w_sussex) # #summary(roads_w_sussex) # # ht = c("motorway", "trunk", "primary", "secondary", "tertiary", "residential", "unclassified", "service") # osm_w_sussex_maj_roads = roads_w_sussex[roads_w_sussex$highway %in% ht, ] # #class(osm_w_sussex_maj_roads) # osm_e_sussex_maj_roads = roads_e_sussex[roads_e_sussex$highway %in% ht, ] # # osm_sussex_maj_roads <- rbind(osm_w_sussex_maj_roads, osm_e_sussex_maj_roads) %>% st_transform(27700) # #summary(osm_sussex_maj_roads) # #summary(brighton_oa) # # #osm_roads_buffered <- st_buffer(osm_sussex_maj_roads, dist = 100) # brighton_oa_buffer <- st_buffer(brighton_oa, dist = 500) # brighton_blob <- st_union(brighton_oa_buffer) # # #tm_shape(brighton_blob) + # # tm_polygons(alpha = 0.5) # # brighton_roads <- st_intersection(osm_sussex_maj_roads, brighton_blob) # # tmap_mode("view") # tm_shape(brighton_roads) + # tm_lines() # # #going to need to write this out as a shapefile and then use JOSM to convert to a PBF file - https://josm.openstreetmap.de/ # #st_write(brighton_roads, "data/brighton_roads.shp", append=FALSE) ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} library(sf) library(raster) library(rayshader) library(ggplot2) bn_postcodes_pop1 <- st_transform(bn_postcodes_pop1, 4326) # Extract the coordinates and the 'Households with Dependent Children' variable coordinates <- st_coordinates(bn_postcodes_pop1) dependent_children <- bn_postcodes_pop1$pcd_dep_ch_hh_count_round idaci <- bn_postcodes_pop1$idaci_decile # Create a data frame for ggplot rayshader_data <- data.frame( x = coordinates[, 1], y = coordinates[, 2], z = dependent_children, w = idaci# Use the number of households with children as the height ) rayshader_data$z <- scales::rescale(rayshader_data$z, to = c(1, 100)) # Normalize heights rayshader_data$w <- scales::rescale(rayshader_data$w, to = c(0, 10)) # Normalize colors ``` ## Rayshader Plot ```{r} # library(basemaps) # # # Define bounding box for your map area # # Calculate bounding box explicitly # bbox <- st_bbox(bn_postcodes_pop1, crs = 4326) # # # set defaults for the basemap # set_defaults(ext = bbox, map_service = "esri", map_type = "world_light_gray_base") # # # Overlay with your height plot or main plot # osm_basemap <- ggplot() + # basemap_gglayer(ext = bbox) + # Pass bbox explicitly here # scale_fill_identity() + # coord_sf() # # # osm_basemap # # # Render the rayshader plot with the basemap # osm_basemap %>% # plot_gg( # ggobj_height = height_plot, # multicore = TRUE, # width = 5, # height = 5, # scale = 50, # Reduce height scaling # windowsize = c(800, 800), # shadow_intensity = 0.5 # Optional: adjust shadow for realism # ) ``` ```{r} height_plot <- ggplot(rayshader_data, aes(x = x, y = y, colour = z)) + geom_point(size = 0.05) + # Raster fills for continuous height mapping scale_colour_gradientn( colors = brewer.pal(9, "YlGnBu"), # Adjust as needed limits = c(min(rayshader_data$z, na.rm = TRUE), max(rayshader_data$z, na.rm = TRUE)), name = "Height" ) + theme_minimal(base_size = 10) + labs( title = "Households with Dependent Children", colour = "IDACI Decile" ) + theme( panel.background = element_rect(fill = "white", color = NA), # White background plot.background = element_rect(fill = "white", color = NA), # White plot background legend.background = element_rect(fill = "white", color = NA), legend.position = "none" # Disable all legends# White legend background )# Minimal theme for height mapping #height_plot main_plot <- ggplot(rayshader_data, aes(x = x, y = y, colour = w)) + geom_point(size = 0.05) + scale_colour_gradientn( colors = brewer.pal(11, "RdYlBu"), limits = c(0, 10), name = "IDACI Decile" ) + theme_minimal(base_size = 10) + labs( title = "IDACI Decile", colour = "IDACI Decile" ) + theme( panel.background = element_rect(fill = "white", color = NA), # White background plot.background = element_rect(fill = "white", color = NA), # White plot background legend.background = element_rect(fill = "white", color = NA), legend.position = "none" # Disable all legends# White legend background ) #main_plot main_plot %>% plot_gg( ggobj_height = height_plot, # Specify height mapping plot multicore = TRUE, width = 5, height = 5, scale = 50, # Adjust scale for height exaggeration windowsize = c(1200, 1200), sunangle=225 ) # Add 3D rendering render_camera(theta = 10, phi = 30, zoom = 0.4) rgl::rglwidget() ``` ```{r eval=FALSE} # # bn_bbox = st_bbox(bn_postcodes_pop1) # # bn_roads <- osmdata::opq(bn_bbox) %>% # osmdata::add_osm_feature("highway") %>% # osmdata::osmdata_sf() # # bn_roads <- bn_roads$osm_lines # # ggplot(bn_roads, aes(color = highway)) + # geom_sf() + # theme(legend.position = "none") + # labs(title = "Open Street Map `highway` attribute in Brighton") # # rayshader_data_sf <- rayshader_data %>% # st_as_sf(coords = c("x", "y"), crs = 4326) # # ggplot(rayshader_data_sf,aes(color=z)) + # geom_sf() + # theme(legend.position = "none") # # map_plot <- ggplot() + # geom_sf(data = rayshader_data_sf, aes(color = z)) + # geom_sf(data = bn_roads, aes(fill = highway), size = 0.5) + # theme_minimal() + # theme( # panel.background = element_rect(fill = "white", color = NA), # White background # plot.background = element_rect(fill = "white", color = NA), # White plot background # legend.background = element_rect(fill = "white", color = NA), # legend.position = "none") # # map_plot %>% # plot_gg( # ggobj_height = height_plot, # Specify height mapping plot # multicore = TRUE, # width = 5, # height = 5, # scale = 50, # Adjust scale for height exaggeration # windowsize = c(1200, 1200), # sunangle=225 # )# Disable all legends# White legend ``` ```{r} library(sf) library(dplyr) library(ggplot2) library(RColorBrewer) # Perform spatial join joined_data_city <- st_join(bn_postcodes_pop1, option0, join = st_within) # Drop spatial attributes before further operations joined_data_city <- st_drop_geometry(joined_data_city) # Aggregate for city-level proportions city_proportions <- joined_data_city %>% group_by(idaci_decile) %>% summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>% mutate(city_prop = total_dep_ch / sum(total_dep_ch)) %>% ungroup() # Aggregate for catchment-level proportions catchment_proportions <- joined_data_city %>% group_by(catchment, idaci_decile) %>% summarise(total_dep_ch = sum(pcd_dep_ch_hh_count_round, na.rm = TRUE)) %>% mutate(catchment_prop = total_dep_ch / sum(total_dep_ch)) %>% ungroup() # Merge city and catchment proportions # comparison_data <- catchment_proportions %>% # left_join(city_proportions %>% select(idaci_decile, city_prop), by = "idaci_decile") # # Plot proportions comparison # ggplot(comparison_data, aes(x = factor(idaci_decile))) + # geom_bar(aes(y = catchment_prop, fill = "Catchment"), stat = "identity", position = "dodge") + # geom_line(aes(y = city_prop, group = 1, color = "City"), size = 1.2) + # scale_fill_manual(values = c("Catchment" = "steelblue")) + # scale_color_manual(values = c("City" = "darkred")) + # labs( # title = "Proportions of Households by IDACI Decile", # subtitle = "Comparison of Catchment Proportions with City-Level Proportions", # x = "IDACI Decile", # y = "Proportion of Households", # fill = "Legend", # color = "Legend" # ) + # facet_wrap(~ catchment, scales = "free_y") + # theme_minimal() + # theme( # legend.position = "bottom", # axis.text.x = element_text(angle = 45, hjust = 1) # ) ``` ```{r} # Initialize an empty data frame to hold the results comparison_data <- data.frame() # Get the unique IDACI deciles from city_proportions unique_deciles <- unique(city_proportions$idaci_decile) # Loop through each catchment in catchment_proportions for (catchment in unique(catchment_proportions$catchment)) { # Subset data for the current catchment catchment_data <- catchment_proportions %>% filter(catchment == !!catchment) # Loop through each IDACI decile for (decile in unique_deciles) { # Check if the current decile exists in the catchment data if (decile %in% catchment_data$idaci_decile) { # If it exists, get the proportion for that decile catchment_prop <- catchment_data$catchment_prop[catchment_data$idaci_decile == decile] } else { # If it does not exist, assign a proportion of 0 catchment_prop <- 0 } # Add the current row to the comparison_data comparison_data <- rbind( comparison_data, data.frame( catchment = catchment, idaci_decile = decile, catchment_prop = catchment_prop ) ) } } comparison_data <- comparison_data %>% left_join(city_proportions, by = "idaci_decile") # Add city-level proportions to the data # comparison_data <- comparison_data %>% # left_join(city_proportions %>% select(idaci_decile, city_prop), by = "idaci_decile") comparison_data <- comparison_data %>% mutate( diff = catchment_prop - city_prop ) ``` ```{r} ggplot(comparison_data, aes(x = factor(idaci_decile), y = diff, fill = idaci_decile)) + geom_bar(stat = "identity") + scale_fill_gradientn( colors = brewer.pal(11, "RdYlBu"), limits = c(0, 10), # Ensure 0 maps to red and 10 maps to blue breaks = 0:10, # Optional: discrete legend ticks for each decile name = "IDACI Decile" ) + facet_wrap(~ catchment, scales = "fixed") + labs( title = "Deviation from City Proportion of Deprivation - Current Catchments", x = "IDACI Decile", y = "Proportion of Households in Catchment" ) + theme_minimal() + theme(legend.position = "bottom") ``` ```{r echo=FALSE, message=FALSE, warning=FALSE, error=FALSE} library(leaflet) optionZ <- st_transform(optionZ, 4326) option0 <- st_transform(option0, 4326) selection <- st_transform(selection, 4326) bn_postcodes_pop1 <- st_transform(bn_postcodes_pop1, 4326) # Create color palette color_palette <- colorNumeric( palette = "Blues", domain = selection$`Total non-white children 0-18` ) color_palette(runif(10, 60, 100)) leaflet() %>% addProviderTiles(providers$CartoDB.Positron) %>% # Add the `optionZ` polygon layer with a red border addPolygons( data = optionZ, color = "red", opacity = 0.5, fillOpacity = 0 ) %>% # Add the `option0` polygon layer addPolygons( data = option0, color = NULL, opacity = 0.5, fillOpacity = 0 ) %>% # Add the `selection` point layer with color mapping for `Total non-white children 0-18` addCircleMarkers( data = selection, radius = ~`Total non-white children 0-18` / 2.5, # Scale radius fillColor = ~color_palette(`Total non-white children 0-18`), # Apply color fillOpacity = 0.5, stroke = FALSE, label = ~paste( Postcode, "Total non-white children 0-18:", `Total non-white children 0-18`, "Total Households:", `Total Households`, "IDACI Decile:", `IDACI Decile` ) ) %>% # Add a legend for the `Total non-white children 0-18` column addLegend( pal = color_palette, values = selection$`Total non-white children 0-18`, title = "Total non-white children 0-18", opacity = 1, position = "bottomright" ) ```