# R/m.surround.R In spqdep: Testing for Spatial Independence of Qualitative Data in Cross Section

#### Documented in m.surround

#'
#' @title A function to generate m-surroundings
#' @usage m.surround(x, m, r = 1, distance = "Euclidean", control = list())
#' @param x input sf object with points/multipolygons geometry or matrix
#'          of spatial coordinates
#' @param m dimension of m-surrounding
#' @param r maximum overlapping between any two m-surroundings.
#' @param distance character. For Cartesian coordinates only: one of Euclidean,
#' Hausdorff or Frechet; for geodetic coordinates, Great Circle distances are
#' computed. (see \code{sf::st_distance()}).
#' Default = "Euclidean".
#' @param control List of additional control arguments.
#'
#' @description This function obtains the m-surroundings
#'   by selecting the *m-1* nearest neighbors
#'   of each observation, allowing for a degree of
#'   overlap of *r*.
#'
#' @return A list of class \code{list} and \code{m_surr} containing the following components:\cr
#'   \tabular{ll}{
#'     \code{ms} \tab a matrix. Each row is a m-surrounding. \cr
#'     \code{R} \tab total number of observations. \cr
#'     \code{rowexcl} \tab index of rows excluded. \cr
#'     \code{mdtms} \tab distances between the observations of
#'       each m-surrounding.\cr
#'     \code{N} \tab total number of symbolized observations.\cr
#'     \code{m} \tab length of the m-surroundings.\cr
#'     \code{r} \tab overlapping degree.\cr
#'     \code{initobs} \tab element to start the generation of m-surroundings.\cr
#'     \code{distance} \tab type of distance.\cr
#'     \code{m} \tab length of the m-surroundings.\cr
#'     \code{x} \tab the input "x" as sf-object.\cr
#'     }
#'
#' @section Control arguments:
#'   \itemize{
#'     \item {\code{initobs}:} {Initial observation to begin the
#'       m-surrounding process. Default = 1.}
#'     \item {\code{dtmaxabs}:} {Threshold of distance (in absolute
#'       value) to prune the m-surroundings. Any m-surrounding
#'       exceeding the threshold is excluded.
#'       If \code{dtmaxabs} = 0 there is no exclusion of m-surroundings.
#'       Default = 0.}
#'     \item {\code{dtmaxpc}:} {Threshold of distance (as a percentage
#'       of the maximum distance between observations) to prune the m-surroundings. Any m-surrounding
#'       exceeding the threshold is excluded.
#'       Example if \code{dtmaxpc} = 0.1 the m-surrounding exceeding
#'       the 10% of the maximum distance between observations are excluded.
#'       If \code{dtmaxpc} = 0 there is no exclusion of m-surroundings.
#'       Default = 0.}
#'     \item {\code{dtmaxknn}:} {Eliminate m-surroundings where some
#'       of the elements are not among the closest
#'       knn (k-nearest-neighbors).
#'       Example, if \code{dtmaxknn} = 4 exclude m-surroundings
#'       where some of the elements are not between the 4 closest.
#'       Default \code{dtmaxknn} = 0 (no exclusion)}
#'    }
#'
#' @author
#'   \tabular{ll}{
#'   Fernando López  \tab \email{fernando.lopez@@upct.es} \cr
#'   Román Mínguez  \tab \email{roman.minguez@@uclm.es} \cr
#'   Antonio Páez \tab \email{paezha@@gmail.com} \cr
#'   Manuel Ruiz \tab \email{manuel.ruiz@@upct.es} \cr
#'   }
#' @references
#'   \itemize{
#'     \item Ruiz M, López FA, A Páez. (2010). Testing for spatial association of qualitative
#'     data using symbolic dynamics. \emph{Journal of Geographical Systems}. 12 (3) 281-309
#'   }
#' @export
#' @examples
#'
#' # Example 1: Obtain m-surroundings with degree of overlapping r
#' N <- 100
#' cx <- runif(N)
#' cy <- runif(N)
#' x <- cbind(cx,cy)
#' m <- 3
#' r <- 1
#' msurr_points <- m.surround(x = x, m = m, r = r)
#' plot(msurr_points, type = 1)
#' plot(msurr_points, type = 2)
#' summary(msurr_points)
#' msurr_points <- m.surround(x = x, m = m, r = r,
#'                 control = list(dtmaxpc = 0.1))
#' plot(msurr_points, type = 1)
#' plot(msurr_points, type = 2)
#' summary(msurr_points)
#' msurr_points <- m.surround(x = x, m = m, r = r,
#'                    control = list(dtmaxknn = 20))
#' plot(msurr_points, type = 1)
#' plot(msurr_points, type = 2)
#' summary(msurr_points)
#'
#' # Example 2:
#' \donttest{
#' data("FastFood.sf")
#' m <- 3
#' r <- 1
#' msurr_points <- m.surround(x = FastFood.sf, m = m, r = r,
#'                            distance = "Euclidean",
#'                            control = list(dtmaxpc = .001))
#' plot(msurr_points, type = 1)
#' plot(msurr_points, type = 2)
#' print(msurr_points)
#' summary(msurr_points)
#' msurr_points <- m.surround(x = FastFood.sf, m = m, r = r,
#'                            distance = "Euclidean",
#'                            control = list(dtmaxknn = 20))
#' plot(msurr_points, type = 1)
#' plot(msurr_points, type = 2)
#' summary(msurr_points)
#' }
#'
#' # Example 3: With isolated areas
#' data(provinces_spain)
#' sf::sf_use_s2(FALSE)
#' plot(sf::st_geometry(provinces_spain))
#' m <- 3
#' r <- 1
#' msurr_points <- m.surround(x = provinces_spain, m = m, r = r,
#'                            distance = "Euclidean",
#'                            control = list(dtmaxknn = 8))
#' plot(msurr_points, type = 1)
#' plot(msurr_points, type = 2)
#'
#' # Example 4: Examples with multipolygons
#' fname <- system.file("shape/nc.shp", package="sf")
#' plot(sf::st_geometry(nc))
#' m <- 3
#' r <- 1
#' msurr_polygonsf <- m.surround(x = nc, m = m, r = r,
#'                    distance = "Great Circle",
#'                    control=list(dtmaxpc = 0.20))
#' plot(msurr_polygonsf, type = 1)
#' plot(msurr_polygonsf, type = 2)
#'
#' # Example 5: With regular lattice
#' sfc = sf::st_sfc(sf::st_polygon(list(rbind(c(0,0), c(1,0), c(1,1), c(0,1), c(0,0)))))
#' hexs <- sf::st_make_grid(sfc, cellsize = 0.1, square = FALSE)
#' hexs.sf <- sf::st_sf(hexs)
#' listw  <- spdep::poly2nb(as(hexs.sf, "Spatial"), queen = FALSE)
#' m <- 3
#' r <- 1
#' msurr_polygonsf <- m.surround(x = hexs.sf, m = m, r = r)
#' plot(msurr_polygonsf, type = 1)
#' plot(msurr_polygonsf, type = 2)
#' summary(msurr_polygonsf)

m.surround <- function(x , m, r = 1, distance = "Euclidean", control = list()) {
con <- list(initobs = 1, dtmaxabs = 0, dtmaxpc = 0, dtmaxknn = 0)
nmsC <- names(con)
con[(namc <- names(control))] <- control
if (length(noNms <- namc[!namc %in% nmsC]))
warning("unknown names in control: ", paste(noNms, collapse = ", "))
if (sum((con$dtmaxabs>0) + (con$dtmaxpc>0) + (con$dtmaxknn>0))>1){ stop("Include only a restriction to create m-surrounding") } if (!is.null(r)) { if (r < 1 || r > (m-1)) warning("overlapping degree must be between 1 and m-1") } # Transform matrix coordinates into SpatialPoints class if (is.matrix(x)) x <- sp::SpatialPoints(x) # Transform Spatial classes into sf class if (inherits(x, "Spatial")) x <- as(x, "sf") # Compute centroids/coordinates/distances from sf objects if (inherits(x, "sf")) { xct <- suppressWarnings(st_centroid(st_geometry(x))) } else stop("object must be either sf, sp or matrix class") N <- length(xct) if (!is.null(r)) R <- trunc((N - m)/(m - r)) + 1 else R <- N if (inherits(st_geometry(x)[1], "sfc_MULTIPOLYGON")) mcoor <- st_coordinates(xct) if (inherits(st_geometry(x)[1], "sfc_POINT")) mcoor <- st_coordinates(x) mcoor <- st_coordinates(xct) rownames(mcoor) <- as.character(1:N) mdtfull <- st_distance(xct, which = distance) # full distance matrix rownames(mdtfull) <- colnames(mdtfull) <- as.character(1:N) ms <- m_surr_no(mdtfull = mdtfull, m = m, r = r, initobs = con$initobs, control_knn = con$dtmaxknn, coor = mcoor) mdtms <- NULL for (i in seq_along(1:nrow(ms))) { rowds <- mdtfull[ms[i,],ms[i,]] rowds <- rowds[1,] mdtms <- rbind(mdtms, rowds) } colnames(mdtms) <- NULL rownames(mdtms) <- ms[, 1] # ms <- matrix(NA, R, m) # mdtms <- matrix(NA, R, m) # m-surrounding distance matrix # rownames(ms) <- rownames(mdtms) <- 1:R # set.seed(con$seedinit)
# Si <- as.character(1:N)
# s0i <- as.character(sample(1:N, 1)) # initial s0
# mdti <- mdtfull
# if (!is.null(r)) {
#   for (i in 1:R) {
#     vdsts0i <- mdti[s0i,] # Vector of distances to s0i
#     vidxnbs0i <- sort(vdsts0i, decreasing = FALSE, index.return = TRUE)
#     # Distances to s0i ordered
#     if (is.list(vidxnbs0i)) vnbs0i <- names(vidxnbs0i$x[2:m]) # if (inherits(vidxnbs0i, "units")) vnbs0i <- names(vidxnbs0i[2:m]) # # Set of (m-1) neighbors to s0i (excluding itself) # ms[i,] <- c(s0i, vnbs0i) # m-surrounding ith # mdtms[i,] <- vdsts0i[c(s0i,vnbs0i)] # distance vector of ms_i # if ((m-r-1) > 0) Ai <- c(s0i,vnbs0i[1:(m-r-1)]) else Ai <- c(s0i) # # if ((m-r-2) > 0) Ai <- c(s0i,vnbs0i[1:(m-r-2)]) else Ai <- c(s0i) # Si <- Si[!(Si %in% Ai)] # mdti <- mdti[Si, Si] # s0i <- vnbs0i[m-r] # } # } else { # Bootstrap resampling # ms <- cbind(1:N, knearneigh(mcoor, k=(m-1))$nn)
#
#
# }
if (is.null(con$dtmaxabs)) dtmaxabs <- 0 else dtmaxabs <- con$dtmaxabs
if (is.null(con$dtmaxpc)) dtmaxpc <- 0 else dtmaxpc <- con$dtmaxpc
if (is.null(con$dtmaxknn)) dtmaxknn <- 0 else dtmaxknn <- con$dtmaxknn

lms <- prunemsdthr(dtmaxabs = dtmaxabs, dtmaxpc = dtmaxpc,
mdtfull = mdtfull, ms = ms,
mdtms = mdtms)
lms$N <- N lms$m <- m
lms$r <- r lms$initobs <- con$initobs lms$distance <- distance
lms$x <- x class(lms) <- c("m_surr","list") return(lms) } # nnlist <- matrix(0, N, m - 1) # Matrix with list of nearest neighbors # nn1 <- cbind(1:N, mdtfull[1,], mcoor) # #nn1 <- cbind(1:N, sqrt((mcoor[1,1] - mcoor[,1])^2 + # # (mcoor[1,2] - mcoor[,2])^2), mcoor) # nn1 <- nn1[order(nn1[,2]), ] # nnlist[1, ] <- nn1[2:m, 1] # ns <- trunc((N - m)/(m - r)) + 1 # list <- rep(0, ns) #zeros(1,ns) # list[1] = 1 # blacklist <- c(1, nnlist[1, 1:(m - (r + 1))]) # t <- 1 # for (v in 2:ns) { # list[v] = nnlist[t, m - r] # h <- list[v] # nn1 <- nn1[!nn1[,1] %in% blacklist,] # #**VIP**: SALE DISTINTO PERO DEBERÍA DAR IGUAL... CONSULTAR CON # # FERNANDO # #nn1[,2] <- mdtfull[nn1[1,1],nn1[,1]] # nn1[,2] <- sqrt((nn1[1,3] - nn1[,3])^2 + (nn1[1,4] - nn1[,4])^2) # nn1 <- nn1[order(nn1[,2]),] # nnlist[h,] <- nn1[2:(m),1] # t = h # blacklist <- c(blacklist, h, nnlist[h, 1:(m - (r + 1))]) # } # nnlist <- cbind(1:N, nnlist) # nnlist1 <- cbind(nnlist[which(!nnlist[,2] == 0),]) # # control debería ser una lista # # Se eliminan aquellas m-historias que contengan vecinos que no estén # # dentro de los k vecinos más próximos # if (!is.null(control)){ # knn <- cbind(1:N, knearneigh(mcoor, control)$nn)
#   int <- numeric()
#   for (i in 1:dim(nnlist1)[1]){
#     int[i] <- length(intersect(nnlist1[i,],knn[nnlist1[i,1],]))
#   }
#   nnlist1 <- nnlist1[int==m,]
# }
#return(nnlist1)
#}


## Try the spqdep package in your browser

Any scripts or data that you put into this service are public.

spqdep documentation built on March 28, 2022, 5:06 p.m.