From adf2be392e4b4685ea778ff690276bbf827cda32 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 00:36:34 +0200 Subject: [PATCH 01/27] Add files via upload Two new variants of the Moran scatterplot. The Moran drop plot visualises indications of p-values and the Moran seismogram allows for visual inspection of potential local spatial weigths misspecifications. --- R/moran.plot.drop.R | 115 ++++++++++++++++++++++++++++++++++++++ R/moran.plot.seismogram.R | 85 ++++++++++++++++++++++++++++ 2 files changed, 200 insertions(+) create mode 100644 R/moran.plot.drop.R create mode 100644 R/moran.plot.seismogram.R diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R new file mode 100644 index 00000000..796e9e93 --- /dev/null +++ b/R/moran.plot.drop.R @@ -0,0 +1,115 @@ +moran.plot.drop <- function(x, listw, nsim = 999, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) { + require(spdep) + if (!inherits(listw, "listw")) + stop(paste(deparse(substitute(listw)), "is not a listw object")) + stopifnot(is.vector(x)) + stopifnot(is.logical(significant)) + stopifnot(is.logical(plain)) + stopifnot(is.logical(return_df)) + if (is.null(zero.policy)) + zero.policy <- get("zeroPolicy", envir = .spdepOptions) + stopifnot(is.logical(zero.policy)) + xname <- deparse(substitute(x)) + if (!is.numeric(x)) + stop(paste(xname, "is not a numeric vector")) + if (any(is.na(x))) + stop("NA in X") + n <- length(listw$neighbours) + if (n != length(x)) + stop("objects of different length") + if (is.null(spChk)) + spChk <- get.spChkOption() + if (spChk && !chkIDs(x, listw)) + stop("Check of data and weights ID integrity failed") + labs <- TRUE + if (is.logical(labels) && !labels) + labs <- FALSE + if (is.null(labels) || length(labels) != n) + labels <- as.character(attr(listw, "region.id")) + if (is.null(xlab)) + xlab <- xname + if (is.null(ylab)) + ylab <- paste("spatially lagged", xname) + Z <- as.vector(scale(x)) + ZLXi <- lag.listw(listw, Z, zero.policy = zero.policy) + locmoran <- localmoran_perm(Z, listw, nsim = nsim) + ZIi <- locmoran[, 4] + b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) + b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) + b[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] <- b2[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] + if(plot) { + if(!plain) { + if(significant) { + x_q1 <- Z[which(Z > 0 & ZLXi > 0 & ZIi >= cv)] + y_q1 <- ZLXi[which(Z > 0 & ZLXi > 0 & ZIi >= cv)] + b_q1 <- b[which(Z > 0 & ZLXi > 0 & ZIi >= cv)] + x_q2 <- Z[which(Z > 0 & ZLXi < 0 & ZIi <= (-1) * cv)] + y_q2 <- ZLXi[which(Z > 0 & ZLXi < 0 & ZIi <= (-1) * cv)] + b_q2 <- b[which(Z > 0 & ZLXi < 0 & ZIi <= (-1) * cv)] + x_q3 <- Z[which(Z < 0 & ZLXi < 0 & ZIi >= cv)] + y_q3 <- ZLXi[which(Z < 0 & ZLXi < 0 & ZIi >= cv)] + b_q3 <- b[which(Z < 0 & ZLXi < 0 & ZIi >= cv)] + x_q4 <- Z[which(Z < 0 & ZLXi > 0 & ZIi <= (-1) * cv)] + y_q4 <- ZLXi[which(Z < 0 & ZLXi > 0 & ZIi <= (-1) * cv)] + b_q4 <- b[which(Z < 0 & ZLXi > 0 & ZIi <= (-1) * cv)] + } else { + x_q1 <- Z[which(Z > 0 & ZLXi > 0 & ZIi < cv)] + y_q1 <- ZLXi[which(Z > 0 & ZLXi > 0 & ZIi < cv)] + b_q1 <- b[which(Z > 0 & ZLXi > 0 & ZIi < cv)] + x_q2 <- Z[which(Z > 0 & ZLXi < 0 & ZIi > (-1) * cv)] + y_q2 <- ZLXi[which(Z > 0 & ZLXi < 0 & ZIi > (-1) * cv)] + b_q2 <- b[which(Z > 0 & ZLXi < 0 & ZIi > (-1) * cv)] + x_q3 <- Z[which(Z < 0 & ZLXi < 0 & ZIi < cv)] + y_q3 <- ZLXi[which(Z < 0 & ZLXi < 0 & ZIi < cv)] + b_q3 <- b[which(Z < 0 & ZLXi < 0 & ZIi < cv)] + x_q4 <- Z[which(Z < 0 & ZLXi > 0 & ZIi > (-1) * cv)] + y_q4 <- ZLXi[which(Z < 0 & ZLXi > 0 & ZIi > (-1) * cv)] + b_q4 <- b[which(Z < 0 & ZLXi > 0 & ZIi > (-1) * cv)] + } + } + lw.lm <- lm(ZLXi ~ Z) + plot(Z, ZLXi, xlab="Z", ylab="WZ", pch = 20, cex = 0.33, col = "lightgrey", xlim = c(min(Z),max(Z)), ylim = c(min(ZLXi, b),max(ZLXi, b))) + abline(h = 0, lty = "dashed", col = "grey30") + abline(v = 0, lty = "dashed", col = "grey30") + abline(lw.lm, lty = "dotted", col = "grey40") + if(!plain) { + for(i in 1:length(x_q1)) { + lines(c(x_q1[i], x_q1[i]), c(y_q1[i], b_q1[i]), type = "l", col = "lightpink", lwd = 1) + points(c(x_q1[i]), c(y_q1[i]), pch = 20, cex = 0.5, col = "firebrick") + if (labs && length(x_q1) > 0) + text(x_q1[i], y_q1[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") + } + + for(i in 1:length(x_q2)) { + lines(c(x_q2[i], x_q2[i]), c(y_q2[i], b_q2[i]), type = "l", col = "lightblue1", lwd = 1) + points(c(x_q2[i]), c(y_q2[i]), pch = 20, cex = 0.5, col = "royalblue") + if (labs && length(x_q2) > 0) + text(x_q2[i], y_q2[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") + } + + for(i in 1:length(x_q3)) { + lines(c(x_q3[i], x_q3[i]), c(y_q3[i], b_q3[i]), type = "l", col = "lightpink", lwd = 1) + points(c(x_q3[i]), c(y_q3[i]), pch = 20, cex = 0.5, col = "firebrick") + if (labs && length(x_q3) > 0) + text(x_q3[i], y_q3[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") + } + + for(i in 1:length(x_q4)) { + lines(c(x_q4[i], x_q4[i]), c(y_q4[i], b_q4[i]), type = "l", col = "lightblue1", lwd = 1) + points(c(x_q4[i]), c(y_q4[i]), pch = 20, cex = 0.5, col = "royalblue") + if (labs && length(x_q4) > 0) + text(x_q4[i], y_q4[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") + } + if (zero.policy) { + n0 <- ZLXi == 0 + if (any(n0)) { + symbols(x[n0], wx[n0], inches = FALSE, circles = rep(diff(range(x))/50, length(which(n0))), bg = "grey", add = TRUE) + } + } + } + } + if(return_df) { + res <- data.frame(z = Z, wz = ZLXi, b = b, line_lengths = abs(ZLXi - b)) + invisible(res) + } +} \ No newline at end of file diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R new file mode 100644 index 00000000..f5023d73 --- /dev/null +++ b/R/moran.plot.seismogram.R @@ -0,0 +1,85 @@ +moran.plot.seismogram <- function(x, listw, nsim = 999, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) { + require(spdep) + if (!inherits(listw, "listw")) + stop(paste(deparse(substitute(listw)), "is not a listw object")) + stopifnot(is.vector(x)) + stopifnot(is.logical(plain)) + stopifnot(is.logical(return_df)) + if (is.null(zero.policy)) + zero.policy <- get("zeroPolicy", envir = .spdepOptions) + stopifnot(is.logical(zero.policy)) + xname <- deparse(substitute(x)) + if (!is.numeric(x)) + stop(paste(xname, "is not a numeric vector")) + if (any(is.na(x))) + stop("NA in X") + n <- length(listw$neighbours) + if (n != length(x)) + stop("objects of different length") + if (is.null(spChk)) + spChk <- get.spChkOption() + if (spChk && !chkIDs(x, listw)) + stop("Check of data and weights ID integrity failed") + if (is.null(xlab)) + xlab <- xname + if (is.null(ylab)) + ylab <- paste("spatially lagged", xname) + Z <- as.vector(scale(x)) + ZLXi <- lag.listw(listw, Z, zero.policy = zero.policy) + locmoran <- localmoran_perm(Z, listw, nsim = nsim) + ZIi <- locmoran[, 4] + b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) + b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) + b[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] <- b2[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] + if(plot) { + if(!plain) { + x_q1 <- Z[which(Z > 0 & ZLXi > 0)] + y_q1 <- ZLXi[which(Z > 0 & ZLXi > 0)] + b_q1 <- b[which(Z > 0 & ZLXi > 0)] + x_q2 <- Z[which(Z > 0 & ZLXi < 0)] + y_q2 <- ZLXi[which(Z > 0 & ZLXi < 0)] + b_q2 <- b[which(Z > 0 & ZLXi < 0)] + x_q3 <- Z[which(Z < 0 & ZLXi < 0)] + y_q3 <- ZLXi[which(Z < 0 & ZLXi < 0)] + b_q3 <- b[which(Z < 0 & ZLXi < 0)] + x_q4 <- Z[which(Z < 0 & ZLXi > 0)] + y_q4 <- ZLXi[which(Z < 0 & ZLXi > 0)] + b_q4 <- b[which(Z < 0 & ZLXi > 0)] + } + lw.lm <- lm(ZLXi ~ Z) + plot(Z, ZLXi, xlab="Z", ylab="WZ", pch = 20, cex = 0.33, col = "lightgrey", xlim = c(min(Z),max(Z)), ylim = c(min(ZLXi, b),max(ZLXi, b))) + abline(h = 0, lty = "dashed", col = "grey30") + abline(v = 0, lty = "dashed", col = "grey30") + abline(lw.lm, lty = "dotted", col = "grey40") + if(!plain) { + df_q1 <- data.frame(x_q1, b_q1) + df_q1 <- df_q1[order(x_q1),] + df_q2 <- data.frame(x_q2, b_q2) + df_q2 <- df_q2[order(x_q2),] + df_q3 <- data.frame(x_q3, b_q3) + df_q3 <- df_q3[order(x_q3),] + df_q4 <- data.frame(x_q4, b_q4) + df_q4 <- df_q4[order(x_q4),] + + for(i in 1:(length(df_q1$x_q1) - 1)) { + lines(c(df_q1$x_q1[i], df_q1$x_q1[i+1]), c(df_q1$b_q1[i], df_q1$b_q1[i+1]), type = "l", col = "firebrick") + } + + for(i in 1:length(df_q2$x_q2)) { + lines(c(df_q2$x_q2[i], df_q2$x_q2[i+1]), c(df_q2$b_q2[i], df_q2$b_q2[i+1]), type = "l", col = "royalblue") + } + + for(i in 1:length(df_q3$x_q3)) { + lines(c(df_q3$x_q3[i], df_q3$x_q3[i+1]), c(df_q3$b_q3[i], df_q3$b_q3[i+1]), type = "l", col = "firebrick") + } + + for(i in 1:length(df_q4$x_q4)) { + lines(c(df_q4$x_q4[i], df_q4$x_q4[i+1]), c(df_q4$b_q4[i], df_q4$b_q4[i+1]), type = "l", col = "royalblue") + } + } + } + if(return_df) { + res <- data.frame(z = Z, wz = ZLXi, b = b) + invisible(res) + } +} \ No newline at end of file From d418f0e762ca54b981429215d039ce6704c0dd71 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 00:37:49 +0200 Subject: [PATCH 02/27] Add files via upload Man files for moran.plot.drop and moran.plot.seismogram. --- man/moran.plot.drop.Rd | 57 ++++++++++++++++++++++++++++++++++++ man/moran.plot.seismogram.Rd | 53 +++++++++++++++++++++++++++++++++ 2 files changed, 110 insertions(+) create mode 100644 man/moran.plot.drop.Rd create mode 100644 man/moran.plot.seismogram.Rd diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd new file mode 100644 index 00000000..5c65fa67 --- /dev/null +++ b/man/moran.plot.drop.Rd @@ -0,0 +1,57 @@ +\name{moran.plot.drop} +\alias{moran.plot.drop} + +\title{Moran drop plot} + +\description{ +A version of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines indicating \emph{p} values for visual inspection of statistical significance. +} +\usage{ +moran.plot.drop(x, listw, nsim = 999, cv = 2.58, + significant = TRUE, plain = FALSE, zero.policy = FALSE, + xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, + spChk = NULL, labels = FALSE) +} + +\arguments{ + \item{x}{a numerical vector holding the attribute of interest} + \item{listw}{a \code{listw} spatial weights object} + \item{nsim}{default 999; number of conditonal permutation simulations} + \item{cv}{default 2.58; the desired critical value assuming approximate normality of standardised local Moran's \emph{I} values} + \item{significant}{default TRUE; a parameter indicating whether to display plot distances of significant (default) or non-significant observations} + \item{plain}{default FALSE; a plain Moran scatterplot is displayed if TRUE} + \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} + \item{xlab}{label for x axis} + \item{ylab}{label for y axis} + \item{plot}{default TRUE; if FALSE, plotting is suppressed} + \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} + \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} + \item{labels}{default FALSE; no labels are plotted by default; character labels for points are assigned to significant observations if provided, region IDs are used if set to TRUE} +} + +\details{ +The Moran drop plot is a version of the Moran scatterplot supplemented by indications of \emph{p} values. The standard Moran scatterplot provides indirect information about the effect size (distance from the trend line), but not about the level of confidence to determine whether the effects shown might be more than just random outcomes. The Moran drop plot marks significant points in red (positive) and blue (negative spatial autocorrelation) and adds so-called drop lines connecting the significant observations to the scatterplot positions of their associated critical values. The coordinates of the latter are determined under the assumption of an approximate standard normality of the z-scores of the local Moran's \emph{I} values. The longer the lines, the higher the associated \emph{p} values. The visualisation thus enables visual inspection of statistical significance while maintaining the relationship to both attribute value ranges and scatterplot quadrants. It is also possible to invert the visualised relationship and display the distances of non-significant observations to their corresponding critical values (if significant is set to FALSE). +} + +\value{ +When return_df is TRUE, a data frame object with the following members is returned: + \item{z}{the standardised attribute values} + \item{wz}{the standardised spatially lagged attribute values} + \item{b}{the y-coordinates of the critical values} + \item{line_lengths}{the absolute distances between b and z} +} + +\references{Westerholt, R. (2024): Extending the Moran scatterplot by indications of critical values and \emph{p}-values: introducing the Moran seismogram and the drop plot. Proceedings of the 32nd Annual GIS Research UK Conference (GISRUK), Leeds, UK. \url{https://doi.org/10.5281/zenodo.10897792}} + +\author{Rene Westerholt \email{rene.westerholt@tu-dortmund.de}} + +\seealso{\code{\link{moran.plot}}} + +\examples{ +# Boston example (CMEDV; owner-occupied housing in USD) +data(boston) +boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) +boston.nb <- poly2nb(boston.tr) +boston.listw <- nb2listw(boston.nb) +moran.plot.drop(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, significant = TRUE, plain = FALSE, labels = FALSE, plot = TRUE) +} \ No newline at end of file diff --git a/man/moran.plot.seismogram.Rd b/man/moran.plot.seismogram.Rd new file mode 100644 index 00000000..0cd40c08 --- /dev/null +++ b/man/moran.plot.seismogram.Rd @@ -0,0 +1,53 @@ +\name{moran.plot.seismogram} +\alias{moran.plot.seismogram} + +\title{Moran seismogram} + +\description{ +A variant of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines connecting location-wise critical values. The plot allows for visual inspection of potential spatial weights misspecifiation. +} +\usage{ +moran.plot.seismogram(x, listw, nsim = 999, cv = 2.58, + plain = FALSE, zero.policy = FALSE, xlab = NULL, + ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) +} + +\arguments{ + \item{x}{a numerical vector holding the attribute of interest} + \item{listw}{a \code{listw} spatial weights object} + \item{nsim}{default 999; number of conditonal permutation simulations} + \item{cv}{default 2.58; the desired critical value assuming approximate normality of standardised local Moran's \emph{I} values} + \item{plain}{default FALSE; a plain Moran scatterplot is displayed if TRUE} + \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} + \item{xlab}{label for x axis} + \item{ylab}{label for y axis} + \item{plot}{default TRUE; if FALSE, plotting is suppressed} + \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} + \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} +} + +\details{ +The Moran seismogram is a version of the Moran scatterplot that is complemented by lines connecting the location-specific critical values of local Moran's \emph{I}. Lines in quadrants with positive spatial autocorrelation are shown in red and lines in quadrants with negative spatial autocorrelation are shown in blue. The representation is similar to a seismogram for detecting earthquakes and thus reveals potentially suspicious local spatial weights configurations by visualising spikes. The latter are displayed in an integrated manner with their positions in the attribute value range and in connection with the types of the associated spatial patterns (by the quadrants of the scatterplot). +} + +\value{ +When return_df is TRUE, a data frame object with the following members is returned: + \item{z}{the standardised attribute values} + \item{wz}{the standardised spatially lagged attribute values} + \item{b}{the y-coordinates of the critical values} +} + +\references{Westerholt, R. (2024): Extending the Moran scatterplot by indications of critical values and \emph{p}-values: introducing the Moran seismogram and the drop plot. Proceedings of the 32nd Annual GIS Research UK Conference (GISRUK), Leeds, UK. \url{https://doi.org/10.5281/zenodo.10897792}} + +\author{Rene Westerholt \email{rene.westerholt@tu-dortmund.de}} + +\seealso{\code{\link{moran.plot}}} + +\examples{ +# Boston example (CMEDV; owner-occupied housing in USD) +data(boston) +boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) +boston.nb <- poly2nb(boston.tr) +boston.listw <- nb2listw(boston.nb) +moran.plot.seismogram(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, plain = FALSE, plot = TRUE) +} \ No newline at end of file From 094805b2d3310310146a82aae36793abf73a60e2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 00:39:46 +0200 Subject: [PATCH 03/27] Update moran.plot.Rd Fixed a typo (y axis instead of x axis for item "ylab"). --- man/moran.plot.Rd | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/man/moran.plot.Rd b/man/moran.plot.Rd index 8175990b..7ee10c22 100644 --- a/man/moran.plot.Rd +++ b/man/moran.plot.Rd @@ -18,7 +18,7 @@ moran.plot(x, listw, zero.policy=attr(listw, "zero.policy"), spChk=NULL, labels= \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} \item{labels}{character labels for points with high influence measures, if set to FALSE, no labels are plotted for points with large influence} \item{xlab}{label for x axis} - \item{ylab}{label for x axis} + \item{ylab}{label for y axis} \item{quiet}{default NULL, use !verbose global option value; if TRUE, output of summary of influence object suppressed} \item{plot}{default TRUE, if false, plotting is suppressed} \item{return_df}{default TRUE, invisibly return a data.frame object; if FALSE invisibly return an influence measures object} From d4e43ea01e327af6a200294cada0ae9f9bcdfab9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 00:52:20 +0200 Subject: [PATCH 04/27] Update moran.plot.drop.R --- R/moran.plot.drop.R | 9 +-------- 1 file changed, 1 insertion(+), 8 deletions(-) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index 796e9e93..e44d9c3a 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -1,5 +1,4 @@ moran.plot.drop <- function(x, listw, nsim = 999, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) { - require(spdep) if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) stopifnot(is.vector(x)) @@ -100,16 +99,10 @@ moran.plot.drop <- function(x, listw, nsim = 999, cv = 2.58, significant = TRUE, if (labs && length(x_q4) > 0) text(x_q4[i], y_q4[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") } - if (zero.policy) { - n0 <- ZLXi == 0 - if (any(n0)) { - symbols(x[n0], wx[n0], inches = FALSE, circles = rep(diff(range(x))/50, length(which(n0))), bg = "grey", add = TRUE) - } - } } } if(return_df) { res <- data.frame(z = Z, wz = ZLXi, b = b, line_lengths = abs(ZLXi - b)) invisible(res) } -} \ No newline at end of file +} From d6d3343d9b93ce831731ed06eb5e8b616039ec77 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 00:52:47 +0200 Subject: [PATCH 05/27] Update moran.plot.seismogram.R --- R/moran.plot.seismogram.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R index f5023d73..3638dc72 100644 --- a/R/moran.plot.seismogram.R +++ b/R/moran.plot.seismogram.R @@ -1,5 +1,4 @@ moran.plot.seismogram <- function(x, listw, nsim = 999, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) { - require(spdep) if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) stopifnot(is.vector(x)) @@ -82,4 +81,4 @@ moran.plot.seismogram <- function(x, listw, nsim = 999, cv = 2.58, plain = FALSE res <- data.frame(z = Z, wz = ZLXi, b = b) invisible(res) } -} \ No newline at end of file +} From 8daea5289de2a4eafe407538084eaf39ac077991 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 01:09:45 +0200 Subject: [PATCH 06/27] Update NAMESPACE --- NAMESPACE | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 0947d07b..9fc5eaf3 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -42,8 +42,9 @@ export(gabrielneigh, geary.test, geary, geary.mc, globalG.test, graph2nb, export(listw2sn, sn2listw, read.gwt2nb, write.sn2gwt, lm.LMtests, lm.RStests, SD.RStests, lm.morantest, localG, localG_perm, localmoran, localmoran_perm, moran, - moran.test, moran.mc, moran.plot, localmoran.sad, lm.morantest.sad, - nb2listw, nb2listwdist, nb2mat, listw2mat, mat2listw, nbdists, nblag, + moran.test, moran.mc, moran.plot, moran.plot.drop, moran.plot.seismogram, + localmoran.sad, lm.morantest.sad, nb2listw, nb2listwdist, nb2mat, + listw2mat, mat2listw, nbdists, nblag, nblag_cumul, poly2nb, read.gal, write.nb.gal, read.geoda, relativeneigh, soi.graph, sp.correlogram, sp.mantel.mc, set.spChkOption, chkIDs, get.spChkOption, spNamedVec, tri2nb, From fb627233945037de530774f880de09204263dc76 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 01:10:32 +0200 Subject: [PATCH 07/27] Update moran.plot.drop.Rd --- man/moran.plot.drop.Rd | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd index 5c65fa67..bb5332ef 100644 --- a/man/moran.plot.drop.Rd +++ b/man/moran.plot.drop.Rd @@ -10,7 +10,7 @@ A version of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented b moran.plot.drop(x, listw, nsim = 999, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, - spChk = NULL, labels = FALSE) + spChk = NULL, labels = NULL) } \arguments{ @@ -26,7 +26,7 @@ moran.plot.drop(x, listw, nsim = 999, cv = 2.58, \item{plot}{default TRUE; if FALSE, plotting is suppressed} \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} - \item{labels}{default FALSE; no labels are plotted by default; character labels for points are assigned to significant observations if provided, region IDs are used if set to TRUE} + \item{labels}{default NULL; no labels are plotted by default; character labels for points are assigned to significant observations if provided, region IDs are used if set to TRUE} } \details{ @@ -53,5 +53,6 @@ data(boston) boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) boston.nb <- poly2nb(boston.tr) boston.listw <- nb2listw(boston.nb) -moran.plot.drop(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, significant = TRUE, plain = FALSE, labels = FALSE, plot = TRUE) -} \ No newline at end of file +moran.plot.drop(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, + significant = TRUE, plain = FALSE, labels = FALSE, plot = TRUE) +} From eee79be060af8a6d5e58b3dd0d7acfde444faf85 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 14 Apr 2024 01:11:06 +0200 Subject: [PATCH 08/27] Update moran.plot.seismogram.Rd --- man/moran.plot.seismogram.Rd | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/man/moran.plot.seismogram.Rd b/man/moran.plot.seismogram.Rd index 0cd40c08..13dce4eb 100644 --- a/man/moran.plot.seismogram.Rd +++ b/man/moran.plot.seismogram.Rd @@ -49,5 +49,6 @@ data(boston) boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) boston.nb <- poly2nb(boston.tr) boston.listw <- nb2listw(boston.nb) -moran.plot.seismogram(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, plain = FALSE, plot = TRUE) -} \ No newline at end of file +moran.plot.seismogram(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, + plain = FALSE, plot = TRUE) +} From 5c104b479b4d30ebab28677f140740793a7ad28a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Tue, 23 Apr 2024 17:56:38 +0200 Subject: [PATCH 09/27] Update moran.plot.drop.R moran.plot.drop is now a function of a fitted object of type locmoran. --- R/moran.plot.drop.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index e44d9c3a..a10cb6c7 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -1,4 +1,4 @@ -moran.plot.drop <- function(x, listw, nsim = 999, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) { +moran.plot.drop <- function(x, listw, locmoran, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) { if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) stopifnot(is.vector(x)) @@ -31,7 +31,6 @@ moran.plot.drop <- function(x, listw, nsim = 999, cv = 2.58, significant = TRUE, ylab <- paste("spatially lagged", xname) Z <- as.vector(scale(x)) ZLXi <- lag.listw(listw, Z, zero.policy = zero.policy) - locmoran <- localmoran_perm(Z, listw, nsim = nsim) ZIi <- locmoran[, 4] b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) From e805e18f96b923bb1c91709521171e742952e531 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Tue, 23 Apr 2024 17:57:32 +0200 Subject: [PATCH 10/27] Update moran.plot.seismogram.R moran.plot.seismogram is now a function of a fitted object of type locmoran. --- R/moran.plot.seismogram.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R index 3638dc72..d95884c3 100644 --- a/R/moran.plot.seismogram.R +++ b/R/moran.plot.seismogram.R @@ -1,4 +1,4 @@ -moran.plot.seismogram <- function(x, listw, nsim = 999, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) { +moran.plot.seismogram <- function(x, listw, locmoran, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) { if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) stopifnot(is.vector(x)) @@ -25,7 +25,6 @@ moran.plot.seismogram <- function(x, listw, nsim = 999, cv = 2.58, plain = FALSE ylab <- paste("spatially lagged", xname) Z <- as.vector(scale(x)) ZLXi <- lag.listw(listw, Z, zero.policy = zero.policy) - locmoran <- localmoran_perm(Z, listw, nsim = nsim) ZIi <- locmoran[, 4] b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) From 60e347fed6f49ead03ab3dee69ea61ebfa8ee39f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Tue, 23 Apr 2024 17:59:11 +0200 Subject: [PATCH 11/27] Update moran.plot.drop.Rd --- man/moran.plot.drop.Rd | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd index bb5332ef..1bed1b75 100644 --- a/man/moran.plot.drop.Rd +++ b/man/moran.plot.drop.Rd @@ -7,7 +7,7 @@ A version of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines indicating \emph{p} values for visual inspection of statistical significance. } \usage{ -moran.plot.drop(x, listw, nsim = 999, cv = 2.58, +moran.plot.drop(x, listw, locmoran, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) @@ -16,7 +16,7 @@ moran.plot.drop(x, listw, nsim = 999, cv = 2.58, \arguments{ \item{x}{a numerical vector holding the attribute of interest} \item{listw}{a \code{listw} spatial weights object} - \item{nsim}{default 999; number of conditonal permutation simulations} + \item{locmoran}{a fitted object of type localmoran} \item{cv}{default 2.58; the desired critical value assuming approximate normality of standardised local Moran's \emph{I} values} \item{significant}{default TRUE; a parameter indicating whether to display plot distances of significant (default) or non-significant observations} \item{plain}{default FALSE; a plain Moran scatterplot is displayed if TRUE} @@ -53,6 +53,6 @@ data(boston) boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) boston.nb <- poly2nb(boston.tr) boston.listw <- nb2listw(boston.nb) -moran.plot.drop(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, - significant = TRUE, plain = FALSE, labels = FALSE, plot = TRUE) +moran.plot.drop(boston.c$CMEDV, boston.listw, localmoran(boston.c$CMEDV, boston.listw), 2.58, + zero.policy = TRUE, significant = TRUE, plain = FALSE, labels = FALSE, plot = TRUE) } From a72808e01e3dd21e1c0f481952dadb996b85b139 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Tue, 23 Apr 2024 18:00:32 +0200 Subject: [PATCH 12/27] Update moran.plot.seismogram.Rd --- man/moran.plot.seismogram.Rd | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/man/moran.plot.seismogram.Rd b/man/moran.plot.seismogram.Rd index 13dce4eb..fb40853b 100644 --- a/man/moran.plot.seismogram.Rd +++ b/man/moran.plot.seismogram.Rd @@ -7,7 +7,7 @@ A variant of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines connecting location-wise critical values. The plot allows for visual inspection of potential spatial weights misspecifiation. } \usage{ -moran.plot.seismogram(x, listw, nsim = 999, cv = 2.58, +moran.plot.seismogram(x, listw, locmoran, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) } @@ -15,7 +15,7 @@ moran.plot.seismogram(x, listw, nsim = 999, cv = 2.58, \arguments{ \item{x}{a numerical vector holding the attribute of interest} \item{listw}{a \code{listw} spatial weights object} - \item{nsim}{default 999; number of conditonal permutation simulations} + \item{locmoran}{a fitted object of type localmoran} \item{cv}{default 2.58; the desired critical value assuming approximate normality of standardised local Moran's \emph{I} values} \item{plain}{default FALSE; a plain Moran scatterplot is displayed if TRUE} \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} @@ -49,6 +49,6 @@ data(boston) boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) boston.nb <- poly2nb(boston.tr) boston.listw <- nb2listw(boston.nb) -moran.plot.seismogram(boston.c$CMEDV, boston.listw, 999, 2.58, zero.policy = TRUE, - plain = FALSE, plot = TRUE) +moran.plot.seismogram(boston.c$CMEDV, boston.listw, localmoran(boston.c$CMEDV, boston.listw), 2.58, +zero.policy = TRUE, plain = FALSE, plot = TRUE) } From 0f36086f8f4fcc115deeeacbed4aeba65ef84373 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Tue, 23 Apr 2024 18:03:29 +0200 Subject: [PATCH 13/27] Update moran.plot.drop.R --- R/moran.plot.drop.R | 2 ++ 1 file changed, 2 insertions(+) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index a10cb6c7..00d729a5 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -1,6 +1,8 @@ moran.plot.drop <- function(x, listw, locmoran, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) { if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) + if (!inherits(locmoran, "localmoran")) + stop(paste(deparse(substitute(locmoran)), "is not a localmoran object")) stopifnot(is.vector(x)) stopifnot(is.logical(significant)) stopifnot(is.logical(plain)) From bb3b688b053e14a48cd98a454b1de4f00071fbd4 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Tue, 23 Apr 2024 18:03:44 +0200 Subject: [PATCH 14/27] Update moran.plot.seismogram.R --- R/moran.plot.seismogram.R | 2 ++ 1 file changed, 2 insertions(+) diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R index d95884c3..225f070f 100644 --- a/R/moran.plot.seismogram.R +++ b/R/moran.plot.seismogram.R @@ -1,6 +1,8 @@ moran.plot.seismogram <- function(x, listw, locmoran, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) { if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) + if (!inherits(locmoran, "localmoran")) + stop(paste(deparse(substitute(locmoran)), "is not a localmoran object")) stopifnot(is.vector(x)) stopifnot(is.logical(plain)) stopifnot(is.logical(return_df)) From 595fad5cf383d49a56a8952a117a54d8b11a5716 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sat, 18 Apr 2026 20:32:25 +0200 Subject: [PATCH 15/27] Update moran.plot.drop.R Updated version of the drop plot implementation. --- R/moran.plot.drop.R | 160 +++++++++++++++++++++++--------------------- 1 file changed, 85 insertions(+), 75 deletions(-) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index 00d729a5..cdd594fc 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -1,23 +1,31 @@ -moran.plot.drop <- function(x, listw, locmoran, cv = 2.58, significant = TRUE, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL, labels = NULL) { +moran.plot.drop <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, significant = TRUE, xlab = NULL, ylab = NULL, return_df = TRUE, spChk = NULL, labels = NULL, zero.policy=attr(listw, "zero.policy")) { if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) if (!inherits(locmoran, "localmoran")) stop(paste(deparse(substitute(locmoran)), "is not a localmoran object")) stopifnot(is.vector(x)) stopifnot(is.logical(significant)) - stopifnot(is.logical(plain)) stopifnot(is.logical(return_df)) - if (is.null(zero.policy)) - zero.policy <- get("zeroPolicy", envir = .spdepOptions) + stopifnot(is.numeric(alpha)) + if (is.null(zero.policy)) + zero.policy <- get.ZeroPolicyOption() stopifnot(is.logical(zero.policy)) xname <- deparse(substitute(x)) if (!is.numeric(x)) stop(paste(xname, "is not a numeric vector")) - if (any(is.na(x))) + if (anyNA(x)) stop("NA in X") n <- length(listw$neighbours) if (n != length(x)) stop("objects of different length") + usePadj <- FALSE + if (!is.null(adjusted_p)) { + if(n != length(adjusted_p)) + stop("objects of different length") + if (anyNA(adjusted_p)) + stop("NA in vector of adjusted p values") + usePadj <- TRUE + } if (is.null(spChk)) spChk <- get.spChkOption() if (spChk && !chkIDs(x, listw)) @@ -27,83 +35,85 @@ moran.plot.drop <- function(x, listw, locmoran, cv = 2.58, significant = TRUE, p labs <- FALSE if (is.null(labels) || length(labels) != n) labels <- as.character(attr(listw, "region.id")) - if (is.null(xlab)) + if (is.null(xlab)) xlab <- xname if (is.null(ylab)) ylab <- paste("spatially lagged", xname) - Z <- as.vector(scale(x)) - ZLXi <- lag.listw(listw, Z, zero.policy = zero.policy) - ZIi <- locmoran[, 4] - b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) - b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) - b[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] <- b2[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] - if(plot) { - if(!plain) { - if(significant) { - x_q1 <- Z[which(Z > 0 & ZLXi > 0 & ZIi >= cv)] - y_q1 <- ZLXi[which(Z > 0 & ZLXi > 0 & ZIi >= cv)] - b_q1 <- b[which(Z > 0 & ZLXi > 0 & ZIi >= cv)] - x_q2 <- Z[which(Z > 0 & ZLXi < 0 & ZIi <= (-1) * cv)] - y_q2 <- ZLXi[which(Z > 0 & ZLXi < 0 & ZIi <= (-1) * cv)] - b_q2 <- b[which(Z > 0 & ZLXi < 0 & ZIi <= (-1) * cv)] - x_q3 <- Z[which(Z < 0 & ZLXi < 0 & ZIi >= cv)] - y_q3 <- ZLXi[which(Z < 0 & ZLXi < 0 & ZIi >= cv)] - b_q3 <- b[which(Z < 0 & ZLXi < 0 & ZIi >= cv)] - x_q4 <- Z[which(Z < 0 & ZLXi > 0 & ZIi <= (-1) * cv)] - y_q4 <- ZLXi[which(Z < 0 & ZLXi > 0 & ZIi <= (-1) * cv)] - b_q4 <- b[which(Z < 0 & ZLXi > 0 & ZIi <= (-1) * cv)] - } else { - x_q1 <- Z[which(Z > 0 & ZLXi > 0 & ZIi < cv)] - y_q1 <- ZLXi[which(Z > 0 & ZLXi > 0 & ZIi < cv)] - b_q1 <- b[which(Z > 0 & ZLXi > 0 & ZIi < cv)] - x_q2 <- Z[which(Z > 0 & ZLXi < 0 & ZIi > (-1) * cv)] - y_q2 <- ZLXi[which(Z > 0 & ZLXi < 0 & ZIi > (-1) * cv)] - b_q2 <- b[which(Z > 0 & ZLXi < 0 & ZIi > (-1) * cv)] - x_q3 <- Z[which(Z < 0 & ZLXi < 0 & ZIi < cv)] - y_q3 <- ZLXi[which(Z < 0 & ZLXi < 0 & ZIi < cv)] - b_q3 <- b[which(Z < 0 & ZLXi < 0 & ZIi < cv)] - x_q4 <- Z[which(Z < 0 & ZLXi > 0 & ZIi > (-1) * cv)] - y_q4 <- ZLXi[which(Z < 0 & ZLXi > 0 & ZIi > (-1) * cv)] - b_q4 <- b[which(Z < 0 & ZLXi > 0 & ZIi > (-1) * cv)] - } - } - lw.lm <- lm(ZLXi ~ Z) - plot(Z, ZLXi, xlab="Z", ylab="WZ", pch = 20, cex = 0.33, col = "lightgrey", xlim = c(min(Z),max(Z)), ylim = c(min(ZLXi, b),max(ZLXi, b))) - abline(h = 0, lty = "dashed", col = "grey30") - abline(v = 0, lty = "dashed", col = "grey30") - abline(lw.lm, lty = "dotted", col = "grey40") - if(!plain) { - for(i in 1:length(x_q1)) { - lines(c(x_q1[i], x_q1[i]), c(y_q1[i], b_q1[i]), type = "l", col = "lightpink", lwd = 1) - points(c(x_q1[i]), c(y_q1[i]), pch = 20, cex = 0.5, col = "firebrick") - if (labs && length(x_q1) > 0) - text(x_q1[i], y_q1[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") - } + Z <- as.vector(scale(x, scale = F)) + WZ <- lag.listw(listw, Z, zero.policy = zero.policy) + if (anyNA(WZ)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") + if(usePadj) + cv <- min(abs(locmoran[(which(adjusted_p <= alpha/2)), 4])) + else + cv <- qnorm(1 - alpha, lower.tail = FALSE) + b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / Z + b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / Z + b[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] <- b2[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] + + if(significant) { + x_q1 <- Z[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] >= cv)] + y_q1 <- WZ[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] >= cv)] + b_q1 <- b[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] >= cv)] + x_q2 <- Z[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] <= (-1) * cv)] + y_q2 <- WZ[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] <= (-1) * cv)] + b_q2 <- b[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] <= (-1) * cv)] + x_q3 <- Z[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] >= cv)] + y_q3 <- WZ[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] >= cv)] + b_q3 <- b[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] >= cv)] + x_q4 <- Z[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] <= (-1) * cv)] + y_q4 <- WZ[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] <= (-1) * cv)] + b_q4 <- b[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] <= (-1) * cv)] + } else { + x_q1 <- Z[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] < cv)] + y_q1 <- WZ[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] < cv)] + b_q1 <- b[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] < cv)] + x_q2 <- Z[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] > (-1) * cv)] + y_q2 <- WZ[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] > (-1) * cv)] + b_q2 <- b[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] > (-1) * cv)] + x_q3 <- Z[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] < cv)] + y_q3 <- WZ[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] < cv)] + b_q3 <- b[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] < cv)] + x_q4 <- Z[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] > (-1) * cv)] + y_q4 <- WZ[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] > (-1) * cv)] + b_q4 <- b[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] > (-1) * cv)] + } + + lw.lm <- lm(WZ ~ Z) + plot(Z, WZ, xlab="X_centred", ylab="WX", pch = 20, cex = 0.33, col = "gray70", xlim = c(min(Z),max(Z)), ylim = c(min(WZ, b),max(WZ, b))) + abline(h = mean(WZ), lty = "dashed", col = "grey30") + abline(v = 0, lty = "dashed", col = "grey30") + abline(lw.lm, lty = "dotted", col = "grey40") + + for(i in 1:length(x_q1)) { + lines(c(x_q1[i], x_q1[i]), c(y_q1[i], b_q1[i]), type = "l", col = "lightpink", lwd = 1) + points(c(x_q1[i]), c(y_q1[i]), pch = 20, cex = 0.5, col = "firebrick") + if (labs && length(x_q1) > 0) + text(x_q1[i], y_q1[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") + } - for(i in 1:length(x_q2)) { - lines(c(x_q2[i], x_q2[i]), c(y_q2[i], b_q2[i]), type = "l", col = "lightblue1", lwd = 1) - points(c(x_q2[i]), c(y_q2[i]), pch = 20, cex = 0.5, col = "royalblue") - if (labs && length(x_q2) > 0) - text(x_q2[i], y_q2[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") - } - - for(i in 1:length(x_q3)) { - lines(c(x_q3[i], x_q3[i]), c(y_q3[i], b_q3[i]), type = "l", col = "lightpink", lwd = 1) - points(c(x_q3[i]), c(y_q3[i]), pch = 20, cex = 0.5, col = "firebrick") - if (labs && length(x_q3) > 0) - text(x_q3[i], y_q3[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") - } + for(i in 1:length(x_q2)) { + lines(c(x_q2[i], x_q2[i]), c(y_q2[i], b_q2[i]), type = "l", col = "lightblue1", lwd = 1) + points(c(x_q2[i]), c(y_q2[i]), pch = 20, cex = 0.5, col = "royalblue") + if (labs && length(x_q2) > 0) + text(x_q2[i], y_q2[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") + } - for(i in 1:length(x_q4)) { - lines(c(x_q4[i], x_q4[i]), c(y_q4[i], b_q4[i]), type = "l", col = "lightblue1", lwd = 1) - points(c(x_q4[i]), c(y_q4[i]), pch = 20, cex = 0.5, col = "royalblue") - if (labs && length(x_q4) > 0) - text(x_q4[i], y_q4[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") - } - } + for(i in 1:length(x_q3)) { + lines(c(x_q3[i], x_q3[i]), c(y_q3[i], b_q3[i]), type = "l", col = "lightpink", lwd = 1) + points(c(x_q3[i]), c(y_q3[i]), pch = 20, cex = 0.5, col = "firebrick") + if (labs && length(x_q3) > 0) + text(x_q3[i], y_q3[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") + } + + for(i in 1:length(x_q4)) { + lines(c(x_q4[i], x_q4[i]), c(y_q4[i], b_q4[i]), type = "l", col = "lightblue1", lwd = 1) + points(c(x_q4[i]), c(y_q4[i]), pch = 20, cex = 0.5, col = "royalblue") + if (labs && length(x_q4) > 0) + text(x_q4[i], y_q4[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") } + if(return_df) { - res <- data.frame(z = Z, wz = ZLXi, b = b, line_lengths = abs(ZLXi - b)) + res <- data.frame(z = Z, wz = WZ, b = b, line_lengths = abs(WZ - b)) #n_sig = length(x_q1) + length(x_q2) + length(x_q3) + length(x_q4) invisible(res) } } From 779056f01f94007a843b586dd19dfb83c5088980 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sat, 18 Apr 2026 20:33:09 +0200 Subject: [PATCH 16/27] Update moran.plot.seismogram.R Updated version of the Moran seismogram. --- R/moran.plot.seismogram.R | 116 ++++++++++++++++++++------------------ 1 file changed, 61 insertions(+), 55 deletions(-) diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R index 225f070f..d58de201 100644 --- a/R/moran.plot.seismogram.R +++ b/R/moran.plot.seismogram.R @@ -1,13 +1,13 @@ -moran.plot.seismogram <- function(x, listw, locmoran, cv = 2.58, plain = FALSE, zero.policy = FALSE, xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) { +moran.plot.seismogram <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, xlab = NULL, ylab = NULL, return_df = TRUE, spChk = NULL, zero.policy=attr(listw, "zero.policy")) { if (!inherits(listw, "listw")) stop(paste(deparse(substitute(listw)), "is not a listw object")) if (!inherits(locmoran, "localmoran")) stop(paste(deparse(substitute(locmoran)), "is not a localmoran object")) stopifnot(is.vector(x)) - stopifnot(is.logical(plain)) stopifnot(is.logical(return_df)) - if (is.null(zero.policy)) - zero.policy <- get("zeroPolicy", envir = .spdepOptions) + stopifnot(is.numeric(alpha)) + if (is.null(zero.policy)) + zero.policy <- get.ZeroPolicyOption() stopifnot(is.logical(zero.policy)) xname <- deparse(substitute(x)) if (!is.numeric(x)) @@ -17,6 +17,14 @@ moran.plot.seismogram <- function(x, listw, locmoran, cv = 2.58, plain = FALSE, n <- length(listw$neighbours) if (n != length(x)) stop("objects of different length") + usePadj <- FALSE + if (!is.null(adjusted_p)) { + if(n != length(adjusted_p)) + stop("objects of different length") + if (any(is.na(adjusted_p))) + stop("NA in vector of adjusted p values") + usePadj <- TRUE + } if (is.null(spChk)) spChk <- get.spChkOption() if (spChk && !chkIDs(x, listw)) @@ -25,61 +33,59 @@ moran.plot.seismogram <- function(x, listw, locmoran, cv = 2.58, plain = FALSE, xlab <- xname if (is.null(ylab)) ylab <- paste("spatially lagged", xname) - Z <- as.vector(scale(x)) - ZLXi <- lag.listw(listw, Z, zero.policy = zero.policy) - ZIi <- locmoran[, 4] + Z <- as.vector(scale(x, scale = F)) + WZ <- lag.listw(listw, Z, zero.policy = zero.policy) + if (anyNA(WZ)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") + if(usePadj) + cv <- min(abs(locmoran[(which(adjusted_p <= alpha/2)),4])) + else + cv <- qnorm(1 - alpha, lower.tail = FALSE) b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) - b[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] <- b2[which((Z < 0 & ZLXi > 0) | (Z > 0 & ZLXi < 0))] - if(plot) { - if(!plain) { - x_q1 <- Z[which(Z > 0 & ZLXi > 0)] - y_q1 <- ZLXi[which(Z > 0 & ZLXi > 0)] - b_q1 <- b[which(Z > 0 & ZLXi > 0)] - x_q2 <- Z[which(Z > 0 & ZLXi < 0)] - y_q2 <- ZLXi[which(Z > 0 & ZLXi < 0)] - b_q2 <- b[which(Z > 0 & ZLXi < 0)] - x_q3 <- Z[which(Z < 0 & ZLXi < 0)] - y_q3 <- ZLXi[which(Z < 0 & ZLXi < 0)] - b_q3 <- b[which(Z < 0 & ZLXi < 0)] - x_q4 <- Z[which(Z < 0 & ZLXi > 0)] - y_q4 <- ZLXi[which(Z < 0 & ZLXi > 0)] - b_q4 <- b[which(Z < 0 & ZLXi > 0)] - } - lw.lm <- lm(ZLXi ~ Z) - plot(Z, ZLXi, xlab="Z", ylab="WZ", pch = 20, cex = 0.33, col = "lightgrey", xlim = c(min(Z),max(Z)), ylim = c(min(ZLXi, b),max(ZLXi, b))) - abline(h = 0, lty = "dashed", col = "grey30") - abline(v = 0, lty = "dashed", col = "grey30") - abline(lw.lm, lty = "dotted", col = "grey40") - if(!plain) { - df_q1 <- data.frame(x_q1, b_q1) - df_q1 <- df_q1[order(x_q1),] - df_q2 <- data.frame(x_q2, b_q2) - df_q2 <- df_q2[order(x_q2),] - df_q3 <- data.frame(x_q3, b_q3) - df_q3 <- df_q3[order(x_q3),] - df_q4 <- data.frame(x_q4, b_q4) - df_q4 <- df_q4[order(x_q4),] - - for(i in 1:(length(df_q1$x_q1) - 1)) { - lines(c(df_q1$x_q1[i], df_q1$x_q1[i+1]), c(df_q1$b_q1[i], df_q1$b_q1[i+1]), type = "l", col = "firebrick") - } - - for(i in 1:length(df_q2$x_q2)) { - lines(c(df_q2$x_q2[i], df_q2$x_q2[i+1]), c(df_q2$b_q2[i], df_q2$b_q2[i+1]), type = "l", col = "royalblue") - } + b[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] <- b2[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] + + x_q1 <- Z[which(Z > 0 & WZ > mean(WZ))] + y_q1 <- WZ[which(Z > 0 & WZ > mean(WZ))] + b_q1 <- b[which(Z > 0 & WZ > mean(WZ))] + x_q2 <- Z[which(Z > 0 & WZ < mean(WZ))] + y_q2 <- WZ[which(Z > 0 & WZ < mean(WZ))] + b_q2 <- b[which(Z > 0 & WZ < mean(WZ))] + x_q3 <- Z[which(Z < 0 & WZ < mean(WZ))] + y_q3 <- WZ[which(Z < 0 & WZ < mean(WZ))] + b_q3 <- b[which(Z < 0 & WZ < mean(WZ))] + x_q4 <- Z[which(Z < 0 & WZ > mean(WZ))] + y_q4 <- WZ[which(Z < 0 & WZ > mean(WZ))] + b_q4 <- b[which(Z < 0 & WZ > mean(WZ))] + + lw.lm <- lm(WZ ~ Z) + plot(Z, WZ, xlab="X_centred", ylab="WX", pch = 20, cex = 0.33, col = "gray70", xlim = c(min(Z),max(Z)), ylim = c(min(WZ, b),max(WZ, b))) + abline(h = mean(WZ), lty = "dashed", col = "grey30") + abline(v = 0, lty = "dashed", col = "grey30") + abline(lw.lm, lty = "dotted", col = "grey40") + + df_q1 <- data.frame(x_q1, b_q1) + df_q1 <- df_q1[order(x_q1),] + df_q2 <- data.frame(x_q2, b_q2) + df_q2 <- df_q2[order(x_q2),] + df_q3 <- data.frame(x_q3, b_q3) + df_q3 <- df_q3[order(x_q3),] + df_q4 <- data.frame(x_q4, b_q4) + df_q4 <- df_q4[order(x_q4),] + + for(i in 1:(length(df_q1$x_q1) - 1)) + lines(c(df_q1$x_q1[i], df_q1$x_q1[i+1]), c(df_q1$b_q1[i], df_q1$b_q1[i+1]), type = "l", col = "firebrick") + + for(i in 1:length(df_q2$x_q2)) + lines(c(df_q2$x_q2[i], df_q2$x_q2[i+1]), c(df_q2$b_q2[i], df_q2$b_q2[i+1]), type = "l", col = "royalblue") - for(i in 1:length(df_q3$x_q3)) { - lines(c(df_q3$x_q3[i], df_q3$x_q3[i+1]), c(df_q3$b_q3[i], df_q3$b_q3[i+1]), type = "l", col = "firebrick") - } - - for(i in 1:length(df_q4$x_q4)) { - lines(c(df_q4$x_q4[i], df_q4$x_q4[i+1]), c(df_q4$b_q4[i], df_q4$b_q4[i+1]), type = "l", col = "royalblue") - } - } - } + for(i in 1:length(df_q3$x_q3)) + lines(c(df_q3$x_q3[i], df_q3$x_q3[i+1]), c(df_q3$b_q3[i], df_q3$b_q3[i+1]), type = "l", col = "firebrick") + + for(i in 1:length(df_q4$x_q4)) + lines(c(df_q4$x_q4[i], df_q4$x_q4[i+1]), c(df_q4$b_q4[i], df_q4$b_q4[i+1]), type = "l", col = "royalblue") + if(return_df) { - res <- data.frame(z = Z, wz = ZLXi, b = b) + res <- data.frame(z = Z, wz = WZ, b = b) invisible(res) } } From 69728a7b749ba0781264b2a9954eb8eec7694079 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sat, 18 Apr 2026 23:42:59 +0200 Subject: [PATCH 17/27] Update moran.plot.drop.Rd Updated documentation. --- man/moran.plot.drop.Rd | 27 +++++++++++++-------------- 1 file changed, 13 insertions(+), 14 deletions(-) diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd index 1bed1b75..c1c16e03 100644 --- a/man/moran.plot.drop.Rd +++ b/man/moran.plot.drop.Rd @@ -7,34 +7,33 @@ A version of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines indicating \emph{p} values for visual inspection of statistical significance. } \usage{ -moran.plot.drop(x, listw, locmoran, cv = 2.58, - significant = TRUE, plain = FALSE, zero.policy = FALSE, - xlab = NULL, ylab = NULL, plot = TRUE, return_df = TRUE, - spChk = NULL, labels = NULL) +moran.plot.drop(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, + significant = TRUE, xlab = NULL, ylab = NULL, return_df = TRUE, + spChk = NULL, labels = NULL), zero.policy=attr(listw, "zero.policy") } \arguments{ \item{x}{a numerical vector holding the attribute of interest} \item{listw}{a \code{listw} spatial weights object} \item{locmoran}{a fitted object of type localmoran} - \item{cv}{default 2.58; the desired critical value assuming approximate normality of standardised local Moran's \emph{I} values} - \item{significant}{default TRUE; a parameter indicating whether to display plot distances of significant (default) or non-significant observations} - \item{plain}{default FALSE; a plain Moran scatterplot is displayed if TRUE} - \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} + \item{alpha}{default 0.05; the desired significance level regarding local Moran's \emph{I} values} + \item{adjusted_p}{default NULL; a vector of \emph{p} values adjusted to account for multiple testing as is returned by \code{\link{p.adjustSP}}; standard normal distribution is used to determine critical values, if NULL} + \item{significant}{default TRUE; a parameter indicating whether to display plot critical value distances of significant (default) or non-significant observations} \item{xlab}{label for x axis} \item{ylab}{label for y axis} - \item{plot}{default TRUE; if FALSE, plotting is suppressed} \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} - \item{labels}{default NULL; no labels are plotted by default; character labels for points are assigned to significant observations if provided, region IDs are used if set to TRUE} + \item{labels}{default NULL; no labels are plotted by default; region IDs are used as labels of significant observations if set to TRUE} + \item{zero.policy}{default option stored in the listw object; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} } \details{ -The Moran drop plot is a version of the Moran scatterplot supplemented by indications of \emph{p} values. The standard Moran scatterplot provides indirect information about the effect size (distance from the trend line), but not about the level of confidence to determine whether the effects shown might be more than just random outcomes. The Moran drop plot marks significant points in red (positive) and blue (negative spatial autocorrelation) and adds so-called drop lines connecting the significant observations to the scatterplot positions of their associated critical values. The coordinates of the latter are determined under the assumption of an approximate standard normality of the z-scores of the local Moran's \emph{I} values. The longer the lines, the higher the associated \emph{p} values. The visualisation thus enables visual inspection of statistical significance while maintaining the relationship to both attribute value ranges and scatterplot quadrants. It is also possible to invert the visualised relationship and display the distances of non-significant observations to their corresponding critical values (if significant is set to FALSE). +The Moran drop plot is a version of the Moran scatterplot supplemented by visual indications of \emph{p} values. The standard Moran scatterplot provides indirect information about the effect size but not about the level of confidence to determine whether the effects shown might be more than just random outcomes. The Moran drop plot marks significant points in red (positive) and blue (negative spatial autocorrelation) and adds so-called drop lines connecting the significant observations to the scatterplot positions of their associated critical values. The coordinates of the latter are determined either under the assumption of approximate standard normality of the z-scores of the local Moran's \emph{I} values, or based on \emph{p}-values provided through \code{adjusted_p}. The longer the lines, the lower the associated \emph{p} values. The visualisation thus enables visual inspection of statistical significance while maintaining the relationship to both attribute value ranges and scatterplot quadrants. It is also possible to invert the visualised relationship and display the distances of non-significant observations to their corresponding critical values (if significant is set to FALSE). } \value{ When return_df is TRUE, a data frame object with the following members is returned: + \item{region.id}{the region ids} \item{z}{the standardised attribute values} \item{wz}{the standardised spatially lagged attribute values} \item{b}{the y-coordinates of the critical values} @@ -50,9 +49,9 @@ When return_df is TRUE, a data frame object with the following members is return \examples{ # Boston example (CMEDV; owner-occupied housing in USD) data(boston) -boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) +boston.tr <- sf::st_read(system.file("shapes/boston_tracts.gpkg", package="spData")[1]) boston.nb <- poly2nb(boston.tr) boston.listw <- nb2listw(boston.nb) -moran.plot.drop(boston.c$CMEDV, boston.listw, localmoran(boston.c$CMEDV, boston.listw), 2.58, - zero.policy = TRUE, significant = TRUE, plain = FALSE, labels = FALSE, plot = TRUE) +moran.plot.drop(boston.c$CMEDV, boston.listw, localmoran(boston.c$CMEDV, boston.listw), 0.01, + significant = TRUE, labels = NULL) } From 62e661c52073de6fb239078ad312b30309c55974 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sat, 18 Apr 2026 23:46:38 +0200 Subject: [PATCH 18/27] Update moran.plot.drop.R --- R/moran.plot.drop.R | 95 ++++++++++++++++++++++++--------------------- 1 file changed, 50 insertions(+), 45 deletions(-) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index cdd594fc..f99ed5a4 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -39,81 +39,86 @@ moran.plot.drop <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, xlab <- xname if (is.null(ylab)) ylab <- paste("spatially lagged", xname) - Z <- as.vector(scale(x, scale = F)) - WZ <- lag.listw(listw, Z, zero.policy = zero.policy) - if (anyNA(WZ)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") + + WX <- lag.listw(listw, scale(x, scale = F), zero.policy = zero.policy) + if (anyNA(WX)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") if(usePadj) cv <- min(abs(locmoran[(which(adjusted_p <= alpha/2)), 4])) else - cv <- qnorm(1 - alpha, lower.tail = FALSE) - b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / Z - b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / Z - b[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] <- b2[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] + cv <- abs(qnorm(1 - alpha, lower.tail = FALSE)) + + b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(x) / (x - mean(x)) + b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(x) / (x - mean(x)) + b[which((x < mean(x) & WX > mean(WX)) | (x > mean(x) & WX < mean(WX)))] <- b2[which((x < mean(x) & WX > mean(WX)) | (x > mean(x) & WX < mean(WX)))] if(significant) { - x_q1 <- Z[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] >= cv)] - y_q1 <- WZ[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] >= cv)] - b_q1 <- b[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] >= cv)] - x_q2 <- Z[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] <= (-1) * cv)] - y_q2 <- WZ[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] <= (-1) * cv)] - b_q2 <- b[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] <= (-1) * cv)] - x_q3 <- Z[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] >= cv)] - y_q3 <- WZ[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] >= cv)] - b_q3 <- b[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] >= cv)] - x_q4 <- Z[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] <= (-1) * cv)] - y_q4 <- WZ[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] <= (-1) * cv)] - b_q4 <- b[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] <= (-1) * cv)] + x_q1 <- x[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] >= cv)] + y_q1 <- WX[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] >= cv)] + b_q1 <- b[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] >= cv)] + labels_q1 <- labels[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] >= cv)] + x_q2 <- x[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] <= (-1) * cv)] + y_q2 <- WX[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] <= (-1) * cv)] + b_q2 <- b[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] <= (-1) * cv)] + labels_q2 <- labels[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] <= (-1) * cv)] + x_q3 <- x[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] >= cv)] + y_q3 <- WX[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] >= cv)] + b_q3 <- b[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] >= cv)] + labels_q3 <- labels[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] >= cv)] + x_q4 <- x[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] <= (-1) * cv)] + y_q4 <- WX[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] <= (-1) * cv)] + b_q4 <- b[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] <= (-1) * cv)] + labels_q4 <- labels[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] <= (-1) * cv)] } else { - x_q1 <- Z[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] < cv)] - y_q1 <- WZ[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] < cv)] - b_q1 <- b[which(Z > 0 & WZ > mean(WZ) & locmoran[, 4] < cv)] - x_q2 <- Z[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] > (-1) * cv)] - y_q2 <- WZ[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] > (-1) * cv)] - b_q2 <- b[which(Z > 0 & WZ < mean(WZ) & locmoran[, 4] > (-1) * cv)] - x_q3 <- Z[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] < cv)] - y_q3 <- WZ[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] < cv)] - b_q3 <- b[which(Z < 0 & WZ < mean(WZ) & locmoran[, 4] < cv)] - x_q4 <- Z[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] > (-1) * cv)] - y_q4 <- WZ[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] > (-1) * cv)] - b_q4 <- b[which(Z < 0 & WZ > mean(WZ) & locmoran[, 4] > (-1) * cv)] + x_q1 <- x[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] < cv)] + y_q1 <- WX[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] < cv)] + b_q1 <- b[which(x > mean(x) & WX > mean(WX) & locmoran[, 4] < cv)] + x_q2 <- x[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] > (-1) * cv)] + y_q2 <- WX[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] > (-1) * cv)] + b_q2 <- b[which(x > mean(x) & WX < mean(WX) & locmoran[, 4] > (-1) * cv)] + x_q3 <- x[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] < cv)] + y_q3 <- WX[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] < cv)] + b_q3 <- b[which(x < mean(x) & WX < mean(WX) & locmoran[, 4] < cv)] + x_q4 <- x[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] > (-1) * cv)] + y_q4 <- WX[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] > (-1) * cv)] + b_q4 <- b[which(x < mean(x) & WX > mean(WX) & locmoran[, 4] > (-1) * cv)] } - lw.lm <- lm(WZ ~ Z) - plot(Z, WZ, xlab="X_centred", ylab="WX", pch = 20, cex = 0.33, col = "gray70", xlim = c(min(Z),max(Z)), ylim = c(min(WZ, b),max(WZ, b))) - abline(h = mean(WZ), lty = "dashed", col = "grey30") - abline(v = 0, lty = "dashed", col = "grey30") + lw.lm <- lm(WX ~ x) + plot(x, WX, xlab="X", ylab="WX", pch = 20, cex = 0.5, col = "gray70", xlim = c(min(x),max(x)), ylim = c(min(WX, b),max(WX, b))) + abline(h = mean(WX), lty = "dashed", col = "grey30") + abline(v = mean(x), lty = "dashed", col = "grey30") abline(lw.lm, lty = "dotted", col = "grey40") for(i in 1:length(x_q1)) { lines(c(x_q1[i], x_q1[i]), c(y_q1[i], b_q1[i]), type = "l", col = "lightpink", lwd = 1) points(c(x_q1[i]), c(y_q1[i]), pch = 20, cex = 0.5, col = "firebrick") - if (labs && length(x_q1) > 0) - text(x_q1[i], y_q1[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") + if (labs) + text(x_q1[i], y_q1[i], labels = labels_q1[i], pos = 2, cex = 0.5, col = "firebrick") } - + for(i in 1:length(x_q2)) { lines(c(x_q2[i], x_q2[i]), c(y_q2[i], b_q2[i]), type = "l", col = "lightblue1", lwd = 1) points(c(x_q2[i]), c(y_q2[i]), pch = 20, cex = 0.5, col = "royalblue") - if (labs && length(x_q2) > 0) - text(x_q2[i], y_q2[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") + if (labs) + text(x_q2[i], y_q2[i], labels = labels_q2[i], pos = 2, cex = 0.5, col = "royalblue") } for(i in 1:length(x_q3)) { lines(c(x_q3[i], x_q3[i]), c(y_q3[i], b_q3[i]), type = "l", col = "lightpink", lwd = 1) points(c(x_q3[i]), c(y_q3[i]), pch = 20, cex = 0.5, col = "firebrick") - if (labs && length(x_q3) > 0) - text(x_q3[i], y_q3[i], labels = labels[i], pos = 2, cex = 0.5, col = "firebrick") + if (labs) + text(x_q3[i], y_q3[i], labels = labels_q3[i], pos = 2, cex = 0.5, col = "firebrick") } for(i in 1:length(x_q4)) { lines(c(x_q4[i], x_q4[i]), c(y_q4[i], b_q4[i]), type = "l", col = "lightblue1", lwd = 1) points(c(x_q4[i]), c(y_q4[i]), pch = 20, cex = 0.5, col = "royalblue") - if (labs && length(x_q4) > 0) - text(x_q4[i], y_q4[i], labels = labels[i], pos = 2, cex = 0.5, col = "royalblue") + if (labs) + text(x_q4[i], y_q4[i], labels = labels_q4[i], pos = 2, cex = 0.5, col = "royalblue") } if(return_df) { - res <- data.frame(z = Z, wz = WZ, b = b, line_lengths = abs(WZ - b)) #n_sig = length(x_q1) + length(x_q2) + length(x_q3) + length(x_q4) + res <- data.frame(region.id = labels, x = x, WX = WX, b = b, line_lengths = abs(WX - b)) invisible(res) } } From 1a67b69bb77d97782e3f483469a49779cbd9e92c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sat, 18 Apr 2026 23:49:47 +0200 Subject: [PATCH 19/27] Update moran.plot.seismogram.R --- R/moran.plot.seismogram.R | 48 +++++++++++++++++++-------------------- 1 file changed, 24 insertions(+), 24 deletions(-) diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R index d58de201..88c722f1 100644 --- a/R/moran.plot.seismogram.R +++ b/R/moran.plot.seismogram.R @@ -33,34 +33,34 @@ moran.plot.seismogram <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = xlab <- xname if (is.null(ylab)) ylab <- paste("spatially lagged", xname) - Z <- as.vector(scale(x, scale = F)) - WZ <- lag.listw(listw, Z, zero.policy = zero.policy) - if (anyNA(WZ)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") + + WX <- lag.listw(listw, scale(x, scale = F), zero.policy = zero.policy) + if (anyNA(WX)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") if(usePadj) cv <- min(abs(locmoran[(which(adjusted_p <= alpha/2)),4])) else - cv <- qnorm(1 - alpha, lower.tail = FALSE) - b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) - b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(Z) / (Z - mean(Z)) - b[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] <- b2[which((Z < 0 & WZ > mean(WZ)) | (Z > 0 & WZ < mean(WZ)))] + cv <- abs(qnorm(1 - alpha, lower.tail = FALSE)) + b <- ((cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(x) / (x - mean(x)) + b2 <- ((-cv * sqrt(locmoran[, 3])) + locmoran[, 2]) * var(x) / (x - mean(x)) + b[which((x < mean(x) & WX > mean(WX)) | (x > mean(x) & WX < mean(WX)))] <- b2[which((x < mean(x) & WX > mean(WX)) | (x > mean(x) & WX < mean(WX)))] - x_q1 <- Z[which(Z > 0 & WZ > mean(WZ))] - y_q1 <- WZ[which(Z > 0 & WZ > mean(WZ))] - b_q1 <- b[which(Z > 0 & WZ > mean(WZ))] - x_q2 <- Z[which(Z > 0 & WZ < mean(WZ))] - y_q2 <- WZ[which(Z > 0 & WZ < mean(WZ))] - b_q2 <- b[which(Z > 0 & WZ < mean(WZ))] - x_q3 <- Z[which(Z < 0 & WZ < mean(WZ))] - y_q3 <- WZ[which(Z < 0 & WZ < mean(WZ))] - b_q3 <- b[which(Z < 0 & WZ < mean(WZ))] - x_q4 <- Z[which(Z < 0 & WZ > mean(WZ))] - y_q4 <- WZ[which(Z < 0 & WZ > mean(WZ))] - b_q4 <- b[which(Z < 0 & WZ > mean(WZ))] + x_q1 <- x[which(x > mean(x) & WX > mean(WX))] + y_q1 <- WX[which(x > mean(x) & WX > mean(WX))] + b_q1 <- b[which(x > mean(x) & WX > mean(WX))] + x_q2 <- x[which(x > mean(x) & WX < mean(WX))] + y_q2 <- WX[which(x > mean(x) & WX < mean(WX))] + b_q2 <- b[which(x > mean(x) & WX < mean(WX))] + x_q3 <- x[which(x < mean(x) & WX < mean(WX))] + y_q3 <- WX[which(x < mean(x) & WX < mean(WX))] + b_q3 <- b[which(x < mean(x) & WX < mean(WX))] + x_q4 <- x[which(x < mean(x) & WX > mean(WX))] + y_q4 <- WX[which(x < mean(x) & WX > mean(WX))] + b_q4 <- b[which(x < mean(x) & WX > mean(WX))] - lw.lm <- lm(WZ ~ Z) - plot(Z, WZ, xlab="X_centred", ylab="WX", pch = 20, cex = 0.33, col = "gray70", xlim = c(min(Z),max(Z)), ylim = c(min(WZ, b),max(WZ, b))) - abline(h = mean(WZ), lty = "dashed", col = "grey30") - abline(v = 0, lty = "dashed", col = "grey30") + lw.lm <- lm(WX ~ x) + plot(x, WX, xlab="X", ylab="WX", pch = 20, cex = 0.5, col = "gray70", xlim = c(min(x),max(x)), ylim = c(min(WX, b),max(WX, b))) + abline(h = mean(WX), lty = "dashed", col = "grey30") + abline(v = mean(x), lty = "dashed", col = "grey30") abline(lw.lm, lty = "dotted", col = "grey40") df_q1 <- data.frame(x_q1, b_q1) @@ -85,7 +85,7 @@ moran.plot.seismogram <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = lines(c(df_q4$x_q4[i], df_q4$x_q4[i+1]), c(df_q4$b_q4[i], df_q4$b_q4[i+1]), type = "l", col = "royalblue") if(return_df) { - res <- data.frame(z = Z, wz = WZ, b = b) + res <- data.frame(x = x, WX = WX, b = b) invisible(res) } } From beba07fa2258eea0077cfb05aee180724b3c603f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 12:41:17 +0200 Subject: [PATCH 20/27] Update moran.plot.drop.R Improved label behaviour. --- R/moran.plot.drop.R | 19 ++++++++++++++----- 1 file changed, 14 insertions(+), 5 deletions(-) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index f99ed5a4..f370e2f3 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -31,14 +31,23 @@ moran.plot.drop <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, if (spChk && !chkIDs(x, listw)) stop("Check of data and weights ID integrity failed") labs <- TRUE - if (is.logical(labels) && !labels) + if (is.logical(labels)) { + if(!labels) + labs <- FALSE + labels <- as.character(attr(listw, "region.id")) + } else if (!is.logical(labels) && !is.null(labels)) { + if(length(labels) != n) { + warning("Length of the labels vector does not match number of regions. region.id from the listw object is used instead.") + labels <- as.character(attr(listw, "region.id")) + } + } else if (is.null(labels)) { labs <- FALSE - if (is.null(labels) || length(labels) != n) labels <- as.character(attr(listw, "region.id")) + } if (is.null(xlab)) xlab <- xname if (is.null(ylab)) - ylab <- paste("spatially lagged", xname) + ylab <- paste("spatially lagged centred", xname) WX <- lag.listw(listw, scale(x, scale = F), zero.policy = zero.policy) if (anyNA(WX)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") @@ -84,7 +93,7 @@ moran.plot.drop <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, } lw.lm <- lm(WX ~ x) - plot(x, WX, xlab="X", ylab="WX", pch = 20, cex = 0.5, col = "gray70", xlim = c(min(x),max(x)), ylim = c(min(WX, b),max(WX, b))) + plot(x, WX, xlab=xlab, ylab=ylab, pch = 20, cex = 0.5, col = "gray70", xlim = c(min(x),max(x)), ylim = c(min(WX, b),max(WX, b))) abline(h = mean(WX), lty = "dashed", col = "grey30") abline(v = mean(x), lty = "dashed", col = "grey30") abline(lw.lm, lty = "dotted", col = "grey40") @@ -118,7 +127,7 @@ moran.plot.drop <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, } if(return_df) { - res <- data.frame(region.id = labels, x = x, WX = WX, b = b, line_lengths = abs(WX - b)) + res <- data.frame(labels = labels, x = x, WX = WX, b = b, line_lengths = abs(WX - b)) invisible(res) } } From f05e022471a9ff9cf05c955978859eec18ba2ce2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 12:42:15 +0200 Subject: [PATCH 21/27] Update moran.plot.drop.Rd --- man/moran.plot.drop.Rd | 28 ++++++++++++++-------------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd index c1c16e03..b8cdbedc 100644 --- a/man/moran.plot.drop.Rd +++ b/man/moran.plot.drop.Rd @@ -15,29 +15,29 @@ moran.plot.drop(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, \arguments{ \item{x}{a numerical vector holding the attribute of interest} \item{listw}{a \code{listw} spatial weights object} - \item{locmoran}{a fitted object of type localmoran} - \item{alpha}{default 0.05; the desired significance level regarding local Moran's \emph{I} values} - \item{adjusted_p}{default NULL; a vector of \emph{p} values adjusted to account for multiple testing as is returned by \code{\link{p.adjustSP}}; standard normal distribution is used to determine critical values, if NULL} - \item{significant}{default TRUE; a parameter indicating whether to display plot critical value distances of significant (default) or non-significant observations} - \item{xlab}{label for x axis} - \item{ylab}{label for y axis} + \item{locmoran}{a fitted object of type \code{localmoran}} + \item{alpha}{default 0.05; the desired significance level regarding local Moran's \emph{I}} + \item{adjusted_p}{default NULL; an optional vector of \emph{p} values adjusted to account for multiple testing as is returned by \code{\link{p.adjustSP}}; if NULL, standard normal approximation is used to determine critical values} + \item{significant}{default TRUE; a parameter indicating whether to display critical value distances of significant (default) or non-significant observations} + \item{xlab}{default NULL; an optional label for the x-axis} + \item{ylab}{default NULL; an optional label for the y-axis} \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} - \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} - \item{labels}{default NULL; no labels are plotted by default; region IDs are used as labels of significant observations if set to TRUE} + \item{spChk}{default NULL to use \code{get.spChkOption()}; should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE} + \item{labels}{default NULL; no labels are plotted by default; region IDs from the \code{listw} object are used as labels of significant observations if set to TRUE; custom labels are used if a character vector is provided} \item{zero.policy}{default option stored in the listw object; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} } \details{ -The Moran drop plot is a version of the Moran scatterplot supplemented by visual indications of \emph{p} values. The standard Moran scatterplot provides indirect information about the effect size but not about the level of confidence to determine whether the effects shown might be more than just random outcomes. The Moran drop plot marks significant points in red (positive) and blue (negative spatial autocorrelation) and adds so-called drop lines connecting the significant observations to the scatterplot positions of their associated critical values. The coordinates of the latter are determined either under the assumption of approximate standard normality of the z-scores of the local Moran's \emph{I} values, or based on \emph{p}-values provided through \code{adjusted_p}. The longer the lines, the lower the associated \emph{p} values. The visualisation thus enables visual inspection of statistical significance while maintaining the relationship to both attribute value ranges and scatterplot quadrants. It is also possible to invert the visualised relationship and display the distances of non-significant observations to their corresponding critical values (if significant is set to FALSE). +The Moran drop plot is a version of the Moran scatterplot supplemented by visual indications of \emph{p} values. The standard Moran scatterplot provides information about the effect size but not about the level of confidence to determine whether the effects shown might be more than just random outcomes. The Moran drop plot marks significant points in red (positive) and blue (negative spatial autocorrelation) and adds so-called drop lines connecting the significant observations to the scatterplot positions of their associated critical values. The coordinates of the latter are determined either under the assumption of approximate standard normality of the z-scores of the local Moran's \emph{I} values, or based on \emph{p} values provided through \code{adjusted_p}. In the latter case, the critical value is approximated based on the observation with the highest adjusted \emph{p} value that still satisfies the selected significance level. The longer the lines, the lower the associated \emph{p} values. The visualisation thus enables visual inspection of statistical significance while maintaining the relationship to both attribute value ranges and scatterplot quadrants. It is also possible to invert the visualised relationship and display the distances of non-significant observations to their corresponding critical values (if significant is set to FALSE). } \value{ When return_df is TRUE, a data frame object with the following members is returned: - \item{region.id}{the region ids} - \item{z}{the standardised attribute values} - \item{wz}{the standardised spatially lagged attribute values} - \item{b}{the y-coordinates of the critical values} - \item{line_lengths}{the absolute distances between b and z} + \item{labels}{either the labels provided or the region ids, if not specified} + \item{x}{the attribute values} + \item{wx}{the spatial lags of the centred attribute values} + \item{b}{the y-coordinates (i.e. hypothetical spatial lags) of the critical values given x} + \item{line_lengths}{the absolute distances between b and x} } \references{Westerholt, R. (2024): Extending the Moran scatterplot by indications of critical values and \emph{p}-values: introducing the Moran seismogram and the drop plot. Proceedings of the 32nd Annual GIS Research UK Conference (GISRUK), Leeds, UK. \url{https://doi.org/10.5281/zenodo.10897792}} From cb7ffd7d7e1f7b7ebe904b4b771c7ed9c0e588ac Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 12:49:24 +0200 Subject: [PATCH 22/27] Update moran.plot.drop.Rd --- man/moran.plot.drop.Rd | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd index b8cdbedc..f28847ea 100644 --- a/man/moran.plot.drop.Rd +++ b/man/moran.plot.drop.Rd @@ -22,7 +22,7 @@ moran.plot.drop(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, \item{xlab}{default NULL; an optional label for the x-axis} \item{ylab}{default NULL; an optional label for the y-axis} \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} - \item{spChk}{default NULL to use \code{get.spChkOption()}; should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE} + \item{spChk}{default NULL to use \code{get.spChkOption()}; should the locmoran names be checked against the listw spatial objects for identity integrity, TRUE, or FALSE} \item{labels}{default NULL; no labels are plotted by default; region IDs from the \code{listw} object are used as labels of significant observations if set to TRUE; custom labels are used if a character vector is provided} \item{zero.policy}{default option stored in the listw object; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} } From 044797f9bf1f755fc3352fd47b4a1e524ef4f519 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 12:50:50 +0200 Subject: [PATCH 23/27] Update moran.plot.drop.R --- R/moran.plot.drop.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/moran.plot.drop.R b/R/moran.plot.drop.R index f370e2f3..57eac58c 100644 --- a/R/moran.plot.drop.R +++ b/R/moran.plot.drop.R @@ -28,7 +28,7 @@ moran.plot.drop <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, } if (is.null(spChk)) spChk <- get.spChkOption() - if (spChk && !chkIDs(x, listw)) + if (spChk && !chkIDs(locmoran, listw)) stop("Check of data and weights ID integrity failed") labs <- TRUE if (is.logical(labels)) { From f792540c91538bcd60db486836b839b4a6c3c0bf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 13:09:55 +0200 Subject: [PATCH 24/27] Update moran.plot.seismogram.R --- R/moran.plot.seismogram.R | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/R/moran.plot.seismogram.R b/R/moran.plot.seismogram.R index 88c722f1..da83256a 100644 --- a/R/moran.plot.seismogram.R +++ b/R/moran.plot.seismogram.R @@ -27,12 +27,12 @@ moran.plot.seismogram <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = } if (is.null(spChk)) spChk <- get.spChkOption() - if (spChk && !chkIDs(x, listw)) - stop("Check of data and weights ID integrity failed") + if (spChk && !chkIDs(locmoran, listw)) + stop("Check of locmoran and weights ID integrity failed") if (is.null(xlab)) xlab <- xname if (is.null(ylab)) - ylab <- paste("spatially lagged", xname) + ylab <- paste("spatially lagged centred", xname) WX <- lag.listw(listw, scale(x, scale = F), zero.policy = zero.policy) if (anyNA(WX)) warning("no-neighbour observation(s) found - use zero.policy=TRUE") @@ -58,7 +58,7 @@ moran.plot.seismogram <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = b_q4 <- b[which(x < mean(x) & WX > mean(WX))] lw.lm <- lm(WX ~ x) - plot(x, WX, xlab="X", ylab="WX", pch = 20, cex = 0.5, col = "gray70", xlim = c(min(x),max(x)), ylim = c(min(WX, b),max(WX, b))) + plot(x, WX, xlab=xlab, ylab=ylab, pch = 20, cex = 0.5, col = "gray70", xlim = c(min(x),max(x)), ylim = c(min(WX, b),max(WX, b))) abline(h = mean(WX), lty = "dashed", col = "grey30") abline(v = mean(x), lty = "dashed", col = "grey30") abline(lw.lm, lty = "dotted", col = "grey40") @@ -85,7 +85,7 @@ moran.plot.seismogram <- function(x, listw, locmoran, alpha = 0.05, adjusted_p = lines(c(df_q4$x_q4[i], df_q4$x_q4[i+1]), c(df_q4$b_q4[i], df_q4$b_q4[i+1]), type = "l", col = "royalblue") if(return_df) { - res <- data.frame(x = x, WX = WX, b = b) + res <- data.frame(labels=as.character(attr(listw, "region.id")), x=x, wx=WX, b=b) invisible(res) } } From 208a99431ae1b0cda90265f84a524773ee11535c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 13:11:38 +0200 Subject: [PATCH 25/27] Update moran.plot.drop.Rd --- man/moran.plot.drop.Rd | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/man/moran.plot.drop.Rd b/man/moran.plot.drop.Rd index f28847ea..527100d8 100644 --- a/man/moran.plot.drop.Rd +++ b/man/moran.plot.drop.Rd @@ -9,7 +9,7 @@ A version of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented b \usage{ moran.plot.drop(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, significant = TRUE, xlab = NULL, ylab = NULL, return_df = TRUE, - spChk = NULL, labels = NULL), zero.policy=attr(listw, "zero.policy") + spChk = NULL, labels = NULL, zero.policy=attr(listw, "zero.policy")) } \arguments{ From 0329de3e3b4cc612f873e2ddf810625f1b82642f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 13:12:46 +0200 Subject: [PATCH 26/27] Update moran.plot.seismogram.Rd --- man/moran.plot.seismogram.Rd | 34 ++++++++++++++++------------------ 1 file changed, 16 insertions(+), 18 deletions(-) diff --git a/man/moran.plot.seismogram.Rd b/man/moran.plot.seismogram.Rd index fb40853b..aa0bda9c 100644 --- a/man/moran.plot.seismogram.Rd +++ b/man/moran.plot.seismogram.Rd @@ -4,37 +4,36 @@ \title{Moran seismogram} \description{ -A variant of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines connecting location-wise critical values. The plot allows for visual inspection of potential spatial weights misspecifiation. +A variant of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented by lines connecting location-wise critical value configurations of attribute values and spatial lags. The plot allows for visual inspection of potential spatial weights misspecifiation. } \usage{ -moran.plot.seismogram(x, listw, locmoran, cv = 2.58, - plain = FALSE, zero.policy = FALSE, xlab = NULL, - ylab = NULL, plot = TRUE, return_df = TRUE, spChk = NULL) +moran.plot.seismogram(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, + xlab = NULL, ylab = NULL, return_df = TRUE, spChk = NULL, zero.policy = FALSE) } \arguments{ \item{x}{a numerical vector holding the attribute of interest} \item{listw}{a \code{listw} spatial weights object} \item{locmoran}{a fitted object of type localmoran} - \item{cv}{default 2.58; the desired critical value assuming approximate normality of standardised local Moran's \emph{I} values} - \item{plain}{default FALSE; a plain Moran scatterplot is displayed if TRUE} - \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} - \item{xlab}{label for x axis} - \item{ylab}{label for y axis} - \item{plot}{default TRUE; if FALSE, plotting is suppressed} + \item{alpha}{default 0.05; the desired significance level regarding local Moran's \emph{I}} + \item{adjusted_p}{default NULL; an optional vector of \emph{p} values adjusted to account for multiple testing as is returned by \code{\link{p.adjustSP}}; if NULL, standard normal approximation is used to determine critical values} + \item{xlab}{default NULL; an optional label for the x-axis} + \item{ylab}{default NULL; an optional label for the y-axis} \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} - \item{spChk}{should the data vector names be checked against the spatial objects for identity integrity, TRUE, or FALSE, default NULL to use \code{get.spChkOption()}} + \item{spChk}{default NULL to use \code{get.spChkOption()}; should the locmoran names be checked against the listw spatial objects for identity integrity, TRUE, or FALSE} + \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} } \details{ -The Moran seismogram is a version of the Moran scatterplot that is complemented by lines connecting the location-specific critical values of local Moran's \emph{I}. Lines in quadrants with positive spatial autocorrelation are shown in red and lines in quadrants with negative spatial autocorrelation are shown in blue. The representation is similar to a seismogram for detecting earthquakes and thus reveals potentially suspicious local spatial weights configurations by visualising spikes. The latter are displayed in an integrated manner with their positions in the attribute value range and in connection with the types of the associated spatial patterns (by the quadrants of the scatterplot). +The Moran seismogram is a version of the Moran scatterplot that is complemented by lines connecting the location-specific critical values of local Moran's \emph{I} in the plot. The y-coordinates associated with the critical value are determined either under the assumption of approximate standard normality of the z-scores of the local Moran's \emph{I} values, or based on \emph{p} values provided through \code{adjusted_p}. In the latter case, the critical value is approximated based on the observation with the highest adjusted \emph{p} value that still satisfies the selected significance level. Lines in quadrants with positive spatial autocorrelation are shown in red and lines in quadrants with negative spatial autocorrelation are shown in blue. The representation is similar to a seismogram for detecting earthquakes and thus reveals potentially suspicious local spatial weights configurations by visualising spikes. The latter are displayed in an integrated manner with their positions in the attribute value range and in connection with the types of the associated spatial patterns (by the quadrants of the scatterplot). } \value{ When return_df is TRUE, a data frame object with the following members is returned: - \item{z}{the standardised attribute values} - \item{wz}{the standardised spatially lagged attribute values} - \item{b}{the y-coordinates of the critical values} + \item{labels}{either the labels provided or the region ids, if not specified} + \item{x}{the attribute values} + \item{wx}{the spatial lags of the centred attribute values} + \item{b}{the y-coordinates (i.e. hypothetical spatial lags) of the critical values given x} } \references{Westerholt, R. (2024): Extending the Moran scatterplot by indications of critical values and \emph{p}-values: introducing the Moran seismogram and the drop plot. Proceedings of the 32nd Annual GIS Research UK Conference (GISRUK), Leeds, UK. \url{https://doi.org/10.5281/zenodo.10897792}} @@ -46,9 +45,8 @@ When return_df is TRUE, a data frame object with the following members is return \examples{ # Boston example (CMEDV; owner-occupied housing in USD) data(boston) -boston.tr <- sf::st_read(system.file("shapes/boston_tracts.shp", package="spData")[1]) +boston.tr <- sf::st_read(system.file("shapes/boston_tracts.gpkg", package="spData")[1]) boston.nb <- poly2nb(boston.tr) boston.listw <- nb2listw(boston.nb) -moran.plot.seismogram(boston.c$CMEDV, boston.listw, localmoran(boston.c$CMEDV, boston.listw), 2.58, -zero.policy = TRUE, plain = FALSE, plot = TRUE) +moran.plot.seismogram(boston.c$CMEDV, boston.listw, localmoran(boston.c$CMEDV, boston.listw), 0.01, zero.policy = TRUE) } From ed1ca9a5532566751619f9da5b365bb54ec727a6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ren=C3=A9=20Westerholt?= Date: Sun, 19 Apr 2026 13:42:26 +0200 Subject: [PATCH 27/27] Update moran.plot.seismogram.Rd --- man/moran.plot.seismogram.Rd | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/man/moran.plot.seismogram.Rd b/man/moran.plot.seismogram.Rd index aa0bda9c..605837e3 100644 --- a/man/moran.plot.seismogram.Rd +++ b/man/moran.plot.seismogram.Rd @@ -8,7 +8,7 @@ A variant of the Moran scatterplot (see \code{\link{moran.plot}}) supplemented b } \usage{ moran.plot.seismogram(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, - xlab = NULL, ylab = NULL, return_df = TRUE, spChk = NULL, zero.policy = FALSE) + xlab = NULL, ylab = NULL, return_df = TRUE, spChk = NULL, zero.policy = attr(listw, "zero.policy")) } \arguments{ @@ -21,7 +21,7 @@ moran.plot.seismogram(x, listw, locmoran, alpha = 0.05, adjusted_p = NULL, \item{ylab}{default NULL; an optional label for the y-axis} \item{return_df}{default TRUE; invisibly return a data.frame object, if FALSE do not return anything} \item{spChk}{default NULL to use \code{get.spChkOption()}; should the locmoran names be checked against the listw spatial objects for identity integrity, TRUE, or FALSE} - \item{zero.policy}{default NULL; use global option value; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} + \item{zero.policy}{default option stored in the listw object; if FALSE stop with error for any empty neighbour sets, if TRUE permit the weights list to be formed with zero-length weights vectors} } \details{