title: Speeding up calculation of beta Pert distribution in Python
body: <p>I am calculating a beta PERT distribution for each iteration of a for loop (amongst other things, but the calculation of the distribution is what takes the most time). 
  Originally coded this in R and it took way too long to do this, thus the attempt to use a faster tool.</p>
  
  <p>Some of my datasets can be quite big, for example I just ran a case with 153413 cases and it still took about 8hrs in Python (better than R but still a bit long). </p>
  
  <p>I am quite new to Python and wondering if there is any way to speed up such a calculation?</p>
  
  <p>Example code:</p>
  
  <pre><code>af = lambda pmu, pmin, pmode, pmax: (pmu-pmin)*(2*pmode-pmin-pmax)/((pmode-pmu)*(pmax-pmin))
  bf = lambda pmu, pmin, pmode, pmax: (pmax-pmu)/(pmu-pmin)*((pmu-pmin)*(2*pmode-pmin-pmax)/((pmode-pmu)*(pmax-pmin)))
  
  e=5.
  shape=4.
  max=10.
  mu_d = np.arange(0, 10, 0.05)                
  d = np.arange(0.025, 60.025, 0.05)
  nlocs=153413  # number of rows in dataset
  
  
  f0_dist = np.zeros(len(mu_d))
  f1_dist = np.zeros(len(mu_d))
  f2_dist = np.zeros(len(mu_d))
  
  f0 = st.norm.cdf(d, 0.9/2., 0.9/6.)
  f1 = st.uniform.cdf(d, 0.001, 0.9)
  
  tic = time.clock()     
      for i in xrange(nlocs):
         for j in xrange(len(mu_d)): # mu_d has 121 values
              Rp_min = mu_d[j] - 1.96*e
              Rp_mode = mu_d[j] - 0.75*e
              Rp_max = max
              Rp_mu=(Rp_min+Rp_max+shape*Rp_mode)/(shape+2)
     dist = st.beta.cdf(d, a=af(Rp_mu, Rp_min, Rp_mode, Rp_max), b=bf(Rp_mu, Rp_min, Rp_mode, Rp_max), loc=Rp_min, scale=1-Rp_min)
  
      f0_dist[j] = 1 - np.sum(dist*f0*0.05)
      f1_dist[j] = 1- np.sum(dist*f1*0.05)
      f2_dist[j] = 1 - np.sum(dist*0.05)
      temp = 0.4*f0_dist + 0.5*f1_dist + 0.1*f1_dist
      aggr_dist = aggr_dist + temp
  
  toc = time.clock() - tic
  print '\nTime elapsed: %.3f seconds\n' % toc
  </code></pre>