FindRoot, a big Jacobian, and bivariate probit model

• To: mathgroup at smc.vnet.net
• Subject: [mg14107] FindRoot, a big Jacobian, and bivariate probit model
• From: "ELLIS, Luci (FS)" <EllisL at rba.gov.au>
• Date: Fri, 25 Sep 1998 03:15:31 -0400
• Sender: owner-wri-mathgroup at wolfram.com

```> Dear MathGroup,
> I am currently trying to write Mathematica code for estimating a bivariate
> probit model by Maximum Likelihood and have run into some problems (to say
> the least).  I am not getting any error messages, but even on a toy
> problem
> (3 explanatory variables, 20 data points), Mathematica fails to give an
> answer in any reasonable amount of time.
>
> (For those familiar with FindRoot but not probit models, the short
> explanation is that they are models where the dependent variable takes the
> values 0 or 1, that is,  y = 1 if b.X + e > 0 and 0 otherwise, where e is
> a
> standard normally-distributed disturbance.  In the bivariate case, you
> have
> two y-variables, and the errors for each of them may be correlated, ie
> distributed bivariate normal)
>
> I have analytic expressions for the first-order conditions and the
> Jacobian
> (second derivatives of the log-likelihood, in that case).  However, the
> Jacobian takes a L-O-N-G time to evaluate because it involves the CDF of a
> bivariate normal distribution, and the parameters are not known, which
> means
> that they have to be substituted into this enormous expression.  I have
> already tried NOT giving Mathematica the Jacobian explicitly, but then it
> asks me for one because it can't get the analytic derivatives of these
> expressions.  How can I speed this up so that the code is usable?
>
> Thanks in advance,
> Luci Ellis
>
> The function currently looks like this: (note I haven't Modulised all the
> local variables yet, so I can look at them.  You'll have to copy and paste
> this into Mathematica to see it properly.)
> (****start from next line ****)
> Cell[BoxData[
>     RowBox[{
>       RowBox[{
>         RowBox[{"EstimateBivariateProbitModel", "[",
>           RowBox[{
>             RowBox[{"y1", ":",
>               RowBox[{"{",
>                 RowBox[{
>                   RowBox[{"(",
>                     RowBox[{"0", "|", "1"}], ")"}], ".."}], "}"}]}], ",",
>             RowBox[{"y2", ":",
>               RowBox[{"{",
>                 RowBox[{
>                   RowBox[{"(",
>                     RowBox[{"0", "|", "1"}], ")"}], ".."}], "}"}]}], ",",
>             RowBox[{"x_", "?", "MatrixQ"}], ",", "opts___Rule"}], "]"}],
>         "/;",
>         RowBox[{"(",
>           RowBox[{
>             RowBox[{"(",
>               RowBox[{
>                 RowBox[{"Length", "[", "y1", "]"}], "==",
>                 RowBox[{"Length", "[", "y2", "]"}], "==",
>                 RowBox[{"Length", "[", "x", "]"}]}], ")"}], " ", "&&", "
> ",
>             RowBox[{"And", " ", "@@",
>               RowBox[{"(", " ",
>                 RowBox[{"NumericQ", "/@",
>                   RowBox[{"Flatten", "[", "x", "]"}]}], ")"}]}]}],
> ")"}]}],
>       ":=", "\n", "\t\t",
>       RowBox[{"Module", "[",
>         RowBox[{
>           RowBox[{"{", "}"}], ",", "\n", "\t",
>           RowBox[{"With", "[",
>             RowBox[{
>               RowBox[{"{",
>                 RowBox[{
>                   RowBox[{"\[CapitalBeta]1", "=",
>                     RowBox[{"ToExpression", "[",
>                       RowBox[{
>                         RowBox[{
>                           RowBox[{"\"\<\[Beta]\>\"", "<>", "#"}], " ",
> "&"}],
>                         " ", "/@", " ",
>                         RowBox[{"(",
>                           RowBox[{"ToString", " ", "/@", " ",
>                             RowBox[{"Range", "[",
>                               RowBox[{
>                                 RowBox[{"Dimensions", "[", "x", "]"}],
>                                 "\[LeftDoubleBracket]", "2",
>                                 "\[RightDoubleBracket]"}], "]"}]}],
> ")"}]}],
>
>                       "]"}]}], ",",
>                   RowBox[{"\[CapitalBeta]2", "=",
>                     RowBox[{"ToExpression", "[",
>                       RowBox[{
>                         RowBox[{
>                           RowBox[{"\"\<\[Gamma]\>\"", "<>", "#"}], " ",
>                           "&"}], " ", "/@", " ",
>                         RowBox[{"(",
>                           RowBox[{"ToString", " ", "/@", " ",
>                             RowBox[{"Range", "[",
>                               RowBox[{
>                                 RowBox[{"Dimensions", "[", "x", "]"}],
>                                 "\[LeftDoubleBracket]", "2",
>                                 "\[RightDoubleBracket]"}], "]"}]}],
> ")"}]}],
>
>                       "]"}]}], ",", " ",
>                   RowBox[{"q1", "=",
>                     RowBox[{
>                       RowBox[{"2", " ", "y1"}], "-", "1"}]}], ",",
>                   RowBox[{"q2", " ", "=",
>                     RowBox[{
>                       RowBox[{"2", " ", "y2"}], "-", "1"}]}]}], "}"}],
> ",",
>               "\n", "\t\t",
>               RowBox[{
>                 RowBox[{"z1", "=",
>                   RowBox[{"x", ".", "\[CapitalBeta]1"}]}], ";", "\t",
>                 RowBox[{"z2", "=",
>                   RowBox[{"x", ".", "\[CapitalBeta]2"}]}], ";",
>                 RowBox[{"\[Rho]star", "=", " ",
>                   RowBox[{"q1", "*", "q2", "*", "\[Rho]"}]}], ";",
>                 RowBox[{"w1", "=",
>                   RowBox[{"q1", "*", "z1"}]}], ";", " ",
>                 RowBox[{"w2", "=",
>                   RowBox[{"q2", "*", "z2"}]}], ";", "\n", "\t\t\t",
>                 RowBox[{"g1", "=",
>                   RowBox[{"Evaluate", "[",
>                     RowBox[{
>                       RowBox[{"\[Phi]", "[", "w1", "]"}], "*",
>                       RowBox[{"\[CapitalPhi]", "[",
>                         FractionBox[
>                           RowBox[{"w2", "-",
>                             RowBox[{"\[Rho]star", " ", "w1"}]}],
>                           SqrtBox[
>                             RowBox[{"(",
>                               RowBox[{"1", "-",
>                                 SuperscriptBox["\[Rho]star", "2"]}],
> ")"}]]],
>                         "]"}]}], "]"}]}], ";", "\n", "\t\t\t",
>                 RowBox[{"g2", "=",
>                   RowBox[{"Evaluate", "[",
>                     RowBox[{
>                       RowBox[{"\[Phi]", "[", "w2", "]"}], "*",
>                       RowBox[{"\[CapitalPhi]", "[",
>                         FractionBox[
>                           RowBox[{"w1", "-",
>                             RowBox[{"\[Rho]star", " ", "w2"}]}],
>                           SqrtBox[
>                             RowBox[{"(",
>                               RowBox[{"1", "-",
>                                 SuperscriptBox["\[Rho]star", "2"]}],
> ")"}]]],
>                         "]"}]}], "]"}]}], ";", "\t\t\t", "\n", "\t\t\t",
>                 RowBox[{"triple", "=",
>                   RowBox[{"Transpose", "[",
>                     RowBox[{"{",
>                       RowBox[{"w1", ",", "w2", ",", "\[Rho]star"}], "}"}],
>
>                     "]"}]}], ";", "\n", "\t\t\t",
>                 RowBox[{"\[Phi]2", "=",
>                   RowBox[{
>                     RowBox[{
>                       RowBox[{"PDF", "[",
>                         RowBox[{
>                           RowBox[{"MultinormalDistribution", "[",
>                             RowBox[{
>                               RowBox[{"{",
>                                 RowBox[{"0", ",", "0"}], "}"}], ",",
>                               TagBox[
>                                 RowBox[{"{",
>                                   RowBox[{
>                                     RowBox[{"{",
>                                       RowBox[{"1", ",",
>                                         RowBox[{
>                                         StyleBox["#",
>                                         StyleBoxAutoDelete->True,
>                                         PrintPrecision->1],
>                                         "\[LeftDoubleBracket]", "3",
>                                         "\[RightDoubleBracket]"}]}],
> "}"}],
>                                     ",",
>                                     RowBox[{"{",
>                                       RowBox[{
>                                         RowBox[{
>                                         StyleBox["#",
>                                         StyleBoxAutoDelete->True,
>                                         PrintPrecision->1],
>                                         "\[LeftDoubleBracket]", "3",
>                                         "\[RightDoubleBracket]"}], ",",
>                                         "1"}], "}"}]}], "}"}],
>                                 Short]}], "]"}], ",",
>                           RowBox[{"{",
>                             RowBox[{
>                               RowBox[{
>                               "#", "\[LeftDoubleBracket]", "1",
>                                 "\[RightDoubleBracket]"}], ",",
>                               RowBox[{
>                               "#", "\[LeftDoubleBracket]", "2",
>                                 "\[RightDoubleBracket]"}]}], "}"}]}],
> "]"}],
>
>                       "&"}], " ", "/@", " ", "triple"}]}], ";", "\n",
>                 "\t\t\t",
>                 RowBox[{"\[CapitalPhi]2", "=",
>                   RowBox[{
>                     RowBox[{
>                       RowBox[{"CDF", "[",
>                         RowBox[{
>                           RowBox[{"MultinormalDistribution", "[",
>                             RowBox[{
>                               RowBox[{"{",
>                                 RowBox[{"0", ",", "0"}], "}"}], ",",
>                               TagBox[
>                                 RowBox[{"{",
>                                   RowBox[{
>                                     RowBox[{"{",
>                                       RowBox[{"1", ",",
>                                         RowBox[{
>                                         StyleBox["#",
>                                         StyleBoxAutoDelete->True,
>                                         PrintPrecision->1],
>                                         "\[LeftDoubleBracket]", "3",
>                                         "\[RightDoubleBracket]"}]}],
> "}"}],
>                                     ",",
>                                     RowBox[{"{",
>                                       RowBox[{
>                                         RowBox[{
>                                         StyleBox["#",
>                                         StyleBoxAutoDelete->True,
>                                         PrintPrecision->1],
>                                         "\[LeftDoubleBracket]", "3",
>                                         "\[RightDoubleBracket]"}], ",",
>                                         "1"}], "}"}]}], "}"}],
>                                 Short]}], "]"}], ",",
>                           RowBox[{"{",
>                             RowBox[{
>                               RowBox[{
>                               "#", "\[LeftDoubleBracket]", "1",
>                                 "\[RightDoubleBracket]"}], ",",
>                               RowBox[{
>                               "#", "\[LeftDoubleBracket]", "2",
>                                 "\[RightDoubleBracket]"}]}], "}"}]}],
> "]"}],
>
>                       "&"}], " ", "/@", " ", "triple"}]}], ";", "\n",
>                 "\t\t\t",
>                 RowBox[{"\[CapitalPhi]2sq", "=",
>                   RowBox[{"\[CapitalPhi]2", "^", "2"}]}], ";", "\n",
>                 "\t\t\t",
>                 RowBox[{"\[Phi]\[CapitalPhi]", "=",
>                   RowBox[{"\[Phi]2", "/", "\[CapitalPhi]2"}]}], ";", "\n",
>
>                 "\t\t\t",
>                 RowBox[{"\[Delta]", "=",
>                   RowBox[{"1", "/",
>                     SqrtBox[
>                       RowBox[{"(",
>                         RowBox[{"1", "-",
>                           RowBox[{"\[Rho]star", "^", "2"}]}], ")"}]]}]}],
>                 ";", "\n", "\t\t\t",
>                 RowBox[{"v1", "=",
>                   RowBox[{"\[Delta]",
>                     RowBox[{"(",
>                       RowBox[{"w2", "-",
>                         RowBox[{"\[Rho]star", " ", "w1"}]}], ")"}]}]}],
> ";",
>
>                 "\n", "\t\t\t",
>                 RowBox[{"v2", "=",
>                   RowBox[{"\[Delta]", " ",
>                     RowBox[{"(",
>                       RowBox[{"w1", "-",
>                         RowBox[{"\[Rho]star", " ", "w2"}]}], ")"}]}]}],
> ";",
>
>                 "\n", "\t\t\t",
>                 RowBox[{"focs", "=",
>                   RowBox[{"Evaluate", "[",
>                     RowBox[{"Thread", "[",
>                       RowBox[{
>                         RowBox[{"Join", "[",
>                           RowBox[{
>                             RowBox[{
>                               RowBox[{"(",
>                                 RowBox[{"q1", " ",
>                                   RowBox[{"g1", "/", "\[CapitalPhi]2"}]}],
>
>                                 ")"}], ".", "x"}], ",",
>                             RowBox[{
>                               RowBox[{"(",
>                                 RowBox[{"q2", " ",
>                                   RowBox[{"g2", "/", "\[CapitalPhi]2"}]}],
>
>                                 ")"}], ".", "x"}], ",",
>                             RowBox[{"{",
>                               RowBox[{
>                                 RowBox[{"(",
>                                   RowBox[{"q1", " ",
>                                     RowBox[{"q2", "/",
> "\[CapitalPhi]2"}]}],
>
>                                   ")"}], ".", "\[Phi]2"}], " ", "}"}]}],
>                           "]"}], "==", "0"}], "]"}], "]"}]}], ";", "\n",
>                 "\t\t\t",
>                 RowBox[{"\[ScriptD]b1b1", "=",
>                   RowBox[{
>                     RowBox[{"Transpose", "[", "x", "]"}], ".",
>                     RowBox[{"(",
>                       RowBox[{"x", "*",
>                         RowBox[{"(", " ",
>                           RowBox[{
>                             RowBox[{"(",
>                               RowBox[{
>                                 RowBox[{"-", "w1"}], " ",
>                                 RowBox[{"g1", "/", "\[CapitalPhi]2"}]}],
>                               ")"}], "-",
>                             RowBox[{"(",
>                               RowBox[{
>                               "\[Rho]star", "  ", "\[Phi]\[CapitalPhi]"}],
>
>                               ")"}], " ", "-",
>                             RowBox[{"(",
>                               RowBox[{
>                                 RowBox[{"g1", "^", "2"}], "/",
>                                 "\[CapitalPhi]2sq"}], ")"}]}], ")"}]}],
>                       ")"}]}]}], ";", "\n", "\t\t\t",
>                 RowBox[{"\[ScriptD]b1b2", "=",
>                   RowBox[{
>                     RowBox[{"Transpose", "[",
>                       RowBox[{"q1", " ", "q2", " ", "x"}], "]"}], ".",
>                     RowBox[{"(",
>                       RowBox[{"x", "*", " ",
>                         RowBox[{"(",
>                           RowBox[{"\[Phi]\[CapitalPhi]", "-",
>                             RowBox[{"(",
>                               RowBox[{"g1", " ",
>                                 RowBox[{
>                                 "g2", " ", "/", "\[CapitalPhi]2sq"}]}],
>                               ")"}]}], ")"}]}], ")"}]}]}], ";", "\n",
>                 "\t\t\t",
>                 RowBox[{"\[ScriptD]b2b2", "=",
>                   RowBox[{
>                     RowBox[{"Transpose", "[", "x", "]"}], ".",
>                     RowBox[{"(",
>                       RowBox[{"x", "*",
>                         RowBox[{"(", " ",
>                           RowBox[{
>                             RowBox[{"(",
>                               RowBox[{
>                                 RowBox[{"-", "w1"}], " ",
>                                 RowBox[{"g1", "/", "\[CapitalPhi]2"}]}],
>                               ")"}], "-",
>                             RowBox[{"(",
>                               RowBox[{
>                               "\[Rho]star", "  ", "\[Phi]\[CapitalPhi]"}],
>
>                               ")"}], " ", "-",
>                             RowBox[{"(",
>                               RowBox[{
>                                 RowBox[{"g1", "^", "2"}], "/",
>                                 "\[CapitalPhi]2sq"}], ")"}]}], ")"}]}],
>                       ")"}]}]}], ";", "\n", "\t\t\t",
>                 RowBox[{"\[ScriptD]b1b\[Rho]", "=",
>                   RowBox[{"{",
>                     RowBox[{
>                       RowBox[{"Transpose", "[",
>                         RowBox[{"q2", "  ", "x"}], "]"}], ".",
>                       RowBox[{"(",
>                         RowBox[{"\[Phi]\[CapitalPhi]",
>                           RowBox[{"(",
>                             RowBox[{
>                               RowBox[{
>                               "\[Rho]star", " ", "\[Delta]", " ", "v1"}],
>                               "-", "w1", " ", "-",
>                               RowBox[{"g1", "/", "\[CapitalPhi]2"}]}],
>                             ")"}]}], ")"}]}], "}"}]}], ";", "\n",
> "\t\t\t",
>                 RowBox[{"\[ScriptD]b2b\[Rho]", "=",
>                   RowBox[{"{",
>                     RowBox[{
>                       RowBox[{"Transpose", "[",
>                         RowBox[{"q1", "  ", "x"}], "]"}], ".",
>                       RowBox[{"(",
>                         RowBox[{"\[Phi]\[CapitalPhi]",
>                           RowBox[{"(",
>                             RowBox[{
>                               RowBox[{
>                               "\[Rho]star", " ", "\[Delta]", " ", "v2"}],
>                               "-", "w2", "-",
>                               RowBox[{"g2", "/", "\[CapitalPhi]2"}]}],
>                             ")"}]}], ")"}]}], "}"}]}], ";", "\n",
> "\t\t\t",
>                 RowBox[{"\[ScriptD]\[Rho]\[Rho]", "=",
>                   RowBox[{"{",
>                     RowBox[{"{",
>                       RowBox[{"\[Phi]\[CapitalPhi]", ".",
>                         RowBox[{"(",
>                           RowBox[{
>                             RowBox[{
>                               RowBox[{"\[Delta]", "^", "2"}], "  ",
>                               "\[Rho]star", " ",
>                               RowBox[{"(",
>                                 RowBox[{"1", "-",
>                                   RowBox[{
>                                     RowBox[{"\[Delta]", "^", "2"}], " ",
>                                     RowBox[{"(",
>                                       RowBox[{
>                                         RowBox[{"w1", "^", "2"}], "+",
>                                         RowBox[{"w2", "^", "2"}], "-",
>                                         RowBox[{
>                                         "2", " ", "\[Rho]star", " ", "w1",
>
>                                         " ", "w2"}]}], ")"}]}]}], ")"}]}],
>
>                             "+",
>                             RowBox[{"\[Delta]", " ", "w1", " ", "w2"}], "
> ",
>
>                             "-", "\[Phi]\[CapitalPhi]"}], ")"}]}], "}"}],
>                     "}"}]}], ";", "\n", "\t\t\t",
>                 RowBox[{"jac", "=",
>                   RowBox[{"BlockMatrix", "[",
>                     RowBox[{"{",
>                       RowBox[{
>                         RowBox[{"{",
>                           RowBox[{
>                           "\[ScriptD]b1b1", ",", "\[ScriptD]b1b2", ",",
>                             RowBox[{
>                             "Transpose", "[", "\[ScriptD]b1b\[Rho]",
> "]"}]}],
>                           "}"}], ",",
>                         RowBox[{"{",
>                           RowBox[{
>                             RowBox[{
>                             "Transpose", "[", "\[ScriptD]b1b2", "]"}],
> ",",
>                             "\[ScriptD]b2b2", ",",
>                             RowBox[{
>                             "Transpose", "[", "\[ScriptD]b2b\[Rho]",
> "]"}]}],
>                           "}"}], ",",
>                         RowBox[{"{",
>                           RowBox[{
>                           "\[ScriptD]b1b\[Rho]", ",",
> "\[ScriptD]b2b\[Rho]",
>
>                             ",", "\[ScriptD]\[Rho]\[Rho]"}], "}"}]}],
> "}"}],
>
>                     "]"}]}], ";", "\n", "\t",
>                 RowBox[{"res", "=", "\t",
>                   RowBox[{"FindRoot", "[",
>                     RowBox[{
>                       RowBox[{"Evaluate", "[", "focs", "]"}], ",",
>                       RowBox[{"Evaluate", "[",
>                         RowBox[{"Sequence", " ", "@@",
>                           RowBox[{"Transpose", "[",
>                             RowBox[{"{",
>                               RowBox[{
>                                 RowBox[{"Join", "[",
>                                   RowBox[{
>                                   "\[CapitalBeta]1", ",",
> "\[CapitalBeta]2",
>
>                                     ",",
>                                     RowBox[{"{", "\[Rho]", "}"}]}], "]"}],
>
>                                 ",",
>                                 RowBox[{"Table", "[",
>                                   RowBox[{"0.01", ",",
>                                     RowBox[{"{",
>                                       RowBox[{
>                                         RowBox[{
>                                         "Length", "[", "\[CapitalBeta]1",
>                                         "]"}], "+",
>                                         RowBox[{
>                                         "Length", "[", "\[CapitalBeta]2",
>                                         "]"}], "+", "1"}], "}"}]}],
> "]"}]}],
>
>                               "}"}], "]"}]}], "]"}], ",",
>                       RowBox[{"Jacobian", "\[Rule]", "jac"}], ",",
>                       RowBox[{"Evaluate", "[",
>                         RowBox[{"FilterOptions", "[",
>                           RowBox[{"FindRoot", ",", "opts", " ", ",",
>                             RowBox[{"Sequence", " ", "@@",
>                               RowBox[{"Options", "[", "FindRoot",
> "]"}]}]}],
>
>                           " ", "]"}], "]"}]}], " ", "]"}]}], ";", " ",
> "\n",
>
>                 "\t\t\t\t",
>                 RowBox[{"{", "res", "}"}]}]}], "\n", "\t", "]"}]}],
> "]"}]}]],
>   "Input"]
>
> ________________________________________________________
> Luci Ellis
> x8786
> Payments Policy Department
> ellisl at rba.gov.au
>

```

• Prev by Date: finite element analysis in mathematica
• Next by Date: Re: Fourier transform
• Previous by thread: finite element analysis in mathematica
• Next by thread: Re: How to import DTA files conatining comment lines starting with '#'