MathGroup Archive 2003

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: Why is Mathematica incredible slow?

  • To: mathgroup at smc.vnet.net
  • Subject: [mg38944] Re: [mg38900] Re: Why is Mathematica incredible slow?
  • From: "Johannes Ludsteck" <johannes.ludsteck at wiwi.uni-regensburg.de>
  • Date: Tue, 21 Jan 2003 07:38:46 -0500 (EST)
  • Organization: Universitaet Regensburg
  • Sender: owner-wri-mathgroup at wolfram.com

Dear Mr. Woll,
I tried to follow your advice to make the jacobian
faster. Below you find the optimized code. (You have
to past it into Mathematica. Otherwise it's unreadable).
Now I would be really surprised if it could be improved
further, since the code is really short and makes extensive
use of the listable attribute and matrix algebra operators.
Nevertheless, evaluation time has increased only by a small
amount.
BUT MY ORIGINAL QUESTION IS OPEN: MATHEMATICA STILL
REMAINS to be slower by a factor of 12 than its competitors.

Best regards,
	Johannes


\!\(\(pdf = 
      Function[{x}, 
        Exp[\(-\(x\^2\/2\)\)]\/\@\(2\ \[Pi]\)];\)\[IndentingNewLine]
  \(cdf = 
      Function[{x}, 0.5 \((1.  + Erf[x/Sqrt[2.0]])\)];\)\[IndentingNewLine]
  SetAttributes[pdf, Listable]; SetAttributes[cdf, Listable];\)

\!\(LogLikJacC = 
    Compile[{{y, _Integer, 1}, {x, _Real, 2}, {b, _Real, 1}}, 
      Module[{p, c, 
          z = Map[If[# \[Equal] 1, 1, \(-1\)] &, y]\ x . 
                b}, \[IndentingNewLine]p = pdf[z]; 
        c = cdf[z]; \(-Transpose[
              x] . \((x \(\( p\ c\ z\)\(+\)\(p\^2\)\(\ \)\)\/c\^2)\)\)]]\)


On 18 Jan 2003 at 0:38, Carl K. Woll wrote:

> Johannes,
> 
> I took a brief look at your problem, and I have a couple comments. I agree
> that there is no reason that FindRoot should be slower in Mathematica than
> in any other program, as the functionality of FindRoot is not complex.
> 
> First, I couldn't get your function Probit to work with your data set. I
> have not tried to figure out the problem here.
> 
> Second, and more importantly, your Probit function uses the nongradient form
> of FindRoot, since it provides two starting points for each variable instead
> of just one. With two starting points for each variable, FindRoot will not
> try to evaluate the gradient, and will instead use a root finding technique
> which does not rely on gradient information. Of course, using FindRoot
> without gradient information is much slower than using FindRoot with
> gradient information, and I believe this is the root of your problem.
> 
> I'm guessing that you used the two starting point, nongradient form of
> FindRoot because Mathematica is unable to find the gradient of your function
> by itself, and so Mathematica returned an error when you tried using
> FindRoot with only one starting point. However, when using the two starting
> point, nongradient form of FindRoot, gradient information is not needed, so
> supplying a gradient to FindRoot does nothing. Try supplying a gradient
> while using only one starting point for each variable.
> 
> Carl Woll
> Physics Dept
> U of Washington
> 
> "Johannes Ludsteck" <johannes.ludsteck at wiwi.uni-regensburg.de> wrote in
> message news:b00ras$plt$1 at smc.vnet.net...
> > Der MathGroup members,
> > I want to find a maximum likelihood estimator using FindRoot.
> > Since the problem seems not be related to specialities of the
> > estimation procedure, I report it here.
> > I know that the maximand is globally concave. Therefore, the
> > maximum can be found be searching the root of the gradient of
> > the maximand.
> >
> > The problem is that Mathematica is slower by a factor of 20-
> > 100 than other comparable packages, though
> > the maximand is compiled and I provide the gradient. An even
> > more amazing finding is that providing the (compiled) Jacobian
> > makes the algorithm even slower that using the default version
> > or QuasiNewton.
> >
> > By the way: Searching for the Maximum of the maximand with
> > FindfMinimum is much slower, even if I provide the gradient
> > function. Find my code below if you are interested in the
> > details or want to play around with it.
> >
> > My questions:
> > [1] Why are FindMinimum and FindRoot so slow. If I compare the
> > Timings of other number crunching functions, e.g.
> > Inverse with the corresponding ones in Gauss and Stata, I find
> > only small differences.
> >
> > [2] Are there any clear inefficiencies in my code below?
> >
> > Here is the Gradient function:
> >
> > LogLikGrad[y_, x_, b_] :=
> >   With[{s = Map[If[# == 1, 1, \(-1\)] &, y],
> >         d = NormalDistribution[0. , 1. ]},
> >     Tr /@ Transpose[s (PDF[d, #]/CDF[d, #]) x &[s x . b]]]
> >
> > I do not add the compiled version of LogLikGrad, since it does
> > not speed up things considerably and is more complex.
> >
> > and here is the estimation function:
> > Probit[y_,x_,bStart_]:=With[
> >   {b0=If[#==0.0,{-0.1,0.1},{0.9,#,1.1,#}]&/@bStart},
> >     #1/.FindRoot[LogLikGrad[y,x,#1],
> >          ##2,Compiled->False]&[
> >           #,Sequence@@Transpose[{#,b0}]]&[
> >              Table[Unique[],{Length[b0]}]]]
> >
> > Probit[y_,x_]:=Probit[y,x,Table[0.0,{Dimensions[x][[2]]}]]
> >
> > In case you want to play around with the stuff, here are data
> > generation commands and a call to Probit:
> >
> > <<Statistics`ContinuousDistributions`;
> > x=Table[2(Random[]-0.5),{10000},{5}];
> > b={0.1,0.2,0.3,0.4,0.5};
> > y=Map[If[#>0,1,0]&, x.b+
> > Table[Random[NormalDistribution[0,1]],{Length[x]}]];
> >
> > Probit[y,x]
> >
> > In case you are interested in having a look at the jacobian,
> > here it is:
> >
> > pdf =  Function[{x},
> >         Exp[-(x\^2/2]/Sqrt[2 \[Pi]];
> > cdf = Function[{x}, 0.5 (1. + Erf[x/2])];
> > SetAttributes[pdf, Listable];
> > SetAttributes[cdf, Listable];
> >
> > \!\(LogLikJacC =
> >     Compile[{{y, _Integer, 1}, {x, _Real, 2}, {b, _Real, 1}},
> >       Module[{s =
> >             Map[If[# \[Equal] 1, 1, \(-1\)] &, y]},
> > \[IndentingNewLine]Apply[
> >           Plus, \(-Map[Outer[Times, #, #] &,
> >                 x]\)*\((\(\(\(pdf[#] cdf[#]\
> > #\)\(+\)\(pdf[#]\^2\)\(\ \
> > \)\)\/cdf[#]\^2\  &\)[s\ x . b])\)]],
> > \[IndentingNewLine]{{pdf[_], _Real,
> >           1}, {cdf[_], _Real, 1}, {s, _Integer, 1}}]\)
> >
> > Best regards,
> > Johannes Ludsteck
> >
> >
> > <><><><><><><><><><><><>
> > Johannes Ludsteck
> > Economics Department
> > University of Regensburg
> > Universitaetsstrasse 31
> > 93053 Regensburg
> > Phone +49/0941/943-2741
> >
> >
> 
> 
> 


<><><><><><><><><><><><>
Johannes Ludsteck
Economics Department
University of Regensburg
Universitaetsstrasse 31
93053 Regensburg
Phone +49/0941/943-2741



  • Prev by Date: Re: Non-mathematician ... please help: Spherical Harmonic Coefficients
  • Next by Date: MatchQ[x^1, x^_Integer]
  • Previous by thread: Re: Why is Mathematica incredible slow?
  • Next by thread: STRAIGHT LINES SEGMENTS !!!