Re: Re: Why is Mathematica incredible slow?
- To: mathgroup at smc.vnet.net
- Subject: [mg38878] Re: [mg38873] Re: [mg38806] Why is Mathematica incredible slow?
- From: "Johannes Ludsteck" <johannes.ludsteck at wiwi.uni-regensburg.de>
- Date: Fri, 17 Jan 2003 05:38:31 -0500 (EST)
- Organization: Universitaet Regensburg
- Sender: owner-wri-mathgroup at wolfram.com
Dear Selwyn,
unfortunately you are not right here. I didn't
put the compiled Gradient function into the mail,
since it is more complicated (and almost unreadable in ASCII
form, but nice if you paste it into Mathematica).
THE SPEED GAINS ASSOCIATED WITH COMPILATION ARE NEGLECTIBLE HERE!
Just compare the evaluation timings
(I append the compiled version below).
The LogLik and LogLikGrad functions are not the problem:
Computation times do not differ much from Gauss and Stata.
It is plausible that compilation increases speed here
only marginally, since the uncompiled function makes extensive
use of the Listable Attribute (the arguments of PDF and CDF
are long Lists!) and the package functions are
optimized for list arguments.
Try out the compiled version of LogLikGradC below and you will
see. Thus the original question remains unanswered. By the
way: The answer is really important for me, since a negative
answer means that I simply cannot use Mathematica for
real world statistical problems (e.g. problems requiring
maximum likelihood or any other minimum estimators) on large
datasets (with 5000 observations and more).
I would be surprised if it were possible to increase
computation time of the following function.
\!\(\(LogLikGradC =
Compile[{{y, _Integer, 1}, {x, _Real, 2}, {b, _Real,
1}}, \[IndentingNewLine]With[{s = \(If[# \[Equal]
1,
1, \(-1\)] &\) /@ y},
\[IndentingNewLine]Tr /@
Transpose[\(\(\(s\ Exp[\(-
\(#\^2\/2\)\)]\)\/\(\(\@\(\(\ \)\(\(\(\
\[Pi]\)\(\ \)\)\/2\)\)\) \((1. + Erf[#\/\@2])\)\)\) x &\)[s\
x . b]]]];\)\)
If you try FindMinimum using LogLikC and LogLikGradC (as
Gradient function) you will see that providing a symbolic
gradient function again does not speed up things.
Best regards and thank you.
By the way: I would be very happy if a Wolfram Research guy
would respond, since Wolfram obtained large amounts of money
for the updates to the "fast numerics" version and the
fast statistics packages, but Mathematica apparently is not
really usable for computationally demanding applications.
Most examples in the statistics packages are toy-examples with
50-200 data points.
My most important problem is that the documentation on
algorithms of functions like FindMinimum and FindRoot is very
thin. THE EVALUATION TIMES OF LOGLIKGRAD DIFFER
ONLY BY A SMALL AMOUNT IN MATHEMATICA AND GAUSS OR STATA, BUT
NOBODY CAN TELL ME WHY FINDMINIMUM IS SO MUCH SLOWER THAN
THE GAUSS AND STATA MINIMIZATION PROCEDURES.
Johannes Ludsteck
On 16 Jan 2003 at 3:21, Selwyn Hollis wrote:
> Johannes,
>
> I'm surprised no one answered already...
> Here are my initial thoughts:
>
> I would not be so quick to blame FindRoot. Your gradient involves both
> PDF and CDF, which are contained in one of the standard packages, not
> implemented in the kernel. I suspect that is where the slowness arises.
>
> "Ordinary Mathematica code can be called from within compiled code",
> but I don't think such code is actually compiled. (??) Perhaps
> creating a compiled, specialized version of PDF/CDF would bring
> Mathematica up to speed.
>
> -----
> Selwyn Hollis
> http://www.appliedsymbols.com
>
>
> On Tuesday, January 14, 2003, at 06:10 AM, Johannes Ludsteck wrote:
>
> > 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
Institut fuer Volkswirtschaftslehre
Lehrstuhl Prof. Dr. Moeller
Universitaet Regensburg
Universitaetsstrasse 31
93053 Regensburg
Tel +49/0941/943-2741