dask / dask/dask-glm

Approximate exponents

Open
#23 18 comments 0 reactions 0 assignees View on GitHub
Dominant language
Python
Stars
78
Forks
47
PR merge metrics
No merged PRs in 30d

Description

Fun fact: the current ADMM implementation can spend almost half of it's time computing `np.exp`

### Script

```python
from dask import persist
import dask.array as da
import numpy as np
from dask_glm.logistic import admm
from dask_glm.utils import make_y

N = 1e7
M = 2
chunks = 1e6
seed = 20009

X = da.random.random((N, M), chunks=(chunks, M))
y = make_y(X, beta=np.array(list(range(M))), chunks=chunks)

X, y = persist(X, y)

%%prun
import dask
with dask.set_options(get=dask.get):
beta = admm(X, y)
```

### Profile results

I use a wrapped version of `np.exp` just so that it shows up in profile results.

```
565331 function calls (532847 primitive calls) in 93.202 seconds

Ordered by: internal time

ncalls tottime percall cumtime percall filename:lineno(function)
2254 41.543 0.018 41.543 0.018 utils.py:25(exp)
1127 30.614 0.027 52.086 0.046 logistic.py:349(logistic_loss)
1127 13.122 0.012 40.358 0.036 logistic.py:364(logistic_gradient)
1127 6.267 0.006 27.191 0.024 logistic.py:343(sigmoid)
1130 0.793 0.001 0.793 0.001 {method 'reduce' of 'numpy.ufunc' objects}
8 0.078 0.010 0.079 0.010 executionengine.py:100(finalize_object)
1127 0.052 0.000 40.409 0.036 logistic.py:373(proximal_logistic_gradient)
1127 0.049 0.000 52.135 0.046 logistic.py:358(proximal_logistic_loss)
8 0.043 0.005 0.043 0.005 passmanagers.py:94(run)
30 0.034 0.001 93.040 3.101 lbfgsb.py:205(_minimize_lbfgsb)
2254 0.031 0.000 41.582 0.018 dispatcher.py:152(__call__)
1127 0.029 0.000 92.998 0.083 lbfgsb.py:277(func_and_grad)
26957/12971 0.027 0.000 0.083 0.000 {method 'format' of 'str' objects}
1127 0.021 0.000 52.417 0.047 optimize.py:290(function_wrapper)
80776 0.019 0.000 0.030 0.000 {built-in method isinstance}
1127 0.013 0.000 0.814 0.001 fromnumeric.py:1743(sum)
```

### Comments

Here are some comments from @stuartarchibald

How many exp() are you doing at a time?
If you are doing many, are the input values of similar magnitude?
Do you care about IEEE754 correctness over NaN and Inf?

When I've deal with this previously, most of the time spent in a loop over exp() was in a) feraiseexcept or similar dealing with inf/nan b) branch misprediction, the split points are usually f(log_e(2)) IIRC.

a) goes away if you don't care
b) goes away if you can guarantee a range or aren't too bothered about inaccuracy out of range. Remez algs are usually used to compute the coefficient table for a polynomial and a bunch of shifts are done to get values into appropriate range.
c) if you have a load of values to compute, perhaps try Intel VML?

cc @seibert @mcg1969

Contributor guide

Open the contributing guide

Assessment

This issue has not been assessed yet.

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.