New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make return_std 20x faster in Gaussian Processes (includes solution) #9234

Closed
andrewww opened this Issue Jun 27, 2017 · 10 comments

Comments

Projects
None yet
7 participants
@andrewww

andrewww commented Jun 27, 2017

Two requests:

(1) Please replace this line:
https://github.com/scikit-learn/scikit-learn/blob/master/sklearn/gaussian_process/gpr.py#L329

from this

    yvar -= np.einsum("ki,kj,ij->k", K_trans, K_trans, K_inv)

to this

    sum1 = np.dot(K_trans,K_inv).T
    yvar -= np.einsum("ki,ik->k", K_trans, sum1)

For an input data set of size 800x1, the time difference is 12.7 seconds to 0.2 seconds. I have validated that the result is the same up to 1e-12 or smaller.

(2) Please cache the result of the K_inv computation. It depends only on the result of training, and can be very costly for repeated calls to the class.

@andrewww

This comment has been minimized.

Show comment
Hide comment
@andrewww

andrewww Jun 27, 2017

Complete solution, starting here https://github.com/scikit-learn/scikit-learn/blob/master/sklearn/gaussian_process/gpr.py#L322

    if not hasattr(self,'K_inv_'):
        L_inv = solve_triangular(self.L_.T, np.eye(self.L_.shape[0]))
        self.K_inv_ = L_inv.dot(L_inv.T)
        
    # Compute variance of predictive distribution
    y_var = self.kernel_.diag(X)
    sum1 = np.dot(K_trans,self.K_inv_).T
    y_var1 = y_var - np.einsum("ki,ik->k", K_trans, sum1)
    # y_var2 = y_var - np.einsum("ki,kj,ij->k", K_trans, K_trans, self.K_inv_)
    # assert np.all(np.abs(y_var1-y_var2)<1e-12)
    y_var = y_var1

andrewww commented Jun 27, 2017

Complete solution, starting here https://github.com/scikit-learn/scikit-learn/blob/master/sklearn/gaussian_process/gpr.py#L322

    if not hasattr(self,'K_inv_'):
        L_inv = solve_triangular(self.L_.T, np.eye(self.L_.shape[0]))
        self.K_inv_ = L_inv.dot(L_inv.T)
        
    # Compute variance of predictive distribution
    y_var = self.kernel_.diag(X)
    sum1 = np.dot(K_trans,self.K_inv_).T
    y_var1 = y_var - np.einsum("ki,ik->k", K_trans, sum1)
    # y_var2 = y_var - np.einsum("ki,kj,ij->k", K_trans, K_trans, self.K_inv_)
    # assert np.all(np.abs(y_var1-y_var2)<1e-12)
    y_var = y_var1
@TomDLT

This comment has been minimized.

Show comment
Hide comment
@TomDLT

TomDLT Jun 27, 2017

Member

Thanks for the suggestion. Can you open a pull-request with this change?

Note that the first part was already solved in #8591, but not the second

Member

TomDLT commented Jun 27, 2017

Thanks for the suggestion. Can you open a pull-request with this change?

Note that the first part was already solved in #8591, but not the second

@andrewww

This comment has been minimized.

Show comment
Hide comment
@andrewww

andrewww Jun 27, 2017

I'm sorry, I can't. I'm posting this from an environment where I can access the website, but none of the other GitHub tools.

andrewww commented Jun 27, 2017

I'm sorry, I can't. I'm posting this from an environment where I can access the website, but none of the other GitHub tools.

@jnothman

This comment has been minimized.

Show comment
Hide comment
@jnothman

jnothman Jun 27, 2017

Member
Member

jnothman commented Jun 27, 2017

@minghui-liu

This comment has been minimized.

Show comment
Hide comment
@minghui-liu

minghui-liu Jun 27, 2017

Contributor

I would like to help and make the changes if that's ok.

Contributor

minghui-liu commented Jun 27, 2017

I would like to help and make the changes if that's ok.

@jmschrei

This comment has been minimized.

Show comment
Hide comment
@jmschrei

jmschrei Jun 27, 2017

Member
Member

jmschrei commented Jun 27, 2017

@ogrisel

This comment has been minimized.

Show comment
Hide comment
@ogrisel

ogrisel Sep 1, 2017

Member

Fixed in 6e01fef from #9236. Thanks @minghui-liu and @andrewww.

Member

ogrisel commented Sep 1, 2017

Fixed in 6e01fef from #9236. Thanks @minghui-liu and @andrewww.

@ogrisel ogrisel closed this Sep 1, 2017

@ogrisel

This comment has been minimized.

Show comment
Hide comment
@ogrisel

ogrisel Sep 1, 2017

Member

BTW I did not observe a 20x speed up a in my tests. The speed stayed approximately the same on a 1000x5 dataset generated with sklearn.datasets.make_regression.

Note: the second call to predict is significantly faster because of the cached K_inv.

@andrewww I would be curious to know more about the kind of data where you observed the initially reported speedup.

Member

ogrisel commented Sep 1, 2017

BTW I did not observe a 20x speed up a in my tests. The speed stayed approximately the same on a 1000x5 dataset generated with sklearn.datasets.make_regression.

Note: the second call to predict is significantly faster because of the cached K_inv.

@andrewww I would be curious to know more about the kind of data where you observed the initially reported speedup.

@andrewww

This comment has been minimized.

Show comment
Hide comment
@andrewww

andrewww Sep 1, 2017

Hmmm.... This was literally a make-or-break change to the code for me, i.e., the code was so slow I could not actually use it without this change (the change to the .einsum() call).

Only thing I can think of is: I'm on Windows 7 x64 / Anaconda 3.1.4. Doesn't numpy sometimes behave differently on different platforms? Maybe the Windows einsum() call is a lot slower than on Linux?

Also, apparently it was already fixed in an earlier issue #8591

andrewww commented Sep 1, 2017

Hmmm.... This was literally a make-or-break change to the code for me, i.e., the code was so slow I could not actually use it without this change (the change to the .einsum() call).

Only thing I can think of is: I'm on Windows 7 x64 / Anaconda 3.1.4. Doesn't numpy sometimes behave differently on different platforms? Maybe the Windows einsum() call is a lot slower than on Linux?

Also, apparently it was already fixed in an earlier issue #8591

@refrigerator

This comment has been minimized.

Show comment
Hide comment
@refrigerator

refrigerator Dec 5, 2017

Is anyone still having problems with this?

I have a model trained on 15,000 data points, and predictions with return_std=True take about 2 minutes for a single prediction, whereas with return_std=False it's basically instant.

refrigerator commented Dec 5, 2017

Is anyone still having problems with this?

I have a model trained on 15,000 data points, and predictions with return_std=True take about 2 minutes for a single prediction, whereas with return_std=False it's basically instant.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment