Top stuff
library(loclm)
show(packageVersion("loclm"))
#> [1] '1.0.0'
show(date())
#> [1] "Thu Jul 30 07:00:27 2026"
Here I’m just setting up a nicer plotting frame.
oldpar <- par(mgp = c(2.5, 0.7, 0), mar = c(4, 4, 2, 0.5), family = "Hershey",
cex.main = 0.8, cex.axis = 0.8, cex.lab = 0.8, las = 1)
op <- par(no.readonly = TRUE)
Create a dataset
Let’s use the multivariate Banana function variant B from
S. Kok and C. Sandrock, 2009, Locating and Characterizing the
Stationary Points of the Extended Rosenbrock Function, Evolutionary
Computing, 17(3), pages 437-453.
As the paper explains, variant B is pathological, especially in
higher dimensions.
#### d-dimensional Banana function
## (d >= 2 is implied by x)
banana <- function(x) {
stopifnot(is.vector(x, mode = "numeric"))
d <- length(x)
if (d == 1) {
x <- c(x, 1)
d <- 2L
}
robj <- 0
for (j in seq_len(d - 1)) {
robj <- robj + 100 * (x[j]^2 - x[j+1])^2 + (x[j] - 1)^2
}
robj
}
Here is a dataset of 101 runs from a 5D banana function, plus a
little noise. Feel free to play with d and
sigma.
## make a dataset
set.seed(1234) # for reproducibility, feel free to change
d <- 3
sigma <- 1
X <- matrix(runif(d * 101), ncol = d)
y <- apply(X, 1, banana) + rnorm(nrow(X), sd = sigma)
## here's where we'll predict
newX <- matrix(runif(d * 23), ncol = d)
newy <- apply(newX, 1, banana) # leave off the noise
LOO cross-validation to set the span
Now use Leave-One-Out (LOO) cross-validation to set the span. We need
a general-purpose LOO function.
## LOO cross-validation
LOO <- function(X, y, reg, ...) {
n <- length(y)
yhat <- vapply(seq.int(n), \(i) {
fit <- reg(X[-i, , drop=FALSE], y[-i], ...)
predict(fit, X[i, , drop=FALSE])
}, FUN.VALUE = 0)
list(yhat = yhat, res = y - yhat)
}
Now call this with loclm as the regression method, and
span as the parameter.
## try different spans
span_vals <- c(seq(0.1, 0.2, 0.01), 0.75)
rmse <- vapply(span_vals, \(span) {
loo <- LOO(X, y, reg = loclm, span = span)
sqrt(mean(loo$res^2))
}, FUN.VALUE = 0)
show(cbind(span = span_vals, rmse = rmse))
#> span rmse
#> [1,] 0.10 8.302191
#> [2,] 0.11 8.136132
#> [3,] 0.12 8.008744
#> [4,] 0.13 7.910766
#> [5,] 0.14 7.944608
#> [6,] 0.15 7.262508
#> [7,] 0.16 7.284905
#> [8,] 0.17 7.398039
#> [9,] 0.18 7.391820
#> [10,] 0.19 7.454959
#> [11,] 0.20 7.617551
#> [12,] 0.75 14.358125
spanhat <- span_vals[which.min(rmse)]
show(c(span = spanhat))
#> span
#> 0.15
Try with the optimal span.
fit <- loclm(X, y, span = spanhat)
pp <- predict(fit, newX)
par(op)
matplot(newX[, 2], cbind(newy, pp), pch = c(1, 2), cex = cex,
col = col,
main = sprintf("span = %s (optimal)", format(spanhat)))
That’s a big improvement.
