Lecture 12 Equality and Inequality Constraints Syllabus Lecture 01 Lecture 02 Lecture 03 Lecture 04 Lecture 05

Lecture 06 Lecture 07 Lecture 08 Lecture 09 Lecture 10 Lecture 11 Lecture 12 Lecture 13 Lecture 14

Lecture 15 Lecture 16 Lecture 17 Lecture 18 Lecture 19 Lecture 20 Lecture 21 Lecture 22 Lecture 23

Lecture 24 Describing Inverse Problems Probability and Measurement Error, Part 1 Probability and Measurement Error, Part 2 The L2 Norm and Simple Least Squares A Priori Information and Weighted Least Squared Resolution and Generalized Inverses Backus-Gilbert Inverse and the Trade Off of Resolution and Variance

The Principle of Maximum Likelihood Inexact Theories Nonuniqueness and Localized Averages Vector Spaces and Singular Value Decomposition Equality and Inequality Constraints L1 , L Norm Problems and Linear Programming Nonlinear Problems: Grid and Monte Carlo Searches Nonlinear Problems: Newtons Method Nonlinear Problems: Simulated Annealing and Bootstrap Confidence Intervals

Factor Analysis Varimax Factors, Empirical Orthogonal Functions Backus-Gilbert Theory for Continuous Problems; Radons Problem Linear Operators and Their Adjoints Frchet Derivatives Exemplary Inverse Problems, incl. Filter Design Exemplary Inverse Problems, incl. Earthquake Location Exemplary Inverse Problems, incl. Vibrational Problems

Purpose of the Lecture Review the Natural Solution and SVD Apply SVD to other types of prior information and to equality constraints Introduce Inequality Constraints and the Notion of Feasibility Develop Solution Methods Solve Exemplary Problems

Part 1 Review the Natural Solution and SVD subspaces model parameters mp can affect data

m0 cannot affect data data dp can be fit by model d0 cannot be fit by any model natural solution determine mp by solving dp-Gmp=0 set m0=0

natural solution determine mp by solving dp-Gmp=0 set m0=0 error reduced to its minimum E=e0Te0 natural solution determine mp by solving dp-Gmp=0 set m0=0

solution length reduced to its minimum L=mpTmp Singular Value Decomposition (SVD) singular value decomposition UTU=I and VTV=I

suppose only p s are non-zero suppose only p s are non-zero only first p columns of U

only first p columns of V UpTUp=I and VpTVp=I since vectors mutually pependicular and of unit length UpUpTI and VpVpTI

since vectors do not span entire space The Natural Solution The Natural Solution natural generalized inverse G-g

resolution and covariance Part 2 Application of SVD to other types of prior information and to equality constraints

general solution to linear inverse problem general minimum-error solution 2 lectures ago general minimum-error solution

natural solution plus amount of null vectors you can adjust to match whatever a priori information you want

for example m= by minimizing L=||m-||2 w.r.t. you can adjust to match whatever a priori information you want for example m=

by minimizing L=||m-||2 w.r.t. get =V0T so m = Vpp-1UpTd + T equality constraints minimize E with constraint Hm=h Step 1 find part of solution constrained by

Hm=h SVD of H (not G) H = VppUpT so m=Vpp-1UpTh + V0 Step 2 convert Gm=d into and equation for

GVpp-1UpTh + GV0 = d and rearrange [GV0] = [d - GVpp-1UpTh] G= d Step 3 solve G= d for using least squares

Step 4 reconstruct m from m=Vpp-1UpTh + V0 Part 3 Inequality Constraints and the Notion of Feasibility

Not all inequality constraints provide new information x>3 x>2 Not all inequality constraints provide new information x>3 x>2

follows from first constraint Some inequality constraints are incompatible x>3 x<2 Some inequality constraints are

incompatible x>3 x<2 nothing can be both bigger than 3 and smaller than 2

every row of the inequality constraint Hm h divides the space of m into two parts one where a solution is feasible one where it is infeasible the boundary is a planar surface when all the constraints are considered together

they either create a feasible volume or they dont if they do, then the solution must be in that volume if they dont, then no solution exists feasible region

(A) m2 (B) m2 m1

m1 now consider the problem of minimizing the error E subject to inequality constraints Hm h if the global minimum is inside the feasible region

then the inequality constraints have no effect on the solution but if the global minimum is outside the feasible region then the solution is on the surface

of the feasible volume but if the global minimum is outside the feasible region then the solution is on the surface of the feasible volume the point on the surface where E is the smallest

Hmh m2 feasible mest -E

infeasible m1 Emin furthermore the feasible-pointing normal to the surface must be parallel to E

else you could slide the point along the surface to reduce the error E Hmh m 2

mest -E Emin m 1 Kuhn Tucker theorem

its possible to find a vector y with yi0 such that its possible to find a vector y with y0 such that feasible-pointing normals to surface its possible to find a vector y with y0 such that feasible-pointing normals to surface

the gradient of the error its possible to find a vector y with y0 such that feasible-pointing normals to surface the gradient

of the error is a non-negative combination of feasible normals its possible to find a vector y with y0 such that feasible-pointing normals to surface

the gradient of the error y specifies the combination is a non-negative combination of feasible normals

its possible to find a vector y with y0 such that for linear case with Gm=d its possible to find a vector y with y0 such that some coefficients yi are

positive its possible to find a vector y with y0 such that some coefficients yi are positive the solution is on the corresponding constraint

surface its possible to find a vector y with y0 such that some coefficients yi are zero its possible to find a vector y with y0 such that

some coefficients yi are zero the solution is on the feasible side of the corresponding constraint surface Part 4

Solution Methods simplest case minimize E subject to mi>0 (H=I and h=0) iterative algorithm with two nested loops Step 1

Start with an initial guess for m The particular initial guess m=0 is feasible It has all its elements in mE constraints satisfied in the equality sense Step 2 Any model parameter mi in mE that has associated with it a negative gradient [E]i can be changed both to decrease the error and to remain feasible.

If there is no such model parameter in mE, the Kuhn Tucker theorem indicates that this m is the solution to the problem. Step 3 If some model parameter mi in mE has a corresponding negative gradient, then the solution can be changed to decrease the prediction error. To change the solution, we select the model parameter corresponding

to the most negative gradient and move it to the set mS. All the model parameters in mS are now recomputed by solving the system GSmS=dS in the least squares sense. The subscript S on the matrix indicates that only the columns multiplying the model parameters in mS have been included in the calculation. All the mEs are still zero. If the new model parameters are all feasible, then we set m = m and return to Step 2. Step 4

If some of the elements of mS are infeasible, however, we cannot use this vector as a new guess for the solution. So, we compute the change in the solution and add as much of this vector as possible to the solution mS without causing the solution to become infeasible. We therefore replace mS with the new guess mS + m, where is the largest choice that can be made without some mS becoming infeasible. At least one of the mSis has its constraint satisfied in the equality sense and must

be moved back to mE. The process then returns to Step 3. In MatLab mest = lsqnonneg(G,dobs); example gravitational field depends upon density

via the inverse square law example gravitational force depends upon density model parameters observations

via the inverse square law theory (B) 200 150

(D) 100 50 0 -2 -1.5

-1 -0.5 0 0.5 1

distance, y yi 1.5 (xi,yi) 2

2.5 3 (E) ij (A)

R gravitydobs(y) force, di (C) ij

(F) (xj,yj) more complicated case minimize ||m||2 subject to Hmh this problem is solved by transformation to the previous

problem solve by non-negative least squares then compute mi as with e=d-Gm In MatLab

(A) m1-m20 infeasible m2 (B) m1-m21 infeasible

m2 (C) m10.2 infeasible m2 m2 (D) Intersection infeasible mest

feasible m1 feasible m1

feasible feasible m1 m1 yet more complicated case

minimize ||d-Gm||2 subject to Hmh this problem is solved by transformation to the previous problem minimize ||m|| subject to Hmh where

and (A) (B) 0 1

0 2 m2 0

(C) 1 2 m2 d 2

infeasible 1.5 1 1 1

0.5 feasible 2 feasible 2

m1 0 m1 z 0

0.5 1 In MatLab [Up, Lp, Vp] = svd(G,0); lambda = diag(Lp); rlambda = 1./lambda; Lpi = diag(rlambda);

% transformation 1 Hp = -H*Vp*Lpi; hp = h + Hp*Up'*dobs; % transformation 2 Gp = [Hp, hp]'; dp = [zeros(1,length(Hp(1,:))), 1]'; mpp = lsqnonneg(Gp,dp); ep = dp - Gp*mpp; mp = -ep(1:end-1)/ep(end);

% take mp back to m mest = Vp*Lpi*(Up'*dobs-mp); dpre = G*mest;