Hi all, Adapting some code from the MEE-AOA repo, I believe I can ca

`aoa()` appears to return incorrect thresholds (different from Meyer & Pebesma 2021) about cast HOT 2 CLOSED

mikemahoney218 commented on June 4, 2024

`aoa()` appears to return incorrect thresholds (different from Meyer & Pebesma 2021)

from cast.

Comments (2)

HannaMeyer commented on June 4, 2024 1

Thanks for finding that! I fixed it according to your suggestion.

from cast.

mikemahoney218 commented on June 4, 2024

Looking at the source of grDevices::boxplot.stats, I see:

> grDevices::boxplot.stats
function (x, coef = 1.5, do.conf = TRUE, do.out = TRUE) 
{
    if (coef < 0) 
        stop("'coef' must not be negative")
    nna <- !is.na(x)
    n <- sum(nna)
    stats <- stats::fivenum(x, na.rm = TRUE)
    iqr <- diff(stats[c(2, 4)])
    if (coef == 0) 
        do.out <- FALSE
    else {
        out <- if (!is.na(iqr)) {
            x < (stats[2L] - coef * iqr) | x > (stats[4L] + coef * 
                iqr)
        }
        else !is.finite(x)
        if (any(out[nna], na.rm = TRUE)) 
            stats[c(1, 5)] <- range(x[!out], na.rm = TRUE)
    }
    conf <- if (do.conf) 
        stats[3L] + c(-1.58, 1.58) * iqr/sqrt(n)
    list(stats = stats, n = n, conf = conf, out = if (do.out) x[out & 
        nna] else numeric())
}

The relevant lines here are:

out <- if (!is.na(iqr)) {
  x < (stats[2L] - coef * iqr) | x > (stats[4L] + coef * iqr)
}
#> [...]
stats[c(1, 5)] <- range(x[!out], na.rm = TRUE)

In context, that means that the AOA threshold winds up equaling:

max(di[!(di > (quantile(di, 0.75) + 1.5 * IQR(di)))])

Which means that the threshold is going to be the value in di closest to, but not more than quantile(di, 0.75) + (1.5 * IQR(di)), which, particularly for smaller data, may be a significantly different value. The returned value will always be lower than (or the same as) the 75th percentile plus 1.5 times the IQR. Is this expected behavior?

from cast.

Recommend Projects