capfit 2.6.1

CapFit: linearly-constrained non-linear least-squares optimization

CapFit: linearly-constrained non-linear least-squares optimization

This CapFit package contains a Python implementation of the linearly-constrained non-linear least-squares optimization method described in Section 3.2 of Cappellari (2023).

Contents

• Attribution

• Installation

• Usage Examples

• CapFit Purpose

• Calling Sequence

• Input Parameters

• Optional Keywords

• Returns

• License

• Changelog

Attribution

If you use this software for your research, please cite Cappellari (2023), where the algorithm was described in Sec. 3.2. The BibTeX entry for the paper is:

@ARTICLE{Cappellari2023,
    author = {{Cappellari}, M.},
    title = "{Full spectrum fitting with photometry in PPXF: stellar population
        versus dynamical masses, non-parametric star formation history and
        metallicity for 3200 LEGA-C galaxies at redshift $z\approx0.8$}",
    journal = {MNRAS},
    eprint = {2208.14974},
    year = 2023,
    volume = 526,
    pages = {3273-3300},
    doi = {10.1093/mnras/stad2597}
}

Installation

install with:

pip install capfit

Without write access to the global site-packages directory, use:

pip install --user capfit

To upgrade CapFit to the latest version use:

pip install --upgrade capfit

Usage Examples

To learn how to use the CapFit package, see the capfit_examples.py file within the capfit/examples directory. It can be found within the main capfit package installation folder inside site-packages. The detailed documentation is contained in the docstring of the file capfit.py, or on PyPi.

---

CapFit Purpose

CapFit solves linearly-constrained least-squares optimization problems.

Linear inequality/equality constraints and bound constraints are supported. Moreover, one can easily tie or fix some parameters without having to modify the fitting function.

CapFit combines two successful ideas:

1. The Sequential Quadratic Programming (SQP) approach, specialized for the case of linear constraints

2. The Levenberg-Marquardt (LM) method

It was designed for the common situations where the user function is not a simple analytic function but is the result of some complex calculations and is more expensive to compute than the small quadratic subproblem.

CapFit performance is generally better than the best unconstrained or bound-constrained least-squares algorithms currently available, in terms of robustness and number of function evaluations. However, in addition CapFit allows for more general constraints.

Given a function of n model parameters x_k returning the m model residuals f_j(x), CapFit finds a local minimum of the cost function:

G(x) = sum[f_j(x)^2]

subject to:

A_ineq @ x <= b_ineq            # Linear Inequality Constraints
A_eq @ x == b_eq                # Linear Equality Constraints
bounds[0] <= x <= bounds[1]     # Bounds
x_k = f(x)                      # Tied Parameters
x_k = a_k                       # Fixed Parameters

Calling Sequence

from capfit.capfit import capfit

res = capfit(func, p0, A_eq=None, A_ineq=None, abs_step=None,
             b_eq=None, b_ineq=None, bounds=(-np.inf, np.inf),
             diff_step=1e-4, fixed=None, ftol=1e-4,
             linear_method='lsq_lin', max_nfev=None, rcond=None,
             tied=None, verbose=0, x_scale='jac', xtol=1e-4, args=(),
             kwargs={})

print(f"solution: {res.x}")

Input Parameters

funcallable

Function which computes the vector of residuals, with the signature fun(x, *args, **kwargs), i.e., the minimization proceeds with respect to its first argument. The argument x passed to this function is an 1-d darray of shape (n,). The function must return a 1-d array of shape (m,).

x0array_like with shape (n,) or float

Initial guess on independent variables. For guaranteed convergence, the initial guess must be feasible (satsfies the constraints) and for this reason an error is returned if it is not the case.

Optional Keywords

abs_stepNone, scalar or array_like, optional

Determines the absolute step size for the finite difference approximation of the Jacobian. If abs_step is given, then the value of diff_step is ignored and the finite differences use this fixed absolute step.

A_eqarray_like with shape (q, n), optional

Defines the linear equality constraints on the fitted parameters:

A_eq @ x == b_eq

The same result can be achieved using the tied keyword, which also allows for non-linear equality constraints.

A_ineqarray_like with shape (p, n), optional

Defines the linear inequality constraints on the fitted parameters:

A_ineq @ x <= b_ineq

b_ineqarray_like with shape (p), optional

See description of A_ineq.

b_eqarray_like with shape (q), optional

See description of A_eq.

bounds2-tuple of array_like, optional

Lower and upper bounds (lb, ub) on independent variables. Defaults to no bounds. Each array must match the size of x0 or be a scalar, in the latter case a bound will be the same for all variables. Use np.inf with an appropriate sign to disable bounds on all or some variables.

diff_stepNone, scalar or array_like, optional

Determines the relative step size for the finite difference approximation of the Jacobian. The actual step is computed as diff_step*maximum(1, abs(x)) (default diff_step=1e-4)

fixed :

Boolean vector set to True where a given parameter has to be held fixed with the value given in x0.

ftolfloat or None, optional

Tolerance for termination by the change of the cost function (default is 1e-4). The optimization process is stopped when both prered < ftol and abs(actred) < ftol and additionally actred <= 2*prered, where actred and prered are the actual and predicted relative changes respectively (as described in More’ et al. 1980).

linear_methodstr, optional

Method used for the solution of the linear least-squares sub-problem. When using linear constraints, this can be either ‘lsq_lin’ or ‘cvxopt’ (default). In the latter case, the cvxopt package must be installed. With only bounds, this keyword is ignored as capfit uses ‘lsq_box’.

max_nfevNone or int, optional

Maximum number of function evaluations before the termination (default is 100*n).

rcondfloat, optional

Cutoff for small singular values used to determine the effective rank of the Jacobian. Singular values smaller than rcond*largest_singular_value are considered zero.

tiedarray_like with shape (n,), optional

A list of string expressions. Each expression “ties” the parameter to other free or fixed parameters. Any expression involving constants and the parameter array p[j] are permitted. Since they are totally constrained, tied parameters are considered to be fixed; no errors are computed for them.

This is a vector with the same dimensions as x0. In practice, for every element of x0 one needs to specify either an empty string '' implying that the parameter is free, or a string expression involving some of the variables p[j], where j represents the index of the vector of parameters. See usage example.

verbose{0, 1, 2}, optional

Level of algorithm’s verbosity:

• 0 (default) : work silently.

• 1 : display a termination report.

• 2 : display progress during iterations.

x_scalearray_like or ‘jac’, optional

Characteristic scale of each variable. Setting x_scale is equivalent to reformulating the problem in scaled variables xs = x/x_scale. An alternative view is that the initial size/2 of the box trust region along j-th dimension is given by x_scale[j] and the box ratains its shape during the optimization. Improved convergence is achieved by setting x_scale such that a step of a given size along any of the scaled variables has a similar effect on the cost function. If set to ‘jac’, the scale is iteratively updated using the inverse norms of the columns of the Jacobian matrix (as described in More’ 1978).

xtolfloat or None, optional

Tolerance for termination by the change dx of the independent variables (default is 1e-4). The condition is norm(dx) < xtol*(xtol + norm(xs)) where xs is the value of x scaled according to x_scale parameter (see below). If None, the termination by this condition is disabled.

args, kwargstuple and dict, optional

Additional arguments passed to fun, empty by default. The calling signature is fun(x, *args, **kwargs).

Returns

The following attributes of the capfit class are defined in output:

.xndarray, shape (n,)

Solution found.

.x_errndarray, shape (n,)

Formal uncertainties of the solution estimated from diagonal elements of the covariance matrix (i.e. ignoring covariance). This is only meaningful if the fun values represent residuals divided by their 1sigma errors. Uncertainties of the non-free parameters are returned as zero.

.costfloat

Value of the cost = 0.5*(fun @ fun) function at the solution.

.covndarray, shape (nfree, nfree)

Covariance matrix for the free (not tied or fixed) parameters.

.funndarray, shape (m,)

Vector of residuals at the solution.

.jacndarray, shape (m, nfree)

Modified Jacobian matrix at the solution, in the sense that J.T @ J is a Gauss-Newton approximation of the Hessian of the cost function. NB: This is the Jacobian of the free (not tied or fixed) parameters.

.gradndarray, shape (m,)

Gradient of the cost function at the solution.

.nfevint

Number of function evaluations done.

.njevint

Number of Jacobian evaluations done.

.statusint

The reason for algorithm termination:

• -1 : improper input parameters status

• 0 : the maximum number of function evaluations is exceeded.

• 2 : ftol termination condition is satisfied.

• 3 : xtol termination condition is satisfied.

• 4 : Both ftol and xtol termination conditions are satisfied.

.messagestr

Verbal description of the termination reason.

.successbool

True if one of the convergence criteria is satisfied (status > 0).

License

Other/Proprietary License

Copyright (c) 2017-2024 Michele Cappellari

This software is provided as is with no warranty. You may use it for non-commercial purposes and modify it for personal or internal use, as long as you include this copyright and disclaimer in all copies. You may not redistribute the code.

---

Changelog

V2.6.1: MC, Oxford, 19 July 2024

• Removed capfit from the ppxf package and made it a separated package.

V2.5.1: MC, Oxford, 22 May 2023

• capfit: Relaxed tolerance when checking initial guess feasibility.

V2.5.0: MC, Oxford, 16 August 2022

• Uses scipy.optimize.linprog to find feasible starting point in lsq_lin.

• Set default linear_method='lsq_lin' in capfit. This eliminates the need to install cvxopt when using general linear constraints.

V2.4.0: MC, Oxford, 04 March 2022

• Remove the non-free variables before the optimization. This reduces the degeneracy of the Jacobian.

V2.3.0: MC, Oxford, 20 December 2020

• New linear_method keyword to select cvxopt or lsq_lin, for cases where the latter stops, when using linear constraints. Thanks to Kyle Westfall (UCO Lick) for a detailed bug report.

V2.2.1: MC, Oxford, 11 September 2020

• Fixed possible infinite loop in lsq_box and lsq_lin. Thanks to Shravan Shetty (pku.edu.cn) for the detailed report.

• Use NumPy rather than SciPy version of linalg.lstsq to avoid a Scipy bug in the default criterion for rank deficiency.

• Pass rcond keyword to cov_err for consistency.

V2.2.0: MC, Oxford, 20 August 2020

• New function lsq_lin implementing a robust linear least-squares linearly-constrained algorithm which works with a rank-deficient matrix and allows for a starting guess. lsq_lin supersedes the former lsqlin.

• Renamed lsqbox to lsq_box and revised its interface.

V2.1.0: MC, Oxford, 09 July 2020

• New function lsqbox implementing a fast linear least-squares box-constrained (bounds) algorithm which allows for a starting guess.

V2.0.2: MC, Oxford, 20 June 2020

• capfit: new keyword cond (passed to lsqlin).

• capfit: Use bvls to solve quadratic subproblem with only bounds.

V2.0.1: MC, Oxford, 24 January 2020

• New keyword cond for lsqlin.

• Relaxed assertion for inconsistent inequalities in lsqlin to avoid false positives. Thanks to Kyle Westfall (UCO Lick) for a detailed bug report.

V2.0.0: MC, Oxford, 10 January 2020

• Use the new general linear least-squares optimization function lsqlin to solve the quadratic sub-problem.

• Allow for linear inequality/equality constraints A_ineq, b_ineq and A_eq, b_eq

V1.0.7: MC, Oxford, 10 October 2019

• Included complete documentation.

• Improved print formatting.

• Return .message attribute.

• Improved xtol convergence test.

• Only accept final move if chi2 decreased.

• Strictly satisfy bounds during Jacobian computation.

V1.0.6: MC, Oxford, 11 June 2019

• Use only free parameters for small-step convergence tests.

• Explain in words convergence status when verbose != 0

• Fixed program stops when abs_step is undefined.

• Fixed capfit ignoring optional max_nfev.

V1.0.5: MC, Oxford, 28 March 2019

• Raise an error if the user function returns non-finite values.

V1.0.4: MC, Oxford, 30 November 2018

• Allow for a scalar abs_step.

V1.0.3: MC, Oxford, 20 September 2018

• Raise an error if one tries to tie parameters to themselves. Thanks to Kyle Westfall (Univ. Santa Cruz) for the feedback.

• Use Python 3.6 f-strings.

V1.0.2: MC, Oxford, 10 May 2018

• Dropped legacy Python 2.7 support.

V1.0.1: MC, Oxford, 13 February 2018

• Make output errors of non-free variables exactly zero.

V1.0.0: MC, Oxford, 15 June 2017

• Written by Michele Cappellari