Skip to content

Commit 31b2121

Browse files
idc9idc9
committed
added cvxpy
1 parent cef8777 commit 31b2121

7 files changed

Lines changed: 392 additions & 0 deletions

File tree

‎ya_glm/backends/cvxpy/__init__.py‎

Whitespace-only changes.

‎ya_glm/backends/cvxpy/glm_solver.py‎

Lines changed: 209 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,209 @@
1+
import cvxpy as cp
2+
from functools import partial
3+
4+
from ya_glm.utils import clip_zero
5+
from ya_glm.cvxpy.penalty import lasso, ridge
6+
from ya_glm.cvxpy.loss_functions import lin_reg_loss, log_reg_loss,\
7+
quantile_reg_loss
8+
from ya_glm.cvxpy.utils import solve_with_backups
9+
from ya_glm.backends.fista.glm_solver import process_param_path
10+
11+
12+
def solve_glm(X, y,
13+
loss_func='lin_reg',
14+
loss_kws={},
15+
fit_intercept=True,
16+
lasso_pen=None,
17+
lasso_weights=None,
18+
groups=None,
19+
L1to2=False,
20+
nuc=False,
21+
ridge_pen=None,
22+
ridge_weights=None,
23+
tikhonov=None,
24+
coef_init=None,
25+
intercept_init=None,
26+
zero_tol=1e-8,
27+
cp_kws={}):
28+
29+
######################
30+
# objective function #
31+
######################
32+
33+
problem, coef, intercept, _, __ = \
34+
setup_problem(X=X, y=y,
35+
loss_func=loss_func,
36+
loss_kws=loss_kws,
37+
fit_intercept=fit_intercept,
38+
lasso_pen=lasso_pen,
39+
lasso_weights=lasso_weights,
40+
groups=groups,
41+
L1to2=L1to2,
42+
nuc=nuc,
43+
ridge_pen=ridge_pen,
44+
ridge_weights=ridge_weights,
45+
coef_init=coef_init,
46+
intercept_init=intercept_init)
47+
48+
solve_with_backups(problem=problem, variable=coef, **cp_kws)
49+
50+
if coef.value is None:
51+
raise RuntimeError("cvxpy solvers failed")
52+
53+
return process_output(problem=problem, coef=coef,
54+
intercept=intercept,
55+
fit_intercept=fit_intercept,
56+
zero_tol=zero_tol)
57+
58+
59+
def solve_glm_path(fit_intercept=True, cp_kws={}, zero_tol=1e-8,
60+
lasso_pen_seq=None, ridge_pen_seq=None,
61+
check_decr=True, **kws):
62+
63+
param_path = process_param_path(lasso_pen_seq=lasso_pen_seq,
64+
ridge_pen_seq=ridge_pen_seq,
65+
check_decr=check_decr)
66+
67+
# make sure we setup the right penalty
68+
if 'lasso_pen' in param_path:
69+
kws['lasso_pen'] = param_path['lasso_pen'][0]
70+
if 'ridge_pen' in param_path:
71+
kws['ridge_pen'] = param_path['ridge_pen'][0]
72+
73+
problem, coef, intercept, lasso_pen, ridge_pen = setup_problem(**kws)
74+
75+
for params in param_path:
76+
if 'lasso_pen' in params:
77+
lasso_pen.value = params['lasso_pen']
78+
79+
if 'ridge_pen' in params:
80+
ridge_pen.value = params['ridge_pen']
81+
82+
solve_with_backups(problem=problem, variable=coef, **cp_kws)
83+
84+
if coef.value is None:
85+
raise RuntimeError("cvxpy solvers failed")
86+
87+
# format output
88+
fit_out = process_output(problem=problem, coef=coef,
89+
intercept=intercept,
90+
fit_intercept=fit_intercept,
91+
zero_tol=zero_tol)
92+
93+
# if generator:
94+
yield fit_out, params
95+
96+
97+
def process_output(problem, coef, intercept, fit_intercept, zero_tol):
98+
99+
opt_data = {**problem.solver_stats.__dict__,
100+
'status': problem.status}
101+
102+
coef_sol = clip_zero(coef.value, zero_tol=zero_tol)
103+
if fit_intercept:
104+
intercept_sol = clip_zero(intercept.value, zero_tol=zero_tol)
105+
else:
106+
intercept_sol = None
107+
108+
return coef_sol, intercept_sol, opt_data
109+
110+
111+
def setup_problem(X, y,
112+
loss_func='lin_reg',
113+
loss_kws={},
114+
fit_intercept=True,
115+
lasso_pen=None,
116+
lasso_weights=None,
117+
groups=None,
118+
L1to2=False,
119+
nuc=False,
120+
ridge_pen=None,
121+
ridge_weights=None,
122+
tikhonov=None,
123+
coef_init=None,
124+
intercept_init=None):
125+
"""
126+
127+
Output
128+
------
129+
problem, coef, intercept, lasso_pen, ridge_pen
130+
131+
"""
132+
glm_loss = get_glm_loss(loss_func)
133+
134+
# TODO: add these
135+
if tikhonov is not None:
136+
raise NotImplementedError
137+
if groups is not None:
138+
raise NotImplementedError
139+
if L1to2:
140+
raise NotImplementedError
141+
if nuc:
142+
raise NotImplementedError
143+
144+
if lasso_pen is None and lasso_weights is not None:
145+
lasso_pen = 1
146+
147+
if ridge_pen is None and ridge_weights is not None:
148+
ridge_pen = 1
149+
150+
if lasso_pen is not None:
151+
lasso_pen = cp.Parameter(nonneg=True, value=lasso_pen)
152+
153+
if ridge_pen is not None:
154+
ridge_pen = cp.Parameter(nonneg=True, value=ridge_pen)
155+
156+
if lasso_pen is not None and ridge_pen is not None:
157+
158+
def objective(coef, intercept):
159+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept) + \
160+
lasso_pen * lasso(coef, weights=lasso_weights) + \
161+
ridge_pen * ridge(coef, weights=ridge_weights)
162+
163+
elif lasso_pen is not None:
164+
def objective(coef, intercept):
165+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept) + \
166+
lasso_pen * lasso(coef, weights=lasso_weights)
167+
168+
elif ridge_pen is not None:
169+
def objective(coef, intercept):
170+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept) + \
171+
ridge_pen * ridge(coef, weights=ridge_weights)
172+
173+
else:
174+
def objective(coef, intercept):
175+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept)
176+
177+
###################
178+
# setup variables #
179+
###################
180+
181+
coef = cp.Variable(shape=X.shape[1], value=coef_init)
182+
if fit_intercept:
183+
intercept = cp.Variable(value=intercept_init)
184+
else:
185+
intercept = None
186+
187+
###########################
188+
# setup and solve problem #
189+
###########################
190+
problem = cp.Problem(cp.Minimize(objective(coef, intercept)))
191+
192+
return problem, coef, intercept, lasso_pen, ridge_pen
193+
194+
195+
def get_glm_loss(loss_func, loss_kws={}):
196+
197+
if loss_func == 'lin_reg':
198+
loss = lin_reg_loss
199+
200+
elif loss_func == 'log_reg':
201+
loss = log_reg_loss
202+
203+
elif loss_func == 'quantile':
204+
loss = quantile_reg_loss
205+
206+
if len(loss_kws) > 0:
207+
loss = partial(loss, **loss_kws)
208+
209+
return loss

‎ya_glm/cvxpy/__init__.py‎

Whitespace-only changes.

‎ya_glm/cvxpy/loss_functions.py‎

Lines changed: 35 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,35 @@
1+
import cvxpy as cp
2+
3+
4+
def lin_reg_loss(X, y, coef, intercept=None):
5+
"""
6+
1/(2 *n_samples) ||X @ ceof + intercept - y||_2^2
7+
"""
8+
pred = X @ coef
9+
if intercept is not None:
10+
pred += intercept
11+
12+
return (0.5 / X.shape[0]) * cp.sum_squares(pred - y)
13+
14+
15+
def log_reg_loss(X, y, coef, intercept=None):
16+
pred = X @ coef
17+
if intercept is not None:
18+
pred += intercept
19+
20+
return (1 / X.shape[0]) * cp.sum(cp.logistic(pred) - cp.multiply(y, pred))
21+
22+
23+
def quantile_reg_loss(X, y, coef, intercept=None, quantile=0.5):
24+
pred = X @ coef
25+
if intercept is not None:
26+
pred += intercept
27+
28+
return (1 / X.shape[0]) * cp.sum(tilted_L1(y - pred, quantile=quantile))
29+
30+
31+
def tilted_L1(u, quantile=0.5):
32+
"""
33+
tilted_L1(u; quant) = quant * [u]_+ + (1 - quant) * [u]_
34+
"""
35+
return 0.5 * cp.abs(u) + (quantile - 0.5) * u

‎ya_glm/cvxpy/penalty.py‎

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
import cvxpy as cp
2+
3+
4+
def lasso(coef, weights=None):
5+
if weights is not None:
6+
return cp.norm1(cp.multiply(coef, weights))
7+
else:
8+
return cp.norm1(coef)
9+
10+
11+
def ridge(coef, weights=None):
12+
if weights is not None:
13+
return 0.5 * sum(cp.multiply(coef ** 2, weights))
14+
else:
15+
return 0.5 * cp.sum_squares(coef)

‎ya_glm/cvxpy/solve_glm.py‎

Lines changed: 101 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,101 @@
1+
import cvxpy as cp
2+
from time import time
3+
from functools import partial
4+
5+
from ya_glm.cvxpy.penalty import lasso, ridge
6+
from ya_glm.cvxpy.loss_functions import lin_reg_loss, log_reg_loss,\
7+
quantile_reg_loss
8+
9+
from ya_glm.cvxpy.utils import solve_with_backups
10+
11+
12+
def solve_glm(X, y,
13+
loss_func='lin_reg',
14+
loss_kws={},
15+
fit_intercept=True,
16+
lasso_pen=None,
17+
lasso_weights=None,
18+
ridge_pen=None,
19+
ridge_weights=None,
20+
coef_init=None,
21+
intercept_init=None,
22+
cp_kws={}):
23+
24+
start_time = time()
25+
######################
26+
# objective function #
27+
######################
28+
glm_loss = get_glm_loss(loss_func)
29+
30+
if lasso_pen is None and lasso_weights is not None:
31+
lasso_pen = 1
32+
33+
if ridge_pen is None and ridge_weights is not None:
34+
ridge_pen = 1
35+
36+
if lasso_pen is not None and ridge_pen is not None:
37+
38+
def objective(coef, intercept):
39+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept) + \
40+
lasso_pen * lasso(coef, weights=lasso_weights) + \
41+
ridge_pen * ridge(coef, weights=ridge_weights)
42+
43+
elif lasso_pen is not None:
44+
def objective(coef, intercept):
45+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept) + \
46+
lasso_pen * lasso(coef, weights=lasso_weights)
47+
48+
elif ridge_pen is not None:
49+
def objective(coef, intercept):
50+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept) + \
51+
ridge_pen * ridge(coef, weights=ridge_weights)
52+
53+
else:
54+
def objective(coef, intercept):
55+
return glm_loss(X=X, y=y, coef=coef, intercept=intercept)
56+
57+
###################
58+
# setup variables #
59+
###################
60+
61+
coef = cp.Variable(shape=X.shape[1], value=coef_init)
62+
if fit_intercept:
63+
intercept = cp.Variable(value=intercept_init)
64+
else:
65+
intercept = None
66+
67+
###########################
68+
# setup and solve problem #
69+
###########################
70+
problem = cp.Problem(cp.Minimize(objective(coef, intercept)))
71+
72+
solve_with_backups(problem=problem, variable=coef, **cp_kws)
73+
74+
if coef.value is None:
75+
raise RuntimeError("cvxpy solvers failed")
76+
77+
opt_data = {'runtime': time() - start_time,
78+
'problem': problem}
79+
80+
if fit_intercept:
81+
return coef.value, intercept.value, opt_data
82+
83+
else:
84+
return coef.value, None, opt_data
85+
86+
87+
def get_glm_loss(loss_func, loss_kws={}):
88+
89+
if loss_func == 'lin_reg':
90+
loss = lin_reg_loss
91+
92+
elif loss_func == 'log_reg':
93+
loss = log_reg_loss
94+
95+
elif loss_func == 'quantile':
96+
loss = quantile_reg_loss
97+
98+
if len(loss_kws) > 0:
99+
loss = partial(loss, **loss_kws)
100+
101+
return loss

‎ya_glm/cvxpy/utils.py‎

Lines changed: 32 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,32 @@
1+
import cvxpy as cp
2+
3+
4+
def solve_with_backups(problem, variable, verbosity=0, **kws):
5+
6+
problem.solve(**kws)
7+
8+
# try back up solvers if the solver did not work
9+
if variable.value is None:
10+
11+
if verbosity >= 1:
12+
print("Solver {} failed".
13+
format(problem.solver_stats.solver_name))
14+
15+
# list available solvers we have not yet tried
16+
avail_solvers = cp.installed_solvers()
17+
avail_solvers.remove(problem.solver_stats.solver_name)
18+
19+
for solver in avail_solvers:
20+
21+
# try each solver
22+
kws['solver'] = solver
23+
problem.solve(**kws)
24+
25+
# if we have succeded then we are done!
26+
if variable.value is not None:
27+
break
28+
29+
else:
30+
if verbosity >= 1:
31+
print("Solver {} failed".
32+
format(problem.solver_stats.solver_name))

0 commit comments

Comments
 (0)