How to perform quadratic optimization

3 views (last 30 days)
Parag on 27 Apr 2016
Answered: Jack on 1 Oct 2019
Hi, I am trying to perform quadratic optimisation for correction the Edge spread function as shown in the figure.Looking at the figure we can visualize optimize solution that the ESF should not have bumps before start and end of slope but want to optimize using optimization methods.
I never used the optimisation tool before. I want to optimize it by minimizing the equation as shown below and constrain.
I tried following thing but it shooting me error fmincon stopped because it exceeded the function evaluation limit
if true
% code
A = zeros(length(ESF),length(ESF));
% ind_c = 0;
for i = 1:length(ESF)-1
A(i,i) = 1;
A(i,i+1)= -1;
X = sym('x',[365,1]);
func = @(y)(sum((y-ESF).^2)); b = zeros(size(A,1),1); y = fmincon(func,Es',A,b); end
i will appreciate your help. Thank you

Accepted Answer

John D'Errico
John D'Errico on 27 Apr 2016
Edited: John D'Errico on 27 Apr 2016
I'll assume that you have the optimization toolbox. If you don't, you could still do this, but I'm not going to teach you how to write the equivalent of lsqlin, and you need lsqlin to do this efficiently. Ok, yes, you can solve the problem using lsqnonneg.
I assume that you have several hundred points from a curve, and you wish to solve the problem where you minimize the difference between the data and an approximate curve, subject to the constraint that the new curve is truly monotone. I don't see any data attached, so I cannot give you an example where I solve your problem.
The absolutely simplest solution is to use my SLM toolbox, to fit a (piecewise linear) spline through your data, subject to the constraints that the curve is monotone.
(I need to run right now for about an hour. I'll return and see if I can cook up an example of how that will work. If you are able to attach some data, that would help.)
More complex is you could solve it using lsqlin. But that will force you to formulate the problem for that tool. Slightly harder than using SLM, but not a lot harder. BRB...
Back, with an example of what I THINK you are asking to do. I'll start with a smooth curve over a few points in x. I'll use erf, to create a basic sigmoidal shape. Then interpolate it using spline. This will incur ringing in the curve.
x = linspace(-15,15,20);
xi = -15:.1:15;
yi = spline(x,erf(x),xi);
So you can see the ringing/edge effects in the curve.
slm = slmengine(xi,yi,'knots',xi,'degree',1,'plot','on','increasing','on');
The plot with the points and the fitted curve is hard to see, but I'll zoom into the shoulder area next.
As you see, the curve follows the original data exactly wherever it was monotone already. It only modifies the shape when it needs to do so, using the requirement of monotonicity.
You can extract the actual fitted y values (for this curve type) directly from the struct it produces as:
yihat = slm.coef';
So this is a very simple way to solve your problem. Simple, because slmengine did all the work for you. Find the SLM toolbox on the file exchange at this link .
I could also have approximated the curve using a nice smooth spline, again using monotonicity as a requirement. But that was not your stated goal.
And yes, I can show you how to solve this problem using other tools, like lsqlin or lsqnonneg. But do I really need to do so? I am being lazy and hoping SLM will suffice. :)
Parag on 28 Apr 2016
Thank you so much for your help and brief explanation. I learned a lot

Sign in to comment.

More Answers (2)

Teja Muppirala
Teja Muppirala on 28 Apr 2016
Here is how you might solve it using QUADPROG. Note that you don't need to do anything symbolic using "sym".
L = length(ESF);
H = eye(L);
f = -ESF;
A = zeros(L-1,L);
for i = 1:L-1
A(i,i) = 1;
A(i,i+1)= -1;
b = zeros(L-1,1);
yopt = quadprog(H,f,A,b);
hold on;
  1 Comment
Parag on 28 Apr 2016
Thank you so much for your help. It works really well by using quadprog

Sign in to comment.

Jack on 1 Oct 2019
The physical Edge Spread Function is monotonic in normal conditions, so the undershoots/overshoots are due to processing (demosaicing or, more likely, sharpening). If one uses raw data, there never are under/overshoots.
This is actual raw data projected onto the normal, from the green channels of a 310x248px edge capture, slanted 13.45 degrees off the vertical:
Esge Profile Projected on Normal GFX50s110mmf28.png
We would like to get an ESF out of that without doing much filtering, if any. Quadprod does a decent job of it with Teja's code (1 and -1 switched in the for loop for this orientation, correct?) compared to a fit by the sum of three sigmoids:
Binned Edge Profile Sigmoids vs Quad Prog Monotonic_.png
There is a little more nuance in quadprod's ESF and that nuance is important for higher frequency accuracy in the OTF.
Thank you Teja.
PS Projected raw data attached for anyone who would like to try it or improve upon this (pr is location on the projected axis, esfr is the raw intensity at pr, mid is an educated guess at the center of the edge.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!