Solving system of non-linear equations in Mathematica

0
votes

I am trying to numerically solve the below system of six equations (g0-g5) for a0-a5 in Mathematica. I am no expert in Mathematica and am not entirely sure how to do this. 

f[x_, y_] := Exp[a0 - 1 + a1*x + a2*y + a3*x*x + a4*y*y + a5*x*y]
g0[x_, y_] := Integrate[f[x, y],{y,-Infinity,Infinity},{x,-Infinity,Infinity}] - 1
g1[x_, y_] := Integrate[x*f[x, y],{y,-Infinity,Infinity},{x,-Infinity,Infinity}]
g2[x_, y_] := Integrate[y*f[x, y],{y,-Infinity,Infinity},{x,-Infinity,Infinity}]
g3[x_, y_] := Integrate[x*x*f[x, y],{y,-Infinity,Infinity},{x,-Infinity,Infinity}] - 1
g4[x_, y_] := Integrate[y*y*f[x, y],{y,-Infinity,Infinity},{x,-Infinity,Infinity}] - 1
g5[x_, y_] := Integrate[x*y*f[x, y],{y,-Infinity,Infinity},{x,-Infinity,Infinity}]

I have, however, spent considerable time trying to get NSolve and FindRoot to yield a solution. Here is that code:

NSolve[{g0[x, y]==0, g1[x, y]==0, g2[x, y]==0, g3[x, y]==0, g4[x, y]==0, 
 g5[x, y]==0}, {a0, a1, a2, a3, a4, a5}, Reals]

FindRoot[{g0[x, y]==0, g1[x, y]==0, g2[x, y]==0, g3[x, y]==0, g4[x, y]==0,
g5[x, y]==0}, {{a0,1}, {a1,1}, {a2,1}, {a3,1}, {a4,1}, {a5,1}}]

One additional piece of information I can offer is that the resulting solution for f(x,y) should be equivalent to the bivariate standard normal density. Any help would be much appreciated. This is my first post on SO, so please let me know if any additional information is necessary.

2
wolfram-mathematica
At least in V9, each of those appear to integrate just fine. They are a little slow. Each returns a ConditionalExpression[ answer, IfThisIsTrue ] and that depends on the sign of something inside a square root in the answer. So I did not Numerically solve those, I just evaluated your code and evaluated things like g0[ x,y ] and got the result. If you want to edit your post to include some numbers for a0-a5 then I could show how to get a numerical result. If you can edit your post to show the NSolve and FindRoot you were using then perhaps we can see the problem. - Bill
you don't show equation to solve, are you after gi[x,y]==0? - agentp
i'll also point out the expressions on the right of all the g definitions are not functions of x,y. You probably ought to be writing g0[a0_,a1_,a2_,a3_,a4_,a5_]:= .. - agentp
I have edited the post to include my NSolve and FindRoot code. Thank you in advance for your feedback. - Heath Henderson

2 Answers

1
votes

I am astonished. I never expected it to finish. But if you subtract off all the time for it to do the integrals then Reduce finishes in the blink of an eye.

f[x_, y_] := Exp[a0 - 1 + a1*x + a2*y + a3*x*x + a4*y*y + a5*x*y];
g0[x_, y_] := Integrate[f[x, y], {y, -Infinity, Infinity}, {x, -Infinity, Infinity}]-1;
g1[x_, y_] := Integrate[x*f[x, y], {y, -Infinity, Infinity}, {x, -Infinity, Infinity}];
g2[x_, y_] := Integrate[y*f[x, y], {y, -Infinity, Infinity}, {x, -Infinity, Infinity}];
g3[x_, y_] := Integrate[x*x*f[x, y], {y, -Infinity, Infinity}, {x, -Infinity, Infinity}]-1;
g4[x_, y_] := Integrate[y*y*f[x, y], {y, -Infinity, Infinity}, {x, -Infinity, Infinity}]-1;
g5[x_, y_] := Integrate[x*y*f[x, y], {y, -Infinity, Infinity}, {x, -Infinity, Infinity}];
Reduce[Simplify[{g0[x,y]==0, g1[x,y]==0, g2[x,y]==0, g3[x,y]==0, g4[x,y]==0, g5[x,y]==0},
 Re[4 a4-a5^2/a3]<0], {a0,a1,a2,a3,a4,a5}]

gives you

C[1] \[Element] Integers && a0==1+2I\[Pi] C[1]-Log[2]-Log[\[Pi]] &&
   a1==0 && a2==0 && a3== -(1/2) && a4== -(1/2) && a5==0

Note: That gives Simplify one assumption which you should verify is justified. That assumption allows it to turn all your ConditionalExpression into presumably valid expressions for your problem. I got that assumption by looking at each of the results returned from Integrate and seeing that they all depended on that for the result to be valid.

0
votes

Here is how to formulate this numerically:

 f[x_, y_, a0_, a1_, a2_, a3_, a4_, a5_] :=
  Exp[a0 - 1 + a1*x + a2*y + a3*x*x + a4*y*y + a5*x*y]
 g0[a0_?NumericQ, a1_, a2_, a3_, a4_, a5_] :=
   Integrate[
   f[x, y, a0, a1, a2, a3, a4, a5], {y, -Infinity, Infinity},
      {x, -Infinity, Infinity}] - 1
 g1[a0_?NumericQ, a1_, a2_, a3_, a4_, a5_] := 
     Integrate[
     x*f[x, y, a0, a1, a2, a3, a4, a5], {y, -Infinity, Infinity},
       {x, -Infinity, Infinity}]
 g2[a0_?NumericQ, a1_, a2_, a3_, a4_, a5_] := 
      Integrate[
        y*f[x, y, a0, a1, a2, a3, a4, a5], {y, -Infinity,  Infinity},  
          {x, -Infinity, Infinity}]
 g3[a0_?NumericQ, a1_, a2_, a3_, a4_, a5_] := 
      Integrate[
       x*x*f[x, y, a0, a1, a2, a3, a4, a5], {y, -Infinity, 
         Infinity}, {x, -Infinity, Infinity}] - 1
 g4[a0_?NumericQ, a1_, a2_, a3_, a4_, a5_] := 
        Integrate[
         y*y*f[x, y, a0, a1, a2, a3, a4, a5], {y, -Infinity, 
           Infinity}, {x, -Infinity, Infinity}] - 1
 g5[a0_?NumericQ, a1_, a2_, a3_, a4_, a5_] := 
      Integrate[
       x*y*f[x, y, a0, a1, a2, a3, a4, a5], {y, -Infinity, 
             Infinity}, {x, -Infinity, Infinity}]

 FindRoot[  {
  g0[a0, a1, a2, a3, a4, a5] == 0, 
  g1[a0, a1, a2, a3, a4, a5] == 0,
  g2[a0, a1, a2, a3, a4, a5] == 0,
  g3[a0, a1, a2, a3, a4, a5] == 0,
  g4[a0, a1, a2, a3, a4, a5] == 0,
  g5[a0, a1, a2, a3, a4, a5] == 0}  ,
      {{a0, -.8379}, {a1, 0}, {a2, 0}, {a3, -.501},
       {a4, -.499}, {a5, 0}}]

Note I've put in an initial guess very close to the known solution (thanks @Bill) and it still takes a very long time to find the answer.

{a0 -> -0.837388 - 1.4099*10^-29 I, a1 -> -6.35273*10^-22 + 7.19577*10^-46 I, a2 -> -1.27815*10^-20 + 6.00264*10^-38 I, a3 -> -0.500489 + 1.41128*10^-29 I, a4 -> -0.5 - 7.13595*10^-44 I, a5 -> -5.55356*10^-28 - 9.23563*10^-47 I}

 Chop@%

{a0 -> -0.837388, a1 -> 0, a2 -> 0, a3 -> -0.500489, a4 -> -0.5, a5 -> 0}