R apply matching speed in Python
I want to effectively apply a complex function to a matrix row in Python (edit: Python 3). In R, this is the apply function and its ilk, which runs quickly.
In Python, I know this can be done in several ways. List understanding, numpy.apply_along_axis, panas.dataframe.apply.
In my coding, these Python methods are very slow. Should I use another method? Or am I not implementing these Python methods correctly?
Here is an example. The math is taken from probabilistic regression models. To be clear, my goal is not to perform probabilistic regression, and I am interested in an effective application method.
In R:
> n = 100000
> p = 7
> x = matrix(rnorm(700000, 0 , 2), ncol = 7)
> beta = rep(1, p)
> start <- Sys.time()
> test <- apply(x, 1, function(t)(dnorm(sum(t*beta))*sum(t*beta)/pnorm(sum(t*beta))) )
> end <- Sys.time()
> print(end - start)
Time difference of 0.6112201 secs
In Python by understanding:
import numpy as np
from scipy.stats import norm
import time
n = 100000
p = 7
x = norm.rvs(0, 2, n * p)
x = x.reshape( (n , p) )
beta = np.ones(p)
start = time.time()
test = [
norm.pdf(sum(x[i,]*beta))*sum(x[i,]*beta)/norm.cdf(sum(x[i,]*beta))
for i in range(100000) ]
end = time.time()
print (end - start)
23.316735982894897
In Python via pandas.dataframe.apply:
frame = DataFrame(x)
f = lambda t: norm.pdf(sum(t))*sum(t)/norm.cdf(sum(t))
start = time.time()
test = frame.apply(f, axis = 1)
end = time.time()
print(end - start)
34.39404106140137
In the answer to this question states apply_along_ AXIS is not for speed. So I don’t include this method.
Again, I’m interested in performing these calculations quickly. Thank you very much for your help!
Solution
List understanding uses python-level loops and is very inefficient. You should modify your code to take advantage of numpy vectorization. If you change the content between time.time()
calls
xbeta = np.sum(x * beta, axis=1)
test = norm.pdf(xbeta) * xbeta / norm.cdf(xbeta)
You’ll see a huge difference. For my machine, it completes in 0.02 seconds.
To give you peace of mind, I’ve tested it against your list understanding and they all give the same results.
xbeta
is something you waste your calculations many times in your calculations. By summing along the second axis, we collapse it into a one-dimensional array, which is your 100000 rows. Now all calculations deal with one-dimensional arrays, so we let numpy handle the rest.