ortho_optim: Orthogonality constrained optimization

Description

A general purpose optimization solver with orthogonality constraint. The orthogonality constrained optimization method is a nearly direct translation from Wen and Yin (2010)'s MATLAB code.

Usage

ortho_optim(
  B,
  fn,
  grad = NULL,
  ...,
  maximize = FALSE,
  control = list(),
  maxitr = 500,
  verbose = FALSE
)

Value

A orthoDr object that consists of a list with named entries of:

B: The optimal B value
fn: The final functional value
itr: The number of iterations
converge: convergence code

Arguments

B: Initial B values. Must be a matrix, and the columns are subject to the orthogonality constrains. Will be processed by Gram-Schmidt if not orthogonal
fn: A function that calculate the objective function value. The first argument should be B. Returns a single value.
grad: A function that calculate the gradient. The first argument should be B. Returns a matrix with the same dimension as B. If not specified, then numerical approximation is used.
...: Arguments passed to fn and grad
maximize: By default, the solver will try to minimize the objective function unless maximize = TRUE
control: A list of tuning variables for optimization. epsilon is the size for numerically approximating the gradient. For others, see Wen and Yin (2013).
maxitr: Maximum number of iterations
verbose: Should information be displayed

References

Wen, Z., & Yin, W. (2013). A feasible method for optimization with orthogonality constraints. Mathematical Programming, 142(1), 397-434. DOI: tools:::Rd_expr_doi("10.1007/s10107-012-0584-1")

Examples

Run this code

# an eigen value problem
library(pracma)
set.seed(1)
n <- 100
k <- 6
A <- matrix(rnorm(n * n), n, n)
A <- t(A) %*% A
B <- gramSchmidt(matrix(rnorm(n * k), n, k))$Q

fx <- function(B, A) -0.5 * sum(diag(t(B) %*% A %*% B))
gx <- function(B, A) -A %*% B
fit <- ortho_optim(B, fx, gx, A = A)
fx(fit$B, A)

# compare with the solution from the eigen function
sol <- eigen(A)$vectors[, 1:k]
fx(sol, A)