EDM simplex calculations test (E = 1)

main.R

library(rEDM)

# generate time series
set.seed(42)
ts_length <- 50
x <- rnorm(ts_length)

# pre-compute distance matrix
x_dists <- as.matrix(dist(x))
x_dists <- x_dists[, -ts_length] # last point can't be a nearest neighbor

# compute simplex predictions "manually"
pred <- rep.int(NaN, ts_length)
for(t in seq(ts_length-1))
{
  # take 2nd and 3rd closest for 2 nearest neighbors
  # (dropping closest because the distance matrix includes itself)
  nn_idx <- order(x_dists[t,])[2:3]
  
  # should be euclidean dist, but abs is fine for E = 1
  nn_wgt <- exp(-abs(x[nn_idx] - x[t])/abs(x[nn_idx[1]]-x[t]))
  nn_wgt <- pmax(nn_wgt, 1e-6) # minimum weighting is 1e-6
  
  pred[t] <- sum(x[nn_idx+1] * nn_wgt) / sum(nn_wgt)
}

# compute simplex predictions using rEDM::simplex()
out_simplex <- simplex(x, E = 1, stats_only = FALSE)
pred_simplex <- out_simplex$model_output[[1]]$pred

all.equal(pred, pred_simplex)

# compute simplex predictions using rEDM::block_lnlp()
out_bl <- block_lnlp(x, columns = 1, method = "simplex", stats_only = FALSE)
pred_bl <- out_bl$model_output[[1]]$pred

all.equal(pred, pred_bl)