Dev ds.morans i

R/ds.lagsarlm.R

@@ -0,0 +1,148 @@
+# library(spatialreg)
+# 
+# model <- lagsarlm(
+#   obesity_qof_2024_25 ~ dist_food_fastfood + dist_food_supermarket + dist_hc_emergency,
+#   data = cells,
+#   listw = W[["band800"]],
+#   zero.policy = TRUE,
+#   method = "Matrix"
+# )
+# 
+#' 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.lagsarlm <- function(formula = NULL, data = NULL, listw = NULL, type='split',
+                       do.checks = TRUE,
+                       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(formula)){
+    stop("Please provide the input formula!", 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)
+  }
+ browser()
+
+  calltext <- paste0("lagsarlmDS(", deparse(formula), ", ", listw, ", ", data, ")", collapse = "")
+
+  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)
+  }
+}
+
+

vignettes/Opal_functioncheck.Rmd

@@ -3,6 +3,7 @@ title: "R Notebook"
 output: html_notebook
 ---
 
+Build package locally.
 
 ```{r, eval=FALSE}
 # r console
@@ -21,9 +22,106 @@ devtools::build()
 # go to administration --> R see if the profiles appear
 # go to administration --> Datashield activate profile and initialize
 ```
-This is an [R Markdown](http://rmarkdown.rstudio.com) Notebook. When you execute code within the notebook, the results appear beneath the code. 
 
-Try executing this chunk by clicking the *Run* button within the chunk or by placing your cursor inside it and pressing *Cmd+Shift+Enter*. 
+Create example data set.
+
+```{r, eval=FALSE}
+# retrieve LSOA boundaries for Cheshire and Merseyside
+# to reduce coverage, remove names from the `within_names`
+cm_lsoa_sf <- boundr::bounds(
+  "lsoa",
+  within_level = "lad",
+  within_names = c(
+    "Cheshire East",
+    "Cheshire West and Chester",
+    "Halton",
+    "Knowsley",
+    "Liverpool",
+    "Sefton",
+    "St. Helens",
+    "Warrington",
+    "Wirral"
+  ),
+  lookup_year = 2011, # can be either 2011 or 2021
+  opts = boundr::boundr_options(resolution = "BFC")
+)
+
+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" |>
+  readr::read_csv() |>
+  # convert to spatial object
+  sf::st_as_sf(coords = c("longitude", "latitude"), crs = 4326) |>
+  # add LSOA details
+  sf::st_join(cm_lsoa_sf)
+
+
+
+qof_asthma_tbl <- "https://pldr.org/download/e6nzv/ng1/QOF_4_03_Asthma_LSOA.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::mutate(
+        dist_gs_dec = dplyr::ntile(distance_local_greenspace, 10)
+      ) |>
+      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)
+    )
+
+    # 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()
+
+
+```
+
+Subset the data for needed coloumns. 
 
 ```{r, eval=FALSE}
 library(opalr)
@@ -35,11 +133,12 @@ dat <- qof_asthma_cohort_dist_gs_sf |>
 
 # Uploading a subset to save upload time
 
-dat <- dat |> dplyr::filter(lsoa11cd %in% sample(levels(factor(dat$lsoa11cd)), 35))
+#dat <- dat |> dplyr::filter(lsoa11cd %in% sample(levels(factor(dat$lsoa11cd)), 35))
 dat$row_id <- seq_len(nrow(dat)) # needs a unique row id
 
 ```
 
+Creating project and uploading table.
 
 ```{r}
 library(opalr)
@@ -76,6 +175,7 @@ opal.table_save(
 
 ```
 
+Login 
 ```{r}
 library(DSOpal)
 library(dsBase)
@@ -122,8 +222,8 @@ ds.summarystat(x = "D$has_asthma", y = "D$lsoa11cd",
 ```{r}
 
 ds.groupcor(x = "D$has_asthma", y = "D$distance_local_greenspace", 
-            index = "D$lsoa11cd", datasources = conns, type = "combine", 
-            names_region = "Liverpool")
+            index = "D$lsoa11cd", datasources = conns, type = "split", 
+            names_region = names_region)
 ```
 
 ```{r}
@@ -143,8 +243,4 @@ fitdat <- ds.glmerSLMA(formula="D$has_asthma ~ D$sex + (1|D$lsoa11cd)",
              family="poisson",datasources=conns
 )
 ```
-Add a new chunk by clicking the *Insert Chunk* button on the toolbar or by pressing *Cmd+Option+I*.
-
-When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the *Preview* button or press *Cmd+Shift+K* to preview the HTML file). 
 
-The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike *Knit*, *Preview* does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.