<br><br><div class="gmail_quote">On Thu, Nov 18, 2010 at 11:43 AM, Skipper Seabold <span dir="ltr">&lt;<a href="mailto:jsseabold@gmail.com">jsseabold@gmail.com</a>&gt;</span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
Hi Brian,<br>
<br>
Thanks for having a look.  I think I figured out that it is an issue<br>
of overhead and the parallel version wins as I scale it up.<br>
<div class="im"><br>
On Tue, Nov 16, 2010 at 5:16 PM, Brian Granger &lt;<a href="mailto:ellisonbg@gmail.com">ellisonbg@gmail.com</a>&gt; wrote:<br>
&gt; Skipper,<br>
&gt; I looked through your code and have a few questions and thoughts:<br>
&gt; * What data are you sending back and forth between the client and the<br>
&gt; engines.  It looks like that is going on here:<br>
<br>
</div>These are different starting values for the optimizer.  In this case<br>
is 4 different float64 arrays of shape (2,)<br>
<div class="im"><br>
&gt;        mec.push(dict(start_params=start_params1), targets=0)<br>
&gt;        mec.push(dict(start_params=start_params2), targets=1)<br>
&gt;        mec.push(dict(start_params=start_params3), targets=2)<br>
&gt;        mec.push(dict(start_params=start_params4), targets=3)<br>
&gt;        mec.execute(&quot;&quot;&quot;<br>
&gt; params = start_params<br>
&gt; ret = optimize.fmin_tnc(obj, params, approx_grad=True,<br>
&gt;                             eta=eta, maxfun=3000, args=(Y,X),<br>
&gt; messages=0)&quot;&quot;&quot;)<br>
&gt;        ret1 = mec.pull(&#39;ret&#39;, targets=0)[0]<br>
&gt;        ret2 = mec.pull(&#39;ret&#39;, targets=1)[0]<br>
&gt;        ret3 = mec.pull(&#39;ret&#39;, targets=2)[0]<br>
&gt;        ret4 = mec.pull(&#39;ret&#39;, targets=3)[0]<br>
<br>
</div>ret is a tuple of len 3 returned by the tnc optimizer.  It contains a<br>
shape (2,) array, a float, and an int.<br>
<div class="im"><br>
&gt; and:<br>
&gt;    mec.push(dict(np=np, optimize=optimize, Y=Y,X=X,eta=eta,dot=dot))<br>
&gt;<br>
<br>
</div>I&#39;m wondering if pushing the whole numpy and scipy.optimize namespace<br>
could cause a performance hit.<br>
<br></blockquote><div><br></div><div>Yes, I would use the MultiEngineClient to run &quot;import numpy as np; from scipy import optimize&quot; on all the engine rather than pushing the entire namespaces.</div><div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
Y and X are arrays of shape (100,1) and (100,2) respectively.<br>
<div class="im"><br>
&gt; One way of getting a sense of the data movement overhead is to run the code<br>
&gt; with just the data movement and no computation.<br>
&gt; * Can you time the parallel and serial code to get a better sense of where<br>
&gt; it is taking a long time.<br>
&gt; * Can you estimate the amount of time spend in the part of the code you are<br>
&gt; parallelizing.  Or, how much time is spent in the part of the code you<br>
&gt; didn&#39;t parallelize.<br>
&gt; Let&#39;s start there..<br>
&gt; Cheers,<br>
&gt; Brian<br>
&gt;<br>
<br>
</div>Comes back later...  It seems that it is indeed a matter of overhead.<br>
I&#39;ve attached a standalone script if anyone is interested that times<br>
only the relevant parts that I want to speed up.  As I make the size<br>
of the arrays bigger, the parallel versions starts to win handily.<br>
<br></blockquote><div><br></div><div>OK great.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
In [18]: run parallel_optimize.py<br>
&lt;snip&gt;<br>
T = 100<br>
0.105112075806 seconds for parallel optimization<br>
0.0321190357208 seconds for serial optimization<br>
<br>
In [19]: run parallel_optimize.py<br>
&lt;snip&gt;<br>
T = 1000<br>
0.117481946945 seconds for parallel optimization<br>
0.062824010849 seconds for serial optimization<br>
<br>
In [20]: run parallel_optimize.py<br>
&lt;snip&gt;<br>
T = 10000<br>
0.262063980103 seconds for parallel optimization<br>
0.464597940445 seconds for serial optimization<br>
<br>
In [21]: run parallel_optimize.py<br>
&lt;snip&gt;<br>
T = 100000<br>
1.49375987053 seconds for parallel optimization<br>
4.25988888741 seconds for serial optimization<br>
<br>
In my case, T will be 10, 50, and 500.  Since what I&#39;m really doing is<br>
running this code as monte carlo simulations many times.  What I could<br>
do is leave the individual fitting in serial and break the many times<br>
over my cores I suppose.  Surely that will give a performance boost.<br>
We&#39;ll see.<br>
<br></blockquote><div><br></div><div>Yes, that is a good alternative.  But I would definitely avoid pushing all of numpy and scipy.optimize.</div><div><br></div><div>Cheers,</div><div><br></div><div>Brian</div><div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
Skipper<br>
<br>
PS.  It seems that the new syntax for ipcluster is<br>
<br>
ipcluster start -n #<br>
<br>
with the most recent git cloned.  Should the nightly build docs reflect this?<br>
<div><div></div><div class="h5"><br>
&gt;<br>
&gt;<br>
&gt; On Sun, Nov 14, 2010 at 9:07 PM, Skipper Seabold &lt;<a href="mailto:jsseabold@gmail.com">jsseabold@gmail.com</a>&gt;<br>
&gt; wrote:<br>
&gt;&gt;<br>
&gt;&gt; I have a class that does some optimizations multiple times using<br>
&gt;&gt; different starting values in order to check for convergence problems<br>
&gt;&gt; due to some &quot;bad&quot; data.  My code was taking about an hour, then I<br>
&gt;&gt; added the multiple calls to scipy.optimize, and it went up to about 4<br>
&gt;&gt; hours.  I am new to this, but I tried to parallelize this code,<br>
&gt;&gt; thinking it could speed things up.  It looks like it&#39;s much actually<br>
&gt;&gt; slower.  Is there anything I can do to optimize this or remove some<br>
&gt;&gt; overhead?  A test case is inlined and attached.<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; In [3]: timeit mod.fit()<br>
&gt;&gt; 10 loops, best of 3: 39.2 ms per loop<br>
&gt;&gt;<br>
&gt;&gt; In [4]: timeit pmod.fit()<br>
&gt;&gt; 1 loops, best of 3: 255 ms per loop<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; import numpy as np<br>
&gt;&gt; from numpy import log, dot, array<br>
&gt;&gt; from scipy import optimize<br>
&gt;&gt; from IPython.kernel import client<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; def pobj(params, Y,X):<br>
&gt;&gt;    &quot;&quot;&quot;<br>
&gt;&gt;    Normal log-likelihood<br>
&gt;&gt;    &quot;&quot;&quot;<br>
&gt;&gt;    nobs2 = len(X)/2.<br>
&gt;&gt;    resid = Y - dot(X,params[:,None])<br>
&gt;&gt;    return -(- nobs2*np.log(2*np.pi)-nobs2*np.log(1/(2*nobs2) *\<br>
&gt;&gt;        dot(resid.T, resid)) - nobs2)<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; class PModel(object):<br>
&gt;&gt;    def __init__(self, Y,X):<br>
&gt;&gt;        self.X = X<br>
&gt;&gt;        self.Y = Y<br>
&gt;&gt;<br>
&gt;&gt;    def fit(self, start_params=[0,0], eta=1e-8):<br>
&gt;&gt;        start_params = np.asarray(start_params)<br>
&gt;&gt;        Y = self.Y<br>
&gt;&gt;        X = self.X<br>
&gt;&gt;<br>
&gt;&gt;        start_params1 = start_params.copy()<br>
&gt;&gt;        start_params2 = start_params + 5<br>
&gt;&gt;        start_params3 = start_params - 5<br>
&gt;&gt;        start_params4 = np.random.randint(-10,10,size=2)<br>
&gt;&gt;        mec.push(dict(start_params=start_params1), targets=0)<br>
&gt;&gt;        mec.push(dict(start_params=start_params2), targets=1)<br>
&gt;&gt;        mec.push(dict(start_params=start_params3), targets=2)<br>
&gt;&gt;        mec.push(dict(start_params=start_params4), targets=3)<br>
&gt;&gt;        mec.execute(&quot;&quot;&quot;<br>
&gt;&gt; params = start_params<br>
&gt;&gt; ret = optimize.fmin_tnc(obj, params, approx_grad=True,<br>
&gt;&gt;                             eta=eta, maxfun=3000, args=(Y,X),<br>
&gt;&gt; messages=0)&quot;&quot;&quot;)<br>
&gt;&gt;        ret1 = mec.pull(&#39;ret&#39;, targets=0)[0]<br>
&gt;&gt;        ret2 = mec.pull(&#39;ret&#39;, targets=1)[0]<br>
&gt;&gt;        ret3 = mec.pull(&#39;ret&#39;, targets=2)[0]<br>
&gt;&gt;        ret4 = mec.pull(&#39;ret&#39;, targets=3)[0]<br>
&gt;&gt;        self.results = ret1<br>
&gt;&gt;        # check results<br>
&gt;&gt;        try:<br>
&gt;&gt;            np.testing.assert_almost_equal(ret1[0], ret2[0], 4)<br>
&gt;&gt;            np.testing.assert_almost_equal(ret1[0], ret3[0], 4)<br>
&gt;&gt;            np.testing.assert_almost_equal(ret1[0], ret4[0], 4)<br>
&gt;&gt;            self.converged = str(ret1[-1])+&quot;:<br>
&gt;&gt; &quot;+optimize.tnc.RCSTRINGS[ret1[-1]]<br>
&gt;&gt;        except:<br>
&gt;&gt;            self.converged = &quot;9: Results sensitive to starting values&quot;<br>
&gt;&gt;<br>
&gt;&gt; class Model(object):<br>
&gt;&gt;    def __init__(self, Y,X):<br>
&gt;&gt;        self.X = X<br>
&gt;&gt;        self.Y = Y<br>
&gt;&gt;<br>
&gt;&gt;    def obj(self, params, Y, X):<br>
&gt;&gt;        &quot;&quot;&quot;<br>
&gt;&gt;        Normal log-likelihood<br>
&gt;&gt;        &quot;&quot;&quot;<br>
&gt;&gt;        X = self.X<br>
&gt;&gt;        Y = self.Y<br>
&gt;&gt;        nobs2 = len(X)/2.<br>
&gt;&gt;        resid = Y - np.dot(X,params[:,None])<br>
&gt;&gt;        return - nobs2*np.log(2*np.pi)-nobs2*np.log(1/(2*nobs2) *\<br>
&gt;&gt;                np.dot(resid.T,resid)) - nobs2<br>
&gt;&gt;<br>
&gt;&gt;    def fit(self, start_params=[0,0], eta=1e-8):<br>
&gt;&gt;        start_params = np.asarray(start_params)<br>
&gt;&gt;        obj = self.obj<br>
&gt;&gt;        Y = self.Y<br>
&gt;&gt;        X = self.X<br>
&gt;&gt;<br>
&gt;&gt;        start_params1 = start_params.copy()<br>
&gt;&gt;        start_params2 = start_params + 5<br>
&gt;&gt;        start_params3 = start_params - 5<br>
&gt;&gt;        start_params4 = np.random.randint(-10,10,size=2)<br>
&gt;&gt;        ret1 = optimize.fmin_tnc(obj, start_params1, approx_grad=True,<br>
&gt;&gt;                             eta=eta, maxfun=3000, args=(Y,X), messages=0)<br>
&gt;&gt;        ret2 = optimize.fmin_tnc(obj, start_params2, approx_grad=True,<br>
&gt;&gt;                             eta=eta, maxfun=3000, args=(Y,X), messages=0)<br>
&gt;&gt;        ret3 = optimize.fmin_tnc(obj, start_params3, approx_grad=True,<br>
&gt;&gt;                             eta=eta, maxfun=3000, args=(Y,X), messages=0)<br>
&gt;&gt;        ret4 = optimize.fmin_tnc(obj, start_params4, approx_grad=True,<br>
&gt;&gt;                             eta=eta, maxfun=3000, args=(Y,X), messages=0)<br>
&gt;&gt;<br>
&gt;&gt;        self.results = ret1<br>
&gt;&gt;        # check results<br>
&gt;&gt;        try:<br>
&gt;&gt;            np.testing.assert_almost_equal(ret1[0], ret2[0], 4)<br>
&gt;&gt;            np.testing.assert_almost_equal(ret1[0], ret3[0], 4)<br>
&gt;&gt;            np.testing.assert_almost_equal(ret1[0], ret4[0], 4)<br>
&gt;&gt;            self.converged = str(ret1[-1])+&quot;:<br>
&gt;&gt; &quot;+optimize.tnc.RCSTRINGS[ret1[-1]]<br>
&gt;&gt;        except:<br>
&gt;&gt;            self.converged = &quot;9: Results sensitive to starting values&quot;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; if __name__ == &quot;__main__&quot;:<br>
&gt;&gt;    np.random.seed(12345)<br>
&gt;&gt;    X = np.random.uniform(0,10,size=(10000,2))<br>
&gt;&gt;    beta = np.array([3.5,-3.5])<br>
&gt;&gt;    Y = np.dot(X,beta[:,None]) + np.random.normal(size=(10000,1))<br>
&gt;&gt;<br>
&gt;&gt;    eta = 1e-8<br>
&gt;&gt;<br>
&gt;&gt;    mec = client.MultiEngineClient()<br>
&gt;&gt;    mec.push_function(dict(obj=pobj))<br>
&gt;&gt;    mec.push(dict(np=np, optimize=optimize, Y=Y,X=X,eta=eta,dot=dot))<br>
&gt;&gt;    pmod = PModel(Y,X)<br>
&gt;&gt;    pmod.fit()<br>
&gt;&gt;<br>
&gt;&gt;    mod = Model(Y,X)<br>
&gt;&gt;    mod.fit()<br>
&gt;&gt;<br>
&gt;&gt; _______________________________________________<br>
&gt;&gt; IPython-User mailing list<br>
&gt;&gt; <a href="mailto:IPython-User@scipy.org">IPython-User@scipy.org</a><br>
&gt;&gt; <a href="http://mail.scipy.org/mailman/listinfo/ipython-user" target="_blank">http://mail.scipy.org/mailman/listinfo/ipython-user</a><br>
&gt;&gt;<br>
&gt;<br>
&gt;<br>
&gt;<br>
&gt; --<br>
&gt; Brian E. Granger, Ph.D.<br>
&gt; Assistant Professor of Physics<br>
&gt; Cal Poly State University, San Luis Obispo<br>
&gt; <a href="mailto:bgranger@calpoly.edu">bgranger@calpoly.edu</a><br>
&gt; <a href="mailto:ellisonbg@gmail.com">ellisonbg@gmail.com</a><br>
&gt;<br>
</div></div></blockquote></div><br><br clear="all"><br>-- <br>Brian E. Granger, Ph.D.<br>Assistant Professor of Physics<br>Cal Poly State University, San Luis Obispo<br><a href="mailto:bgranger@calpoly.edu">bgranger@calpoly.edu</a><br>
<a href="mailto:ellisonbg@gmail.com">ellisonbg@gmail.com</a><br>