1
$\begingroup$

I need to solve, numerically, a very large system of equations (Logit-type). The system has one parameter $d$ which is numerically defined in order to find the solution.

My code works, but it requires 11 minutes to run. I want:

  1. To increase the efficiency of my code (I have to estimate also a different "pay" function which will generate around 10.000 equations instead of the actual 250).
  2. I want Mathematica to run my program across many values of $d$ (user-defined).
  3. I want Mathematica to save solutions for each $d$ in a dataset which can be exported to Stata for a maximum likelihood estimation (my task will be to find the parameter $d$ which better fit my data). Or, maybe, I can import my observations in Mathematica for the maximum likelihood estimation? What's best?

Here is the code:

ClearAll[qmax, qmax2, pay]
qmax[h_] := Floor[8 Sqrt[2] Sqrt[h]];
qmax2[h_, q_] := Min[qmax[h], (q - 1)];
pay[h_, q_] := (32 - q^2/(4 h))*(1/10)*Sum[
   Sum[Symbol["t" <> ToString[t] <> "q" <> ToString[s]], {s, 0, qmax2[t, q]}]+ 
   If[q <= qmax[t], 0.5*Symbol["t" <> ToString[t] <> "q" <> ToString[q]], 0],
   {t, 1, 10}];

equations = Flatten@Table[
  Table[
  Symbol["t" <> ToString[h] <> "q" <> ToString[q]] == 
  Exp[d*pay[h,q]]/Sum[Exp[d*pay[h, s]], {s, 0, qmax[h]}],
  {q, 0, qmax[h]}],
  {h, 1, 10}];

vars = Flatten[Table[
  {Symbol["t" <> ToString[h] <> "q" <> ToString[q]], 0.01}, 
  {h,1, 10}, {q, 0, qmax[h]}], 1];

d = 1;
FindRoot[equations, vars]
$\endgroup$
3
  • $\begingroup$ What have you tried so far to work towards your goals? If you encounter specific problems, people here will be able to help. However, as it is currently written, yours is not a request for help, but unfortunately it rather reads like a job advertisement. The issue is that any contributions towards your problem would likely only help this specific case, and not be of general value to this community. $\endgroup$
    – MarcoB
    Commented Jun 27, 2018 at 15:59
  • 1
    $\begingroup$ A general way to speed up FindRoot is to start with a good guess. Do you know more about the solution than 0.01 for each component? When you do get to varying the parameter d, try using the previous d's answer as an initial guess for the next root -- or get more sophisticated with linear or quadratic extrapolation. $\endgroup$
    – Chris K
    Commented Jun 27, 2018 at 16:40
  • $\begingroup$ This equations = Total[Abs[Flatten@Table[Table[Symbol["t"<>ToString[h]<>"q" <>ToString[q]]*Sum[Exp[d*pay[h, s]], {s,0,qmax[h]}]-Exp[d*pay[h,q]], {q,0,qmax[h]}], {h,1,10}]]] //. {0.5*xxx_+0.5*yyy_->0.5*(xxx+yyy), 2*xxx_+2*yyy_->2*(xxx+yyy)}; vars=Flatten[Table[{Symbol["t"<>ToString[h]<>"q"<>ToString[q]]}, {h,1,10}, {q,0,qmax[h]}]]; d=1; AbsoluteTiming[NMinimize[equations, vars]] takes about 30 seconds to find an approximate solution. That only times the solve. Can you confirm that this works? Can you show a problem with 1000 equations? Can you give a list of example d values? $\endgroup$
    – Bill
    Commented Sep 4, 2018 at 2:09

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.