Extending Occam's inversion with lasso fusion, overcomplete dictionaries, and isotropic total variation regularisation

Abstract

Occam's inversion is a robust algorithm to perform nonlinear geophysical inversion. It provides the smoothest model within observation noise, thereby discouraging geological overinterpretation. While Occam originally penalised l2 model roughness, l1 can be used to provide models that are visually sharp. However, l1 regularised geophysical inversion has diverged from the larger body of statistics and imaging literature. For example, l1 regularisation with a difference operator (i.e., total variation) in multiple dimensions has two distinct forms, only one of which is invariant to edge orientation -- a distinction often overlooked in geophysics. In one dimension this reduces to the fused lasso problem in statistics. Within the Occam framework, for l2 data norm and l1 model regularisation, we show that lasso fusion through either synthesis or analysis leads to the same 1D problem. We extend the framework to multiple dimensions and to operators such as wavelet transforms in a dictionary, unravelling the mathematics behind three commonly used solvers for l1 regularised problems. These are coordinate descent, iteratively reweighted least squares (IRLS), and split Bregman. We use them to solve a sequence of problems that are linear (1D regression, 2D deblurring) and nonlinear (1D airborne transient electromagnetics) including a field data example. We recommend coordinate descent for 1D problems, and IRLS over split Bregman for 2D problems, with IRLS requiring up to an order of magnitude fewer least squares solves. We hope this work will further encourage geophysicists to adopt l1 regularised inversion, by presenting its various forms as a familiar Occam's inversion.