userFcnLeastSquares¶

Computes a least-squares fit using user-supplied functions.

Synopsis¶

userFcnLeastSquares (fcn, nbasis, xdata, ydata)

Required Arguments¶

float fcn (k, x) (Input): User-supplied function that defines the subspace from which the least-squares fit is to be performed. The k-th basis function evaluated at x is f(k, x) where k = 1, 2, …, nbasis.
int ndata (Input): Number of data points.
float xdata[] (Input): Array with ndata components containing the abscissas of the least-squares problem.
float ydata[] (Input): Array with ndata components containing the ordinates of the least-squares problem.

Return Value¶

The vector containing the coefficients of the basis functions. If a fit cannot be computed, then None is returned.

Optional Arguments¶

intercept (Output)

This option adds an intercept to the model. Thus, the least-squares fit is computed using the user-supplied basis functions augmented by the constant function. The coefficient of the constant function is stored in intercept.

sse (Output)

This option returns the error sum of squares.

weights, float weights[] (Input)

This option requires the user to provide the weights.

Default: all weights equal one

Description¶

The function userFcnLeastSquares computes a best least-squares approximation to given univariate data of the form

\[\left\{\left(x_i,f_i\right)\right\}_{i=0}^{n-1}\]

by M basis functions

\[\left\{F_j\right\}_{j=1}^{M}\]

(where M = nbasis). In particular, the default for this function returns the coefficients a which minimize

\[\sum_{i=0}^{n-1} w_i \left[ f_i - \sum_{j=1}^{M} a_{j-1}F_j\left(x_i\right) \right]^2\]

where w = weights, n = ndata, x = xdata, and f = ydata.

If the optional argument intercept is chosen, then an intercept is placed in the model, and the coefficients a, returned by userFcnLeastSquares, minimize the error sum of squares as indicated below.

\[\sum_{i=0}^{n-1} w_i \left[ f_i - \mathit{intercept} - \sum_{j=1}^{M} a_{j-1}F_j\left(x_i\right) \right]^2\]

Remarks¶

Note that although the system is linear, for very large problems the Chapter 8 function, nonlinLeastSquares, might be better suited. This is because userFcnLeastSquares will gather the matrix entries one at a time by calls to the user‑supplied function. By using the nonlinear solver and supplying the Jacobian the user can sidestep this issue and likely achieve accurate results since the nonlinear solver utilizes regularization and iterative refinement.

Examples¶

Example 1¶

This example fits the following two functions (indexed by δ):

\[1 + sinx + 7 sin3x + δε\]

where ɛ is a random uniform deviate over the range [‑1, 1] and δ is 0 for the first function and 1 for the second. These functions are evaluated at 90 equally spaced points on the interval [0, 6]. Four basis functions are used: 1, sinx, sin2x, sin3x.

from __future__ import print_function
from numpy import *
from pyimsl.math.randomSeedSet import randomSeedSet
from pyimsl.math.randomUniform import randomUniform
from pyimsl.math.userFcnLeastSquares import userFcnLeastSquares
from pyimsl.math.writeMatrix import writeMatrix

# Define functions


def F(x):
    return (1. + sin(x) + 7. * sin(3.0 * x))


def fcn(n, x):
    if (n == 1):
        return 1.0
    else:
        return sin((n - 1.0) * x)


# Generate random numbers
ndata = 90
nbasis = 4
randomSeedSet(1234567)
random = randomUniform(ndata)

# Set up data
xdata = empty(ndata)
ydata = empty(ndata)
for delta in range(0, 2):
    for i in range(0, ndata):
        xdata[i] = 6. * i / (ndata - 1)
        ydata[i] = F(xdata[i]) + delta * 2. * (random[i] - .5)
    coef = userFcnLeastSquares(fcn, nbasis, xdata, ydata)
    print("\nFor delta = ", delta)
    writeMatrix("the computed coefficients are", coef)

Output¶

For delta =  0

For delta =  1
 
           the computed coefficients are
          1            2            3            4
          1            1            0            7
 
           the computed coefficients are
          1            2            3            4
      0.979        0.998        0.096        6.839

Example 2¶

Recall that the first example fitted the following two functions (indexed by δ):

1 + sinx + 7 sin3x + δε

where ε is a random uniform deviate over the range[−1, 1], and δ is 0 for the first function and 1 for the second. These functions are evaluated at 90 equally spaced points on the interval [0, 6]. Previously, the four basis functions were used: 1, sinx, sin2x, sin3x. This example uses the four basis functions: sinx, sin2x, sin3x, sin4x, combined with the intercept option.

from __future__ import print_function
from numpy import *
from pyimsl.math.randomSeedSet import randomSeedSet
from pyimsl.math.randomUniform import randomUniform
from pyimsl.math.userFcnLeastSquares import userFcnLeastSquares
from pyimsl.math.writeMatrix import writeMatrix

# Define functions


def F(x):
    return (1. + sin(x) + 7. * sin(3.0 * x))


def fcn(n, x):
    return sin(n * x)


# Generate random numbers
ndata = 90
nbasis = 4
randomSeedSet(1234567)
random = randomUniform(ndata)

# Set up data
xdata = empty(ndata)
ydata = empty(ndata)
for delta in range(0, 2):
    for i in range(0, ndata):
        xdata[i] = 6. * i / (ndata - 1)
        ydata[i] = F(xdata[i]) + delta * 2. * (random[i] - .5)
    intercept = []
    coef = userFcnLeastSquares(fcn, nbasis, xdata, ydata,
                               intercept=intercept)
    print("\nFor delta = ", delta)
    print("The predicted intercept value is %10.3f"
          % (intercept[0]))
    writeMatrix("the computed coefficients are", coef)

Output¶

For delta =  0
The predicted intercept value is      1.000

For delta =  1
The predicted intercept value is      0.978
 
           the computed coefficients are
          1            2            3            4
          1           -0            7            0
 
           the computed coefficients are
          1            2            3            4
      0.998        0.097        6.841        0.075

Warning Errors¶

`IMSL_LINEAR_DEPENDENCE`	Linear dependence of the basis functions exists. One or more components of coef are set to zero.
`IMSL_LINEAR_DEPENDENCE_CONST`	Linear dependence of the constant function and basis functions exists. One or more components of coef are set to zero.

Fatal Errors¶

IMSL_NEGATIVE_WEIGHTS_2 All weights must be greater than or equal to zero.

IMSL_STOP_USER_FCN

Request from user supplied function to stop algorithm.

User flag = “#”.