Commit 2e77c2fb authored by davidkep's avatar davidkep

Initial commit

parent fbfbe0fb
Package: nsoptim
Type: Package
Title: Utility Library to Implement Non-Smooth Optimization in C++
Version: 0.1.0
Date: 2019-01-30
Authors@R: person("David", "Kepplinger", , "david.kepplinger@gmail.com",
role = c("aut", "cre"))
Copyright: See the file COPYRIGHTS for copyright details on some of the
functions and algorithms used.
Encoding: UTF-8
SystemRequirements: C++11
URL: https://github.com/dakep/nsoptim-rpkg
BugReports: https://github.com/dakep/nsoptim-rpkg/issues
Description: A C++ library to help solve non-smooth optimziation problems
in R packages.
Depends:
R (>= 3.4.0)
Imports:
Rcpp
LinkingTo:
Rcpp,
RcppArmadillo (>= 0.9.100)
License: GPL (>= 2)
RoxygenNote: 6.0.1
# Empty namespace
inlineCxxPlugin <- function(...) {
plugin <- Rcpp::Rcpp.plugin.maker(package = 'nsoptim',
include.before = '#include <nsoptim.hpp>')
settings <- plugin()
settings$env$PKG_CPPFLAGS <- '-I../inst/include'
settings$env$USE_CXX11 <- 'yes'
settings
}
# nsoptim R Package
## Interfaces
### Loss Function
Loss functions which operate on the data-type `Data` must derive from the base class `LossFunction<Data>`.
Loss functions must implement a copy constructor which should be designed for fast and shallow copying of the loss.
Loss functions must support the following operations for one or more `Coefficients` types and the loss function's
data type `Data`.
#### `LossFunction(const LossFunction& other)` and `LossFunction& operator=(const LossFunction& other)`
Copy constructor/assignment operator, potentially creating a shallow copy of the `other` loss function.
This is called several times often and should be designed with speed in mind.
#### `LossFunction Clone()`
Create a deep copy of the loss function, i.e., should not share any data with this loss function, except for the
constant, underlying data.
#### `Coefficients ZeroCoefficients() const`
Create an appropriate `Coefficients` object which can be evaluated by this loss function, with all values set to zero.
#### `const Data& data() const`
Return a constant reference to the data the loss object operates on.
#### `double RelativeDifference(const Coefficients& x, const Coefficients& y) const`
Compute the relative difference between the two `Coefficients` objects `x` and `y`.
#### `double operator()(const Coefficients& where) const`
Evaluate the loss function at the coefficient values `where`.
### Penalty Function
Penalty functions must support evaluation of the penalty at one or more `Coefficients` types.
Loss functions must implement a copy constructor which should be designed for fast and shallow copying of the loss.
#### `double operator()(const Coefficients& where) const`
Evaluate the penalty at the coefficient values `where` and return the numeric value.
### Convex Functions or Convex Surrogates
Any loss or penalty function can define a method to create a convex surrogate at a particual value of the coefficients.
If a function `Function` derives from the type `ConvexFunction<Function>`, the convex surrogate is automatically
determined to be the function itself.
A non-convex function which wants to provide a convex surrogate must implement the following method for all supported
coefficient types `Coefficients`.
#### `const Function& ConvexSurrogate(const Coefficients& where) const`
Return a reference (or constant reference) to the convex surrogate of this function at the coefficients values `where`.
A convex surrogate function must be a convex function which is above this function everywhere and equals this
function at `where`.
In other words, if this function is denoted by _f_ and the convex surrogate is denoted by _g_, _g_ must be a convex
function such that _f(x) <= g(x)_ for all coefficients values `x` and _f(where) = g(where)_ at the given coefficients
values `where`.
### Optimizer
Optimizer operate on the additive objective function _f(x; d) = l(x; d) + p(x)_ where _l_ is the loss function operating
on data _d_, and _p_ is the penalty function. These have the following types:
* _l_ is of type `LossFunction`,
* _p_ is of type `PenaltyFunction`,
* and _x_ is of type `Coefficients`.
The type of _d_ is determined by the loss function _l_, nameley `LossFunction::DataType`.
An optimizer must implement the following methods.
#### `Optimizer::Optimum Optimize()` and `Optimizer::Optimum Optimize(const Coefficients& start)`
Optimize the current objective function, optionally starting at the coefficients value `start`.
#### `LossFunction& loss()`
Give access to the loss function currently in use. This allows the change of variable hyper-parameters which
don't change the loss in a significant way.
#### `void loss(const LossFunction& loss)`
Replace the current loss function with a new one.
#### `PenaltyFunction& penalty()`
Give access to the penalty function currently in use. This allows the change of variable hyper-parameters which
don't change the penalty in a significant way.
#### `void penalty(const Penalty& loss)`
Replace the current penalty function with a new one.
#### `Optimizer Clone()`
Clone the optimizer such that it operates on it's own clones of the loss and penalty functions and does not share
any state with this optimizer.
#### `Optimizer(const Optimizer& other)` and `Optimizer& (const Optimizer& other)`
Copy constructor/assignment operator to create a new optimizer which inherits the state of the `other` optimizer and
maybe share some internal data.
### Iterative Optimizer
Iterative optimizer should additionally provide a method to set the convergence tolerance to different levels.
Convergence is determined by the optimizer in different ways, but the optimizer should ensure that the convergence
criteria is comparable to the relative difference between two values of the coefficients as determined by the
loss function in use.
#### `void convergence_tolerance(const double tolerance)`
Set the convergence tolerance to `tolerance`.
//
// nsoptim.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_HPP_
#define NSOPTIM_HPP_
#if defined(Rcpp_hpp)
#error "The file 'Rcpp.h' should not be included. Please correct to include only 'nsoptim.hpp'."
#endif
#if defined(RcppArmadillo__RcppArmadillo__h)
#error "The file 'RcppArmadillo.h' should not be included. Please correct to include only 'nsoptim.hpp'."
#endif
#include <nsoptim_forward.hpp>
#include <nsoptim/rcpp_integration.hpp>
#endif // NSOPTIM_HPP_
//
// armadillo.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_ARMADILLO_HPP_
#define NSOPTIM_ARMADILLO_HPP_
#define ARMA_USE_CXX11 1
#define ARMA_DONT_USE_OPENMP 1
#ifdef __clang__
# pragma clang diagnostic push
# pragma clang diagnostic ignored "-Weverything"
#endif
#ifdef HAVE_RCPP
// For RCPP
# include <RcppArmadilloForward.h>
#else
// For stand-alone
# include <armadillo>
#endif // HAVE_RCPP
#ifdef __clang__
# pragma clang diagnostic pop
#endif
#endif // NSOPTIM_ARMADILLO_HPP_
//
// armadillo.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_CONFIG_HPP_
#define NSOPTIM_CONFIG_HPP_
#ifndef NSOPTIM_METRICS_ENABLED
#define NSOPTIM_METRICS_ENABLED_BOOL false
#else
#define NSOPTIM_METRICS_ENABLED_BOOL true
#endif
#endif // NSOPTIM_CONFIG_HPP_
//
// container.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_CONTAINER_HPP_
#define NSOPTIM_CONTAINER_HPP_
#include <nsoptim/container/data.hpp>
#include <nsoptim/container/metrics.hpp>
#include <nsoptim/container/regression_coefficients.hpp>
#endif // NSOPTIM_CONTAINER_HPP_
//
// data.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_CONTAINER_DATA_HPP_
#define NSOPTIM_CONTAINER_DATA_HPP_
#include <iostream>
#include <nsoptim/armadillo.hpp>
namespace nsoptim {
//! Simple structure holding a matrix with predictor data *x* (with `n_obs` rows and `n_pred` colums)
//! and the response vector *y* (with `n_obs` rows)
class PredictorResponseData {
public:
//! Initialize the predictor-response data with empty x and y
PredictorResponseData() noexcept : x_(), y_(), n_obs_(0), n_pred_(0) {}
//! Initialize predictor-response data with the given x and y.
//! @note the given data will be copied!
//!
//! @param other_x predictor matrix to copy.
//! @param other_y response vector to copy.
PredictorResponseData(const arma::mat& other_x, const arma::vec& other_y) noexcept
: x_(other_x), y_(other_y), n_obs_(other_x.n_rows), n_pred_(other_x.n_cols) {}
//! Initialize predictor-response data with the given x and y.
//! @note the given data will be moved to this container!
//!
//! @param other_x predictor matrix to move.
//! @param other_y response vector to move.
PredictorResponseData(arma::mat&& other_x, arma::vec&& other_y) noexcept
: x_(std::move(other_x)), y_(std::move(other_y)), n_obs_(x_.n_rows), n_pred_(x_.n_cols) {}
//! Copy the given predictor-response data, but pointing to the same underlying data!
//!
//! @param other predictor-response data to copy.
PredictorResponseData(const PredictorResponseData& other) = default;
PredictorResponseData& operator=(const PredictorResponseData& other) = default;
//! Get a data set with the observations at the requested indices.
//!
//! @param indices the indicies of the observations to get.
//! @return the subset of the data with the requested observations.
inline PredictorResponseData observations(const arma::uvec& indices) const {
return PredictorResponseData(x_.rows(indices), y_.rows(indices));
}
//! Get a data set with the first `n_obs` observations of the data.
//!
//! @param n_obs number of observations to extract.
//! @return the subset of the data with the requested observations.
inline PredictorResponseData HeadRows(const arma::uword n_obs) const {
return PredictorResponseData(x_.head_rows(n_obs), y_.head_rows(n_obs));
}
//! Get a data set with the last `n_obs` observations of the data.
//!
//! @param n_obs number of observations to extract.
//! @return the subset of the data with the requested rows.
inline PredictorResponseData TailRows(const arma::uword n_obs) const {
return PredictorResponseData(x_.tail_rows(n_obs), y_.tail_rows(n_obs));
}
//! Get a constant reference to the predictor matrix.
//! Only valid as long as the PredictorResponseData object is in scope.
//!
//! @return constant reference to the predictor matrix
inline const arma::mat& cx() const noexcept {
return x_;
}
//! Get a constant reference to the response vector.
//! Only valid as long as the PredictorResponseData object is in scope.
//!
//! @return constant reference to the response vector
inline const arma::vec& cy() const noexcept {
return y_;
}
//! Get non-const references to the data
//! Get a reference to the predictor matrix.
//! Only valid as long as the PredictorResponseData object is in scope.
//!
//! @return reference to the predictor matrix
inline arma::mat& x() noexcept {
return x_;
}
//! Get a reference to the response vector.
//! Only valid as long as the PredictorResponseData object is in scope.
//!
//! @return reference to the response vector
inline arma::vec& y() noexcept {
return y_;
}
//! Get the number of observations in this data set.
inline arma::uword n_obs() const noexcept {
return n_obs_;
}
//! Get the number of observations in this data set.
inline arma::uword n_pred() const noexcept {
return n_pred_;
}
private:
arma::mat x_;
arma::vec y_;
arma::uword n_obs_; //< The number of observations in the data.
arma::uword n_pred_; //< The number of variables in the data.
};
} // namespace nsoptim
#endif // NSOPTIM_CONTAINER_DATA_HPP_
//
// metrics.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_CONTAINER_METRICS_HPP_
#define NSOPTIM_CONTAINER_METRICS_HPP_
#include <memory>
#include <string>
#include <forward_list>
#include <nsoptim/config.hpp>
namespace nsoptim {
//! Simple wrapper around a named metric.
template<typename T>
struct Metric {
Metric(const std::string& _name, const T _value) : name(_name), value(_value) {}
Metric(Metric&&) = default;
Metric& operator=(Metric&&) = default;
Metric(const Metric&) = default;
Metric& operator=(const Metric&) = default;
const std::string name;
const T value;
};
#ifdef NSOPTIM_METRICS_ENABLED
#else
#endif
namespace _metrics_internal {
template<bool>
class Metrics {
public:
//! Create a metrics collection with the given namespace
//!
//! @param name the name(space) for the metrics collection.
explicit Metrics(const std::string&) noexcept {}
Metrics(Metrics&&) = default;
Metrics& operator=(Metrics&&) = default;
Metrics(const Metrics&) = default;
Metrics& operator=(const Metrics&) = default;
//! Add a metric to the collection.
//!
//! @param name the name of the metric.
//! @param value the value for the metric.
void AddMetric(const std::string&, const int) noexcept {}
//! Add a metric to the collection.
//!
//! @param name the name of the metric.
//! @param value the value for the metric.
void AddMetric(const std::string&, const double) noexcept {}
//! Add a sub-collection of metrics to the collection.
//!
//! @param metrics a collection of metrics to be added to this metrics collection.
void AddSubMetrics(const Metrics&) noexcept {}
//! Get the name(space) of this Metrics object.
//!
//! @return the name(space).
std::string name() const noexcept {
return "";
}
//! Get all integer metrics in this collection.
//!
//! @return a reference to the named integer metrics.
const std::forward_list<const Metric<int>>& IntegerMetrics() const {
static const std::forward_list<const Metric<int>> tmp;
return tmp;
}
//! Get all double (numeric) metrics in this collection.
//!
//! @return a reference to the named double (numeric) metrics.
const std::forward_list<const Metric<double>>& DoubleMetrics() const {
static const std::forward_list<const Metric<double>> tmp;
return tmp;
}
//! Get all sub-metrics in this collection.
//!
//! @return a reference to the sub-metrics.
const std::forward_list<const Metrics>& SubMetrics() const {
static const std::forward_list<const Metrics> tmp;
return tmp;
}
};
template<>
class Metrics<true> {
public:
//! Create a metrics collection with the given namespace
//!
//! @param name the name(space) for the metrics collection.
explicit Metrics(const std::string& name) : name_(name) {}
Metrics(Metrics&&) = default;
Metrics& operator=(Metrics&&) = default;
Metrics(const Metrics&) = default;
Metrics& operator=(const Metrics&) = default;
//! Add a metric to the collection.
//!
//! @param name the name of the metric.
//! @param value the value for the metric.
void AddMetric(const std::string& name, const int value) {
int_metrics_.emplace_front(name, value);
}
//! Add a metric to the collection.
//!
//! @param name the name of the metric.
//! @param value the value for the metric.
void AddMetric(const std::string& name, const double value) {
dbl_metrics_.emplace_front(name, value);
}
//! Add a sub-collection of metrics to the collection.
//!
//! @param metrics a collection of metrics to be added to this metrics collection.
void AddSubMetrics(const Metrics& metrics) {
sub_metrics_.push_front(metrics);
}
//! Get the name(space) of this Metrics object.
//!
//! @return the name(space).
std::string name() const noexcept {
return name_;
}
//! Get all integer metrics in this collection.
//!
//! @return a reference to the named integer metrics.
const std::forward_list<const Metric<int>>& IntegerMetrics() const {
return int_metrics_;
}
//! Get all double (numeric) metrics in this collection.
//!
//! @return a reference to the named double (numeric) metrics.
const std::forward_list<const Metric<double>>& DoubleMetrics() const {
return dbl_metrics_;
}
//! Get all sub-metrics in this collection.
//!
//! @return a reference to the sub-metrics.
const std::forward_list<const Metrics>& SubMetrics() const {
return sub_metrics_;
}
private:
const std::string name_;
std::forward_list<const Metric<int>> int_metrics_;
std::forward_list<const Metric<double>> dbl_metrics_;
std::forward_list<const Metrics> sub_metrics_;
};
} // namespace _metrics_internal
//! Export the default templated Metrics collection as `OptionalMetrics`.
using Metrics = _metrics_internal::Metrics<NSOPTIM_METRICS_ENABLED_BOOL>;
} // namespace nsoptim
#endif // NSOPTIM_CONTAINER_METRICS_HPP_
//
// regression_coefficients.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_CONTAINER_REGRESSION_COEFFICIENTS_HPP_
#define NSOPTIM_CONTAINER_REGRESSION_COEFFICIENTS_HPP_
#include <memory>
#include <type_traits>
#include <nsoptim/armadillo.hpp>
namespace nsoptim {
//! Regression coefficients in the form of a scalar intercept and a (sparse) vector of slope coefficients.
//! The slope coefficients must be either of type `arma::vec` or `arma::sp_vec`.
template <class T>
class RegressionCoefficients {
static_assert(std::is_same<T, arma::vec>() || std::is_same<T, arma::sp_vec>(),
"T must be a (sparse) vector.");
public:
//! The dimensions of regression coefficients is the number if *slope* parameters.
//! The intercept is not included in this number!
using Dimensions = int;
//! The type of the regression slope vector.
using SlopeCoefficient = T;
//! Initialize an empty coefficient vector and an intercept of 0.
RegressionCoefficients() noexcept : intercept(0), beta() {}
//! Initialize the coefficients as 0-vector
//!
//! @param n_pred the number of predictors, i.e., number of slope coefficients.
explicit RegressionCoefficients(const arma::uword n_pred) noexcept
: RegressionCoefficients(n_pred, typename std::is_same<T, arma::sp_vec>::type{}) {}
//! Copy the slope coefficients from the given vector. The intercept will be set to 0.
//!
//! @param beta the value for the slope coefficients.
explicit RegressionCoefficients(const SlopeCoefficient& beta) noexcept : intercept(0), beta(beta) {}
//! Copy the coefficients from the given intercept and slope.
//!
//! @param intercept the value for the intercept.
//! @param beta the value for the slope coefficients.
RegressionCoefficients(const double intercept, const SlopeCoefficient& beta) noexcept
: intercept(intercept), beta(beta) {}
//! Copy constructor
//!
//! @param coef the RegressionCoefficients object to copy.
RegressionCoefficients(const RegressionCoefficients& coef) = default;
//! Move constructor.
//!
//! @param coef the RegressionCoefficients object to move.
RegressionCoefficients(RegressionCoefficients&& coef) = default;
//! Copy assignment
//!
//! @param coef the RegressionCoefficients object to copy.
RegressionCoefficients& operator=(const RegressionCoefficients& other) = default;
//! Move assignment
//!
//! @param coef the RegressionCoefficients object to move.
RegressionCoefficients& operator=(RegressionCoefficients&& other) = default;
//! Reset the coefficients to an empty vector and an intercept of 0.
inline void Reset() noexcept {
intercept = 0;
beta.reset();
}
//! Print the coefficients to the given output stream
//!
//! @param user_stream the stream to print the coefficients to.
inline void Print(std::ostream &user_stream) const {
user_stream << "Intercept: " << intercept << "\n";
beta.t().print(user_stream, "Beta:");
}
//! Print the coefficients to stdout
inline void Print() const {
std::cout << "Intercept: " << intercept << "\n";
beta.t().print("Beta:");
}
double intercept; //< The intercept coefficient
SlopeCoefficient beta; //< The slope coefficients
private:
// Initialize the regression coefficients with all zeros if the slope is a dense vector.
RegressionCoefficients(const arma::uword n_pred, std::false_type) noexcept
: intercept(0), beta(n_pred, arma::fill::zeros) {}
// Initialize the regression coefficients with all zeros if the slope is a sparse vector.
RegressionCoefficients(const arma::uword n_pred, std::true_type) noexcept
: intercept(0), beta(n_pred) {}
};
} // namespace nsoptim
#endif // NSOPTIM_CONTAINER_REGRESSION_COEFFICIENTS_HPP_
//
// objective.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_OBJECTIVE_HPP_
#define NSOPTIM_OBJECTIVE_HPP_
// Loss functions
#include <nsoptim/objective/ls_loss.hpp>
// Penalty functions
#include <nsoptim/objective/en_penalty.hpp>
#include <nsoptim/objective/adaptive_en_penalty.hpp>
#endif // NSOPTIM_OBJECTIVE_HPP_
//
// adaptive_en_penalty.hpp
// nsoptim
//
// Created by David Kepplinger on 2018-11-30.
// Copyright © 2018 David Kepplinger. All rights reserved.
//
#ifndef NSOPTIM_OBJECTIVE_ADAPTIVE_EN_PENALTY_HPP_
#define NSOPTIM_OBJECTIVE_ADAPTIVE_EN_PENALTY_HPP_
#include <nsoptim/armadillo.hpp>
#include <nsoptim/objective/convex.hpp>
#include <nsoptim/container/regression_coefficients.hpp>
#include <nsoptim/objective/penalty.hpp>
namespace nsoptim {
//! The adaptive elastic net penalty function defined as
//! $P(\beta; \lambda_1, \lambda_2, w) = \sum_{j=1}^p \lambda_1 w_j |\beta_j| + \lambda_2 \beta_j^2$
//! with a p-dimensional vector of penalty loadings, $w$, and hyper-parameters $\lambda_1$ and $\lambda_2$.
class AdaptiveEnPenalty : public PenaltyFunction, public ConvexFunction<AdaptiveEnPenalty> {
public:
//! Declare this penalty function as an EN penalty.
struct is_en_penalty_tag {};
template<typename T>
using GradientType = typename T::SlopeCoefficient;
explicit AdaptiveEnPenalty(const arma::vec& loadings) noexcept : loadings_(loadings), lambda_1_(0), lambda_2_(0) {}
AdaptiveEnPenalty(const arma::vec& loadings, const double lambda_1, const double lambda_2) noexcept
: loadings_(loadings), lambda_1_(lambda_1), lambda_2_(lambda_2) {}
AdaptiveEnPenalty(const AdaptiveEnPenalty& other) = default;
AdaptiveEnPenalty& operator=(const AdaptiveEnPenalty& other) = default;
AdaptiveEnPenalty(AdaptiveEnPenalty&& other) = default;
AdaptiveEnPenalty& operator=(AdaptiveEnPenalty&& other) = default;
~AdaptiveEnPenalty() = default;
AdaptiveEnPenalty Copy() {
return AdaptiveEnPenalty(*this);
}
AdaptiveEnPenalty Clone() const {
return AdaptiveEnPenalty(*this);
}
void lambda_1(const double lambda_1) noexcept {
lambda_1_ = lambda_1;
}
double lambda_1() const noexcept {
return lambda_1_;
}
void lambda_2(const double lambda_2) noexcept {
lambda_2_ = lambda_2;
}
double lambda_2() const noexcept {
return lambda_2_;
}