diff --git a/.Rbuildignore b/.Rbuildignore index f2ca8ef..5539a12 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -2,3 +2,4 @@ ^\.Rproj\.user$ ^README\.Rmd$ ^\.github$ +^LICENSE\.md$ diff --git a/.Rhistory b/.Rhistory index 46f0bdf..af94a48 100644 --- a/.Rhistory +++ b/.Rhistory @@ -1,327 +1,6 @@ -results, -function(x) x$`Correlation Matrix`[1, 2] -), -server = if (type == "split") names(results) else rep("combine", length(results)), -row.names = NULL -) -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -source("dsExampleClient/R/ds.groupcor.R") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -results_data %>% head() -if(type == 'split' && length(names_region) != length(stdnames)){ -names_region <- rep(names_region, times = length(stdnames)) -} -shape_list <- lapply(names_region, function(region) { -boundr::bounds( -"lsoa", -within_level = "lad", -within_names = region, -lookup_year = 2011, -opts = boundr::boundr_options(resolution = "BFC") -) |> -dplyr::select(lsoa11cd, geometry) -}) -head(shape_list[[1]]) -shape_list[[1]] %>% head() -shape_list[[2]] %>% head() -shape_list[[3]] %>% head() -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -shape_list %>% names() -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -names(shape_list) -shape_list[[1]] -shape_list[[1]] %>% head() -shape_list[[2]] %>% head() -stdnames -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -head(results_data) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -correlation_data %>% head() -correlation_data %>% tail() -shape_list[[1]] %>% head() -server_name -shape_list[[server_name]] %>% head() -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -head(results_data) -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -head(results_data) -levels(factor(results_data$server)) -ggplot2::ggplot(results_data) + -ggplot2::geom_sf(ggplot2::aes(fill = correlation.coef)) + -ggplot2::scale_fill_gradientn( -colours = grDevices::colorRampPalette(c( -"#440154", -"#414487", -"#2A788E", -"#22A884", -"#7AD151", -"#FDE725" -))(100), -na.value = "grey90" -) + -ggplot2::labs(fill = "Correlation Coefficient") + -ggplot2::facet_wrap(~server) + -ggplot2::theme_minimal() -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "split") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -shape_list[[1]] -shape_list[[2]] -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/groupcorDS.R") -source("dsExample/R/groupcorDS.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -names(result) -source("dsExample/R/groupcorDS.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExample/R/groupcorDS.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -result[["Nfilter.tab"]] -names(result) %>% tail() -source("dsExample/R/groupcorDS.R") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -names(output) -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -!is.null(common_lsoas) -is.null(common_lsoas) -length(common_lsoas) -length(unique(common_lsoas)) -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -pwd -getwd() -pdf("Personal/images/combined_correlation_map.pdf", width = 8, height = 6) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -png("Personal/images/combined_correlation_map.png", width = 8, height = 6) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -png("Personal/images/combined_correlation_map.png", width = 8, height = 6) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -is.null(common_lsoas) || length(common_lsoas) < output[[1]][["Nfilter.tab"]] -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -reference <- output[[1]][[lsoa]] -reference -combined.sums.of.products <- matrix( -0, -nrow = nrow(reference[[1]]), -ncol = ncol(reference[[1]]), -dimnames = dimnames(reference[[1]]) -) -combined.sums.of.products -combined.sums.of.products <- matrix( -0, -nrow = nrow(reference[[1]]), -ncol = ncol(reference[[1]]), -dimnames = dimnames(reference[[1]]) -) -combined.sums <- matrix( -0, -nrow = nrow(reference[[2]]), -ncol = ncol(reference[[2]]), -dimnames = dimnames(reference[[2]]) -) -combined.complete.cases <- matrix( -0, -nrow = nrow(reference[[3]]), -ncol = ncol(reference[[3]]), -dimnames = dimnames(reference[[3]]) -) -combined.missing.cases.vector <- matrix( -0, -nrow = nrow(reference[[4]][[1]]), -ncol = ncol(reference[[4]][[1]]), -dimnames = dimnames(reference[[4]][[1]]) -) -combined.missing.cases.matrix <- matrix( -0, -nrow = nrow(reference[[4]][[2]]), -ncol = ncol(reference[[4]][[2]]), -dimnames = dimnames(reference[[4]][[2]]) -) -combined.sums.of.squares <- matrix( -0, -nrow = nrow(reference[[5]]), -ncol = ncol(reference[[5]]), -dimnames = dimnames(reference[[5]]) -) -reference[[4]][[2]] -dimnames(reference[[4]][[2]]) -nrow(reference[[4]][[2]]) -ncol(reference[[4]][[2]]) -reference[[6]] -reference[[5]] -tail(common_lsoas) -names(output[[1]]) %>% tail() -common_lsoas <- sort(Reduce(intersect, lapply(head(output, -1), names))) -tail(common_lsoas) -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -png("Personal/images/combined_correlation_map.png", width = 8, height = 6) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -dev.off() -source("dsExampleClient/R/ds.groupcor.R") -png("Personal/images/combined_correlation_map.png", width = 8, height = 6) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -dev.off() -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -png("Personal/images/combined_correlation_map.png", width = 12, height = 10) -ds.groupcor("D$has_asthma", "D$distance_local_greenspace", "D$lsoa11cd", datasources = conns, type = "combine") -dev.off() -png( -filename = "Personal/images/combined_correlation_map.png", -width = 12, -height = 10, -units = "in", -res = 300 -) -p <- ds.groupcor( -"D$has_asthma", -"D$distance_local_greenspace", -"D$lsoa11cd", -datasources = conns, -type = "combine" -) -print(p) -dev.off() -png( -filename = "Personal/images/split_correlation_map.png", -width = 12, -height = 10, -units = "in", -res = 300 -) -p <- ds.groupcor( -"D$has_asthma", -"D$distance_local_greenspace", -"D$lsoa11cd", -datasources = conns, -type = "split" -) -source("dsExampleClient/R/ds.groupcor.R") -source("dsExampleClient/R/ds.groupcor.R") -png( -filename = "Personal/images/split_correlation_map.png", -width = 12, -height = 10, -units = "in", -res = 300 -) -p <- ds.groupcor( -"D$has_asthma", -"D$distance_local_greenspace", -"D$lsoa11cd", -datasources = conns, -type = "split" -) -source("dsExampleClient/R/ds.groupcor.R") -png( -filename = "Personal/images/split_correlation_map.png", -width = 12, -height = 10, -units = "in", -res = 300 -) -p <- ds.groupcor( -"D$has_asthma", -"D$distance_local_greenspace", -"D$lsoa11cd", -datasources = conns, -type = "split" -) -source("dsExampleClient/R/ds.groupcor.R") -ds.groupcor( -"D$has_asthma", -"D$distance_local_greenspace", -"D$lsoa11cd", -datasources = conns, -type = "split" -) -length(output$server_1) -output[[1]][["Nfilter.tab"]] -nfilter.tab <- output[[1]][["Nfilter.tab"]] -output <- lapply(output, function(x) { -x[-length(x)] -}) -length(output) -length(output[[1]]) -length(output[[2]]) -source("dsExampleClient/R/ds.groupcor.R") -png( -filename = "Personal/images/split_correlation_map.png", -width = 12, -height = 10, -units = "in", -res = 300 -) -p <- ds.groupcor( -"D$has_asthma", -"D$distance_local_greenspace", -"D$lsoa11cd", -datasources = conns, -type = "split" -) -print(p) -dev.off() -getwd() -dat <- dd -getwd() -setwd("dsGeospatial/") -usethis::create_package(".", open = FALSE) -usethis::use_roxygen_md() -usethis::use_readme_rmd() -usethis::use_news_md() -usethis::use_vignette() -usethis::use_vignette("getting_started") -usethis::use_package_doc() -usethis::use_spell_check() -getwd() -library(devtools) -devtools::document() -devtools::document() -?dsGeospatial -rm(list = ls()) -knitr::opts_chunk$set(echo = TRUE) -# retrieve LSOA boundaries for Cheshire and Merseyside -# to reduce coverage, remove names from the `within_names` +library(opalr) +library(readr) +library(magrittr) cm_lsoa_sf <- boundr::bounds( "lsoa", within_level = "lad", @@ -336,177 +15,309 @@ within_names = c( "Warrington", "Wirral" ), -lookup_year = 2011, # can be either 2011 or 2021 +lookup_year = 2021, +# can be either 2011 or 2021 opts = boundr::boundr_options(resolution = "BFC") ) +# Get Liverpool LSOAs liverpool_sf <- cm_lsoa_sf |> -dplyr::filter(lad22nm %in% "Liverpool") -uprn_gs_sf <- "https://pldr.org/download/2k6r3/n12/UPRN_2_1_greenspace_distances_with_coords.csv" |> +dplyr::filter(lad25nm == "Liverpool") +satellite_sf <- "https://pldr.org/download/2o64n/jqv/UPRN_5_2_satellite_measures_cm_multiple_years_with_coords.csv" |> readr::read_csv() |> -# convert to spatial object +sf::st_as_sf( +coords = c("longitude", "latitude"), +crs = 4326, +remove = FALSE +) +# Match CRS with Liverpool boundary +satellite_sf <- st_transform(satellite_sf, st_crs(liverpool_sf)) +library(sf) +satellite_sf <- st_transform(satellite_sf, st_crs(liverpool_sf)) +satellite_liverpool_sf <- satellite_sf[liverpool_sf, ] +# convert sf to dataframe +satellite_liverpool <- satellite_liverpool_sf |> +st_drop_geometry() +ggplot() + +geom_sf(data = satellite_liverpool_sf, +fill = "lightgrey", +colour = "white", +linewidth = 0.1) + +geom_sf(data = pts_sf, +colour = "red", +size = 0.3) + +theme_minimal() +library(ggplot2) +ggplot() + +geom_sf(data = satellite_liverpool_sf, +fill = "lightgrey", +colour = "white", +linewidth = 0.1) + +geom_sf(data = pts_sf, +colour = "red", +size = 0.3) + +theme_minimal() +ggplot() + +geom_sf(data = cm_lsoa_sf, +fill = "lightgrey", +colour = "white", +linewidth = 0.1) + +geom_sf(data = satellite_liverpool_sf, +colour = "red", +size = 0.3) + +theme_minimal() +o <- opalr::opal.login(username = "administrator", +password = "password", +url = "https://opal-demo.obiba.org") +?opal.project_create +opal.project_create(o, "UPRN", database = NULL, title = NULL) +opal.table_save( +o, +satellite_liverpool, +"UPRN", +"satellite_metrics", +overwrite = TRUE, +force = TRUE, +id.name = "UPRN" +) +opal.project_create(o, "UPRN", database = "DataSHIELD", title = NULL) +opal.project_delete(o, "UPRN") +opal.project_create(o, "UPRN", database = "DataSHIELD", title = NULL) +opal.project_create(o, "UPRN", database = "mongodb", title = NULL) +opal.table_save( +o, +satellite_liverpool, +"UPRN", +"satellite_metrics", +overwrite = TRUE, +force = TRUE, +id.name = "UPRN" +) +uprn_dist_gs <- "https://pldr.org/download/2k6r3/n12/UPRN_2_1_greenspace_distances_with_coords.csv" |> +readr::read_csv() +uprn_dist_bs <- "https://pldr.org/download/2g6k0/ksp/UPRN_4_1_bluespace_distances_with_coords.csv" |> +readr::read_csv() +uprn_evi <- "https://pldr.org/download/2o64n/0l9/UPRN_5_2_EVI_multiple_years_recategorised_with_coords.csv"|> +readr::read_csv() +uprn_330300 <- "https://pldr.org/download/2z71q/1me/UPRN_6_1_3_30_300_measures_cm_with_coords.csv" |> +readr::read_csv() |> +dplyr::select(UPRN, composite_3_30_300) +head(uprn_dist_gs) +head(satellite_sf) +uprn_dist_gs |> +dplyr::select(UPRN, distance_local_greenspace, latitude, longitude) |> head() +uprn_dist_gs_v2 <- uprn_dist_gs |> +dplyr::select(UPRN, distance_local_greenspace, latitude, longitude) |> +dplyr::mutate( +distance_local_greenspace = round(distance_local_greenspace / 5) * 5, # round to nearest 5th +distance_local_greenspace = ifelse(distance_local_greenspace > 1500, 1500, distance_local_greenspace), +distance_local_greenspace_100cat = distance_local_greenspace |> +cut(breaks = c(-Inf, seq(0, 1500, 100))) +) +uprn_dist_gs_v2 <- uprn_dist_gs |> +dplyr::select(UPRN, distance_local_greenspace, latitude, longitude) |> +dplyr::mutate( +distance_local_greenspace = round(distance_local_greenspace / 5) * 5, # round to nearest 5th +distance_local_greenspace = ifelse(distance_local_greenspace > 1500, 1500, distance_local_greenspace), +distance_local_greenspace_100cat = distance_local_greenspace |> +cut(breaks = c(-Inf, seq(0, 1500, 100))) +) +head(uprn_dist_gs_v2) +levels(factor(uprn_dist_gs_v2$distance_local_greenspace_100cat)) +uprn_dist_gs_v2[which(uprn_dist_gs_v2$distance_local_greenspace == "(-Inf,0]"),] +uprn_dist_gs_v2[which(uprn_dist_gs_v2$distance_local_greenspace == "(600,700]"),] +uprn_dist_gs_v2[which(uprn_dist_gs_v2$distance_local_greenspace == "(100,200]"),] +uprn_dist_gs_v2[which(uprn_dist_gs_v2$distance_local_greenspace_100cat == "(100,200]"),] +uprn_dist_gs_v2[which(uprn_dist_gs_v2$distance_local_greenspace_100cat == "(-Inf,0]"),] +uprn_dist_bs_v2 <- uprn_dist_bs |> +dplyr::filter(!is.na(distance_any_bluespace)) |> +dplyr::mutate( +dist_bs_100cat = distance_any_bluespace |> +cut(breaks = c(-Inf, seq(0, 4200, 100), Inf)) |> +as.factor() +) |> +dplyr::arrange(UPRN) +uprn_dist_bs_sf <- uprn_dist_bs_v2 |> sf::st_as_sf(coords = c("longitude", "latitude"), crs = 4326) |> -# add LSOA details +# attach LSOA information sf::st_join(cm_lsoa_sf) -qof_asthma_tbl <- "https://pldr.org/download/e6nzv/ng1/QOF_4_03_Asthma_LSOA.csv" |> +head(uprn_evi) +head(uprn_330300) +uprn_dist_bs |> +dplyr::filter(!is.na(distance_any_bluespace)) |> head() +dist_bs_evi_330300_dist_gs_v2 %>% head() +uprn_dist_gs <- "https://pldr.org/download/2k6r3/n12/UPRN_2_1_greenspace_distances_with_coords.csv" |> +readr::read_csv() +uprn_dist_bs <- "https://pldr.org/download/2g6k0/ksp/UPRN_4_1_bluespace_distances_with_coords.csv" |> +readr::read_csv() +uprn_evi <- "https://pldr.org/download/2o64n/0l9/UPRN_5_2_EVI_multiple_years_recategorised_with_coords.csv"|> +readr::read_csv() +uprn_330300 <- "https://pldr.org/download/2z71q/1me/UPRN_6_1_3_30_300_measures_cm_with_coords.csv" |> readr::read_csv() |> -# filter out CM LSOAs -dplyr::filter(lsoa11 %in% cm_lsoa_sf$lsoa11cd) |> -# filter out last year -dplyr::filter(year == 2024) -# set seed for reproducibility -set.seed(20250603) -qof_asthma_cohort_tbl <- qof_asthma_tbl |> -purrr::pmap(function(lsoa11, den, num, ...) { -# recalculate pop to be a multiple of 2 -pop <- ifelse(den %% 2 != 0, den + 1, den) -# extract UPRNs for this LSOA -uprn_tbl <- uprn_gs_sf |> -sf::st_drop_geometry() |> -dplyr::filter(lsoa11cd == lsoa11) |> +dplyr::select(UPRN, composite_3_30_300) +# re-categorise the distance to local greenspace +uprn_dist_gs_v2 <- uprn_dist_gs |> +dplyr::select(UPRN, distance_local_greenspace, latitude, longitude) |> dplyr::mutate( -dist_gs_dec = dplyr::ntile(distance_local_greenspace, 10) +distance_local_greenspace = round(distance_local_greenspace / 5) * 5, # round to nearest 5th +distance_local_greenspace = ifelse(distance_local_greenspace > 1500, 1500, distance_local_greenspace), +distance_local_greenspace_100cat = distance_local_greenspace |> +cut(breaks = c(-Inf, seq(0, 1500, 100))) +) +# re-categorise the distance to bluespace +uprn_dist_bs_v2 <- uprn_dist_bs |> +dplyr::filter(!is.na(distance_any_bluespace)) |> +dplyr::mutate( +dist_bs_100cat = distance_any_bluespace |> +cut(breaks = c(-Inf, seq(0, 4200, 100), Inf)) |> +as.factor() +) |> +dplyr::arrange(UPRN) +# WRANGLING --- +# create spatial version of the UPRN metric +uprn_dist_bs_sf <- uprn_dist_bs_v2 |> +sf::st_as_sf(coords = c("longitude", "latitude"), crs = 4326) |> +# attach LSOA information +sf::st_join(cm_lsoa_sf) +# COMBINING DATASETS ---- +# inspect that there are not categories in an LSOA for which the counts are smaller than 10 +# combining all 4 datasets +dist_bs_evi_330300_dist_gs <- uprn_dist_bs_sf |> +dplyr::inner_join(uprn_evi, by = "UPRN") |> +dplyr::inner_join(uprn_330300, by = "UPRN") |> +dplyr::inner_join(uprn_dist_gs_v2, by = "UPRN") +distance_bs_100_levels <- +c("<0m", "0_100m", "101_200m", "201_300m", "301_400m","401_500m", "501_600m", +"601_700m", "701_800m", "801_900m", "901_1000m","1km_1.1km", "1.1km_1.2km", +"1.2km_1.3km","1.3km_1.4km", "1.4km_1.5km", "1.5km_1.6km", "1.6km_1.7km", +"1.7km_1.8km", "1.8km_1.9km", "1.9km_2km", "2km_2.1km", "2.1km_2.2km", +"2.2km_2.3km", "2.3km_2.4km", "2.4km_2.5km", "2.5km_2.6km", "2.6km_2.7km", +"2.7km_2.8km", "2.8km_2.9km","2.9km_3km", "3km_3.1km", "3.1km_3.2km", +"3.2km_3.3km", "3.3km_3.4km", "3.4km_3.5km","3.5km_3.6km", "3.6km_3.7km", +"3.7km_3.8km", "3.8km_3.9km", "3.9km_4km","4km_4.1km", "4.1km_4.2km", +">4.2km" +) +head(dist_bs_evi_330300_dist_gs_v2 ) +head(dist_bs_evi_330300_dist_gs) +dist_bs_evi_330300_dist_gs %>% dplyr::select("UPRN","distance_any_bluespace", "dist_bs_100cat", "distance_local_greenspace", "distance_local_greenspace_100cat") +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) ) |> -dplyr::arrange(dist_gs_dec) -# create cohort table -cohort_tbl <- tibble::tibble( -lsoa11cd = lsoa11, -sex = rep(c(0, 1), each = pop / 2), -has_asthma = FALSE, -uprn = sample(uprn_tbl$UPRN, size = pop, replace = TRUE) +dplyr::select(-lad22cd, -distance_any_bluespace, +-latitude.x, -longitude.x,-latitude.y, -longitude.y, +-lad22nm, -distance_local_greenspace) |> +dplyr::relocate(lsoa11cd, lsoa11nm, UPRN, distance_local_greenspace_100cat, .before = 1) +dist_bs_evi_330300_dist_gs <- uprn_dist_bs_sf |> +dplyr::inner_join(uprn_evi, by = "UPRN") |> +dplyr::inner_join(uprn_330300, by = "UPRN") |> +dplyr::inner_join(uprn_dist_gs_v2, by = "UPRN") +distance_bs_100_levels <- +c("<0m", "0_100m", "101_200m", "201_300m", "301_400m","401_500m", "501_600m", +"601_700m", "701_800m", "801_900m", "901_1000m","1km_1.1km", "1.1km_1.2km", +"1.2km_1.3km","1.3km_1.4km", "1.4km_1.5km", "1.5km_1.6km", "1.6km_1.7km", +"1.7km_1.8km", "1.8km_1.9km", "1.9km_2km", "2km_2.1km", "2.1km_2.2km", +"2.2km_2.3km", "2.3km_2.4km", "2.4km_2.5km", "2.5km_2.6km", "2.6km_2.7km", +"2.7km_2.8km", "2.8km_2.9km","2.9km_3km", "3km_3.1km", "3.1km_3.2km", +"3.2km_3.3km", "3.3km_3.4km", "3.4km_3.5km","3.5km_3.6km", "3.6km_3.7km", +"3.7km_3.8km", "3.8km_3.9km", "3.9km_4km","4km_4.1km", "4.1km_4.2km", +">4.2km" +) +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) +) |> +dplyr::select(-lad22cd, -distance_any_bluespace, +-latitude.x, -longitude.x,-latitude.y, -longitude.y, +-lad22nm, -distance_local_greenspace) |> +dplyr::relocate(lsoa11cd, lsoa11nm, UPRN, distance_local_greenspace_100cat, .before = 1) +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) ) -# randomly assigned the outcome to `num` patients -idx <- sample(seq_len(pop), ceiling(num), replace = FALSE) -cohort_tbl$has_asthma[idx] <- TRUE -# randomly re-assigned UPRNs with higher distances to patients with outcome -## create vector of probabilities, so further distances are favoured -p <- c(0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.1, 0.2, 0.2, 0.2) -idx_uprn_tile <- sample(1:10, ceiling(num), replace = TRUE, prob = p) -## subset UPRNs with match tile as per `idx_uprn_tile` -cohort_tbl$uprn[idx] <- purrr::map(idx_uprn_tile, function(x) { -# subset based on tile decile -aux <- dplyr::filter(uprn_tbl, dist_gs_dec == x) -# extract UPRN -aux$UPRN[sample(seq_len(nrow(aux)), 1)] -}) |> -purrr::list_c() -# return final cohort table -return(cohort_tbl) -}) |> -purrr::list_c() -qof_asthma_cohort_dist_gs_sf <- qof_asthma_cohort_tbl |> -dplyr::left_join(uprn_gs_sf, by = c("uprn" = "UPRN", "lsoa11cd")) |> -dplyr::arrange(lsoa11cd, uprn) |> -sf::st_as_sf() -library(dsGeospatial) -knitr::opts_chunk$set(echo = TRUE) -library(DSLite) -library(devtools) +colnames(dist_bs_evi_330300_dist_gs_v2) +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) +) |> +dplyr::select(-lad25cd, -distance_any_bluespace, +-latitude.x, -longitude.x,-latitude.y, -longitude.y, +-lad25nm, -distance_local_greenspace) |> +dplyr::relocate(lsoa11cd, lsoa11nm, UPRN, distance_local_greenspace_100cat, .before = 1) +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) +) |> +dplyr::select(-lad25cd, -distance_any_bluespace, +-latitude.x, -longitude.x,-latitude.y, -longitude.y, +-lad25nm, -distance_local_greenspace) |> +dplyr::relocate(lsoa21cd, lsoa21nm, UPRN, distance_local_greenspace_100cat, .before = 1) +head(dist_bs_evi_330300_dist_gs_v2) +?dplyr::relocate +df <- tibble(a = 1, b = 1, c = 1, d = "a", e = "a", f = "a") +library(dplyr) +df <- tibble(a = 1, b = 1, c = 1, d = "a", e = "a", f = "a") +df |> relocate(f) +df |> relocate(a, .after = c) +df |> relocate(f, .before = b) +df |> relocate(f, .before = 1) +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) +) |> +dplyr::select(-lad25cd, -distance_any_bluespace, +-latitude.x, -longitude.x,-latitude.y, -longitude.y, +-lad25nm, -distance_local_greenspace) +head(dist_bs_evi_330300_dist_gs_v2) +dist_bs_evi_330300_dist_gs_v2 <- dist_bs_evi_330300_dist_gs |> +dplyr::mutate( +distance_bluespace_description = +ifelse(is.na(distance_any_bluespace), NA, distance_bs_100_levels[dist_bs_100cat]) +) |> +dplyr::select(-lad25cd, -distance_any_bluespace, +-latitude.x, -longitude.x,-latitude.y, -longitude.y, +-lad25nm, -distance_local_greenspace) |> +dplyr::relocate(lsoa21cd, lsoa21nm, UPRN, distance_local_greenspace_100cat, .before = 1) +colna, +head(dist_bs_evi_330300_dist_gs_v2 ) +levels(factor(dist_bs_evi_330300_dist_gs_v2$distance_bluespace_description)) +levels(factor(dist_bs_evi_330300_dist_gs_v2$distance_bluespace_description))[1:11] +levels(factor(dist_bs_evi_330300_dist_gs_v2$dist_bs_100cat)) +getwd() +setwd("/Users/tkariya/Documents/DataSHIELD/dsGeospatialClient/") +library(opalr) +# creating and uploading table +# creating project +o <- opalr::opal.login(username = "administrator", +password = "password", +url = "http://localhost:8880") +library(DSOpal) library(dsBase) library(dsBaseClient) -library(dsGeospatial) -dat <- qof_asthma_cohort_dist_gs_sf |> -dplyr::select("lsoa11cd", "has_asthma","uprn","distance_local_greenspace") |> -sf::st_drop_geometry() -devtools::load_all("/Users/tkariya/Documents/DataSHIELD/dsExample") -#devtools::load_all("/Users/tkariya/Documents/DataSHIELD/dsExampleClient") -dslite.server1 <- newDSLiteServer( -tables = list( -uprn_asthma_green1 = dat -) -) -dslite.server2 <- newDSLiteServer( -tables = list( -uprn_asthma_green2 = dat -) -) -dslite.server1$config(defaultDSConfiguration(include=c("dsBase", "dsExample"))) -dslite.server1$aggregateMethod("MeanSdGpDS", "MeanSdGpDS") -dslite.server2$config(defaultDSConfiguration(include=c("dsBase", "dsExample"))) -dslite.server2$aggregateMethod("MeanSdGpDS", "MeanSdGpDS") -options( -datashield.privacyControlLevel = "banana", -nfilter.tab = 3, -nfilter.subset = 3, -nfilter.glm = 0.33, -nfilter.string = 80, -nfilter.stringShort = 20, -nfilter.kNN = 3, -nfilter.levels.density = 0.33, -nfilter.levels.max = 40, -nfilter.noise = 0.25, -nfilter.privacy.old = 5 -) builder <- DSI::newDSLoginBuilder() -builder$append( -server = "server_1", -url = "dslite.server1", -table = "uprn_asthma_green1", -driver = "DSLiteDriver" -) -builder$append( -server = "server_2", -url = "dslite.server2", -table = "uprn_asthma_green2", -driver = "DSLiteDriver" -) +builder$append(server = "study1", +url = "http://localhost:8880", +user = "administrator", +password = "password", +table = "UPRN.greenspace_distance_metrics", +profile = "ds-geospatial") +builder$append(server = "study2", +url = "http://localhost:8880", +user = "administrator", +password = "password", +table = "UPRN.greenspace_distance_metrics", +profile = "ds-geospatial") logindata <- builder$build() -conns <- DSI::datashield.login(logins = logindata, assign = FALSE) -datashield.assign.table(conns, "D", logindata) -nfilter.tab <- 3 -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns) -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns, -type='split') -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns) -p <- ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns) -print(p) -source("R/ds.geoheatmapPlot.R") -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns) -source("R/ds.geoheatmapPlot.R") -dev.off() -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns) -library(devtools) -devtools::document() -library(dsGeospatial) -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns) -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns, -type='split') -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns, -type='split', names_region = c("Liverpool", "Sefton")) -ds.geoheatmapPlot(x="D$distance_local_greenspace", y = "D$lsoa11cd", -datasources = conns) -ds.geoheatmapPlot(x="D$has_asthma", y = "D$lsoa11cd", datasources = conns, -type = 'split', names_region = c("Liverpool", "Sefton")) -install.packages("dsGeospatial") -rm(list = ls()) -pak::pak("FederatedMethods/dsGeospatial") -pak::pak("git::git@github.com:FederatedMethods/dsGeospatial.git") -pak::pak("git::git@github.com:FederatedMethods/dsGeospatial.git") -pak::pak( -"git::ssh://git@github.com/FederatedMethods/dsGeospatial.git" -) -pak::pkg_install( -"git::ssh://git@github.com/FederatedMethods/dsGeospatial.git" -) -pak::local_install( -"/Users/tkariya/Documents/DataSHIELD/dsGeospatial/" -) -pak::local_install( -"/Users/tkariya/Documents/DataSHIELD/dsGeospatial/" -) -pak::local_install( -"/Users/tkariya/Documents/DataSHIELD/dsGeospatial/" -) -pak::local_install( -"/Users/tkariya/Documents/DataSHIELD/dsGeospatial/" -) -getwd() -setwd("/Users/tkariya/Documents/DataSHIELD/") -setwd("/Users/tkariya/Documents/DataSHIELD/dsGeospatialClient/") -usethis::create_package() -usethis::create_package(".", open = FALSE) -usethis::use_roxygen_md() -usethis::use_readme_rmd() -getwd() -usethis::use_readme_rmd() -usethis::use_news_md() -usethis::use_vignette("getting-started") +conns <- DSI::datashield.login(logins = logindata, assign = TRUE, symbol = "D") +ds.colnames("D", conns) +library(opalr) +# creating and uploading table +# creating project +o <- opalr::opal.login(username = "administrator", +password = "password", +url = "http://localhost:8880") diff --git a/.github/.gitignore b/.github/.gitignore index cc7b82f..926e968 100644 --- a/.github/.gitignore +++ b/.github/.gitignore @@ -1,3 +1,4 @@ .Rproj.user inst/doc *.html +.Rhistory \ No newline at end of file diff --git a/DESCRIPTION b/DESCRIPTION index 00de032..33601b8 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -9,10 +9,9 @@ Description: Provides privacy-preserving methods for the visualisation and remain within the data-holding organisations, while only disclosure- controlled aggregate results are returned to the analyst. The package is being developed as part of the GROVE project within the DARE UK programme. -License: `use_mit_license()`, `use_gpl3_license()` or friends to pick a - license -Depends:R (>= 4.0.0) -Imports:cli, +License: MIT + file LICENSE +Depends: R (>= 4.0.0) +Imports: cli, methods Encoding: UTF-8 Roxygen: list(markdown = TRUE) diff --git a/LICENSE b/LICENSE new file mode 100644 index 0000000..d8f6df1 --- /dev/null +++ b/LICENSE @@ -0,0 +1,2 @@ +YEAR: 2026 +COPYRIGHT HOLDER: University of Liverpool diff --git a/LICENSE.md b/LICENSE.md new file mode 100644 index 0000000..632fa02 --- /dev/null +++ b/LICENSE.md @@ -0,0 +1,21 @@ +# MIT License + +Copyright (c) 2026 University of Liverpool + +Permission is hereby granted, free of charge, to any person obtaining a copy +of this software and associated documentation files (the "Software"), to deal +in the Software without restriction, including without limitation the rights +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the Software is +furnished to do so, subject to the following conditions: + +The above copyright notice and this permission notice shall be included in all +copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE +SOFTWARE. diff --git a/R/ds.localmoran.R b/R/ds.localmoran.R new file mode 100644 index 0000000..7653466 --- /dev/null +++ b/R/ds.localmoran.R @@ -0,0 +1,138 @@ +#' Calculate Local Moran's I +#' +#' Calculates Local Moran's I for a numeric variable using a spatial weights +#' object across one or more DataSHIELD studies. +#' +#' @param x Character string specifying the name of the numeric vector on the +#' server to be analysed. +#' @param listw Character string specifying the name of the spatial weights +#' object of class \code{listw} on the server. +#' @param type Character string specifying the output format. Currently only +#' \code{"split"} is supported. +#' @param vars Optional variable specification. Currently unused. +#' @param do.checks Logical. If \code{TRUE}, checks that the input objects are +#' defined and have consistent classes across the studies. +#' @param nsim Numeric value specifying the number of simulations used for +#' calculating the Monte Carlo p-value. Default is 199. +#' @param datasources A list of \code{DSConnection} objects. If \code{NULL}, +#' active DataSHIELD connections are used. +#' +#' @return A named list with one element per study. Each study contains a +#' named numeric vector with: +#' \itemize{ +#' \item \code{moran_I}: Moran's I statistic. +#' \item \code{z}: standardised Moran's I statistic. +#' } +#' +#' @examples +#' \dontrun{ +#' DSOpal::Opal() +#' builder <- DSI::newDSLoginBuilder() +#' builder$append( +#' server = "mersey-demo", +#' url = "http://localhost:8880", +#' user = "administrator", +#' password = "password", +#' driver = "OpalDriver", +#' profile = "ds-geospatial" +#' ) +#' +#' logindata <- builder$build() +#' conns <- DSI::datashield.login(logins = logindata) +#' DSI::datashield.assign.resource( +#' conns, +#' symbol = "cells", +#' resource = paste0(project, ".", "cells") +#' ) +#' DSI::datashield.assign.resource( +#' conns, +#' symbol = "W", +#' resource = paste0(project, ".", "W_D_5km") +#' ) +#' ## 2.8 Materialise arbitrary R objects ---- +#' DSI::datashield.assign.resource( +#' conns, +#' symbol = "cells_resource", +#' resource = paste0(project, ".cells") +#' ) +#' +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "cells_object", +#' expr = quote(as.resource.object(cells_resource)) +#' ) +#' ## 2.9 Map specific columns ---- +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "obesity_qof_2024_25", +#' expr = quote(cells_object$obesity_qof_2024_25) +#' ) +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "W_object", +#' expr = quote(as.resource.object(W)) +#' ) +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "W_band800", +#' expr = quote(W_object$band800) +#' ) +#' +#' ds.moransI( +#' x = "obesity_qof_2024_25", +#' listw = "W_band800", +#' nsim = 199, +#' datasources = conns +#' ) +#' +#' +#' @export + +ds.localmoran <- function(x = NULL, listw = NULL, type='split', + vars = NULL, do.checks = TRUE, + nsim = 199, datasources=NULL){ + + # look for DS connections + if(is.null(datasources)){ + datasources <- datashield.connections_find() + } + + # ensure datasources is a list of DSConnection-class + if(!(is.list(datasources) && all(unlist(lapply(datasources, function(d) {methods::is(d,"DSConnection")}))))){ + stop("The 'datasources' were expected to be a list of DSConnection-class objects", call.=FALSE) + } + + if(is.null(x)){ + stop("Please provide the input vector 'x', a numeric!", call.=FALSE) + } + + if(is.null(listw)){ + stop("Please provide the input vector 'listw'!, a list object", call.=FALSE) + } + + if(do.checks){ + + # check if the input objects are defined in all the studies + isDefined(datasources, x) + isDefined(datasources, listw) + + # call the internal function that checks the input object is of the same class in all studies. + typ1 <- checkClass(datasources, x) + typ2 <- checkClass(datasources, listw) + } + + if(!is.numeric(nsim) & length(nsim != 1)){ + stop("nsim should be numeric of length 1") + } + calltext <- paste0("Localmoran(", x, ", ", listw, ")") + + output <- DSI::datashield.aggregate(datasources, as.symbol(calltext)) + + ## Only type split supported now + if(type == "split"){ + return(output) + } + +} + + diff --git a/R/ds.moransI.R b/R/ds.moransI.R new file mode 100644 index 0000000..0306257 --- /dev/null +++ b/R/ds.moransI.R @@ -0,0 +1,140 @@ +#' Calculate Moran's I +#' +#' Calculates Moran's I for a numeric variable using a spatial weights +#' object across one or more DataSHIELD studies. +#' +#' @param x Character string specifying the name of the numeric vector on the +#' server to be analysed. +#' @param listw Character string specifying the name of the spatial weights +#' object of class \code{listw} on the server. +#' @param type Character string specifying the output format. Currently only +#' \code{"split"} is supported. +#' @param vars Optional variable specification. Currently unused. +#' @param do.checks Logical. If \code{TRUE}, checks that the input objects are +#' defined and have consistent classes across the studies. +#' @param nsim Numeric value specifying the number of simulations used for +#' calculating the Monte Carlo p-value. Default is 199. +#' @param datasources A list of \code{DSConnection} objects. If \code{NULL}, +#' active DataSHIELD connections are used. +#' +#' @return A named list with one element per study. Each study contains a +#' named numeric vector with: +#' \itemize{ +#' \item \code{moran_I}: Moran's I statistic. +#' \item \code{z}: standardised Moran's I statistic. +#' \item \code{p_mc}: Monte Carlo p-value. +#' } +#' +#' @examples +#' \dontrun{ +#' DSOpal::Opal() +#' builder <- DSI::newDSLoginBuilder() +#' builder$append( +#' server = "mersey-demo", +#' url = "http://localhost:8880", +#' user = "administrator", +#' password = "password", +#' driver = "OpalDriver", +#' profile = "ds-geospatial" +#' ) +#' +#' logindata <- builder$build() +#' conns <- DSI::datashield.login(logins = logindata) +#' DSI::datashield.assign.resource( +#' conns, +#' symbol = "cells", +#' resource = paste0(project, ".", "cells") +#' ) +#' DSI::datashield.assign.resource( +#' conns, +#' symbol = "W", +#' resource = paste0(project, ".", "W_D_5km") +#' ) +#' ## 2.8 Materialise arbitrary R objects ---- +#' DSI::datashield.assign.resource( +#' conns, +#' symbol = "cells_resource", +#' resource = paste0(project, ".cells") +#' ) +#' +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "cells_object", +#' expr = quote(as.resource.object(cells_resource)) +#' ) +#' ## 2.9 Map specific columns ---- +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "obesity_qof_2024_25", +#' expr = quote(cells_object$obesity_qof_2024_25) +#' ) +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "W_object", +#' expr = quote(as.resource.object(W)) +#' ) +#' DSI::datashield.assign.expr( +#' conns, +#' symbol = "W_band800", +#' expr = quote(W_object$band800) +#' ) +#' +#' ds.moransI( +#' x = "obesity_qof_2024_25", +#' listw = "W_band800", +#' nsim = 199, +#' datasources = conns +#' ) +#' +#' +#' @export + +ds.moransI <- function(x = NULL, listw = NULL, type='split', + vars = NULL, do.checks = TRUE, + nsim = 199, datasources=NULL){ + + # look for DS connections + if(is.null(datasources)){ + datasources <- datashield.connections_find() + } + + # ensure datasources is a list of DSConnection-class + if(!(is.list(datasources) && all(unlist(lapply(datasources, function(d) {methods::is(d,"DSConnection")}))))){ + stop("The 'datasources' were expected to be a list of DSConnection-class objects", call.=FALSE) + } + + if(is.null(x)){ + stop("Please provide the input vector 'x', a numeric!", call.=FALSE) + } + + if(is.null(listw)){ + stop("Please provide the input vector 'listw'!, a list object", call.=FALSE) + } + + if(do.checks){ + + # check if the input objects are defined in all the studies + isDefined(datasources, x) + isDefined(datasources, listw) + + # call the internal function that checks the input object is of the same class in all studies. + typ1 <- checkClass(datasources, x) + typ2 <- checkClass(datasources, listw) + } + + if(!is.numeric(nsim) & length(nsim != 1)){ + stop("nsim should be numeric of length 1") + } + calltext <- paste0("MoransI(", x, ", ", listw, ", ", nsim, ")") + + output <- DSI::datashield.aggregate(datasources, as.symbol(calltext)) + + ## Only type split supported now + if(type == "split"){ + + result <- lapply(output, function(x) unlist(x)) + return(result) + } +} + + diff --git a/R/ds.summarystat.R b/R/ds.summarystat.R index c1e0380..c70df63 100644 --- a/R/ds.summarystat.R +++ b/R/ds.summarystat.R @@ -198,18 +198,71 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, plot_matrix <- list() split <- type == 'split' - if(!is.null(names_region) && all(is.character(names_region))){ + ## names_region check + # check if list and all elements are character + + if (is.list(names_region) && + length(names_region) > 0 && + all(is.character(unlist(names_region)))) { names_region.valid = TRUE - if(!split && length(names_region) > 1){ - warning(paste0("more than one region selected for type = ", type, - ", using only first names_region to return")) - names_region = names_region[1] - } - if(split && length(names_region) != numsources){ - names_region <- rep(names_region, times = numsources) - } + names_region <- lapply(names_region, function(i){ + if("CheshireMercyside" %in% i){ + indx <- which(unlist(i) == "CheshireMercyside") + i <- i[-indx] + i <- c(i, "Cheshire East","Cheshire West and Chester", + "Halton","Knowsley","Liverpool","Sefton", + "St. Helens","Warrington","Wirral") + }else{ + i <- i + } + }) + + # names_region is valid check type combination + + if (split) { + if (numsources == 1) { + warning( + "type 'split', expects length of datasources to be more than 1, + defaulting to type 'combine'" + ) + split = FALSE + } else{ + # check list length = numsources + # check is all entries match when length is >1 + if (length(names_region) != numsources) { + warning( + paste0( + "length of names_region should match length of datasources, + using only first elements in names_region to return" + ) + ) + names_region <- lapply(1:numsources, function(i) { + names_region[[i]] <- names_region[[1]] + }) + names_region <- lapply(names_region, unique) + } else{ + names_region <- lapply(names_region, unique) + } + } + } else if(!split) { + # when combine, list must have length 1 + # ifnot repeat list[[1]] for all numsources + if (length(names_region) != 1) { + warning( + "length of names_region must be 1 for type 'combine', + using only first elements in names_region to return" + ) + names_region = list(names_region[[1]]) + names_region <- lapply(names_region, unique) + + } else{ + names_region = list(names_region[[1]]) + names_region <- lapply(names_region, unique) + } + } + shape_list <- lapply(names_region, function(region) { boundr::bounds( @@ -221,16 +274,20 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, ) |> dplyr::select(lsoa11cd, geometry) }) + q.names <- c("5th_Quantile", "10th_Quantile", "25th_Quantile", "50th_Quantile", + "75th_Quantile", "90th_Quantile", "95th_Quantile") - if(is.null(metric) || !metric %in% c("Mean", "SD","SEM")){ + if(is.null(metric) || !metric %in% c("Mean", "SD","SEM", q.names)){ warning("Invalid metric type, returning available metric table") draw.plot = FALSE names_region.valid =FALSE } else { + if(!split){ metric_map <- list("Mean" = "mean.gp", "SD" = "SD.gp", "SEM" = "SEM.gp") + } else { metric_map <- list("Mean" = "mean.gp.study", "SD" = "SD.gp.study", @@ -239,11 +296,11 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, } }else{ - warning("invalid region names, returning the whole table") + warning("invalid region names/format, returning the whole table") draw.plot = FALSE names_region.valid = FALSE } - + ########## @@ -253,6 +310,7 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, numsources <- length(output) mean.matrix <- NULL sd.matrix <- NULL + qq.matrix <- NULL n.matrix <- NULL Nvalid <- 0 Nmissing <- 0 @@ -261,10 +319,11 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, for(j in 1:numsources){ mean.matrix <- rbind(mean.matrix,as.numeric(unlist(output[[j]][2]))) sd.matrix <- rbind(sd.matrix,as.numeric(unlist(output[[j]][3]))) - n.matrix <- rbind(n.matrix,as.numeric(unlist(output[[j]][4]))) - Nvalid <- Nvalid+as.numeric(unlist(output[[j]][5])) - Nmissing <- Nmissing+as.numeric(unlist(output[[j]][6])) - Ntotal <- Ntotal+as.numeric(unlist(output[[j]][7])) + qq.matrix <- rbind(qq.matrix,as.numeric(unlist(output[[j]][4]))) + n.matrix <- rbind(n.matrix,as.numeric(unlist(output[[j]][5]))) + Nvalid <- Nvalid+as.numeric(unlist(output[[j]][6])) + Nmissing <- Nmissing+as.numeric(unlist(output[[j]][7])) + Ntotal <- Ntotal+as.numeric(unlist(output[[j]][8])) } @@ -274,6 +333,29 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, # Calculate weighted means across studies in each group mean.gp <- (diag(t(mean.matrix)%*%n.matrix))/(t(n.matrix)%*%nsum.vector) + # Calculate weighted quantiles across studies in each group + # Repeat each LSOA sample size for its 7 quantiles + qq.n.matrix <- t( + apply(n.matrix, 1, function(x) rep(x, each = 7)) + ) + + # Weighted quantiles across studies + qq.gp <- (diag(t(qq.matrix) %*% qq.n.matrix)) / + (t(qq.n.matrix) %*% nsum.vector) + + qq.gp <- matrix( + qq.gp, + ncol = 7, + byrow = TRUE + ) + + # Add names + rownames(qq.gp) <- names(output[[1]][[4]]) + + colnames(qq.gp) <- c("5%_gp", "10%_gp", "25%_gp", "50%_gp", "75%_gp", "90%_gp", "95%_gp") + + qq.gp <- as.data.frame(qq.gp) + # Calculate weighted SDs across studies in each group var.gp <- (diag(t(var.matrix)%*%n.matrix))/(t(n.matrix)%*%nsum.vector) SD.gp <- sqrt(var.gp) @@ -292,7 +374,6 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, dimnames(N.gp) <- c(list(names.gp),list("Nvalid_gp")) dimnames(SEM.gp) <- c(list(names.gp),list("SEM_gp")) - } # SPLIT @@ -300,6 +381,7 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, numsources <- length(output) mean.matrix <- NULL sd.matrix <- NULL + qq.matrix <- list() n.matrix <- NULL Nvalid <- 0 Nmissing <- 0 @@ -308,23 +390,32 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, for(j in 1:numsources){ mean.matrix <- rbind(mean.matrix,as.numeric(unlist(output[[j]][2]))) sd.matrix <- rbind(sd.matrix,as.numeric(unlist(output[[j]][3]))) - n.matrix <- rbind(n.matrix,as.numeric(unlist(output[[j]][4]))) - Nvalid <- Nvalid+as.numeric(unlist(output[[j]][5])) - Nmissing <- Nmissing+as.numeric(unlist(output[[j]][6])) - Ntotal <- Ntotal+as.numeric(unlist(output[[j]][7])) + qq.matrix[[j]] <- rbind(qq.matrix,as.numeric(unlist(output[[j]][4]))) + n.matrix <- rbind(n.matrix,as.numeric(unlist(output[[j]][5]))) + Nvalid <- Nvalid+as.numeric(unlist(output[[j]][6])) + Nmissing <- Nmissing+as.numeric(unlist(output[[j]][7])) + Ntotal <- Ntotal+as.numeric(unlist(output[[j]][8])) } var.matrix <- sd.matrix^2 - mean.gp.study <- t(mean.matrix) SD.gp.study <- t(sd.matrix) N.gp.study <- t(n.matrix) SEM.gp.study <- SD.gp.study/sqrt(N.gp.study) - + qq.gp.study <- lapply(qq.matrix, function(i){ + matrix(i, ncol = 7, + byrow = TRUE)}) - # create names + qq.gp.study <- lapply(qq.gp.study, function(i) { + colnames(i) <- c("5%_gp", "10%_gp", "25%_gp", "50%_gp", + "75%_gp", "90%_gp", "95%_gp") + i + }) + + + # create names names.gp <- rep(NA,dim(mean.gp.study)[1]) for(k in 1:dim(mean.gp.study)[1]){ @@ -338,7 +429,6 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, dimnames(N.gp.study) <- c(list(names.gp),list(names.study)) dimnames(SEM.gp.study) <- c(list(names.gp),list(names.study)) - } if(type=="combine"){ @@ -360,20 +450,28 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, if(draw.plot){ ## If plot is TRUE please select which metric to plot else table is returned - sel.metric <- get(metric_map[[metric]]) - sel.metric.data <- data.frame(lsoa11cd = lsoanames, + if(metric %in% q.names){ + + sel.metric.data <- data.frame(lsoa11cd = rownames(qq.gp)) + sel.metric.data$value <- qq.gp[ , which(metric == q.names)] + sel.metric.data$server <- 'combine' + + } else { + sel.metric <- get(metric_map[[metric]]) + + sel.metric.data <- data.frame(lsoa11cd = lsoanames, value = as.numeric(sel.metric), server = 'combine') - + } shape_sf <- shape_list[[1]] |> dplyr::left_join(sel.metric.data, by = "lsoa11cd") |> sf::st_as_sf() } else { - result <- list(mean.gp,SD.gp,N.gp,SEM.gp,Nvalid,Nmissing,Ntotal, lsoa_names[[1]]) - names(result) <- list("Mean_gp","StDev_gp","Nvalid_gp","SEM_gp","Total_Nvalid", + result <- list(mean.gp,SD.gp,N.gp,SEM.gp,qq.gp,Nvalid,Nmissing,Ntotal, lsoa_names[[1]]) + names(result) <- list("Mean_gp","StDev_gp","Nvalid_gp","SEM_gp", "Q_gp","Total_Nvalid", "Total_Nmissing","Total_Ntotal", "LSOAnames") return(result) } @@ -383,6 +481,20 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, if(type=="split"){ if(draw.plot){ + browser() + if(metric %in% q.names){ + + sel.metric.data <- + lapply(1:length(qq.gp.study), function(i) { + data.frame( + server = names(datasources)[i], + value = qq.gp.study[[i]][ ,which(metric == q.names)], + lsoa11cd = lsoa_names[[i]], + row.names = NULL + ) + }) + + }else{ sel.metric <- get(metric_map[[metric]]) sel.metric.data <- @@ -394,7 +506,7 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, row.names = NULL ) }) - + } shape_sf <- do.call( rbind, lapply(1:numsources, function(i){ @@ -409,10 +521,10 @@ ds.summarystat <- function(x=NULL, y=NULL, type='combine', do.checks=FALSE, lsoanames <- matrix(unlist(lsoa_names), ncol = numsources) dimnames(lsoanames) <- dimnames(mean.gp.study) - result <- list(mean.gp.study,SD.gp.study,N.gp.study,SEM.gp.study,Nvalid, + result <- list(mean.gp.study,SD.gp.study,N.gp.study,SEM.gp.study,qq.gp.study,Nvalid, Nmissing,Ntotal, lsoanames) names(result) <- list("Mean_gp_study","StDev_gp_study","Nvalid_gp_study", - "SEM_gp_study","Total_Nvalid","Total_Nmissing", + "SEM_gp_study","Q_gp","Total_Nvalid","Total_Nmissing", "Total_Ntotal", "LSOAnames") return(result) } diff --git a/vignettes/Demin's_work_functions.Rmd b/vignettes/Demin's_work_functions.Rmd new file mode 100644 index 0000000..e7f2d9f --- /dev/null +++ b/vignettes/Demin's_work_functions.Rmd @@ -0,0 +1,115 @@ +--- +title: "R Notebook" +output: html_notebook +--- + + +```{r} +DSOpal::Opal() +builder <- DSI::newDSLoginBuilder() +builder$append( + server = "mersey-demo", + url = "http://localhost:8880", + user = "administrator", + password = "password", + driver = "OpalDriver", + profile = "ds-geospatial" +) +``` + +```{r} +logindata <- builder$build() +conns <- DSI::datashield.login(logins = logindata) + +``` + + +```{r} +DSI::datashield.assign.resource( +conns, +symbol = "cells", +resource = paste0(project, ".", "cells") +) + +DSI::datashield.assign.resource( +conns, +symbol = "W", +resource = paste0(project, ".", "W_D_5km") +) +``` + +## Materialise arbitrary R objects ---- + + +```{r} + +DSI::datashield.assign.resource( + conns, + symbol = "cells_resource", + resource = paste0(project, ".cells") +) + +DSI::datashield.assign.expr( + conns, + symbol = "cells_object", + expr = quote(as.resource.object(cells_resource)) +) + +``` + +## Map specific columns ---- + + +```{r} +DSI::datashield.assign.expr( + conns, + symbol = "obesity_qof_2024_25", + expr = quote(cells_object$obesity_qof_2024_25) +) +DSI::datashield.assign.expr( + conns, + symbol = "W_object", + expr = quote(as.resource.object(W)) +) +DSI::datashield.assign.expr( + conns, + symbol = "W_band800", + expr = quote(W_object$band800) +) + +``` + + + +```{r} +ds.moransI( + x = "obesity_qof_2024_25", + listw = "W_band800", + nsim = 199, + datasources = conns +) + +``` + +```{r} +out <- ds.localmoran( + x = "obesity_qof_2024_25", + listw = "W_band800", + datasources = conns +) + +quadr_mean <- out$`mersey-demo`$quadr$mean + +counts <- table(quadr_mean) + +barplot( + counts, + xlab = "Local Moran's I quadrant", + ylab = "Count", + main = "Local Moran's I Quadrants" +) + + +``` + + diff --git a/vignettes/Opal_functioncheck.Rmd b/vignettes/Opal_functioncheck.Rmd index bcb8269..533984e 100644 --- a/vignettes/Opal_functioncheck.Rmd +++ b/vignettes/Opal_functioncheck.Rmd @@ -110,9 +110,11 @@ ds.colnames("D", conns) ```{r} library(dsGeospatialClient) +names_region <- list(c("Liverpool", "Sefton"), c("Liverpool")) + ds.summarystat(x = "D$has_asthma", y = "D$lsoa11cd", - type='combine', do.checks=FALSE, - draw.plot = TRUE, names_region = "Liverpool", + type='combe', do.checks=FALSE, + draw.plot = TRUE, names_region = names_region, metric = "Mean", datasources=conns) ```