## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

## ----setup--------------------------------------------------------------------
library(covercorr)

## -----------------------------------------------------------------------------
set.seed(1)
n <- 100
X <- rnorm(n)
Y <- rnorm(n)
result <- coverage_correlation(Y, X, visualise = TRUE)
result

## -----------------------------------------------------------------------------
set.seed(2)
n <- 100
X <- rnorm(n)
Z <- rnorm(n)
rho <- 0.9
Y <- rho * X + sqrt(1 - rho^2) * Z
result <- coverage_correlation(Y, X, visualise = TRUE)
result

## -----------------------------------------------------------------------------
set.seed(3)
n <- 100
p <- 2
X <- matrix(rnorm(p * n), ncol = p)
Y <- matrix(0, nrow = n, ncol = p)
Y[, 1] <- X[, 1]^2
Y[, 2] <- X[, 1] * X[, 2]
result <- coverage_correlation(Y, X)
result

## -----------------------------------------------------------------------------
set.seed(4)
n <- 50
p <- 2
X <- matrix(rnorm(p * n), ncol = p)
Y <- matrix(rnorm(p * n), ncol = p)
result <- coverage_correlation(Y, X, method = "approx")
result

## -----------------------------------------------------------------------------
set.seed(5)
n <- 100
X <- rnorm(n)
Y <- sin(3 * X) + rnorm(n) * 0.05

# Deterministic: repeated calls give identical results
r1 <- coverage_correlation_grid(Y, X)
r2 <- coverage_correlation_grid(Y, X)
c(r1$stat, r2$stat)

## -----------------------------------------------------------------------------
set.seed(6)
n <- 100
X <- rnorm(n)
Y <- sin(3 * X) + rnorm(n) * 0.05

# Pairwise, random reference points
covercorr(X, Y)

# Pairwise, deterministic grid
covercorr(X, Y, reference = "deterministic")

## -----------------------------------------------------------------------------
set.seed(7)
n <- 100
X <- rnorm(n)
Y <- sin(3 * X) + rnorm(n) * 0.05
covercorr(cbind(X, Y))

## -----------------------------------------------------------------------------
set.seed(8)
n <- 100
X1 <- rnorm(n)
X2 <- sin(3 * X1) + rnorm(n) * 0.05
X3 <- X1 * X2 + rnorm(n) * 0.1
coverage_correlation_K(list(X1, X2, X3))

## -----------------------------------------------------------------------------
set.seed(9)
n <- 100
X1 <- rnorm(n)
X2 <- sin(3 * X1) + rnorm(n) * 0.05
X3 <- X1 * X2 + rnorm(n) * 0.1
covercorr(cbind(X1, X2, X3))

## -----------------------------------------------------------------------------
set.seed(10)
n <- 100
X1 <- rnorm(n)
X2 <- sin(3 * X1) + rnorm(n) * 0.05
X3 <- X1 * X2 + rnorm(n) * 0.1
res <- coverage_correlation_K_grid(list(X1, X2, X3))
res$stat
res$pval
res$pval_available

## ----fig.width = 7.2, fig.height = 3.6, out.width = "100%"--------------------
data(CD8T)
n <- nrow(CD8T)

i <- 1
j <- 2
x <- CD8T[, i]
y <- CD8T[, j]

# The coefficient and its p-value for this pair
covercorr(x, y)

# Monge-Kantorovich ranks, using random reference points
set.seed(1)
x_rank <- MK_rank(as.matrix(x), runif(n))
y_rank <- MK_rank(as.matrix(y), runif(n))

op <- par(mfrow = c(1, 2))
plot(x, y, pch = ".", xlab = colnames(CD8T)[i], ylab = colnames(CD8T)[j])
visualise_density(x_rank, y_rank)
par(op)

