Hessian Update Strategies
This module provides various generic Hessian approximation strategies that can be employed when the calculating the exact Hessian or an approximation is computationally too demanding.
Classes Summary
|
Broydens "bad" method as introduced in [Broyden 1965](https://doi.org/10.1090%2FS0025-5718-1965-0198670-6). |
|
Broyden-Fletcher-Goldfarb-Shanno update strategy. |
|
Broydens "good" method as introduced in [Broyden 1965](https://doi.org/10.1090%2FS0025-5718-1965-0198670-6). |
|
BroydenClass Update scheme as described in [Nocedal & Wright]( http://dx.doi.org/10.1007/b98874) Chapter 6.3. |
|
Davidon-Fletcher-Powell update strategy. |
|
|
|
|
|
Abstract class from which Hessian update strategies should subclass |
|
|
|
|
|
|
|
|
|
Iterative update schemes that only use s and y values for update. |
|
Symmetric Rank 1 update strategy as described in [Nocedal & Wright](http://dx.doi.org/10.1007/b98874) Chapter 6.2. |
|
Structured Secant Method as introduced by [Dennis et al 1989](https://doi.org/10.1007/BF00962795), which is compatible with BFGS, DFP update schemes. |
|
|
|
Totally Structured Secant Method as introduced by [Huschens 1994](https://doi.org/10.1137/0804005), which uses a self-adjusting update method for the second order term. |
Functions Summary
|
Scale free implementation of the broyden class update scheme. |
Classes
- class fides.hessian_approximation.BB(init_with_hess=False)[source]
Broydens “bad” method as introduced in [Broyden 1965](https://doi.org/10.1090%2FS0025-5718-1965-0198670-6). This is a rank 1 update strategy that does not preserve symmetry or positive definiteness.
This scheme only works with a function that returns (fval, grad)
- update(s, y)[source]
Update the Hessian approximation
- Parameters
s (
numpy.ndarray
) – step in optimization variablesy (
numpy.ndarray
) – step in gradient
- Return type
- class fides.hessian_approximation.BFGS(init_with_hess=False, enforce_curv_cond=True)[source]
Broyden-Fletcher-Goldfarb-Shanno update strategy. This is a rank 2 update strategy that preserves symmetry and positive-semidefiniteness.
This scheme only works with a function that returns (fval, grad)
- __init__(init_with_hess=False, enforce_curv_cond=True)[source]
Create a Hessian update strategy instance
- Parameters
init_with_hess (
typing.Optional
[bool
]) – Whether the hybrid update strategy should be initialized according to the user-provided objective function
- class fides.hessian_approximation.BG(init_with_hess=False)[source]
Broydens “good” method as introduced in [Broyden 1965](https://doi.org/10.1090%2FS0025-5718-1965-0198670-6). This is a rank 1 update strategy that does not preserve symmetry or positive definiteness.
This scheme only works with a function that returns (fval, grad)
- class fides.hessian_approximation.Broyden(phi, init_with_hess=False, enforce_curv_cond=True)[source]
BroydenClass Update scheme as described in [Nocedal & Wright]( http://dx.doi.org/10.1007/b98874) Chapter 6.3. This is a generalization of BFGS/DFP methods where the parameter \(phi\) controls the convex combination between the two. This is a rank 2 update strategy that preserves positive-semidefiniteness and symmetry (if \(\phi \in [0,1]\)).
This scheme only works with a function that returns (fval, grad)
- Parameters
phi (
float
) – convex combination parameter interpolating between BFGS (phi==0) and DFP (phi==1).enforce_curv_cond (
typing.Optional
[bool
]) – boolean that controls whether the employed broyden class update should attempt to preserve positive definiteness. If set to True, updates from steps that violate the curvature condition will be discarded.
- __init__(phi, init_with_hess=False, enforce_curv_cond=True)[source]
Create a Hessian update strategy instance
- Parameters
init_with_hess (
typing.Optional
[bool
]) – Whether the hybrid update strategy should be initialized according to the user-provided objective function
- class fides.hessian_approximation.DFP(init_with_hess=False, enforce_curv_cond=True)[source]
Davidon-Fletcher-Powell update strategy. This is a rank 2 update strategy that preserves symmetry and positive-semidefiniteness.
This scheme only works with a function that returns (fval, grad)
- __init__(init_with_hess=False, enforce_curv_cond=True)[source]
Create a Hessian update strategy instance
- Parameters
init_with_hess (
typing.Optional
[bool
]) – Whether the hybrid update strategy should be initialized according to the user-provided objective function
- class fides.hessian_approximation.FX(happ=<fides.hessian_approximation.BFGS object>, hybrid_tol=0.2)[source]
- __init__(happ=<fides.hessian_approximation.BFGS object>, hybrid_tol=0.2)[source]
Hybrid method HY2 as introduced by [Fletcher & Xu 1986](https://doi.org/10.1093/imanum/7.3.371). This approximation scheme employs a dynamic approximation as long as function values satisfy \(\frac{f_k - f_{k+1}}{f_k} < \epsilon\) and employs the iterative scheme applied to the last dynamic approximation if not.
This scheme only works with a function that returns (fval, grad, hess)
- Parameters
hybrid_tol (
typing.Optional
[float
]) – switch tolerance \(\epsilon\)
- update(delta, gamma, r, rprev, hess)[source]
Update the Hessian approximation
- Parameters
delta (
numpy.ndarray
) – step in optimization variablesgamma (
numpy.ndarray
) – step in gradientr (
numpy.ndarray
) – residuals after current steprprev (
numpy.ndarray
) – residuals befor current stephess (
numpy.ndarray
) – user-provided (Gauss-Newton) Hessian approximation
- Return type
- class fides.hessian_approximation.GNSBFGS(hybrid_tol=1e-06, enforce_curv_cond=True)[source]
- __init__(hybrid_tol=1e-06, enforce_curv_cond=True)[source]
Hybrid Gauss-Newton Structured BFGS method as introduced by [Zhou & Chen 2010](https://doi.org/10.1137/090748470), which combines ideas of hybrid switching methods and structured secant methods.
This scheme only works with a function that returns (res, sres)
- Parameters
hybrid_tol (
float
) – switching tolerance that controls switching between update methods
- update(s, y, r, hess, yb)[source]
Update the structured approximation
- Parameters
s (
numpy.ndarray
) – step in optimization parametersy (
numpy.ndarray
) – step in gradient parametersr (
numpy.ndarray
) – residual vectorhess (
numpy.ndarray
) – user-provided (Gauss-Newton) Hessian approximationyb (
numpy.ndarray
) – approximation to A*s, where A is structured approximation matrix
- Return type
- class fides.hessian_approximation.HessianApproximation(init_with_hess=False)[source]
Abstract class from which Hessian update strategies should subclass
- __init__(init_with_hess=False)[source]
Create a Hessian update strategy instance
- Parameters
init_with_hess (
typing.Optional
[bool
]) – Whether the hybrid update strategy should be initialized according to the user-provided objective function
- get_diff()[source]
Getter for the Hessian approximation update
- Return type
- Returns
Hessian approximation update
- init_mat(dim, hess=None)[source]
Initializes this approximation instance and checks the dimensionality
- Parameters
dim (
int
) – dimension of optimization variableshess (
typing.Optional
[numpy.ndarray
]) – user provided initialization
- Return type
- set_mat(mat)[source]
Setter for the Hessian approximation
- Parameters
mat (
numpy.ndarray
) – Hessian approximation- Return type
- class fides.hessian_approximation.HybridApproximation(happ=<fides.hessian_approximation.BFGS object>)[source]
- __init__(happ=<fides.hessian_approximation.BFGS object>)[source]
Create a Hybrid Hessian update strategy that switches between an iterative approximation and a dynamic approximation
- Parameters
happ (
fides.hessian_approximation.IterativeHessianApproximation
) – Iterative Hessian Approximation
- init_mat(dim, hess=None)[source]
Initializes this approximation instance and checks the dimensionality
- Parameters
dim (
int
) – dimension of optimization variableshess (
typing.Optional
[numpy.ndarray
]) – user provided initialization
- class fides.hessian_approximation.HybridFixed(happ=<fides.hessian_approximation.BFGS object>, switch_iteration=20)[source]
- __init__(happ=<fides.hessian_approximation.BFGS object>, switch_iteration=20)[source]
Switch from a dynamic approximation to the user provided iterative scheme after a fixed number of successive iterations without trust-region update. The switching is non-reversible. The iterative scheme is initialized and updated rom the beginning, but only employed after the specified number of iterations.
This scheme only works with a function that returns (fval, grad, hess)
- Parameters
switch_iteration (
typing.Optional
[int
]) – Number of iterations without trust region update after which switch occurs.
- class fides.hessian_approximation.HybridFraction(happ=<fides.hessian_approximation.BFGS object>, switch_threshold=0.8)[source]
- __init__(happ=<fides.hessian_approximation.BFGS object>, switch_threshold=0.8)[source]
Switch from a dynamic approximation to the user provided iterative scheme as soon as the fraction of iterations where the step is accepted but the trust region is not update exceeds the user provided threshold.Threshold check is only performed after 25 iterations. The switching is non-reversible. The iterative scheme is initialized and updated rom the beginning, but only employed after the switching.
This scheme only works with a function that returns (fval, grad, hess)
- Parameters
switch_threshold (
typing.Optional
[float
]) – Threshold for fraction of iterations where step is accepted but trust region is not updated, which when exceeded triggers switch of approximation.
- class fides.hessian_approximation.HybridSwitchApproximation(happ=<fides.hessian_approximation.BFGS object>)[source]
- class fides.hessian_approximation.IterativeHessianApproximation(init_with_hess=False)[source]
Iterative update schemes that only use s and y values for update.
- update(s, y)[source]
Update the Hessian approximation
- Parameters
s (
numpy.ndarray
) – step in optimization variablesy (
numpy.ndarray
) – step in gradient
- Return type
- class fides.hessian_approximation.SR1(init_with_hess=False)[source]
Symmetric Rank 1 update strategy as described in [Nocedal & Wright](http://dx.doi.org/10.1007/b98874) Chapter 6.2. This is a rank 1 update strategy that preserves symmetry but does not preserve positive-semidefiniteness.
This scheme only works with a function that returns (fval, grad)
- class fides.hessian_approximation.SSM(phi=0.0, enforce_curv_cond=True)[source]
Structured Secant Method as introduced by [Dennis et al 1989](https://doi.org/10.1007/BF00962795), which is compatible with BFGS, DFP update schemes.
This scheme only works with a function that returns (res, sres)
- update(s, y, r, hess, yb)[source]
Update the structured approximation
- Parameters
s (
numpy.ndarray
) – step in optimization parametersy (
numpy.ndarray
) – step in gradient parametersr (
numpy.ndarray
) – residual vectorhess (
numpy.ndarray
) – user-provided (Gauss-Newton) Hessian approximationyb (
numpy.ndarray
) – approximation to A*s, where A is structured approximation matrix
- Return type
- class fides.hessian_approximation.StructuredApproximation(phi=0.0, enforce_curv_cond=True)[source]
- __init__(phi=0.0, enforce_curv_cond=True)[source]
This is the base class for structured secant methods (SSM). SSMs approximate the hessian by combining the Gauss-Newton component C(x) and an iteratively updated component that approximates the difference S to the true Hessian.
- Parameters
phi (
typing.Optional
[float
]) – convex combination parameter interpolating between BFGS (phi==0) and DFP (phi==1) update schemes.enforce_curv_cond (
typing.Optional
[bool
]) – boolean that controls whether the employed broyden class update should attempt to preserve positive definiteness. If set to True, updates from steps that violate the curvature condition will be discarded.
- init_mat(dim, hess=None)[source]
Initializes this approximation instance and checks the dimensionality
- Parameters
dim (
int
) – dimension of optimization variableshess (
typing.Optional
[numpy.ndarray
]) – user provided initialization
- update(s, y, r, hess, yb)[source]
Update the structured approximation
- Parameters
s (
numpy.ndarray
) – step in optimization parametersy (
numpy.ndarray
) – step in gradient parametersr (
numpy.ndarray
) – residual vectorhess (
numpy.ndarray
) – user-provided (Gauss-Newton) Hessian approximationyb (
numpy.ndarray
) – approximation to A*s, where A is structured approximation matrix
- Return type
- class fides.hessian_approximation.TSSM(phi=0.0, enforce_curv_cond=True)[source]
Totally Structured Secant Method as introduced by [Huschens 1994](https://doi.org/10.1137/0804005), which uses a self-adjusting update method for the second order term.
This scheme only works with a function that returns (res, sres)
- update(s, y, r, hess, yb)[source]
Update the structured approximation
- Parameters
s (
numpy.ndarray
) – step in optimization parametersy (
numpy.ndarray
) – step in gradient parametersr (
numpy.ndarray
) – residual vectorhess (
numpy.ndarray
) – user-provided (Gauss-Newton) Hessian approximationyb (
numpy.ndarray
) – approximation to A*s, where A is structured approximation matrix
- Return type
Functions
- fides.hessian_approximation.broyden_class_update(y, s, mat, phi=0.0, enforce_curv_cond=True)[source]
Scale free implementation of the broyden class update scheme.
- Parameters
y (
numpy.ndarray
) – difference in gradients (
numpy.ndarray
) – search direction in previous stepmat (
numpy.ndarray
) – current hessian approximationphi (
float
) – convex combination parameter, interpolates between BFGS (phi=0) and DFP (phi=1).enforce_curv_cond (
bool
) – boolean that controls whether the employed broyden class update should attempt to preserve positive definiteness. If set to True, updates from steps that violate the curvature condition will be discarded.
- Return type