Package {VNDesign}

Build a covariance matrix \Sigma

Description

Generates a covariance matrix \Sigma for a set of locations coords, using different covariance kernel families (exponential, Gaussian, Matérn, spherical, triangular, Brownian, AR(1), or a user-defined covariance).

Usage

build_Sigma(
  coords,
  type = c("exp", "gauss", "matern32", "matern52", "spherical", "triangular",
    "triangular_trunc", "brownian", "ar1", "custom"),
  sigma2 = 1,
  range = 1,
  nugget = 0,
  phi = NULL,
  cov_fun = NULL,
  return_dist = FALSE
)

Arguments

Details

The implemented covariance structures include exponential, Gaussian, Matérn, spherical, triangular, Brownian and AR(1) covariance models.

For distance-based covariance structures, the pairwise distance matrix D is computed as:

• abs(outer(coords, coords, "-")) when coords is one-dimensional,

• as.matrix(dist(coords, method = "euclidean")) when coords is multidimensional.

For type = "brownian", the Brownian covariance kernel K(u,v)=\min(u,v)/range is used, and therefore coords must be one-dimensional.

For type = "custom", the covariance matrix is constructed by evaluating cov_fun(u,v) for each pair of candidate points.

Value

If return_dist = FALSE (default), returns the covariance matrix Sigma of dimension N \times N.

If return_dist = TRUE, returns a list containing:

• Sigma: covariance matrix.

• D: pairwise distance matrix.

Build \chi using an analytical gradient

Description

Constructs the candidate regressor matrix Xchi by evaluating the analytical gradient of the predictor \eta(x,\theta) with respect to the parameter vector theta at each candidate point in chi.

Usage

build_Xchi_grad(chi, theta, grad_fun)

Arguments

Value

Numeric matrix Xchi of dimension N \times p, where row i contains the gradient evaluated at the candidate point x_i.

Examples

## ---------------------------------------------
## Example 1D 
## ---------------------------------------------

chi <- seq(-1, 1, length.out = 5)
theta <- c(0, 1)

# eta(x) = theta1 + theta2 * x
grad_fun <- function(x, theta) {
  c(1, x)
}

build_Xchi_grad(chi, theta, grad_fun)

## ---------------------------------------------
## Example 2D 
## ---------------------------------------------

chi2 <- cbind(x1 = c(0, 1), x2 = c(0.5, 0.8))
theta <- c(0, 1, -1)

grad_fun2 <- function(x, theta) {
  c(1, x[1], x[2])
}

build_Xchi_grad(chi2, theta, grad_fun2)

Build \chi via centered finite differences

Description

Constructs the candidate regressor matrix Xchi by numerically approximating the gradient of the predictor \eta(x,\theta) with respect to the parameter vector theta using centered finite differences.

Usage

build_Xchi_numeric(chi, theta, eta_fun, step = NULL, rel_step = 1e-06)

Arguments

Details

Uses the centered scheme

\partial \eta(x, \theta) / \partial \theta_j \approx \{\eta(x, \theta + h_j e_j) - \eta(x, \theta - h_j e_j)\} / (2 h_j)

which has error order O(h^2).

Value

Numeric matrix N \times p with the finite-difference gradient approximation.

Compute an efficiency lower bound for a given design

Description

Computes a lower bound efficiency of a design \xi under the D-optimality criterion. The bound is based on the formulation proposed in the literature for correlated observations, and is used as a stopping criterion for the Virtual Noise algorithm.

Usage

efficiency_bound(xi_vec, Xchi, Sigma_chi, M_tilde, kappa, n_target)

Arguments

Details

The efficiency is expressed in terms of the information matrix M(\xi) and a deviation term \delta derived from the sensitivity function h(x, \xi).

The efficiency diagnostic is computed as

\mathrm{eff}_{\Phi}(\xi) \ge \frac{n \, \Phi(M(\xi))}{n \, \Phi(M(\xi)) + \kappa \, \delta},

where \Phi(M) = |M|^{1/p} and \delta = \max_x h(x, \xi) - \sum_x \xi(x) h(x, \xi).

The sensitivity function h(x, \xi) is constructed using the transformed matrix T(\xi) and the gradient of \Phi.

This quantity provides a practical measure of convergence and can be used as a stopping rule in iterative design algorithms.

Value

A list with the following components:

• eff_lb: lower bound for the efficiency of the design.

• delta: deviation term measuring how far the design is from optimality.

• phi: value of the D-optimality criterion in its multiplicative form |M(\xi)|^{1/p}.

• h_max: maximum value of the sensitivity function.

• h_xi_mean: weighted average of the sensitivity function under \xi.

Note

This function currently supports only the D-optimality criterion.

Convert a Virtual Noise design measure into an exact design

Description

Converts an approximate design measure obtained from vndesign_Dopt or virtual_noise_Dopt into an exact design with a fixed number of design points.

Usage

exact_design(
  res,
  Xchi,
  Sigma_chi,
  chi = NULL,
  n,
  method = c("random_best", "quantile", "quantile_endpoints"),
  B = 100,
  seed = NULL
)

Arguments

Value

A list with components:

• idx: selected candidate indices.

• points: selected candidate points.

• n: number of selected points.

• criterion_value: D-criterion value of the exact design.

Examples

chi <- seq(-1, 1, length.out = 10)
Xchi <- cbind(1, chi)
Sigma_chi <- build_Sigma(chi, type = "exp", range = 0.4)
res <- list(xi = rep(1 / length(chi), length(chi)))

exact_design(res, Xchi, Sigma_chi, chi = chi, n = 3, B = 5, seed = 1)

Plot convergence of the D-optimality criterion

Description

Plots the value of the D-optimality criterion at each iteration of the Virtual Noise algorithm in order to visualize the convergence of the approximate design sequence. The initial value at iteration 0 is omitted from the plot.

Usage

plot_convergence(res, main = NULL, xlab = "Iteration", ylab = NULL)

Arguments

Value

Produces a base R convergence plot and returns NULL invisibly.

Examples

## Run the Virtual Noise algorithm
chi <- seq(-1, 1, length.out = 50)
res <- vndesign_Dopt(
  chi = chi,
  theta = c(0, 1),
  grad_fun = function(x, theta) c(1, x),
  kappa = 0.99 / length(chi),
  max_iter = 200
)

## Plot convergence of the D-criterion
plot_convergence(res)

Plot design weights over the candidate set

Description

Plots the design weights returned by vndesign_Dopt or virtual_noise_Dopt. One-dimensional designs are shown as vertical lines; two-dimensional designs are shown as points with size proportional to the square root of the weight.

Usage

plot_design_weights(res, chi = NULL, main = NULL)

Arguments

Value

Produces a base R visualization of the design weights over the candidate set and returns NULL invisibly.

Examples

## Run the Virtual Noise algorithm
chi <- seq(-1, 1, length.out = 50)
res <- vndesign_Dopt(
  chi = chi,
  theta = c(0, 1),
  grad_fun = function(x, theta) c(1, x),
  kappa = 0.99 / length(chi),
  max_iter = 200
)

## Plot design weights
plot_design_weights(res, chi = chi)

Plot efficiency diagnostic

Description

Plots the efficiency diagnostic values stored in an object returned by vndesign_Dopt or virtual_noise_Dopt with return_eff_bound = TRUE. The initial value at iteration 0 is omitted from the plot.

Usage

plot_efficiency(res, main = NULL)

Arguments

Value

Produces a base R visualization of the efficiency diagnostic over the iterations of the Virtual Noise algorithm and returns NULL invisibly.

Examples

## Run the Virtual Noise algorithm with efficiency tracking
chi <- seq(-1, 1, length.out = 50)
res <- vndesign_Dopt(
  chi = chi,
  theta = c(0, 1),
  grad_fun = function(x, theta) c(1, x),
  kappa = 0.99 / length(chi),
  max_iter = 200,
  return_eff_bound = TRUE,
  n_target = 3
)

## Plot efficiency diagnostic
plot_efficiency(res)

Virtual Noise algorithm for D-optimal design

Description

Implements the Virtual Noise (VN) algorithm for constructing approximate D-optimal experimental designs under correlated observations.

Usage

virtual_noise_Dopt(
  Xchi,
  Sigma_chi,
  chi = NULL,
  kappa,
  delta = 0.001,
  max_iter = NULL,
  xi_init = NULL,
  verbose = TRUE,
  vn_alt = 3,
  vn_lambda = 40,
  vn_delta = 1e-06,
  stop_rule = c("max_iter", "efficiency"),
  return_eff_bound = FALSE,
  n_target = NULL,
  eff_tol = NULL
)

Arguments

Details

The algorithm iteratively updates a design measure over a candidate set chi by selecting candidate points according to D-optimality sensitivity rules and incorporating an additional virtual-noise penalization mechanism.

At each iteration, the method:

1. Constructs a modified covariance matrix including the virtual-noise term,

2. Computes the corresponding information matrix,

3. Evaluates D-optimality selection quantities over the candidate set,

4. Selects a candidate point according to the VN decision rules,

5. Updates the design measure exactly.

The algorithm updates the weight vector xi over the candidate set by selecting, at each iteration, a candidate point x_s and performing the exact update

\xi^{(s)} = \frac{s-1}{s}\xi^{(s-1)} + \frac{1}{s}e_{x_s},

where e_{x_s} denotes the unit vector associated with the selected candidate point.

The VN component modifies the covariance structure through an additional penalization term controlled by kappa and the selected VN alternative vn_alt.

The D-criterion is implemented as -log det(M_tilde) where M_tilde = X^T Sigma_tilde^{-1} X and Sigma_tilde = Sigma_chi + diag(sigma2_xi).

This implementation follows the classical VN heuristic. The optional efficiency values returned by return_eff_bound = TRUE are useful convergence diagnostics; they should not be read as certificates from the convex linear-programming formulation of Pázman et al. (2022).

Value

A list with components:

• xi: final weights (length N).

• crit_vals: D-criterion values over iterations.

• iterations: number of executed iterations.

• chosen_indices: indices selected at each iteration.

• chosen_points: selected points from chi.

• support_idx: indices with positive weight at the end.

• support_points: support points from chi.

• support_weights: normalized support weights.

If return_eff_bound = TRUE, the list also contains eff_lb_vals, delta_vals, phi_vals, h_max_vals, and h_xi_mean_vals.

References

Pázman, A. (1999). Further developments on virtual noise methods for exact optimal designs.

Pázman, A., Hainy, M., and Müller, W. G. (2022). A convex approach to optimum design of experiments with correlated observations.

López-Fidalgo, J., and Wong, W. K. (2025). Optimal Designs for Correlated Data. Annual Review of Statistics and Its Application.

Main interface for constructing D-optimal designs using the Virtual Noise algorithm

Description

This is the main user-level function of the package. It provides a complete workflow for constructing approximate D-optimal designs under correlated observations using the Virtual Noise (VN) algorithm implemented in virtual_noise_Dopt.

Usage

vndesign_Dopt(
  chi,
  Xchi = NULL,
  theta = NULL,
  grad_fun = NULL,
  eta_fun = NULL,
  kappa,
  max_iter = NULL,
  xi_init = NULL,
  delta = 0.001,
  verbose = TRUE,
  stop_rule = c("max_iter", "efficiency"),
  cov_type = c("exp", "gauss", "matern32", "matern52", "spherical", "triangular",
    "triangular_trunc", "brownian", "ar1", "custom"),
  sigma2 = 1,
  range = 1,
  nugget = 0,
  phi = NULL,
  cov_fun = NULL,
  vn_alt = 3,
  vn_lambda = 40,
  vn_delta = 1e-06,
  num_step = NULL,
  num_rel_step = 1e-06,
  return_components = FALSE,
  return_eff_bound = FALSE,
  n_target = NULL,
  eff_tol = NULL
)

Arguments

Details

The function acts as a high-level wrapper integrating most of the package functionality, including:

• construction or validation of the candidate regressor matrix Xchi,

• analytical or numerical gradient evaluation,

• covariance matrix construction via build_Sigma,

• execution of the Virtual Noise algorithm,

• computation of optional efficiency diagnostics,

• and generation of visualization-ready outputs.

The candidate set chi may be:

• a numeric vector (one-dimensional case),

• or a matrix / data.frame of dimension N \times d, where each row represents a candidate point.

Two modes of operation are supported:

Direct mode: If Xchi is provided, it is used directly and no model specification through theta, grad_fun, or eta_fun is required.

Model-based mode: If Xchi is NULL, the function constructs it from theta using either:

• an analytical gradient function grad_fun,

• or numerical finite differences via eta_fun.

This function provides a high-level interface to the Virtual Noise methodology implemented in virtual_noise_Dopt.

For the D-criterion, the implementation minimizes -\log |M|, which is equivalent to maximizing the determinant of the information matrix.

When return_eff_bound = TRUE, the function also computes a Pázman-type D-efficiency diagnostic along the iterations. If stop_rule = "efficiency", this diagnostic is additionally used as a stopping criterion, while max_iter acts as a safety limit. Since this package implements the classical VN heuristic rather than the convex linear-programming formulation of Pázman et al. (2022), these values should be interpreted as convergence diagnostics.

Value

A list containing the output of virtual_noise_Dopt. Main components include:

• xi: final design weights over the candidate set.

• crit_vals: D-criterion values over the iterations.

• iterations: number of executed iterations.

• support_idx: indices with positive final weights.

• support_points: support points from chi.

• support_weights: normalized support weights.

If return_eff_bound = TRUE, the output also includes the efficiency diagnostic values and related quantities. If return_components = TRUE, the output additionally includes Xchi, Sigma_chi, and D.

See Also

virtual_noise_Dopt, efficiency_bound, build_Sigma, build_Xchi_grad, build_Xchi_numeric

Examples

chi <- seq(-1, 1, length.out = 15)
theta <- c(0, 1)
grad_fun <- function(x, theta) c(1, x)
kappa <- 0.99 / length(chi)

res <- vndesign_Dopt(
  chi = chi,
  theta = theta,
  grad_fun = grad_fun,
  kappa = kappa,
  max_iter = 5,
  verbose = FALSE,
  cov_type = "exp",
  range = 0.4,
  vn_alt = 3,
  vn_lambda = 2
)

head(res$xi)