# disable possible multithreading from the
# OPENBLAS and MKL linear algebra backends
import os
os.environ['OPENBLAS_NUM_THREADS'] = '1'
os.environ['MKL_NUM_THREADS'] = '1'
import numpy as np
import scipy as sp
from scipy import stats
mu, sigma = 5, 2
x = np.random.normal(loc=mu, scale=sigma, size=np.int(2e7))
def negll(par, x):
return -np.sum(sp.stats.norm.logpdf(x=x, loc=par[0], scale=par[1]))
y=negll(par=np.array([1,1]), x=x)
'''
o1 = sp.optimize.minimize(fun=negll,
x0=np.array([1,1]),
args=x,
method="L-BFGS-B",
bounds=((-np.inf, np.inf), (0.0001, np.inf)))
print(o1.x)
'''
from optimparallel import minimize_parallel
x=np.linspace(0,100,100000)
y1=x*2+1
print(x)
def f(param,x):
return np.sum((param[0]**x+param[1])**2)
if __name__ == '__main__':
o2 = minimize_parallel(fun=f,
x0=np.array([1,1]),
args=x,
bounds=((-np.inf, np.inf), (0.0001, np.inf)))
#np.all(o1.x == o2.x)
print(o2.x)
'''
o2 = minimize_parallel(fun=negll,
x0=np.array([1,1]),
args=x,
bounds=((-np.inf, np.inf), (0.0001, np.inf)))
np.all(o1.x == o2.x)
optimparallel: A Python Package Providing a Parallel Version of the L-BFGS-B Optimization Method
By Florian Geber, Lewis Blake
The Python package optimparallel provides a parallel version of the L-BFGS-B optimization method of scipy.optimize.minimize(). For an objective function with an execution time of more than 0.1 seconds and p
p
parameters the optimization speed increases by up to factor 1+p
1+p
when no analytic gradient is specified and 1+p
1+p
processor cores with sufficient memory are available.
The purpose of this jupyter notebook is to illustrate the usage of optimparallel. Note that optimparallel is the Python version of the R package optimParallel. The following examples are similar to the examples from this R Journal article.
minimize_parallel() by examples
The main function of the optimparallel package is minimize_parallel(), which has the same usage and output as scipy.optimize.minimize(), but evaluates the objective function fun()
fun()
and its gradient jac()
jac()
in parallel. For illustration, consider 2×10 7
2×107
samples from a normal distribution with mean μ=5
μ=5
and standard deviation σ=2
σ=2
.
'''
"""
Example of `minimize_parallel()` on a Windows OS.
On a Windows OS it might be necessary to run `minimize_parallel()`
in the main scope. Here is an example.
"""
from optimparallel import minimize_parallel
from scipy.optimize import minimize
import numpy as np
import time
## objective function, has to be global
def f(x, sleep_secs=.5):
print('fn')
time.sleep(sleep_secs)
return sum((x-14)**2)
def main():
"""Function to be called in the main scope."""
## start value
x0 = np.array([10,20])
## minimize with parallel evaluation of 'fun'
## and its approximate gradient
o1 = minimize_parallel(fun=f, x0=x0, args=.5, parallel={'loginfo':True})
print(o1)
## test against scipy.optimize.minimize(method='L-BFGS-B')
o2 = minimize(fun=f, x0=x0, args=.5, method='L-BFGS-B')
print(o2)
print(all(np.isclose(o1.x, o2.x, atol=1e-10)),
np.isclose(o1.fun, o2.fun, atol=1e-10),
all(np.isclose(o1.jac, o2.jac, atol=1e-10)))
## timing results
o1_start = time.time()
_ = minimize_parallel(fun=f, x0=x0, args=.5)
o1_end = time.time()
o2_start = time.time()
_ = minimize(fun=f, x0=x0, args=.5, method='L-BFGS-B')
o2_end = time.time()
print("Time parallel {:2.2} Time standard {:2.2} ".
format(o1_end - o1_start, o2_end - o2_start))
if __name__ == '__main__':
main()
"""
Example of `minimize_parallel()`.
On a Windows OS it might be necessary to run `minimize_parallel()`
in the main scope. See `example_windows_os.py`.
"""
from optimparallel import minimize_parallel
from scipy.optimize import minimize
import numpy as np
import time
import matplotlib.pyplot as plt
## objective function
def f(x, sleep_secs=.5):
print('fn')
time.sleep(sleep_secs)
return sum((x-14)**2)
## start value
x0 = np.array([10,20])
## minimize with parallel evaluation of 'fun'
## and its approximate gradient
o1 = minimize_parallel(fun=f, x0=x0, args=.5)
print(o1)
## test against scipy.optimize.minimize(method='L-BFGS-B')
o2 = minimize(fun=f, x0=x0, args=.5, method='L-BFGS-B')
print(all(np.isclose(o1.x, o2.x, atol=1e-10)),
np.isclose(o1.fun, o2.fun, atol=1e-10),
all(np.isclose(o1.jac, o2.jac, atol=1e-10)))
## timing results
o1_start = time.time()
_ = minimize_parallel(fun=f, x0=x0, args=.5)
o1_end = time.time()
o2_start = time.time()
_ = minimize(fun=f, x0=x0, args=.5, method='L-BFGS-B')
o2_end = time.time()
print("Time parallel {:2.2} Time standard {:2.2} ".
format(o1_end - o1_start, o2_end - o2_start))
## loginfo -------------------------------------
o1 = minimize_parallel(fun=f, x0=x0, args=.5, parallel={'loginfo': True})
print(o1.loginfo)
x1, x2 = o1.loginfo['x'][:,0], o1.loginfo['x'][:,1]
plt.plot(x1, x2, '-o')
for i, _ in enumerate(x1):
plt.text(x1[i]+.2, x2[i], 'f = {a:3.3f}'.format(a=o1.loginfo['fun'][i,0]))
plt.xlabel('x[0]')
plt.ylabel('x[1]')
plt.xlim(right=x1[-1]+1)
plt.show()
## example with gradient -----------------------
def g(x, sleep_secs=.5):
print('gr')
time.sleep(sleep_secs)
return 2*(x-14)
o3 = minimize_parallel(fun=f, x0=x0, jac=g, args=.5)
o4 = minimize(fun=f, x0=x0, jac=g, args=.5, method='L-BFGS-B')
print(all(np.isclose(o3.x, o4.x, atol=1e-10)),
np.isclose(o3.fun, o4.fun, atol=1e-10),
all(np.isclose(o3.jac, o4.jac, atol=1e-10)))
o3 = minimize_parallel(fun=negll,
x0=np.array([1,1]),
args=x,
bounds=((-np.inf, np.inf), (0.0001, np.inf)),
parallel={'loginfo': True, 'time': True})
o3.time
# show the three last evaluated parameters (columns 1 and 2)
# and the corresponding 'fun' values (column 3)
np.c_[o3.loginfo['x'][-3:], o3.loginfo['fun'][-3:]]
import matplotlib.pyplot as plt
x1, x2 = o3.loginfo['x'][:,0], o3.loginfo['x'][:,1]
plt.plot(x1, x2, '-o')
plt.plot(x1[-1:], x2[-1:], '-o', color='red', markersize=50, marker='+')
plt.text(1.2, 1, 'start')
plt.xlabel('$mu$'); plt.ylabel('$sigma$')
plt.show()