Weighted regression without priority zone, but with finding the minimum possible value for all the differences

Dear colleagues,

I need to perform a regression. The input data is shown in Figure 1 and is also attached as an MS Excel file.

Figure 1 (Figures are shown in the end of the letter.

The equation is z = a*x^m/(y+c)^n.

The parameters of the fit are shown in Figure 2 and are also attached.

Figure 2

The output data is shown in Figure 3.

Figure 3

The individual differences in percentage are shown in Figure 4.

Figure 4

I know that by using a weighted matrix, it is possible to decrease some differences (priority zones) at the expense of increasing others. However, in this case, my goal is to estimate the minimal possible error (one border value) for all of the values. For example, in the first fit, the maximum difference is 93%. It is easy to decrease it by applying individual weight coefficients, but of course, somewhere, the difference will increase. Therefore, the question is: what method can be used to find the minimal possible difference? If this for example is 18%, it means that everywhere the difference will be smaller, and it is likely that there is no 0 % anymore. Example by hand (not real) is given in Figure 5.

Figure 5

And if we try to decrease the difference at the point with an 18% difference, this will lead to an increase in the difference somewhere above 18%. Also, if in the first fit I have a 0% or 1% difference at a given point, in the second fit, the difference at this point is likely to be close (but smaller) to 18%. And one more thing – if there is a way for doing the whole procedure, is it possible to use curve fitter app or need to use a code?

Answer 1

[Torsten]

x = [1 2 3 4 5 10 20 50].';
y = [5 10 15 20 25 30 40 50 60 90];
z = [265.110074796797	195.008035083193	156.258229529605	129.742584500194	107.062285017337	91.0540850617739	68.5080545479447	50.5038878827341	35.8057303140135	23.9109370743307	
    321.446939682946	249.023136173216	204.538587797874	172.392832029864	146.102087493907	125.302681110766	98.3496300763363	79.6443293949692	60.4259292015808	28.967085981117	
    358.948929784307	282.185871826562	232.988776882345	197.578293515688	169.230643083598	145.977264293244	115.850743148756	96.2843426347109	74.9936831153838	40.5561051536658	
    387.580526618466	306.229191510467	253.083199667196	215.390669679886	185.62217968629	160.805889976744	128.174360996678	107.804496871799	85.3098905504615	48.5261569465696	
    411.022720927616	325.178097515035	268.617658401558	229.174534346536	198.325915732664	172.4000748844	137.680141245939	116.580164412001	93.3004757115783	54.5677932289588	
    491.330306710826	385.749717456981	316.583948069254	271.810779459274	237.728742469111	208.93517622606	166.914170433146	142.97295566989	118.058749869664	72.5800929043424	
    584.345341008974	448.998217493589	364.117335543095	314.175511964579	277.043092535072	246.270752777746	195.70913612502	168.102854357484	142.722367271741	89.5041059266947	
    730.266495461484	536.918175592217	426.295310861942	369.764885508864	328.878387857907	296.884055308305	233.111367401215	199.486893792689	175.18124008861	110.317421963392];
f = @(p) abs(p(1)*x.^p(2)./(y+p(3)).^p(4) - z);
F = @(p) reshape(f(p),[],1);
p0 = [1 1 1 1];
F(p0)

ans = 80×1

264.9434 321.1136 358.4489 386.9139 410.1894 489.6636 581.0120 721.9332 194.9171 248.8413

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

max(F(p0))

ans = 721.9332

format long
[p,fval] = fminimax(F,p0)

Local minimum possible. Constraints satisfied.

fminimax stopped because the size of the current search direction is less than
twice the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.

p = 1×4

1.0e+04 * 1.012399111966462 0.000025177933233 0.001641135844188 0.000118976592444

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

fval = 80×1

0.749551595760067 6.679548910905453 10.350141109652270 12.794953138220137 14.577849760131755 19.293422491732827 22.302958149803658 22.382430607261313 10.938069716167234 3.808316417156561

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

max(fval)

ans = 
  22.382430608024322

[Boyan]

Dear Torsten,

Do you think that in Curve fitter app there is some option that makes same like fminimax? This possibility of minimizing maximum error often is very usefull as a goal in regression analyses and in my opinion it is strange that fminimax is not embeded in Curve fitter (if it is not). What is you opinion?

[Torsten]

Do you think that in Curve fitter app there is some option that makes same like fminimax?

The Curve Fitter App always tries to minimize the 2-norm of the error, not the Inf-norm. The reason is that a theoretical background for statistical analysis of the results is built on minimizing the sum of errors squared, not on minimizing the maximum error. But I agree: depending on the application, minimizing the maximum error can be useful.

[Boyan]

Dear Torsten,

When I repeat the example above, I have had no problem. But when I used a little bit different data, some problems has appeared, and I can't find my mistake.

the data is:

x = [1 2 3 4 5 10 20 50 100].';

y = [5 10 15 20 25 30 40 50 60 90];

z = [213.10	157.80	127.60	107.10	91.40	74.40	56.70	38.10	28.80	25.60
    258.70	199.90	164.70	138.90	120.10	103.50	81.70	60.60	51.70	32.20
    288.40	225.30	186.20	157.00	136.30	119.30	95.40	73.40	63.90	45.40
    310.60	243.60	201.10	169.50	147.40	129.90	104.60	82.30	72.00	53.10
    328.60	257.90	212.50	179.00	155.80	137.80	111.60	88.90	78.00	58.20
    389.00	302.80	247.20	207.40	180.70	160.60	131.70	109.00	94.90	70.70
    456.70	348.90	280.60	234.20	204.10	181.00	149.90	127.80	109.60	79.10
    559.00	411.30	323.00	267.50	232.70	204.70	171.30	151.20	126.20	86.30
    647.40	459.80	353.80	291.10	252.70	220.40	185.80	167.60	136.90	89.70
    ];

f = @(p) abs(p(1)*x.^p(2)./(y+p(3)).^p(4) - z);

F = @(p) reshape(f(p),[],1);

p0 = [1 1 1 1];

F(p0)

ans = 90×1

212.9333 258.3667 287.9000 309.9333 327.7667 387.3333 453.3667 550.6667 630.7333 157.7091

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

ans =

1.0e+02 *

2.129333333333333

2.583666666666667

2.879000000000000

and so on - the values are very good.

The problem appears when input this:

format long
[p,fval] = fminimax(F,p0)

Local minimum possible. Constraints satisfied.

fminimax stopped because the size of the current search direction is less than
twice the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.

p = 1×4

5.380638654412424 0.914518590478233 -4.999977679457625 0.087235693178237

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

fval = 90×1

1.994040300913148 2.328839186763963 2.509950524156043 2.619382360394732 2.689220554897128 2.765106561291370 2.446643377744690 0.688461177790053 2.765106517365292 1.531241720214055

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

And it is that wrong coefficients "p" appeared:

p =

5.380638654412424 0.914518590478233 -4.999977679457625 0.087235693178237

Could you show me my mistake, please

Best regards,

Boyan

[Torsten]

What is "wrong" with the coefficients ?

If you want them to be positive, e.g., set

lb = [0 0 0 0];
ub = [Inf Inf Inf Inf];

and call fminimax as

[p,fval] = fminimax(F,p0,[],[],[],[],lb,ub)

[Boyan]

Thank you for giving me a way for controling the sign of the coefficients.

I think that coefficients is wrong because when I substitute them, I don't get the good results for the z values which is shown in the begining of the code, i.e.

good results for z values :"....

212.9333

258.3667

287.9000

309.9333 and so on"

Furthermore:

The fval is above 100 % (I forgot to copy previous time), that the error calculated by the procedure itself is also inacceptable, above 100 %:

"fval =

1.0e+02 *

1.994040300913148

2.328839186763963

2.509950524156043

[Torsten]

good results for z values :"....

212.9333

258.3667

287.9000

309.9333 and so on"

??

F(p0) does not give you the fitted z-values, but the error for each z-value. So the above coefficients p0 are really bad because the error is almost as large as the z-values themselves.

The fval is above 100 % (I forgot to copy previous time), that the error calculated by the procedure itself is also inacceptable, above 100 %:

The maximum error max(fval) is 19.7173. Why do you say it's above 100 % ? The fval must be as small as possible - ideally 0 - to signal a good fit.

Below, I first tried to determine optimum coefficients in the 2-norm and use them as initial values for the call to fminimax.

x = [1 2 3 4 5 10 20 50 100].';
y = [5 10 15 20 25 30 40 50 60 90];
z = [213.10	157.80	127.60	107.10	91.40	74.40	56.70	38.10	28.80	25.60
    258.70	199.90	164.70	138.90	120.10	103.50	81.70	60.60	51.70	32.20
    288.40	225.30	186.20	157.00	136.30	119.30	95.40	73.40	63.90	45.40
    310.60	243.60	201.10	169.50	147.40	129.90	104.60	82.30	72.00	53.10
    328.60	257.90	212.50	179.00	155.80	137.80	111.60	88.90	78.00	58.20
    389.00	302.80	247.20	207.40	180.70	160.60	131.70	109.00	94.90	70.70
    456.70	348.90	280.60	234.20	204.10	181.00	149.90	127.80	109.60	79.10
    559.00	411.30	323.00	267.50	232.70	204.70	171.30	151.20	126.20	86.30
    647.40	459.80	353.80	291.10	252.70	220.40	185.80	167.60	136.90	89.70
    ];
f = @(p) p(1)*x.^p(2)./(y+p(3)).^p(4) - z;
F = @(p) reshape(f(p),[],1);
p0 = [1 1 1 1];
p = lsqnonlin(F,p0)

Local minimum possible.

lsqnonlin stopped because the final change in the sum of squares relative to 
its initial value is less than the value of the function tolerance.

p = 1×4

1.0e+03 * 2.5206 0.0002 0.0083 0.0009

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

f = @(p)abs(p(1)*x.^p(2)./(y+p(3)).^p(4) - z)

f = function_handle with value:
    @(p)abs(p(1)*x.^p(2)./(y+p(3)).^p(4)-z)

F = @(p)reshape(f(p),[],1);
p0 = p;
[p,fval] = fminimax(F,p0)

Local minimum possible. Constraints satisfied.

fminimax stopped because the size of the current search direction is less than
twice the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.

p = 1×4

1.0e+03 * 3.7505 0.0002 0.0101 0.0010

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

fval = 90×1

19.7173 11.5999 6.5642 3.2169 0.6662 6.7234 12.8786 18.3582 19.7173 16.0419

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

max(fval)

ans = 19.7173

[Boyan]

Dear Torsten,

Until now I thought that F(p0) is fitted "z" values. Because they are similar to input values. Of course second look gives the answer that this is not so, because F(p0) is before fminimax to be used. Sorry for the question.

But .. this situation arise new questions, sorry also for that.

1. Is it possible to calculate fitted values automatically, to see at first look how close they are to the input values?
2. For the initial example (in the begining of our conversation) the p_i values are, right?:

p = 1×4

1.0e+04 *

1.012399111966462 0.000025177933233 0.001641135844188 0.000118976592444

So, when I substitute them with MS Excel in the formula f = @(p) abs(p(1)*x.^p(2)./(y+p(3)).^p(4) - z); the fitted z values are very very different than input z values.

I attached the MS Excel fitted data file, sheet fminimax. Also I put it here like a picture.

[Torsten]

Is it possible to calculate fitted values automatically, to see at first look how close they are to the input values?

Just evaluate p(1)*x.^p(2)./((y+p(3)).^p(4)) for the optimal parameters to get the fitted z-values:

x = [1 2 3 4 5 10 20 50].';
y = [5 10 15 20 25 30 40 50 60 90];
z = [265.110074796797	195.008035083193	156.258229529605	129.742584500194	107.062285017337	91.0540850617739	68.5080545479447	50.5038878827341	35.8057303140135	23.9109370743307	
    321.446939682946	249.023136173216	204.538587797874	172.392832029864	146.102087493907	125.302681110766	98.3496300763363	79.6443293949692	60.4259292015808	28.967085981117	
    358.948929784307	282.185871826562	232.988776882345	197.578293515688	169.230643083598	145.977264293244	115.850743148756	96.2843426347109	74.9936831153838	40.5561051536658	
    387.580526618466	306.229191510467	253.083199667196	215.390669679886	185.62217968629	160.805889976744	128.174360996678	107.804496871799	85.3098905504615	48.5261569465696	
    411.022720927616	325.178097515035	268.617658401558	229.174534346536	198.325915732664	172.4000748844	137.680141245939	116.580164412001	93.3004757115783	54.5677932289588	
    491.330306710826	385.749717456981	316.583948069254	271.810779459274	237.728742469111	208.93517622606	166.914170433146	142.97295566989	118.058749869664	72.5800929043424	
    584.345341008974	448.998217493589	364.117335543095	314.175511964579	277.043092535072	246.270752777746	195.70913612502	168.102854357484	142.722367271741	89.5041059266947	
    730.266495461484	536.918175592217	426.295310861942	369.764885508864	328.878387857907	296.884055308305	233.111367401215	199.486893792689	175.18124008861	110.317421963392];
f = @(p) abs(p(1)*x.^p(2)./(y+p(3)).^p(4) - z);
F = @(p) reshape(f(p),[],1);
p0 = [1 1 1 1];
%format long
[p,fval] = fminimax(F,p0)

Local minimum possible. Constraints satisfied.

fminimax stopped because the size of the current search direction is less than
twice the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.

p = 1×4

1.0e+04 * 1.0124 0.0000 0.0016 0.0001

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

fval = 80×1

0.7496 6.6795 10.3501 12.7950 14.5778 19.2934 22.3030 22.3824 10.9381 3.8083

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

z_opt = p(1)*x.^p(2)./((y+p(3)).^p(4))

z_opt = 8×10

264.3605 205.9461 167.5595 140.5547 120.6030 105.3074 83.4902 68.7559 58.1882 39.2383 314.7674 245.2148 199.5088 167.3549 143.5989 125.3869 99.4097 81.8659 69.2832 46.7201 348.5988 271.5707 220.9522 185.3423 159.0330 138.8635 110.0943 90.6649 76.7298 51.7416 374.7856 291.9711 237.5501 199.2653 170.9796 149.2950 118.3646 97.4757 82.4937 55.6284 396.4449 308.8444 251.2784 210.7810 180.8607 157.9229 125.2050 103.1089 87.2611 58.8433 472.0369 367.7333 299.1909 250.9716 215.3463 188.0348 149.0784 122.7691 103.8996 70.0632 562.0424 437.8507 356.2390 298.8256 256.4074 223.8883 177.5039 146.1781 123.7107 83.4225 707.8841 551.4665 448.6777 376.3664 322.9413 281.9840 223.5636 184.1092 155.8117 105.0694

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

So, when I substitute them with MS Excel in the formula f = @(p) abs(p(1)*x.^p(2)./(y+p(3)).^p(4) - z); the fitted z values are very different than input z values.

The value you use for p(4) in your Excel-Sheet is wrong. And the values you use should include as many digits as possible - don't round after the second digit !

[Boyan]

Dear Torsten,

Everytning is clear now. Thank you again.

Now I have to solve a slighty different task. I have to minimize not the maximum absolute error |z_input - z_fitted|, but the maximal relative error |(z_input - z_fitted)/z_input|.

Could you help with the appropriate change of the code, which I have already have given, in order to minimize the parameter "relative error".

I have tried to ask this as a new question, with link to the current discussion as a base, but the system didn't accept the question. That is why I'm using current question, which is solved. Sorry for that.

Best regards,

Boyan

[Torsten]

x = [1 2 3 4 5 10 20 50].';
y = [5 10 15 20 25 30 40 50 60 90];
z = [265.110074796797	195.008035083193	156.258229529605	129.742584500194	107.062285017337	91.0540850617739	68.5080545479447	50.5038878827341	35.8057303140135	23.9109370743307	
    321.446939682946	249.023136173216	204.538587797874	172.392832029864	146.102087493907	125.302681110766	98.3496300763363	79.6443293949692	60.4259292015808	28.967085981117	
    358.948929784307	282.185871826562	232.988776882345	197.578293515688	169.230643083598	145.977264293244	115.850743148756	96.2843426347109	74.9936831153838	40.5561051536658	
    387.580526618466	306.229191510467	253.083199667196	215.390669679886	185.62217968629	160.805889976744	128.174360996678	107.804496871799	85.3098905504615	48.5261569465696	
    411.022720927616	325.178097515035	268.617658401558	229.174534346536	198.325915732664	172.4000748844	137.680141245939	116.580164412001	93.3004757115783	54.5677932289588	
    491.330306710826	385.749717456981	316.583948069254	271.810779459274	237.728742469111	208.93517622606	166.914170433146	142.97295566989	118.058749869664	72.5800929043424	
    584.345341008974	448.998217493589	364.117335543095	314.175511964579	277.043092535072	246.270752777746	195.70913612502	168.102854357484	142.722367271741	89.5041059266947	
    730.266495461484	536.918175592217	426.295310861942	369.764885508864	328.878387857907	296.884055308305	233.111367401215	199.486893792689	175.18124008861	110.317421963392];
% Relative error in percent
f = @(p) 100*abs((p(1)*x.^p(2)./(y+p(3)).^p(4))./z - 1);
F = @(p) reshape(f(p),[],1);
p0 = [1 1 1 1];
F(p0)

ans = 80×1

99.9371 99.8963 99.8607 99.8280 99.7973 99.6608 99.4296 98.8589 99.9534 99.9270

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

max(F(p0))

ans = 99.9646

[p,fval] = fminimax(F,p0,[],[],[],[],[],[],[],optimset('MaxFunEvals',100000,'MaxIter',100000))

Local minimum possible. Constraints satisfied.

fminimax stopped because the size of the current search direction is less than
twice the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.

p = 1×4

1.0e+05 * 2.1974 0.0000 0.0004 0.0000

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

fval = 80×1

21.5325 18.3483 16.2280 14.5580 13.1702 8.3527 2.7741 5.7864 13.5687 14.6030

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

max(fval)

ans = 21.5325

z_opt = p(1)*x.^p(2)./((y+p(3)).^p(4))

z_opt = 8×10

208.0253 168.5480 139.5108 117.5053 100.4143 86.8646 66.9815 53.3207 43.5156 26.3106 262.4670 212.6582 176.0217 148.2572 126.6934 109.5977 84.5111 67.2752 54.9039 33.1963 300.6987 243.6347 201.6616 169.8528 145.1479 125.5621 96.8212 77.0747 62.9014 38.0317 331.1564 268.3124 222.0878 187.0572 159.8499 138.2802 106.6282 84.8815 69.2727 41.8839 356.8901 289.1625 239.3459 201.5931 172.2716 149.0258 114.9142 91.4776 74.6557 45.1387 450.2908 364.8384 301.9844 254.3514 217.3563 188.0269 144.9880 115.4179 94.1937 56.9518 568.1350 460.3191 381.0158 320.9170 274.2400 237.2348 182.9324 145.6235 118.8448 71.8565 772.5230 625.9201 518.0873 436.3677 372.8985 322.5806 248.7427 198.0120 161.5995 97.7070

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

[Boyan]

Thank you, Torsten!

[Boyan]

Hi, Torsten,

I have another issue, but I have posted like a new question in the link below.

Because in this question too much subquestions appeared.

If you are still interested in, so here is the link:

Where I make a mistake with fminimax? - MATLAB Answers - MATLAB Central

Answer 2

[Torsten]

Use "fminimax". It minimizes the maximum error and is the nonlinear equivalent to minimizing the Inf-norm in linear regression problems.

https://uk.mathworks.com/help/optim/ug/fminimax.html

[Boyan]

Dear colleague,

Could you give me more information about your advice?