Tidy and geospatial kernel density based local two-sample comparison tests

Tidy and geospatial versions of kernel density based local two-sample comparison tests for 1- and 2-dimensional data.

Usage

tidy_kde_local_test(data1, data2, labels, ...)
st_kde_local_test(x1, x2, labels, ...)

Arguments

data1, data2: data frames/tibbles of data values
x1,x2: sf objects with point geometry
labels: flag or vector of strings for legend labels
...: other parameters in ks::kde.local.test function

Value

The output has the same structure as the kernel density estimate from *_kde, except that estimate is the difference between the density values data1-data2 rather than the density values, and label becomes an indicator factor of the local comparison test result: "f1<f2" = data1 < data2, 0 = data1 = data2, "f2>f1" = data1 > data2.

The output from st_kde_local_test has two contours, with contlabel=-50 (for f1<f2) and contlabel=50 (for f1>f2), as multipolygons which delimit the significant difference regions.

Details

A kernel local density based two-sample comparison is a modification of the standard kernel density estimate where the two data samples are compared. A Hochberg procedure is employed to control the significance level for multiple comparison tests.

For details of the computation of the kernel local density based two-sample comparison test and the bandwidth selector procedure, see ?ks::kde.local.test. The bandwidth matrix of smoothing parameters is computed as in ks::kde per data sample.

If labels is missing, then the first sample label is taken from x1, and the second sample label from x2. If labels="default" then these are "f1" and "f2". Otherwise, they are assigned to the values of the input vector of strings.

Examples

## tidy local test between unsuccessful and successful grafts
library(ggplot2)
data(hsct, package="ks")
hsct <- dplyr::as_tibble(hsct)
hsct <- dplyr::filter(hsct, PE.Ly65Mac1 >0 & APC.CD45.2>0)
hsct6 <- dplyr::filter(hsct, subject==6)   ## unsuccessful graft 
hsct6 <- dplyr::select(hsct6, PE.Ly65Mac1, APC.CD45.2)
hsct12 <- dplyr::filter(hsct, subject==12) ## successful graft 
hsct12 <- dplyr::select(hsct12, PE.Ly65Mac1, APC.CD45.2)
t1 <- tidy_kde_local_test(data1=hsct6, data2=hsct12)
gt <- ggplot(t1, aes(x=PE.Ly65Mac1, y=APC.CD45.2)) + coord_fixed() 
gt + geom_contour_filled_ks(aes(fill=after_stat(contregion)))


t2 <- tidy_kfs(hsct12)
gt + geom_contour_ks(data=t1, aes(colour=after_stat(contregion))) +
    geom_contour_filled_ks(data=t2, aes(fill=after_stat(contregion)), colour=1) + 
    scale_fill_manual(values=7)
#> Scale for fill is already present.
#> Adding another scale for fill, which will replace the existing scale.


## geospatial local test between Grevillea species
data(wa)
data(grevilleasf)
hakeoides <- dplyr::filter(grevilleasf, species=="hakeoides")
paradoxa <- dplyr::filter(grevilleasf, species=="paradoxa")
s1 <- st_kde_local_test(x1=hakeoides, x2=paradoxa)
s2 <- st_kfs(hakeoides)

## base R plot
xlim <- c(1.2e5, 1.1e6); ylim <- c(6.1e6, 7.2e6)
plot(wa, xlim=xlim, ylim=ylim)
plot(s1, add=TRUE)


plot(wa, xlim=xlim, ylim=ylim)
plot(s1, col=NA, add=TRUE,
     border=colorspace::qualitative_hcl(n=2, palette="Dark2", rev=TRUE))
plot(s2, col=7, add=TRUE, pos="bottom")


## geom_sf plot
gs <- ggplot(s1) + geom_sf(data=wa, fill=NA) + theme_sf() 
gs + geom_sf(data=st_get_contour(s1), aes(fill=contregion)) +
    coord_sf(xlim=xlim, ylim=ylim)

gs + geom_sf(data=st_get_contour(s1), aes(colour=contregion), fill=NA) +
    geom_sf(data=st_get_contour(s2), aes(fill=contregion), colour=1) + 
    scale_fill_manual(values=7) + coord_sf(xlim=xlim, ylim=ylim)
#> Scale for fill is already present.
#> Adding another scale for fill, which will replace the existing scale.