CurveFit.CurveFitModule

Simple Least Squares fitting

The CurveFit module provides functions that implement a few least squares approximations.

It is simple in an engineering sense. It simply returns the coefficients and does not do any error analysis.

It does, however, provide a simple and common interface to the routines.

The package also includes nonlinear least squares fitting using a Newton type algorithm

The fitting algorithms include

• Straight lines
• Polynomial fitting
• Power laws
• Log laws
• Exp laws
• Rational polynomial fitting
• A generic non-linear fitting algorithm
• King's law (used in hotwire anemometry)
CurveFit.ExpSumFitMethod
(sol::ExpSumFit)(x)

Calculate the sum of exponentials using solution sol at points x.

CurveFit.RationalPolyType

Type defining a rational polynomial

A rational polynomial is the ratio of two polynomials and it is often useful in approximating functions.

CurveFit.apply_fitMethod

#Uses the object created by curve_fit to estimate values

The call method is overloaded so that the fit object can be used as a function:

Example:

x = 1.0:10.0 @. y = 2*x + 1 + randn()

fit = curve_fit(LinearFit, x, y)

y1 = fit(5.1) y2 = apply_fit(fit, 5.1)

CurveFit.cumints!Method
cumints!(init::ExpSumFitInit, x::AbstractVector{T}, y::AbstractVector{T}) where {T <: Real}

Calculates in-place init.n cumulative integrals using coeffients init.coeff.

CurveFit.curve_fitMethod

Generic interface for curve fitting.

The same function curve_fit can be used to fit the data depending on fit type, shich is specified in the first parameter. This function returns an object that can be used to estimate the value of the fitting model using function apply_fit.

A few examples:

• f = curve_fit(LinearFit, x, y)
• f = curve_fit(Polynomial, x, y, n)
CurveFit.expsum_fill_X!Method
expsum_fill_X!(init::ExpSumFitInit, x, λ, X)

Fills matrix init.X

$$$[e^{x\lambda_1} e^{x\lambda_2} \ldots e^{x\lambda_n}]$$$

or if fitting is done with constant

$$$[e^{x\lambda_1} e^{x\lambda_2} \ldots e^{x\lambda_n} 1]$$$
CurveFit.expsum_fill_Y!Method
expsum_fill_Y!(init::ExpSumFitInit, x, y)

Fills matrix Y with

$$$[\int^1 \int^2 \ldots \int^n x^0 x^1 \ldots x^{m-n-1}]$$$

where ∫ⁱ means the ith cumulative integral of y vs x, and m is the number of columns in matrix Y.

CurveFit.expsum_fitMethod
    expsum_fit(x::Union{T,AbstractVector{T}}, y::AbstractVector{T}, n::Int; m::Int = 1, withconst::Bool = true, init::Union{Nothing,ExpSumFitInit} = nothing) where {T <: Real}

Fits the sum of n exponentials and a constant.

$$$y = k + p_1 e^{\lambda_1 t} + p_2 e^{\lambda_2 t} + \ldots + p_n e^{\lambda_n t}$$$

If the keyword withconst is set to false, the constant is not fitted but set k=0.

Uses numerical integration with m strips, where the default m=1 uses linear interpolation. m=2 and higher require uniform interval and usually lead to better accuracy.

Passing init preallocates most needed memory, and can be initialized with expsum_init.

Returns a ExpSumFit struct containing a constant k and vectors p, λ, τ, where τ = -1/λ.

The algorithm is from Matlab code of Juan Gonzales Burgos.

CurveFit.fit_linear_modelMethod
Calculates the coefficients of a linear model using Least Squares

Given a Vandermonde matrix, this function uses the QR factorization to compute the Least Squares fit of a linear model. The code probably could be optimized.

CurveFit.gauss_newton_fitMethod

a = gaussnewtonfit(x, y, fun, ∇fun!, a0[[, eps,] maxiter])

Gauss-Newton nonlinear least squares. Given vectors x and y, the tries to fit parameters a to a function f using least squares approximation:

$y = f(x, a₁, ..., aₙ)$

For more general approximations, see gauss_newton_fit.

Arguments:

• x Vector with x values
• y Vector with y values
• fun a function that is called as fun(x, a) where a is a vector of parameters.
• ∇fun! A function that calculares the derivatives with respect to parameters a
• a0 Vector with the initial guesses of the parameters
• eps Maximum residual for convergence
• maxiter Maximum number of iterations for convergence

Return value

A vector with the convrged array. If no convergence is achieved, the function throws an error.

Specification of the fitting function

The function that should be fitted shoud be specified by Julia funcion with the following signature:

fun(x::T, a::AbstractVector{T}) where {T<:Number}

The derivatives with respect to each fitting parameter a[i] should have the following signature:

∇fun!(x::T, a::AbstractVector{T}, df::AbstractVector{T}) where {T<:Number}

No return value is expected and the derivatives are returned in argument df.

Initial approximation (guess)

If the initial approximation is not good enough, divergence is possible.

Careful with parameters close to 0. The initial guess should never be 0.0 because the initial value of the parameter is used as reference value for computing resiuduals.

Convergence criteria

The argumento maxiter specifies the maximum number of iterations that should be carried out. At each iteration,

$aₖⁿ⁺¹ = aₖⁿ + δₖ$

Convergence is achieved when

$|δᵏ / aₖ⁰| < ε$

Example

x = 1.0:10.0
a = [3.0, 2.0, 1.0]
y = a[1] + a[2]*x + a[3]*x^2
fun(x, a) = a[1] + a[2]*x + a[3]*x^2

function ∇fun!(x, a, df)
df[1] = 1.0
df[2] = x
df[3] = x^2
end

a = gauss_newton_fit(x, y, fun, ∇fun!, [0.5, 0.5, 0.5], 1e-8, 30)
CurveFit.gauss_newton_generic_fitMethod

a = gaussnewtongeneric_fit(x, y, fun, ∇fun!, a0[[, eps,] maxiter])

Gauss-Newton nonlinear least squares. Given matrix x the function tries to fit parameters a to an implicit function f using least squares approximation:

$f(x₁,..., xₘ, a₁, ..., aₙ) = 0$

This function is very similar to gauss_newton_fit but more generic. It doesn't limit the number of independent variables and doesn't require a y LHS.

Arguments:

• x Matrix where each column is a variable
• 
• fun a function that is called as fun(x, a) where a is a vector of parameters.
• ∇fun! A function that calculares the derivatives with respect to parameters a
• a0 Vector with the initial guesses of the parameters
• eps Maximum residual for convergence
• maxiter Maximum number of iterations for convergence

Return value

A vector with the convrged array. If no convergence is achieved, the function throws an error.

Specification of the fitting function

The function that should be fitted shoud be specified by Julia funcion with the following signature:

fun(x::T, a::AbstractVector{T}) where {T<:Number}

The derivatives with respect to each fitting parameter a[i] should have the following signature:

∇fun!(x::T, a::AbstractVector{T}, df::AbstractVector{T}) where {T<:Number}

No return value is expected and the derivatives are returned in argument df.

Initial approximation (guess)

If the initial approximation is not good enough, divergence is possible.

Careful with parameters close to 0. The initial guess should never be 0.0 because the initial value of the parameter is used as reference value for computing resiuduals.

Convergence criteria

The argumento maxiter specifies the maximum number of iterations that should be carried out. At each iteration,

$aₖⁿ⁺¹ = aₖⁿ + δₖ$

Convergence is achieved when

$|δᵏ / aₖ⁰| < ε$

Example

x = 1.0:10.0
a = [3.0, 2.0, 1.0]
y = a[1] + a[2]*x + a[3]*x^2
fun(x, a) = a[1] + a[2]*x + a[3]*x^2

function ∇fun!(x, a, df)
df[1] = 1.0
df[2] = x
df[3] = x^2
end

a = gauss_newton_fit(x, y, fun, ∇fun!, [0.5, 0.5, 0.5], 1e-8, 30)
CurveFit.king_fitFunction

Uses nonlinear least squares to fit the modified King's law:

E^2 = A + B*U^n

The Original (linear) King's law is used to estimate A and B when n=1/2. This initial value is used as an initial guess for fitting the nonlinear modified King's law using the function nonlinear_fit.

CurveFit.linear_king_fitMethod

Original Kings law (1910) represents the relation between voltage and velocity in a hotwire anemometer. The law is given by:

E^2 = A + B*U^0.5

This function estimates A and B.

CurveFit.linear_rational_fitMethod

Linear Rational LeastSquares

The following curvit is done:

y = p(x) / q(x)

where p(x) and q(x) are polynomials.

The linear case is solved by doing a least square fit on

y*q(x) = p(x)

where the zero order term o q(x) is assumed to be 1.

CurveFit.nonlinear_fitFunction

Nonlinear multi-variable least squares

Uses a Newton like algorithm to compute the least squares fit of a model

Parameters:

• x an array where each colum is a variable
• fun The function that should be fitted (more later on)
• a0 An initial guess for each fitting parameter
• eps Convergence criterion
• maxiter Maximum number of iterations

Specifying the model

The model is specified in argument fun that should be callable according to the following signature:

r = fun(x, a)

where both x and a are vectors

The way the algorithm is implemented allows for the function to be implicit. This means that we are trying to fit the relationship fun(x, a) = 0.

Initial guess

The initial guess is important and should be as close as possible to the expeted values. It should, at least, give an order of magnitude of each parameter. Zero values are not recommended since in this case no order of magnitude is available and it is assumed. The initial step is assumed to be a0/10.

Convergence criterion

The initial guess a0 will probably not fit the data very well. The largest value of maxerr = maxabs(fun(x, a0)) (maximum error) is used as a reference and every time a new approximation is computed, if the change in paramaters is compared to this reference error. If it is small, maxabs(da) < eps*maxerr, the algorithm converged.

Return values

A tuple containing:

• The estimated coefficients
• A Bool stating whether the algorithm converged
• Number of iterations it took to converge.
CurveFit.poly_fitMethod

Fits a polynomial of degree n through a set of points.

Simple algorithm, doesn't use orthogonal polynomials or any such thing and therefore unconditioned matrices are possible. Use it only for low degree polynomials.

This function returns a the coefficients of the polynomial.

CurveFit.rational_fitFunction

Carry out a nonlinear least squares of rational polynomials

Find the polynomial coefficients that best approximate the points given by x and y`.