Problem with weighted center mean

[JoseLastra]

https://github.com/JosiahParry/sfdep/blob/7e5cddbf9be70911644d31cb0367a69d6a077a74/R/point-pattern-centers.R#L38

I've been testing this functions and considering this formula: https://desktop.arcgis.com/es/arcmap/10.6/tools/spatial-statistics-toolbox/GUID-6844F6F5-B0F9-4B4A-B0EA-04AC6AA8C1BC-web.png

When you put a vetcor of weights this line should be like this:

res <- colSums(coords * weights)/sum(weights, na.rm = TRUE)

[JosiahParry]

Thanks! Can you please provide a reproducible example? Or known output for a dataset you can share?

On Thu, Oct 6, 2022 at 11:30 JoseLastra @.***> wrote:

https://github.com/JosiahParry/sfdep/blob/7e5cddbf9be70911644d31cb0367a69d6a077a74/R/point-pattern-centers.R#L38

I've been testing this functions and considering this formula: https://desktop.arcgis.com/es/arcmap/10.6/tools/spatial-statistics-toolbox/GUID-6844F6F5-B0F9-4B4A-B0EA-04AC6AA8C1BC-web.png

When you put a vetcor of weights this line should be like this:

res <- colSums(coords weights)/sum(weights, na.rm = TRUE)*

— Reply to this email directly, view it on GitHub https://github.com/JosiahParry/sfdep/issues/35, or unsubscribe https://github.com/notifications/unsubscribe-auth/ADHIKLBBEDAMUGFOXAVWW2LWB3WCHANCNFSM6AAAAAAQ6X26HQ . You are receiving this because you are subscribed to this thread.Message ID: @.***>

[JoseLastra]

I was using the test points from PySal http://pysal.org/notebooks/explore/pointpats/centrography.html

library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.4.3, PROJ 7.2.1; sf_use_s2() is TRUE
library(sfdep)
library(tidyverse,warn.conflicts = F)

#points creation 
points <- rbind(c(66.22, 32.54), c(22.52, 22.39), c(31.01, 81.21),
          c(9.47, 31.02),  c(30.78, 60.10), c(75.21, 58.93),
          c(79.26,  7.68), c(8.23, 39.93),  c(98.73, 77.17),
          c(89.78, 42.53), c(65.19, 92.08), c(54.46, 8.48)) %>% 
  as.data.frame() %>% rename(x = V1, y = V2) %>% mutate(weights = seq(0,11,1))

# to sf 
points_sf <- points %>% st_as_sf(coords = c('x','y'))

# Test center mean no W
mean_c <- center_mean(points_sf)

#plot results
plot(points_sf$geometry, axes = T)
plot(mean_c, col = 'red', add = T)

# Test center with W, default sfdep function
mean_cw <- center_mean(points_sf, weights = points_sf$weights)
#plot results
plot(points_sf$geometry, axes = T)
plot(mean_cw, col = 'blue', add = T) #result outside of bbox

## changing the line and possible solution
 center_mean_ps <- function (geometry, weights = NULL) {
  #geometry <- check_polygon(geometry) # removed just for this example
  crs <- sf::st_crs(geometry)
  coords <- sf::st_coordinates(geometry)
  n <- nrow(coords)
  if (!is.null(weights)) {
    res <- colSums(coords * weights)/sum(weights, na.rm = T)
  }
  else {
    res <- colSums(coords)/n
  }
  sf::st_sfc(sf::st_point(res), crs = crs)
}

#testing 
 mean_cw <- center_mean_ps(points_sf, weights = points_sf$weights)
 #plot results
 plot(points_sf$geometry, axes = T)
 plot(mean_cw, col = 'blue', add = T) #result outside of bbox

 ## check both mean center
 #plot results
 plot(points_sf$geometry, axes = T)
 plot(mean_c, col = 'red', add = T, pch = 13) 
 plot(mean_cw, col = 'blue', add = T,pch = 15) #result outside of bbox

^{Created on 2022-10-06 with reprex v2.0.2}