sporco.fista.fista

Base classes for FISTA algorithms

Classes

FISTA(Nx, xshape, dtype[, opt]) Base class for Fast Iterative Shrinkage/Thresholding algorithm (FISTA) algorithms [4].
FISTADFT(xshape, dtype[, opt]) Base class for FISTA algorithms with gradients and updates computed in the frequency domain.

Class Descriptions

class sporco.fista.fista.FISTA(Nx, xshape, dtype, opt=None)[source]

Bases: sporco.common.IterativeSolver

Base class for Fast Iterative Shrinkage/Thresholding algorithm (FISTA) algorithms [4].

Solve optimisation problems of the form

\[\mathrm{argmin}_{\mathbf{x}} \; f(\mathbf{x}) + g(\mathbf{x}) \;\;,\]

where \(f, g\) are convex functions and \(f\) is smooth.

This class is intended to be a base class of other classes that specialise to specific optimisation problems.

After termination of the solve method, attribute itstat is a list of tuples representing statistics of each iteration. The default fields of the named tuple IterationStats are:

Iter : Iteration number

ObjFun : Objective function value

FVal : Value of smooth objective function component \(f\)

GVal : Value of objective function component \(g\)

F_Btrack : Value of objective function \(f + g\) (see Sec. 2.2 of [4])

Q_Btrack : Value of Quadratic approximation \(Q_L\) (see Sec. 2.3 of [4])

IterBtrack : Number of iterations in backtracking

Rsdl : Residual

L : Inverse of gradient step parameter

Time : Cumulative run time

Parameters:
Nx : int

Size of variable \(\mathbf{x}\) in objective function

xshape : tuple of ints

Shape of working variable X

dtype : data-type

Data type for working variables (overridden by ‘DataType’ option)

opt : FISTA.Options object

Algorithm options

class Options(opt=None)[source]

Bases: sporco.cdict.ConstrainedDict

ADMM algorithm options.

Options:

FastSolve : Flag determining whether non-essential computation is skipped. When FastSolve is True and Verbose is False, the functional value and related iteration statistics are not computed. If FastSolve is True residuals are also not calculated, in which case the residual-based stopping method is also disabled, with the number of iterations determined only by MaxMainIter.

Verbose : Flag determining whether iteration status is displayed.

StatusHeader : Flag determining whether status header and separator are displayed.

DataType : Specify data type for solution variables, e.g. np.float32.

X0 : Initial value for X variable.

Callback : Callback function to be called at the end of every iteration.

MaxMainIter : Maximum main iterations.

IterTimer : Label of the timer to use for iteration times.

RelStopTol : Relative convergence tolerance for fixed point residual (see Sec. 4.3 of [21]).

L : Inverse of gradient step parameter \(L\).

AutoStop : Options for adaptive stoping strategy (fixed point residual, see Sec. 4.3 of [21]).

Enabled : Flag determining whether the adaptive stopping relative parameter strategy is enabled.

Tau0 : numerator in adaptive criterion (\(\tau_0\) in [21]).

BackTrack : Options for adaptive L strategy (backtracking, see Sec. 4 of [4]).

Enabled : Flag determining whether adaptive inverse step size parameter strategy is enabled.

Eta : Multiplier applied to L when updated (\(L\) in [4]).

MaxIter : Maximum iterations of updating L when backtracking.

Parameters:
opt : dict or None, optional (default None)

FISTA algorithm options

fwiter = 4

Field width for iteration count display column

fpothr = 2

Field precision for other display columns

itstat_fields_objfn = ('ObjFun', 'FVal', 'GVal')

Fields in IterationStats associated with the objective function; see eval_objfun

itstat_fields_alg = ('Rsdl', 'F_Btrack', 'Q_Btrack', 'IterBTrack', 'L')

Fields in IterationStats associated with the specific solver algorithm

itstat_fields_extra = ()

Non-standard fields in IterationStats; see itstat_extra

hdrtxt_objfn = ('Fnc', 'f', 'g')

Display column headers associated with the objective function; see eval_objfun

hdrval_objfun = {'Fnc': 'ObjFun', 'f': 'FVal', 'g': 'GVal'}

Dictionary mapping display column headers in hdrtxt_objfn to IterationStats entries

xinit(xshape)[source]

Return initialiser for working variable X.

solve()[source]

Start (or re-start) optimisation. This method implements the framework for the iterations of a FISTA algorithm. There is sufficient flexibility in overriding the component methods that it calls that it is usually not necessary to override this method in derived clases.

If option Verbose is True, the progress of the optimisation is displayed at every iteration. At termination of this method, attribute itstat is a list of tuples representing statistics of each iteration, unless option FastSolve is True and option Verbose is False.

Attribute timer is an instance of util.Timer that provides the following labelled timers:

init: Time taken for object initialisation by __init__

solve: Total time taken by call(s) to solve

solve_wo_func: Total time taken by call(s) to solve, excluding time taken to compute functional value and related iteration statistics

solve_wo_rsdl : Total time taken by call(s) to solve, excluding time taken to compute functional value and related iteration statistics as well as time take to compute residuals

solve_wo_btrack : Total time taken by call(s) to solve, excluding time taken to compute functional value and related iteration statistics as well as time take to compute residuals and implemented BackTrack mechanism

getmin()[source]

Get minimiser after optimisation.

proximal_step(grad=None)[source]

Compute proximal update (gradient descent + regularization).

combination_step()[source]

Build next update by a smart combination of previous updates.

compute_backtracking()[source]

Estimate step size L by computing a linesearch that guarantees that F <= Q.

eval_grad()[source]

Compute gradient.

Overriding this method is required.

eval_proxop(V)[source]

Compute proximal operator of \(g\).

Overriding this method is required.

eval_R(V)[source]

Evaluate smooth function \(f\) in V.

Overriding this method is required.

eval_Rx()[source]

Evaluate smooth function \(f\) in X.

store_prev()[source]

Store previous X state.

eval_Dxy()[source]

Evaluate difference of state and auxiliary state updates.

compute_residuals()[source]

Compute residuals and stopping thresholds.

classmethod hdrtxt()[source]

Construct tuple of status display column title.

classmethod hdrval()[source]

Construct dictionary mapping display column title to IterationStats entries.

iteration_stats(k, frcxd)[source]

Construct iteration stats record tuple.

eval_objfn()[source]

Compute components of objective function as well as total contribution to objective function.

itstat_extra()[source]

Non-standard entries for the iteration stats record tuple.

getitstat()[source]

Get iteration stats as named tuple of arrays instead of array of named tuples.

display_start()[source]

Set up status display if option selected. NB: this method assumes that the first entry is the iteration count and the last is the L value.

display_status(fmtstr, itst)[source]

Display current iteration status as selection of fields from iteration stats tuple.

display_end(nsep)[source]

Terminate status display if option selected.

var_x()[source]

Get \(\mathbf{x}\) variable.

obfn_f(X)[source]

Compute \(f(\mathbf{x})\) component of FISTA objective function.

Overriding this method is required if eval_objfun is not overridden.

obfn_g(X)[source]

Compute \(g(\mathbf{x})\) component of FISTA objective function.

Overriding this method is required if eval_objfun is not overridden.

rsdl()[source]

Compute fixed point residual.

Overriding this method is required.

Lchange()[source]

Action to be taken, if any, when L parameter is changed.

Overriding this method is optional.

class sporco.fista.fista.FISTADFT(xshape, dtype, opt=None)[source]

Bases: sporco.fista.fista.FISTA

Base class for FISTA algorithms with gradients and updates computed in the frequency domain.


Inheritance diagram of FISTADFT


Solve optimisation problems of the form

\[\mathrm{argmin}_{\mathbf{x}} \; f(\mathbf{x}) + g(\mathbf{x}) \;\;,\]

where \(f, g\) are convex functions and \(f\) is smooth.

This class specialises class FISTA, but remains a base class for other classes that specialise to specific optimisation problems.

Parameters:
xshape : tuple of ints

Shape of working variable X (the primary variable)

dtype : data-type

Data type for working variables

opt : FISTADFT.Options object

Algorithm options

class Options(opt=None)[source]

Bases: sporco.fista.fista.Options

FISTADFT algorithm options.

Options include all of those defined in FISTA.Options.

Parameters:
opt : dict or None, optional (default None)

FISTADFT algorithm options

proximal_step(gradf=None)[source]

Compute proximal update (gradient descent + constraint).

combination_step()[source]

Update auxiliary state by a smart combination of previous updates in the frequency domain.

compute_backtracking()[source]

Estimate step size L by computing a linesearch that guarantees that F <= Q.

store_prev()[source]

Store previous X in frequency domain.

eval_Dxyf()[source]

Evaluate difference of state and auxiliary state.

eval_Rxf()[source]

Evaluate smooth term in X.

eval_gradf()[source]

Compute gradient in Fourier domain. Also computes argument of smooth function \(f\).

Overriding this method is required.

eval_proxop(V)[source]

Compute proximal operator of \(g\).

Overriding this method is required.