--- title: "Using splines2 with Rcpp" author: Wenjie Wang date: "`r Sys.Date()`" bibliography: - ../inst/bib/splines2.bib vignette: > %\VignetteIndexEntry{Using splines2 with Rcpp} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} output: rmarkdown::html_vignette: number_sections: yes toc: yes --- # Introduction In this vignette, we introduce how to use the C++ header-only library that **splines2** contains with the **Rcpp** package [@eddelbuettel2013springer] for constructing spline basis functions directly in C++. The introduction is intended for package developers who would like to use **splines2** package in C++ by adding **splines2** to the `LinkingTo` field of the package `DESCRIPTION` file. # Header Files and Namespace Different from the procedure-based functions in the R interface, the C++ interface follows the commonly-used object-oriented design in C++ for ease of usage and maintenance. The implementations use the **Armadillo** [@sanderson2016armadillo] library with the help of **RcppArmadillo** [@eddelbuettel2014csda] and require C++11. We assume that C++11 is enabled and the header file named `splines2Armadillo.h` is included for access to all the classes and implementations in the namespace `splines2` henceforth. ```cpp #include #include // include header files from splines2 // [[Rcpp::plugins(cpp11)]] ``` To use `Rcpp::sourceCpp()`, one may need to add `[[Rcpp::depends()]]` as follows: ```cpp // [[Rcpp::depends(RcppArmadillo)]] // [[Rcpp::depends(splines2)]] ``` For ease of demonstration, we assume the following *using-directives*: ```cpp using namespace arma using namespace splines2 ``` # Classes for Spline Basis Functions A virtual base class named `SplineBase` is implemented to support a variety of classes for spline basis functions including - `BSpline` for B-splines; - `MSpline` for M-splines; - `ISpline` for I-splines; - `CSpline` for C-splines; - `NaturalSpline` and `NaturalSplineK` for natural cubic splines; - `PeriodicMSpline` for periodic M-splines; - `PeriodicBSpline` for periodic B-splines; ## Constructors of `BSpline`, `MSpline`, `ISpline`, and `CSpline` The `BSpline`, `MSpline`, `ISpline`, and `CSpline` classes share the same constructors inherited from the `SplineBase` class. There are four constructors in addition to the default constructor. The first non-default constructor is invoked when internal knots are explicitly specified as the second argument. Taking B-splines as an example, the first non-default constructor of a `BSpline` object is ```cpp // 1. specify x, internal knots, degree, and boundary knots BSpline(const vec& x, const vec& internal_knots, const unsigned int degree = 3, const vec& boundary_knots = vec()); ``` The second non-default constructor is called when an unsigned integer is specified as the second argument, which represents the degree of freedom (DF) of the *complete spline basis functions* (different from the `df` argument in the R interface) is specified. Then the number of internal knots is computed as `spline_df - degree - 1` and the placement of internal knots uses quantiles of specified `x` within the boundary. ```cpp // 2. specify x, spline DF, degree, and boundary knots BSpline(const vec& x, const unsigned int spline_df, const unsigned int degree = 3, const vec& boundary_knots = vec()); ``` The third non-default constructor is intended for the basis functions with an extended knot sequence, where the multiplicities of the knots can be more than one. ```cpp // 3. specify x, degree, and (extended) knot sequence BSpline(const vec& x, const unsigned int degree, const vec& knot_sequence); ``` The fourth non-default constructor is explicit and takes a pointer to a base class object, which can be useful when we want to create a new object using the same specification (`x`, `degree`, `internal_knots`, and `boundary_knots`) of an existing object (not necessarily a `BSpline` object). ```cpp // 4. create a new object from a base class pointer BSpline(const SplineBase* pSplineBase); ``` This constructor also allows us to easily switch between different types of splines. For example, we can create a `BSpline` object named `bsp_obj` from an existing `MSpline` object named `msp_obj` with the same specification as follows: ```cpp BSpline bsp_obj { &msp_obj }; ``` ## Constructors of `PeriodicMSpline` and `PeriodicBSpline` The `PeriodicMSpline` and `PeriodicBSpline` classes are intended for constructing the periodic M-splines and periodic B-splines, respectively, which provide the same set of non-default constructors with `BSpline`. The only difference is that the knot sequence specified for the third non-default constructor must be a *simple knot sequence*. ## Constructors of `NaturalSpline` and `NaturalSplineK` The classes `NaturalSpline` and `NaturalSplineK` are intended for natural cubic splines. The former corresponds to the function `splines2::naturalSpline()` (or `splines2::nsp()`) in R, while the latter is the engine of the function `splines2::nsk()`. They have the same constructors that do not allow the specification of the `degree`. Taking `NaturalSpline` as an example, the first non-default constructor is called when internal knots are explicitly specified. ```cpp // 1. specify x, internal knots, and boundary knots NaturalSpline(const vec& x, const vec& internal_knots, const vec& boundary_knots = vec()); ``` The second non-default constructor is called when an unsigned integer representing the degree of freedom of the *complete spline basis functions* (different from the `df` argument in the R interface) is specified. Then the number of internal knots is computed as `spline_df - 2` and the placement of internal knots uses quantiles of specified `x`. ```cpp // 2. specify x, spline DF, and boundary knots NaturalSpline(const vec& x, const unsigned int spline_df, const vec& boundary_knots = vec()); ``` The third non-default constructor is explicit and takes a pointer to a base class object. It can be useful when we want to create a new object using the same specification (`x`, `internal_knots`, `boundary_knots`, etc.) of an existing object. ```cpp // 3. create a new object from a base class pointer NaturalSpline(const SplineBase* pSplineBase); ``` ## Function Members The main methods are - `basis()` for spline basis matrix - `derivative()` for derivatives of spline basis - `integral()` for integrals of spline basis (except for the `CSpline` class) The specific function signatures are as follows: ```cpp mat basis(const bool complete_basis = true); mat derivative(const unsigned int derivs = 1, const bool complete_basis = true); mat integral(const bool complete_basis = true); ``` We can set and get the spline specifications through the following *setter* and *getter* functions, respectively. ```cpp // setter functions SplineBase* set_x(const vec&); SplineBase* set_x(const double); SplineBase* set_internal_knots(const vec&); SplineBase* set_boundary_knots(const vec&); SplineBase* set_knot_sequence(const vec&); SplineBase* set_degree(const unsigned int); SplineBase* set_order(const unsigned int); // getter functions vec get_x(); vec get_internal_knots(); vec get_boundary_knots(); vec get_knot_sequence(); unsigned int get_degree(); unsigned int get_order(); unsigned int get_spline_df(); ``` The *setter* function returns a pointer to the current object so that the specification can be chained for convenience. For example, ```cpp vec x { arma::regspace(0, 0.1, 1) }; // 0, 0.1, ..., 1 BSpline obj { x, 5 }; // df = 5 (and degree = 3, by default) // change degree to 2 and get basis mat basis_mat { obj.set_degree(2)->basis() }; ``` The corresponding first derivatives and integrals of the basis functions can be obtained as follows: ```cpp mat derivative_mat { bs.derivative() }; mat integral_mat { bs.integral() }; ``` Notice that there is no available `integral()` method for `CSpline` and no meaningful `degree` related methods for `NaturalSpline`. # Generalized Bernstein Polynomials The `BernsteinPoly` class is provided for the generalized Bernstein polynomials. ## Constructors The main non-default constructor is as follows: ```cpp BernsteinPoly(const vec& x, const unsigned int degree, const vec& boundary_knots = vec()); ``` In addition, two explicit constructors are provided for `BernsteinPoly*` and `SplineBase*`, which set `x`, `degree`, and `boundary_knots` from the objects that the pointers point to. ## Function Members The main methods are - `basis()` for the basis functions - `derivative()` for the derivatives of basis functions - `integral()` for the integrals of basis functions The specific function signatures are as follows: ```cpp mat basis(const bool complete_basis = true); mat derivative(const unsigned int derivs = 1, const bool complete_basis = true); mat integral(const bool complete_basis = true); ``` In addition, we may *set* and *get* the specifications through the following *setter* and *getter* functions, respectively. ```cpp // setter functions BernsteinPoly* set_x(const vec&); BernsteinPoly* set_x(const double); BernsteinPoly* set_degree(const unsigned int); BernsteinPoly* set_order(const unsigned int); BernsteinPoly* set_internal_knots(const vec&); // placeholder, does nothing BernsteinPoly* set_boundary_knots(const vec&); // getter functions vec get_x(); unsigned int get_degree(); unsigned int get_order(); vec get_boundary_knots(); ``` The *setter* function returns a pointer to the current object. # Reference {-}