Merge pull request #25 from lixun910/spatial_validation

[new feature] Spatial validation; join count ratio; make spatial; bivariate local moran

.vscode/settings.json

@@ -100,6 +100,8 @@
         "scoped_allocator": "cpp",
         "valarray": "cpp",
         "geodaweight.h": "c",
-        "filesystem": "cpp"
+        "filesystem": "cpp",
+        "numbers": "cpp",
+        "semaphore": "cpp"
     }
 }

DESCRIPTION

@@ -1,8 +1,8 @@
 Package: rgeoda
 Type: Package
 Title: R Library for Spatial Data Analysis
-Version: 0.0.8-4
-Date: 2021-08-04
+Version: 0.0.8-6
+Date: 2021-09-07
 Authors@R: 
     c(person(given = "Xun", family = "Li", email="[email protected]", role=c("aut","cre")),
       person(given = "Luc", family = "Anselin", email="[email protected]", role="aut"))

NAMESPACE

@@ -24,6 +24,7 @@ export(has_isolates)
 export(hinge15_breaks)
 export(hinge30_breaks)
 export(is_symmetric)
+export(join_count_ratio)
 export(kernel_knn_weights)
 export(kernel_weights)
 export(knn_weights)
@@ -36,6 +37,7 @@ export(lisa_num_nbrs)
 export(lisa_pvalues)
 export(lisa_values)
 export(local_bijoincount)
+export(local_bimoran)
 export(local_g)
 export(local_geary)
 export(local_gstar)
@@ -46,6 +48,7 @@ export(local_multigeary)
 export(local_multijoincount)
 export(local_multiquantilelisa)
 export(local_quantilelisa)
+export(make_spatial)
 export(max_neighbors)
 export(maxp_greedy)
 export(maxp_sa)
@@ -72,6 +75,7 @@ export(sf_to_geoda)
 export(skater)
 export(sp_to_geoda)
 export(spatial_lag)
+export(spatial_validation)
 export(stddev_breaks)
 export(update_weights)
 export(weights_sparsity)

R/RcppExports.R

@@ -37,6 +37,18 @@ p_azp_tabu <- function(p, xp_w, data, n_vars, tabu_length, conv_tabu, bound_vals
     .Call('_rgeoda_p_azp_tabu', PACKAGE = 'rgeoda', p, xp_w, data, n_vars, tabu_length, conv_tabu, bound_vals, min_bound, inits, init_regions, scale_method, distance_method, seed)
 }
 
+p_spatialvalidation <- function(xp_geoda, clusters, xp_w) {
+    .Call('_rgeoda_p_spatialvalidation', PACKAGE = 'rgeoda', xp_geoda, clusters, xp_w)
+}
+
+p_joincount_ratio <- function(clusters, xp_w) {
+    .Call('_rgeoda_p_joincount_ratio', PACKAGE = 'rgeoda', clusters, xp_w)
+}
+
+p_make_spatial <- function(clusters, xp_w) {
+    .Call('_rgeoda_p_make_spatial', PACKAGE = 'rgeoda', clusters, xp_w)
+}
+
 p_LISA__Run <- function(xp) {
     invisible(.Call('_rgeoda_p_LISA__Run', PACKAGE = 'rgeoda', xp))
 }
@@ -89,6 +101,10 @@ p_localmoran <- function(xp_w, data, permutations, permutation_method, significa
     .Call('_rgeoda_p_localmoran', PACKAGE = 'rgeoda', xp_w, data, permutations, permutation_method, significance_cutoff, cpu_threads, seed)
 }
 
+p_bi_localmoran <- function(xp_w, data1, data2, permutations, permutation_method, significance_cutoff, cpu_threads, seed) {
+    .Call('_rgeoda_p_bi_localmoran', PACKAGE = 'rgeoda', xp_w, data1, data2, permutations, permutation_method, significance_cutoff, cpu_threads, seed)
+}
+
 p_eb_rate <- function(event_data, base_data) {
     .Call('_rgeoda_p_eb_rate', PACKAGE = 'rgeoda', event_data, base_data)
 }

R/clustering.R

@@ -622,3 +622,98 @@ azp_tabu<- function(p, w, df, tabu_length=10, conv_tabu=10, bound_variable=data.
 
   return(p_azp_tabu(p, w$GetPointer(), df, n_vars, tabu_length, conv_tabu, bound_values, min_bound, inits, initial_regions, scale_method, distance_method, random_seed))
 }
+
+############################################################
+#' @title Spatial Validation
+#' @description Spatial validation provides a collection of validation measures including
+#' 1. fragmentations (entropy, simpson), 2. join count ratio, 3. compactness (isoperimeter quotient) 
+#' and 4. diameter.
+#' @param sf_obj An sf (simple feature) object 
+#' @param clusters A cluster classification variable (categorical values from a dataframe or values returned from cluster functions)
+#' @param w An instance of Weight class
+#' @return A list with names "Is Spatially Constrained", "Fragmentation", "Join Count Ratio",
+#' "Compactness", and "Diameter".
+#' @examples
+#' \dontrun{
+#' library(sf)
+#' guerry_path <- system.file("extdata", "Guerry.shp", package = "rgeoda")
+#' guerry <- st_read(guerry_path)
+#' queen_w <- queen_weights(guerry)
+#' data <- guerry[c('Crm_prs','Crm_prp','Litercy','Donatns','Infants','Suicids')]
+#' clusters <- skater(5, queen_w, data)
+#' results <- spatial_validation(guerry, clusters, queen_w)
+#' results
+#' }
+#' @export
+spatial_validation <- function(sf_obj, clusters, w) {
+  if (w$num_obs < 1) {
+    stop("The weights is not valid.")
+  }
+
+  geoda_obj <- getGeoDaObj(sf_obj)
+
+  if (geoda_obj$GetNumObs() <=0) {
+    stop("gda object is not valid.")
+  }
+
+  return (p_spatialvalidation(geoda_obj$GetPointer(), clusters, w$GetPointer()));
+}
+
+############################################################
+#' @title Join Count Ratio
+#' @description Join count ratio is the join counts, the number of times a category is surrounded 
+#' by neighbors of the same category, over the total number of neighbors after converting
+#' each category to a dummy variable.
+#' @param clusters A cluster classification variable (categorical values from a dataframe or values returned from cluster functions)
+#' @param w An instance of Weight class
+#' @return A data.frame with names "Cluster", "N", "Neighbors", "Join Count", "Ratio"
+#' @examples
+#' \dontrun{
+#' library(sf)
+#' guerry_path <- system.file("extdata", "Guerry.shp", package = "rgeoda")
+#' guerry <- st_read(guerry_path)
+#' queen_w <- queen_weights(guerry)
+#' data <- guerry[c('Crm_prs','Crm_prp','Litercy','Donatns','Infants','Suicids')]
+#' clusters <- skater(5, queen_w, data)
+#' results <- join_count_ratio(clusters, queen_w)
+#' results
+#' }
+#' @export
+join_count_ratio<- function(clusters, w) {
+  if (w$num_obs < 1) {
+    stop("The weights is not valid.")
+  }
+
+  return (p_joincount_ratio(clusters, w$GetPointer()));
+}
+
+############################################################
+#' @title Make Spatial
+#' @description Make spatially constrained clusters from spatially non-constrained clusters
+#' using the contiguity information from the input weights
+#' @param clusters A cluster classification variable (categorical values from a dataframe or values returned from cluster functions)
+#' @param w An instance of Weight class
+#' @return A vector of categorical values (cluster classification)
+#' @examples
+#' \dontrun{
+#' library(sf)
+#' guerry_path <- system.file("extdata", "Guerry.shp", package = "rgeoda")
+#' guerry <- st_read(guerry_path)
+#' data <- guerry[c('Crm_prs','Crm_prp','Litercy','Donatns','Infants','Suicids')]
+#' clusters <- kmeans(5, data)
+#' queen_w <- queen_weights(guerry)
+#' results <- make_spatial(clusters, queen_w)
+#' results
+#' }
+#' @export
+make_spatial <- function(clusters, w) {
+  if (w$num_obs < 1) {
+    stop("The weights is not valid.")
+  }
+
+  if (w$num_obs != length(clusters)) {
+    stop("The weights doesn not match with the size of input clusters.") 
+  }
+
+  return (p_make_spatial(clusters, w$GetPointer()));
+}

R/lisa.R

@@ -316,6 +316,60 @@ local_moran <- function(w, df, permutations=999, permutation_method="complete",
   return(LISA$new(p_LISA(lisa_obj)))
 }
 
+#################################################################
+#' @title  Bivariate Local Moran Statistics
+#' @description The function to apply bivariate local Moran statistics
+#' @param w An instance of Weight object
+#' @param df A data frame with two selected variable. E.g. guerry[c('Crm_prs','Litercy')]
+#' @param permutations (optional) The number of permutations for the LISA
+#' computation
+#' @param permutation_method (optional) The permutation method used for the
+#' LISA computation. Options are {'complete', 'lookup'}. Default is 'complete'.
+#' @param significance_cutoff  (optional) A cutoff value for significance
+#' p-values to filter not-significant clusters
+#' @param cpu_threads (optional) The number of cpu threads used for parallel
+#' LISA computation
+#' @param seed (optional) The seed for random number generator
+#' @return An instance of LISA-class
+#' @examples
+#' library(sf)
+#' guerry_path <- system.file("extdata", "Guerry.shp", package = "rgeoda")
+#' guerry <- st_read(guerry_path)
+#' queen_w <- queen_weights(guerry)
+#' lisa <- local_bimoran(queen_w, guerry[c('Crm_prs','Litercy')])
+#' lms <- lisa_values(lisa)
+#' lms
+#' @export
+local_bimoran <- function(w, df, permutations=999, permutation_method="complete",
+                        significance_cutoff=0.05, cpu_threads=6,
+                        seed=123456789) {
+  if (w$num_obs <= 0) {
+    stop("Weights object is not valid.")
+  }
+
+  if (inherits(df, "data.frame") == FALSE) {
+    stop("The input data needs to be a data.frame.")
+  }
+
+  num_vars <- length(df)
+
+  if (inherits(df, "sf")) {
+    num_vars <- num_vars - 1
+  }
+
+  if (num_vars != 2) {
+    stop("Two variables are required for bivariate local moran.")
+  }
+
+  data1 <- df[[1]]
+  data2 <- df[[2]]
+
+  lisa_obj <- p_bi_localmoran(w$GetPointer(), data1, data2, permutations,
+                           permutation_method, significance_cutoff,
+                           cpu_threads, seed)
+  return(LISA$new(p_LISA(lisa_obj)))
+}
+
 #################################################################
 #' @title  Local Moran with Empirical Bayes(EB) Rate
 #' @description The function to apply local Moran with EB Rate statistics. The

README.md

@@ -61,6 +61,7 @@ https://geodacenter.github.io/rgeoda/reference/
 
 * Spatial Autocorrelation
     * Local Moran
+    * Bivariate Local Moran
     * Local Moran EB Rates
     * Local Geary
     * Local Getis-Ord 
@@ -93,6 +94,13 @@ https://geodacenter.github.io/rgeoda/reference/
       * greedy
       * Tabu Search
       * Simulated Annealing
+    * Join Count Ratio
+    * Spatial Validation
+      * Fragmentation
+      * Join Count Ratio
+      * Compactness
+      * Diameter
+    * Make Spatial
 
 ## Build and install from source code