c18-4 (779611)
Текст из файла
804Chapter 18.Integral Equations and Inverse Theory18.4 Inverse Problems and the Use of A PrioriInformationminimize:A[u]orminimize:B[u](18.4.1)(Of course these will generally give different answers for u.) As another possibility,now suppose that we want to minimize A[u] subject to the constraint that B[u] havesome particular value, say b. The method of Lagrange multipliers gives the variationδδ{A[u] + λ1 (B[u] − b)} =(A[u] + λ1 B[u]) = 0δuδu(18.4.2)where λ1 is a Lagrange multiplier.
Notice that b is absent in the second equality,since it doesn’t depend on u.Next, suppose that we change our minds and decide to minimize B[u] subjectto the constraint that A[u] have a particular value, a. Instead of equation (18.4.2)we haveδδ{B[u] + λ2 (A[u] − a)} =(B[u] + λ2 A[u]) = 0δuδu(18.4.3)with, this time, λ2 the Lagrange multiplier. Multiplying equation (18.4.3) by theconstant 1/λ2 , and identifying 1/λ2 with λ1 , we see that the actual variations areexactly the same in the two cases. Both cases will yield the same one-parameterfamily of solutions, say, u(λ1 ).
As λ1 varies from 0 to ∞, the solution u(λ1 )varies along a so-called trade-off curve between the problem of minimizing A andthe problem of minimizing B. Any solution along this curve can equally wellbe thought of as either (i) a minimization of A for some constrained value of B,or (ii) a minimization of B for some constrained value of A, or (iii) a weightedminimization of the sum A + λ1 B.The second preliminary point has to do with degenerate minimization principles.In the example above, now suppose that A[u] has the particular formA[u] = |A · u − c|2(18.4.4)for some matrix A and vector c.
If A has fewer rows than columns, or if A is squarebut degenerate (has a nontrivial nullspace, see §2.6, especially Figure 2.6.1), thenminimizing A[u] will not give a unique solution for u. (To see why, review §15.4,and note that for a “design matrix” A with fewer rows than columns, the matrixAT · A in the normal equations 15.4.10 is degenerate.) However, if we add anymultiple λ times a nondegenerate quadratic form B[u], for example u · H · u with Ha positive definite matrix, then minimization of A[u] + λB[u] will lead to a uniquesolution for u. (The sum of two quadratic forms is itself a quadratic form, with thesecond piece guaranteeing nondegeneracy.)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).Later discussion will be facilitated by some preliminary mention of a coupleof mathematical points. Suppose that u is an “unknown” vector that we plan todetermine by some minimization principle. Let A[u] > 0 and B[u] > 0 be twopositive functionals of u, so that we can try to determine u by either80518.4 Inverse Problems and the Use of A Priori InformationWe can combine these two points, for this conclusion: When a quadraticminimization principle is combined with a quadratic constraint, and both arepositive, only one of the two need be nondegenerate for the overall problem to bewell-posed. We are now equipped to face the subject of inverse problems.Suppose that u(x) is some unknown or underlying (u stands for both unknownand underlying!) physical process, which we hope to determine by a set of Nmeasurements ci , i = 1, 2, .
. . , N . The relation between u(x) and the ci ’s is thateach ci measures a (hopefully distinct) aspect of u(x) through its own linear responsekernel ri , and with its own measurement error ni . In other words,Zci ≡ si + ni = ri (x)u(x)dx + ni(18.4.5)(compare this to equations 13.3.1 and 13.3.2). Within the assumption of linearity,this is quite a general formulation. The ci ’s might approximate values of u(x) atcertain locations xi , in which case ri (x) would have the form of a more or lessnarrow instrumental response centered around x = xi . Or, the ci ’s might “live” in anentirely different function space from u(x), measuring different Fourier componentsof u(x) for example.The inverse problem is, given the ci ’s, the ri (x)’s, and perhaps some informationabout the errors ni such as their covariance matrixSij ≡ Covar[ni , nj ](18.4.6)how do we find a good statistical estimator of u(x), call it ub(x)?It should be obvious that this is an ill-posed problem. After all, how can wereconstruct a whole function ub(x) from only a finite number of discrete values ci ?Yet, whether formally or informally, we do this all the time in science.
We routinelymeasure “enough points” and then “draw a curve through them.” In doing so, weare making some assumptions, either about the underlying function u(x), or aboutthe nature of the response functions ri (x), or both. Our purpose now is to formalizethese assumptions, and to extend our abilities to cases where the measurements andunderlying function live in quite different function spaces. (How do you “draw acurve” through a scattering of Fourier coefficients?)We can’t really want every point x of the function ub(x).
We do want somelarge number M of discrete points xµ , µ = 1, 2, . . ., M , where M is sufficientlylarge, and the xµ ’s are sufficiently evenly spaced, that neither u(x) nor ri (x) variesmuch between any xµ and xµ+1 . (Here and following we will use Greek letters likeµ to denote values in the space of the underlying process, and Roman letters like ito denote values of immediate observables.) For such a dense set of xµ ’s, we canreplace equation (18.4.5) by a quadrature likeXRiµu(xµ ) + ni(18.4.7)ci =µwhere the N × M matrix R has componentsRiµ ≡ ri (xµ )(xµ+1 − xµ−1 )/2(18.4.8)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).The Inverse Problem with Zeroth-Order Regularization806Chapter 18.Integral Equations and Inverse Theory(or any other simple quadrature — it rarely matters which). We will view equations(18.4.5) and (18.4.7) as being equivalent for practical purposes.How do you solve a set of equations like equation (18.4.7) for the unknownu(xµ )’s? Here is a bad way, but one that contains the germ of some correct ideas:Form a χ2 measure of how well a model ub(x) agrees with the measured data,χ =NN XXi=1 j=1"ci −MX"#Riµ ub(xµ )µ=1"#2PMNXci − µ=1 Riµ ub(xµ )≈σi−1Sijcj −MXµ=1#Rjµub(xµ )(18.4.9)i=1(compare with equation 15.1.5).
Here S−1 is the inverse of the covariance matrix,and the approximate equality holds if you can neglect the off-diagonal covariances,with σi ≡ (Covar[i, i])1/2 .Now you can use the method of singular value decomposition (SVD) in §15.4to find the vector bu that minimizes equation (18.4.9). Don’t try to use the methodof normal equations; since M is greater than N they will be singular, as we alreadydiscussed.
The SVD process will thus surely find a large number of zero singularvalues, indicative of a highly non-unique solution. Among the infinity of degeneratesolutions (most of them badly behaved with arbitrarily large ub(xµ )’s) SVD willselect the one with smallest |bu| in the sense ofX[bu(xµ )]2a minimum(18.4.10)µ(look at Figure 2.6.1).
This solution is often called the principal solution. Itis a limiting case of what is called zeroth-order regularization, corresponding tominimizing the sum of the two positive functionalsminimize:u] + λ(bu·bu)χ2 [b(18.4.11)in the limit of small λ. Below, we will learn how to do such minimizations, as wellas more general ones, without the ad hoc use of SVD.What happens if we determine bu by equation (18.4.11) with a non-infinitesimalvalue of λ? First, note that if M N (many more unknowns than equations), thenu will often have enough freedom to be able to make χ2 (equation 18.4.9) quiteunrealistically small, if not zero. In the language of §15.1, the number of degrees offreedom ν = N − M , which is approximately the expected value of χ2 when ν islarge, is being driven down to zero (and, not meaningfully, beyond). Yet, we knowthat for the true underlying function u(x), which has no adjustable parameters, thenumber of degrees of freedom and the expected value of χ2 should be about ν ≈ N .Increasing λ pulls the solution away from minimizing χ2 in favor of minimizingbu·bu.
From the preliminary discussion above, we can view this as minimizing bu·busubject to the constraint that χ2 have some constant nonzero value. A popularchoice, in fact, is to find that value of λ which yields χ2 = N , that is, to get about asmuch extra regularization as a plausible value of χ2 dictates. The resulting ub(x) iscalled the solution of the inverse problem with zeroth-order regularization.Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).2best smoothness(independent of agreement)807achievable solutionsstbenstiolusobest agreement(independent of smoothness)Better SmoothnessFigure 18.4.1.
Almost all inverse problem methods involve a trade-off between two optimizations:agreement between data and solution, or “sharpness”of mapping between true and estimated solution (heredenoted A), and smoothness or stability of the solution (here denoted B). Among all possible solutions,shown here schematically as the shaded region, those on the boundary connecting the unconstrainedminimum of A and the unconstrained minimum of B are the “best” solutions, in the sense that everyother solution is dominated by at least one solution on the curve.The value N is actually a surrogate for any value drawn from a Gaussiandistribution with mean N and standard deviation (2N )1/2 (the asymptotic χ2distribution).
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.















