Menke W. Geophysical Data Analysis: Discrete Inverse Theory

Подождите немного. Документ загружается.

Linear,

Gaussian Inverse Problem, Viewpoint

spread(N)

IIN

Ill:

C, C,

[N,

Z,]’

i-I

j-l

(4.1 1)

spread(R)

IIR

Ill:

C, C,

[R,

Z,],

1-1

J=I

These measures of the goodness of the resolution spread are based on

the

norm

the difference between the resolution matrix and an

identity matrix. They are sometimes called the

Dirichlet

spread

func-

tions.

When

spread(R)

Since the unit standard deviation of the model parameters is a

measure of the amount of error amplification mapped from data to

model parameters, this quantity can be used to estimate the size of the

unit covariance matrix as

size([cov,

m])

Il[var,

m]L/211$

[cov,

m],

(4.12)

where the square root is interpreted component by component. Note

that this measure of covariance size does not take into account the size

the off-diagonal elements in the unit covariance matrix.

4.7

Generalized Inverses with Good

Resolution and Covariance

Having found a way to measure quantitatively the goodness of the

resolution and covariance of a generalized inverse, we now consider

whether it is possible to use these measures as guiding principles for

deriving generalized inverses. This procedure is analogous to that of

chapter

which involves first defining measures of solution prediction

error and simplicity and then using those measures to derive the least

squares and minimum length estimates of the model parameters.

4.7.1

OVERDETERMINED CASE

We first consider a purely overdetermined problem of the form

We postulate that this problem has a solution of the form

mest

G-gd

and try to determine

G-g

by minimizing some combina-

tion of the above measures of goodness. Since we previously noted that

4.7

Generalized Inverses with Good Resolution and Covariance

the overdetermined least squares solution had perfect model resolu-

tion, we shall try to determine

G-g

by minimizing only the spread

the data resolution. We begin by examining the spread of the

row

Since each of the

Jk)s

is positive, we can minimize

spread(N)

by minimizing each individual

Jk.

(4.13)

the total

therefore

insert the definition of the data resolution matrix

GG-g

into the

formula for

and minimize it with respect to the elements of the

generalized inverse matrix:

aJk/aG;:

(4.14)

We shall perform the differentiation separately for each of the three

terms of

Jk.

The first term is given by

The second term is given by

NM NM

(4.16)

The third term is zero, since its is not a function of the generalized

inverse. Writing the complete equation in matrix form yields

GTGG-B

(4.17)

70 4

Linear, Gaussian Inverse Problem, Viewpoint

Since

GTG

is square, we can premultiply by its inverse to solve for the

generalized inverse,

G-g

[GTG]-'GT,

which is precisely the same as

the formula for the least squares generalized inverse. The least squares

generalized inverse can be interpreted either as the inverse that mini-

mizes the

norm of the prediction error or as the inverse that

minimizes the Dirichlet spread of the data resolution.

4.7.2

UNDERDETERMINED CASE

The data can be satisfied exactly in a purely underdetermined

problem. The data resolution matrix is, therefore, precisely an identity

matrix and its spread is zero. We might therefore try to derive a

generalized inverse for this problem by minimizing the spread of the

model resolution matrix with respect to the elements of the general-

ized inverse. It is perhaps not particularly surprising that the general-

ized inverse obtained by this method is exactly the minimum length

generalized inverse

G-g

GT[GGT]-'.

The minimum length solution

can be interpreted either as the inverse that minimizes the

norm

the solution length or as the inverse that minimizes the Dirichlet

spread of the model resolution. This is another aspect

the symmetri-

cal relationship between the least squares and minimum length

solu-

tions.

4.7.3

THE GENERAL CASE WITH DIRICHLET

SPREAD FUNCTIONS

We seek the generalized inverse

G-B

that minimizes the weighted

sum of Dirichlet measures of resolution spread and covariance size.

Minimize:

spread(N)

spread(R)

size([cov,

m])

(4.18)

where the

a's

are arbitrary weighting factors. This problem is done in

exactly the same fashion as the one in Section

4.7.1

except that there is

now three times as much algebra. The result is an equation for the

generalized inverse:

a,[GT'G]G-g

G-g[~,GGT

CY~[COV,

d]]

[a,

a2]GT

(4.19)

This equation has no explicit solution for

G-g

in terms of an algebraic

function of the various matrices. Explicit solutions can be written,

however, for a variety of special choices of the weighting factors. The

least squares solution is recovered if

and

and the

4.8

Sidelobes

and

the

Backus-Gilbert

Spread

Function

minimum length solution is recovered if

and

more interest is the case in which

cyI

equals some

constant (say,

c2)

and [cov,

The generalized inverse is then

given by

G-B

[GTG

cZ21]-IGT

(4.20)

This formula is precisely the damped least squares inverse, which we

derived in the previous chapter by minimizing a combination of

prediction error and solution length. The damped least squares solu-

tion can also be interpreted as the inverse that minimizes a weighted

combination of data resolution spread and covariance size.

Note that it is quite possible for these generalized inverses to possess

resolution matrices containing

negative

off-diagonal elements. For

interpreting the rows of the matrices as localized averages, this is an

unfortunate property. Physically, an average would make more sense

it contained only positive weighting factors. In principle, it is

possible to include nonnegativity as a constraint when choosing the

generalized inverse by minimizing the spread functions. However, in

practice this constraint

never implemented, because it makes the

calculation of the generalized inverse very difficult.

4.8

Sidelobes and the Backus-Gilbert

Spread Function

When there is a natural ordering of data or model parameters, the

Dirichlet spread function is not a particularly appropriate measure of

the goodness of resolution, because the off-diagonal elements of the

resolution matrix are all weighted equally, regardless of whether they

are close or far from the main diagonal. If there is a natural ordering,

we would much prefer that any large elements be near the main

diagonal (Fig.

4.4).

The rows of the resolution matrix are then local-

ized'averagmg functions.

If one uses the Dirichlet spread function to compute a generalized

inverse, it will often have

sidelobes,

that is, large amplitude regions in

the resolution matrices far from the main diagonal. We would prefer to

find a generalized inverse without sidelobes, even at the expense of

widening the band of nonzero elements near the main diagonal, since

Linear, Gaussian Inverse Problem, Viewpoint

index

Fig.

4.4.

(a and b) Resolution matrices have the same spread, when measured by the

Dirichlet spread function. Nevertheless, if the model parameters possess a natural

ordering, then (a)

better resolved. The Backus-Gilbert spread function

designed to

measure (a) as having a smaller spread than (b).

a solution with such a resolution matrix is then interpretable as a

localized average of physically adjacent model parameters.

We therefore add a weighting factor

w(i,

to the measure of spread

that weights the

(i,

element

R according to its physical distance

from the diagonal element. This weighting preferentially selects reso-

lution matrices that are “spiky,” or “deltalike.” If the natural ordering

were a simple linear one, then the choice

w(i,

-j)2

would be

reasonable. If the ordering is multidimensional, a more complicated

weighting factor is needed. It is usually convenient to choose the

spread function

that the diagonal elements have no weight, i.e.,

w(i,

and

that

w(i,

is always nonnegative and symmetric in

and

The new spread function, often called the Backus-Gilbert

spread function, is then given by

MM MM

spread(R)

w(i,

j)[R,

ZijI2

w(i,

j)Ri

(4.21)

I-I

j-1

I-lj-I

[Refs.

1,2].

similar expression holds for the spread

the data resolution. One

can now use this measure of spread to derive new generalized inverses.

Their sidelobes will be smaller than those based on the Dirichlet spread

functions. On the other hand, they are sometimes worse when judged

by other criteria.

we shall see, the Backus-Gilbert generalized

inverse for the completely underdetermined problem does not exactly

satisfy the data, even though the analogous minimum-length general-

ized inverse does. These facts demonstrate that there are unavoidable

trade-offs inherent in finding solutions to inverse problems.

4.9

The Backus-Gilbert Generalized Inverse

4.9

The Backus

Gilbert Generalized

Inverse for the Underdetermined

Problem

This problem is analogous to deriving the minimum length solution

by minimizing the Dirichlet spread of model resolution. Since it is very

easy to satisfy the data when the problem is underdetermined

(so

that

the data resolution has small spread), we shall find a generalized

inverse that minimizes the spread of the model resolution alone.

We seek the generalized inverse

G-g

that minimizes the Backus

Gilbert spread of model resolution. Since the diagonal elements of the

model resolution matrix are given

no weight, we also require that the

resulting model resolution matrix satisfy the equation

RtJ=

’;

(4.22)

This constraint ensures that the diagonal of the resolution matrix is

finite and that the rows are unit averaging functions acting on the true

model parameters. Writing the spread of one row of the resolution

matrix as

and inserting the expression for the resolution matrix, we

have

I-

W(1,

k)Rk/R,

where the quantity

[s;,],

is defined as

(4.23)

(4.24)

The left-hand side of the constraint equation

can also be

written in terms

the generalized inverse

Linear, Gaussian Inverse

Problem,

Viewpoint

Here the quantity

is defined as

(4.26)

The problem of minimizing

with respect to the elements of the

generalized inverse (under the given constraints) can be solved

through the use of Lagrange multipliers. We first define a Lagrange

function

such that

(4.27)

i-I

j-1

where

-22

is the Lagrange multiplier. We then differentiate

with

respect to the elements of the generalized inverse and set the result

equal to zero as

a@,/dGij

[Spr]kGi:

2AUP

(4.28)

(Note that one can solve for each row of the generalized inverse

separately,

that it is only necessary to take derivatives with respect to

the elements in the kth row.) The above equation must be solved along

with the original constraint equation. Treating the kth row of

G-8

as a

vector and the quantity

[S,]!,

as a matrix in indices

andj, we can write

these equations as the matnx equation

r=l

(4.29)

This is a square

system oflinear equations that must

be solved for the

elements of the kth row of the generalized inverse

and for the one Lagrange multiplier

The matrix equation can be solved explicitly using a variant of the

“bordering method” of linear algebra, which is used to construct the

inverse of a matrix by partitioning it into submatrices with simple

properties. Suppose that the inverse

the matrix in

Eq.

(4.29) exists

and that we partition it into an

square matrix

vector

and

scalar

By assumption, premultiplication by the inverse yields the

identity matrix

4.10

Including the

Covariance

Size

A[S,J]k

buT Au]

(4.30)

[

bT[s,J]k

CUT

bTu

The unknown submatrices

and c can now be determined by

equating the submatrices

A[S,],

buT

that

[S,];l[I

buT]

[Sl,IklU

that

[S,J;lu

buT[S,J;'u

and

UTISIJlklU

(4.31)

-1

bT[SjJ],

cuT

0 so

that

S,],

Once the submatrices

A, b,

and c are known, it is easy to solve the

equation

for

the generalized inverse. The result is written as

(4.32)

4.10

Including the Covariance Size

The measure

goodness that was used to determine the Backus-

Gilbert inverse can be modified to include a measure

the covariance

size of the model parameters [Ref.

31.

We shall use the same measure

as we did when considering the Dihchlet spread

goodness is measured

functions,

that

spread(R)

(

size([cov,

m])

w(i,

j)Rt

[cov,

mIli

(4.33)

i-I

j-1

where

is a weighting factor that determines the relative

contribution of model resolution and covariance to the measure of the

goodness

the generalized inverse. The goodness

of the kth row is

Linear, Gaussian Inverse Problem, Viewpoint

then

w(k,

I)R,&

a)[cov,

mIkk

I=I

(4.34)

where the quantity

[S&

is defined by the equation

[S:,Ik

a[S,,Ik

a"ov,

41,

(4.35)

Since the function

has exactly the same form as

had in the

previous section, the generalized inverse is just the previous result with

[S,Ik

replaced by

[S:,],:

(4.36)

4.1

The Trade-off

Resolution and

Variance

Suppose that one is attempting to determine a set of model parame-

ters that represents a discretized version of a continuous function, such

as x-ray opacity in the medical tomography problem (Fig.

4.5).

If the

discretization is made very fine, then the x rays will not sample every

box; the problem will be underdetermined. If we try to determine the

opacity of each box individually, then estimates of opacity will tend to

have rather large variance. Few boxes will have several x rays passing

through them,

that little averaging out of the errors will take place.

On the other hand, the boxes are very small -and very small features

can be detected (the resolution is very good). The large variance can be

reduced by increasing the box size (or alternatively, averaging several

neighboring boxes). Each of these larger regions will then contain

several x rays, and noise will tend to be averaged out. But because the

regions are now larger, small features can

longer be detected and the

resolution of the x-ray opacity has become poorer.

4.1

The Trade-off

Resolution and Variance

’

:p=’

Fig.

4.5.

(a) X-ray paths (straight lines) through body in the tomography problem. (b)

Coarse parameterizations

the opacity have good variance, since several rays pass

through each box. But since they are coarse, they only poorly resolve the structure ofthe

body. (c, Fine parameterizations have poorer variance but better resolution.

This scenario illustrates an important trade-off between model

resolution spread and variance size. One can be decreased only at the

expense

increasing the other. We can study this trade-off by choos-

ing a generalized inverse that minimizes a weighted sum

resolution

spread and covariance size:

spread(R)

size([cov,

m])

(4.37)

spread

model resolution

Fig.

4.6.

Trade-off curve

resolution and variance

for

a given discretization of a

continuous function. The larger the

the more weight resolution is given (relative to

variance) when forming the generalized inverse. The details

the trade-off curve

depend on the parameterization. The resolution can be no better than the smallest

element in the parameterization and no worse than the sum of all the elements.