Skip to content
Browse files

Performance improvements for LARS precomputed Gram matrix.

Greatly improves performance in the case of more samples than
features. In the case of more features than samples, results are
roughly the same.
  • Loading branch information...
1 parent 043cd96 commit 4c0d2702a756e2f27b8e906f8dfb30a989b0f729 @fabianp fabianp committed Sep 27, 2010
Showing with 19 additions and 19 deletions.
  1. +3 −2 scikits/learn/glm/benchmarks/bench_lasso.py
  2. +16 −17 scikits/learn/glm/lars.py
View
5 scikits/learn/glm/benchmarks/bench_lasso.py
@@ -15,7 +15,7 @@
from time import time
import numpy as np
-from bench_glm import make_data
+from bench_glmnet import make_data
def bench(clf, X_train, Y_train, X_test, Y_test, Gram=None):
gc.collect()
@@ -69,8 +69,9 @@ def LassoLARSFactory(alpha):
larslasso_results.append(bench(LassoLARSFactory(alpha),
X, Y, X_test, Y_test))
print "benching LassoLARS (precomp. Gram): "
+ Gram = np.dot(X.T, X)
larslasso_gram_results.append(bench(LassoLARSFactory(alpha),
- X, Y, X_test, Y_test, Gram=np.dot(X.T, X)))
+ X, Y, X_test, Y_test, Gram=Gram))
return lasso_results, larslasso_results, larslasso_gram_results
View
33 scikits/learn/glm/lars.py
@@ -15,7 +15,6 @@
from .base import LinearModel
from ..utils import arrayfuncs
-
def lars_path(X, y, Gram=None, max_features=None, alpha_min=0,
method="lar", verbose=False):
""" Compute Least Angle Regression and LASSO path
@@ -96,7 +95,11 @@ def lars_path(X, y, Gram=None, max_features=None, alpha_min=0,
# whole array, which can cause problems.
L = np.zeros((max_features, max_features), dtype=np.float64)
- X = np.asfortranarray(X) # make sure data are contiguous in memory
+ if Gram is None:
+ # setting the array to be fortran-ordered speeds up the
+ # calculation of the (partial) Gram matrix, but not very
+ # useful if the Gram matrix is already precomputed.
+ X = np.asfortranarray(X)
Xt = X.T
if Gram is not None:
@@ -120,7 +123,6 @@ def lars_path(X, y, Gram=None, max_features=None, alpha_min=0,
# To get the most correlated variable not already in the active set
arrayfuncs.dot_over(Xt, res, active_mask, np.False_, Cov)
else:
- # could use dot_over
arrayfuncs.dot_over(Gram, coefs[n_iter], active_mask, np.False_, a)
Cov = Xty[inactive_mask] - a[:n_inactive]
@@ -180,26 +182,23 @@ def lars_path(X, y, Gram=None, max_features=None, alpha_min=0,
# (Golub & Van Loan, 1996)
# compute eqiangular vector
+ b = linalg.cho_solve((L[:n_active, :n_active], True),
+ sign_active[:n_active])
+ AA = 1. / sqrt(np.sum(b * sign_active[:n_active]))
+ b *= AA
+
if Gram is None:
- b = linalg.cho_solve((L[:n_active, :n_active], True),
- sign_active[:n_active])
- AA = 1. / sqrt(np.sum(b * sign_active[:n_active]))
- b *= AA
+ eqdir = np.dot(X[:,active], b) # equiangular direction (unit vector)
+ # correlation between active variables and eqiangular vector
+ arrayfuncs.dot_over(Xt, eqdir, active_mask, np.False_, a)
else:
- S = sign_active[:n_active][:,None] * sign_active[:n_active][None,:]
- b = linalg.inv(Gram[active][:,active] * S)
- b = np.sum(b, axis=1)
- AA = 1. / sqrt(b.sum())
- b *= sign_active[:n_active]
- b *= AA
-
- eqdir = np.dot(X[:,active], b) # equiangular direction (unit vector)
+ # in the case of huge number of features, this takes 50% of time
+ # to avoid the copies, we could reorder the rows of Gram ...
+ arrayfuncs.dot_over (Gram[active].T, b, active_mask, np.False_, a)
if n_active >= n_features:
gamma_ = C / AA
else:
- # correlation between active variables and eqiangular vector
- arrayfuncs.dot_over(Xt, eqdir, active_mask, np.False_, a)
# equation 2.13
g1 = (C - Cov[:n_inactive]) / (AA - a[:n_inactive])
g2 = (C + Cov[:n_inactive]) / (AA + a[:n_inactive])

0 comments on commit 4c0d270

Please sign in to comment.
Something went wrong with that request. Please try again.