Minimizer#
Once we defined the objective function for our inversion, we need to find a model vector that would minimize it. To do so we need to make use of a minimizer.
In the new framework, a minimizer can be a function that takes the objective function (and other arguments, like an initial model) and returns the inverted model or recovered model, i.e. the model that minimize that objective function.
Linear minimizers#
If our inversion is linear (e.g. potential field inversions), we can make use of linear minimizers, like inversion_ideas.minimize.conjugate_gradient() to minimize the objective function.
The conjugate_gradient() function will ask the objective function and an initial model from which it’ll start the conjugate gradient algorithm, and it will return the recovered model after convergence has been achieved.
For example, let’s consider the LinearRegressor simulation and the data_misfit object we built when we introduced the objective function:
data_misfit = ii.DataMisfit(data, uncertainties, simulation)
data_misfit
where the true_model was:
true_model
array([-5., 2., -3., 4., 1.])
Let’s define an objective function by adding a SimpleSmallness regularization:
beta = 30.0
phi = data_misfit + beta * ii.SimpleSmallness(n_params)
phi
And use the conjugate_gradient() to minimize it by setting an initial model full of zeros:
initial_model = np.zeros(n_params, dtype=np.float64)
model = ii.minimize.conjugate_gradient(phi, initial_model)
model
array([-5.20251337, 1.72959343, -2.26550802, 4.04188586, 1.56209423])
We were able to recover a model that is close to the true model used to generate the synthetic data.
We can make use of the inversion_ideas.DataMisfit.chi_factor() method to see if we are over- or under-fitting the data.
\(\chi \approx 1.0\) would mean that we are recovering the data within the known uncertainties.
\(\chi > 1.0\) would mean that we are under-fitting the data.
\(\chi < 1.0\) would mean that we are over-fitting the data.
data_misfit.chi_factor(model)
np.float64(0.8665542696537072)
Non-linear minimizers#
If our inversion is non-linear (DC-IP, EM inversions, etc.), we’ll need to use a suitable minimizer for non-linear inversion problems, like the Gauss-Newton method.
The inversion_ideas.minimize.GaussNewtonConjugateGradient object implements a Gauss-Newton algorithm using Conjugate Gradient to find search directions, and performing a backtracking line search.
We can use it to minimize the non-linear objective function through the inversion_ideas.minimize.GaussNewtonConjugateGradient.run() method:
minimizer = ii.minimize.GaussNewtonConjugateGradient()
inverted_model = minimizer.run(phi, initial_model)
Alternatively we can iterate over it after calling the minimizer:
minimizer = ii.minimize.GaussNewtonConjugateGradient()
for model in minimizer(phi, initial_model):
# The model is the one obtained after each Gauss-Newton iteration
...
# We can perform stuff here, like plot the model, print information,
# save stuff to disk, etc.
...
# The inverted model is the last value of the `model` variable
inverted_model = model.copy()
Third-party minimizers#
It’s also possible to plug our Objective function objects into third-party minimizers, likes the ones available in scipy.
For example, we could use scipy.optimize.minimize() to minimize our objective function using a conjugate gradient algorithm:
from scipy.optimize import minimize
# Pass the phi.gradient method as the `jac` argument
result = minimize(phi, x0=initial_model, method="CG", jac=phi.gradient)
inverted_model = result.x
inverted_model
array([-5.20251337, 1.72959343, -2.26550802, 4.04188586, 1.56209423])