Can't find explicit solution

Hi Guys,

MATLAB says it can't find an explicit solution for following equation. What i want is to estimate parameters, by doing partial derivations, they should be equal to zero and then solve them.

I tried it for equation after equation or all at once with solve and vpasolve. Nothing worked, but i know, that there must be a solution (from another programm).

Is there a formal error in my code? Or do you have any idea what i could try?

syms A B b y0 real

x11=-1;x21=-1;x31=1;x41=1;
x1=[x11 x21 x31 x41];

x12=-1;x22=1;x32=-1;x42=1;
x2=[x12 x22 x32 x42];

t1=27;             t2=50;             t3=25;             t4=55;

Versuchsdaten=[t1 t2 t3 t4];

%log. Likelihood-Funktion%
logL=0;
for i=1:length(Versuchsdaten)
    logL=logL+log(b/(exp((y0+A*x1(i)+B*x2(i)))))+(b-1)*log((Versuchsdaten(i))/(exp((y0+A*x1(i)+B*x2(i)))))-((Versuchsdaten(i))/exp((y0+A*x1(i)+B*x2(i))))^(b);
end

logL

%differentielle Ableitungen%

dLdb=diff(logL,b);
dLdy0=diff(logL,y0);
dLdA=diff(logL,A);
dLdB=diff(logL,B);

Try 1:

eqn=dLdb==0;
b_hat=solve(eqn,b)

When MATLAB would find a solution for this equation i would youse subs to integrate this solution for b into dLby0 and so one till i get a numerical solution.

Try 2:

[b_hat, y0_hat, A_hat, B_hat] = vpasolve([dLdb==0, dLdy0==0, dLdA==0, dLdB==0], [b, y0, A, B])

[Torsten]

Why don't you simply use an optimizer like "fminsearch" or "fmincon" to determine your parameters ?

Best wishes

Torsten.

[Ronald Singer]

I will try some ways to optimize my function. Just were so focused about solving this problem directly that I forgot about optimazation. Thanks once again, Torsten.

[Walter Roberson]

I tried in Maple; when I set the calculation precision to be low, it suggests

A = -7.39041359316*10^7, B = 4.30445800579, b = 0., y0 = 7.39041458029*10^7

However, it fails when increasing the calculation precision. I went back and substituted b = 0 into the equation, but that generates a division by 0 error. That suggests that possibly there might not be a direct solution in reals, only a solution in limits.

Answer 1

[Karan Gill]

You can check Walter's suggestion that the solution is in limits as

>> limit(dLdb,b,0)
ans =
NaN
>> limit(dLdy0,b,0)
ans =
0
>> limit(dLdA,b,0)
ans =
0
>> limit(dLdB,b,0)
ans =
0

Continuing with Walter's suggestion, I reduced the precision of vpasolve using digits. Then I specified the range in vpasolve to limit b to about 0. Read the doc for digits and vpasolve.

>> digitsOld = digits(4); % change precision
eqns = [dLdb==0, dLdy0==0, dLdA==0, dLdB==0];
vars = [b, y0, A, B];
ranges = [-10^(-5) 10^(-5); NaN NaN; NaN NaN; NaN NaN];
[b_hat, y0_hat, A_hat, B_hat] = vpasolve(eqns, vars, ranges)

b_hat =
-7.358e-6
y0_hat =
98.4
A_hat =
61.14
B_hat =
23.11

vpasolve fails at b = 0. Walter, thanks again for the insight :) And ake sure you change digits back!

digits(digitsOld)

Answer 2

[Walter Roberson]

There are at least two solutions, in two different clusters.

There is at least one solution near

A = 0.00458728466696763, B = 0.351160874946993, b = -27.8555787371441, y0 = 3.58721372891737

The exact position would take more work. The minimum residue appears to be along a line relating A and B, and along a line relating b and y0, but there is no obvious correlation between A and b or y0 or between B and b or y0. Hmmmm, now I see planes of relationships when I view the variables 3 at a time. I am not sure what to make of it all; the situation is as if there is a hyperplane of solutions. It might be the case that there is a point solution in theory but that in practice the round-off errors are causing what appear to be random fuzz around it.

There is a second solution (or cluster of solutions) near

A = 0.00458728466706846433, B = 0.351160874946970347, b = 27.8555787443984215, y0 = 3.62982071185302235

This is pretty much the same A and B, and a b that might plausibly be the negative of the other b, and a y0 that might plausibly be the same as the other one.

When I examine mathematically, it is not obvious to me that substituting -b for b in the equations would produce the same solutions, so at the moment I am not sure what is happening in that regards: the above solutions were produced experimentally by residue minimization, not based on theory.

[Ronald Singer]

Hey Walter, thank you for your answers. This is pretty much the exact solution i need; or in other words, the other Programm has calculated (the second one, because b needs to be >0). You're saying that you found this solution experimentally. So not in MATLAB? And do you see any chance to get this solutions directly without optimization? For my thesis it would be great, if i can calculate it in both ways.

[Walter Roberson]

I have been doing more testing, and finding a number of potential solution clusters. There is something suspicious in the area near -81003023.5509659648, -83630380.8575624228, -9.50685400681583273e-09, -99921437.7476787567

As I search, it keeps finding points with fairly low residue, and the optimal keeps moving towards lower A and lower B.

At the moment the search range is

A = -82060697.1417102218 .. -79988127.5028114915, B = -83751616.8654009104 .. -80982369.0845304728, b = -9.63048143368570503e-09 .. -9.49956191432576237e-09, y0 = -99990723.4130454659 .. -99425994.6932563037

and I have found quite a number of points with low residue within that range -- but there are also many points with fairly high residue as well. I have not found any good pattern as to which points are good and which are not. I tested against the possibility that the limit of A approaching -inf or B approaching -inf would be 0, but if it is then it is not obviously the case.

[Walter Roberson]

That B value turns out to be (1/4)*ln(110/27) which can be written other ways such as (1/4)*ln(2)+(1/4)*ln(5)+(1/4)*ln(11)-3*ln(3)*(1/4)

As to why it turns out to be that, I am still trying to reproduce that. It popped up for me earlier and now when I go back to try to re-create it, I am having reproducibility difficulties, so I am not sure how to prove it quite yet. But the value works out. The line of investigation was approximately: solve dLdb for A, substitute that A into dldy0, solve that for B, and poke at that result until your symbolic engine realizes that the various b cancel out giving you a constant result.