Functions
AMA_MltvLstsqr.c File Reference
#include <AMA.h>

Functions

long int AMA_MltvLstsqr (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *wht, long int *degree, long int *mlamda, double **lamda, double theta, AMA_SPLINE **spline)
 Least Squares Approximation of Multivariate Data More...
 

Function Documentation

long int AMA_MltvLstsqr ( AMA_OPTIONS options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  wht,
long int *  degree,
long int *  mlamda,
double **  lamda,
double  theta,
AMA_SPLINE **  spline 
)

Least Squares Approximation of Multivariate Data

This function employs cnspla to compute a least squares spline approximation of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$. For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, dependent variable data $z_\ell$ and weights $w_\ell$, for $\ell=1,\ldots,N$, this function computes the spline

\[ s({\bf X}) = \sum_{j_n=1}^{m_n}\cdots\sum_{j_1=1}^{m_1}\alpha_{j_1,\ldots,j_n}B_{d_1,j_1}(x_1|{\bf\Lambda}_1)\cdots B_{d_n,j_n}(x_n|{\bf\Lambda_n}) \]

that minimizes

\[ (1.0-\theta)\sum_\ell^N w_\ell ( s({\bf X}_\ell) - z_\ell )^2 + \theta F^{(2)}\left(s({\bf X})\right) \]

for $0.0\le\theta< 1.0$. The function $F^{(2)}\left(s({\bf X})\right)$ imposes a penalty term on the least squares approximation.

In the above definition of $s({\bf X})$ the $\alpha_{j_1,\ldots,j_n}$, for $j_k=1,\ldots,m_k$ and $1\le k\le n$, are the coefficients of the tensor product B-splines

\[ \Phi_{j_1,\ldots,j_n}({\bf X}) = B_{d_1,j_1}(x_1|{\bf\Lambda}_1)\cdots B_{d_n,j_n}(x_n|{\bf\Lambda}_n) \]

where $B_{d_k,j_k}(x_k|{\bf\Lambda}_k)$ are the $m_k$ univariate B-splines of degree $1 \le d_k\le 5$ defined by the knot vector ${\bf\Lambda}_k\in{\rm R}^{M_k}$. The knot vector ${\bf\Lambda}_k$ is given as

\[ {\bf\Lambda}_k = ( \lambda_{k,1},\cdots,\lambda_{k,d_k+1}, \lambda_{k,d_k+2}, \cdots, \lambda_{k,m_k}, \lambda_{k,m_k+1},\cdots,\lambda_{k,m_k+d_k+1} )^T \]

where the knots $\lambda_{k,i}$ must satisfy the conditions

  • $\lambda_{k,1} = \ldots = \lambda_{k,d_k+1}$;
  • $\lambda_{k,i} \le \lambda_{k,i+1}$, for $1\le i\le m_k+d_k$;
  • $\lambda_{k,i} < \lambda_{i+d_k+1}$, for $1\le i\le m_k$;
  • $\lambda_{k,m_k+1} = \ldots = \lambda_{k,m_k+d_k+1}$.

If the knot vector ${\bf\Lambda}_k$ is represented by its distinct knot vector

\[ {\bf\rm K}^{{\bf\Lambda}_k} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\bf\rm H}^{{\bf\Lambda}_k} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T \]

then the aforementioned conditions are

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = d_k+1$.

Also the first and last distinct knots must satisfy the conditions $\kappa_{k,1}\le \min_\ell\lbrace x_{k,\ell}\rbrace$ and $\kappa_{k,m^d_k}\ge \max_\ell\lbrace x_{k,\ell}\rbrace$, respectively. In this case the number of knots is $M_k = m_k + d_k + 1$.

For convenience this function does not require the definition of the $d_k+1$ order knots at the left and right hand endpoints and also accepts the knot vector ${\bf\Lambda}^o_k\in{\rm R}^{M_k}$ given as

\[ {\bf\Lambda}^o_k = ( \lambda_{k,d_k+1}, \lambda_{k,d_k+2}, \cdots, \lambda_{k,m_k}, \lambda_{k,m_k+1} )^T \]

where the knots $\lambda_{k,i}$ satisfy the conditions

  • $\lambda_{k,i} \le \lambda_{k,i+1}$, for $d_k+2\le i\le m_k-1$;
  • $\lambda_{k,i} < \lambda_{i+d_k+1}$, for $d_k+1\le i\le m_k-d_k$;
  • $\lambda_{k,d_k+1} < \lambda_{k,d_k+2}$ and $\lambda_{k,m_k} < \lambda_{k,m_k+1}$.

If the knot vector ${\bf\Lambda}^o$ is represented by its distinct knot vector

\[ {\bf\rm K}^{{\bf\Lambda}^o} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\bf\rm H}^{{\bf\Lambda}^o} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T \]

then the aforementioned conditions are

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = 1$.

As before, the first and last distinct knots must satisfy the conditions $\kappa_{k,1}\le \min_\ell\lbrace x_{k,\ell}\rbrace$ and $\kappa_{k,m^d_k}\ge \max_\ell\lbrace x_{k,\ell}\rbrace$, respectively. In this case the number of knots is $M_k = m_k - d_k + 1$.

The function $F^{(2)}(s({\bf X}))$ is used to impose a penalty term on the least squares approximation and is defined as

\[ F^{(2)}(s({\bf X})) = \sum_{k=1}^n\left[\int_R\left(\partial^2 s({\bf X})\over\partial^2 x_k\right)^2dR + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right] \]

where $R = [\lambda_{1,1},\lambda_{1,m_1+d_1+1}] \times [\lambda_{2,1},\lambda_{2,m_2+d_2+1}] \times\cdots\times [\lambda_{n,1},\lambda_{n,m_n+d_n+1}]$. The cross partial terms can be excluded from the penalty term $F^{(2)}(s({\bf X}))$ with AMA_OptionsSetPenaltyTerm().

Additionally, the spline $s({\bf X})$ is subject to the bounds

\[ \alpha_l \le s({\bf X}) \le \alpha_u \]

where $\alpha_l=-\alpha_\infty$, $\alpha_u=\alpha_\infty$ and $\alpha_\infty$ equals AMA_SplineInfbnd(). By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().

Depending on the number and distribution of the knots relative to the independent variable data it is possible for the least squares approximation problem to have a non-unique solution. Hence, this function provides the capability of performing a minimum norm optimization subsequent to computing the least squares approximation. The minimum norm optimization computes an alternate spline $\bar s({\bf X})$ by minimizing $F^{(2)}(\bar s({\bf X}))$ subject to the approximation constraints

\[ z_l - \epsilon_\ell \le \bar s({\bf X}_\ell) \le z_l + \epsilon_\ell, \]

for $\ell=1,\ldots,N$, where $\epsilon_l = s({\bf X}_\ell)-z_\ell$ and $s({\bf X}_\ell)$ is the least squares approximation. These constraints insure the condition

\[ \sum_\ell^N w_\ell ( \bar s({\bf X}_\ell) - z_\ell )^2 \le \sum_\ell^N w_\ell ( s({\bf X}_\ell) - z_\ell )^2, \]

that is, the minimum norm optimization does not increase the least squares error.

This function does the following:

  • Checks input parameters for validity, see AMA_MltvInputCheck().
  • Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$ or ${\bf\Lambda}^o_k$, for $1\le k\le n$.
  • Invokes cnspla to compute a least squares spline approximation.
  • Invokes cnspla to perform a minimum norm optimization.
  • Stores approximation into spline.
Note
By default the spline coefficients are initialized to zero but a different initial value for the spline coefficients can be set with AMA_OptionsSetCoefficients().
By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().
By default the cross partial terms are included in the penalty term but they can be excluded with AMA_OptionsSetPenaltyTerm().
By default the minimum norm optimization is not performed but it can be enabled with AMA_OptionsSetMinimumNorm().
The penalty term and minimum norm optimization options are mutually exclusive; that is, the minimum norm optimization is not performed if $\theta > 0.0$.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvLstsqr().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
wht[in] Array of size n containing the weights $w_\ell$, for $\ell=1,\ldots,N$. Must satisfy $w_\ell\ge 0.0$ for all $\ell=1,\ldots,N$. If wht = NULL, then this function uses $w_\ell = 1.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\le$ degree[k] $\le 5$.
mlamda[in] Array of size nind containing the number of knots $M_k$ where mlamda[k] $= M_k$. Must satisfy mlamda[k] $\ge 2$.
lamda[in] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knot vector ${\bf\Lambda}_k$ or ${\bf\Lambda}^o_k$.
theta[in] The penalty term weight $\theta$. Must satisfy $0.0\le$ theta $< 1.0$.
spline[out] Pointer to AMA_SPLINE pointer containing the least squares spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 101715