How to add L1 norm as a constraint in PCA
Answered
I am trying to solve the PCA problem adding an extra \(L_1\) constraint into it. Given a data matrix \(X\), my objective function to solve is
\[max_w{w^tX^tXw} \quad s.t. \quad \|w\|_2=1, \quad \|w\|_1 \leq t\]
I have already worked out how to solve the objective function + the first constraint but I have no idea on how to include the \(L_1\) constraint.
So far this is what I have,
import gurobipy as gp
from gurobipy import GRB
import numpy as np
from sklearn.datasets import load_boston
from itertools import product
# Load data
boston = load_boston()
x = boston.data
x = x[:, [0, 2, 4, 5, 6, 7, 10, 11, 12]] # select non-categorical variables
x_center = x - np.mean(x, 0)
n, p = x_center.shape
Quad = np.dot(x_center.T, x_center)
# Define model
pca_gurobi = gp.Model()
weights = pca_gurobi.addVars(p, lb=-GRB.INFINITY, name="weights") # pca loading
obj_pca = (1/n)*sum(Quad[i, j] * weights[i] * weights[j] for i, j in product(range(p), repeat=2))
pca_gurobi.addConstr(sum(weights[i]**2 for i in range(p)) == 1)
# Add L1 constraint
# lasso_threshold = 5.0
# First try (Error due to sum over GenExpr)
# pca_gurobi.addConstr(sum(gp.abs_(weights[i]) for i in range(p)) <= lasso_threshold)
# Second try: Error
# abs_weights = pca_gurobi.addVars(p, name="abs_weights")
# for v, absv in zip(weights, abs_weights):
# pca_gurobi.addConstr(gp.abs_(v)==absv)
# pca_gurobi.addConstr(sum(abs_w for abs_w in abs_weights.tolist()) <= lasso_threshold)
# Solve model
pca_gurobi.setObjective(obj_pca, GRB.MAXIMIZE)
pca_gurobi.params.NonConvex = 2 # Let gurobi know that constraint is non convex
pca_gurobi.optimize()
weights_pca_gurobi = np.array([weights[i].X for i in range(p)])
My second approach on adding the L1 norm was based on this post here:
And I saw that the answer was based on using MVar instead of normal variables, so I converted my problem to use MVar, but it somehow resulted into an infeasible problem. This is the code for my second approach:
spca_gurobi2 = gp.Model()
s_weights2 = spca_gurobi2.addMVar(p, lb=-GRB.INFINITY, name="sparse_weights") # pca loading
lasso_threshold = 5.0
obj_spca2 = 0
for i, j in product(range(p), repeat=2):
list_w = s_weights2.tolist()
obj_spca2 += Quad[i, j] * list_w[i] * list_w[j]
obj_spca2 = (1/n)*obj_spca2
spca_gurobi2.addConstr(sum(w**2 for w in s_weights2.tolist()) == 1)
# First try
abs_s_weights2 = spca_gurobi2.addMVar(p, name="abs_sparse_weights")
for v, absv in zip(s_weights2.tolist(), abs_s_weights2.tolist()):
spca_gurobi2.addConstr(absv == gp.abs_(v))
spca_gurobi2.addConstr(sum(abs_w for abs_w in abs_s_weights2.tolist()) <= lasso_threshold)
spca_gurobi2.setObjective(obj_spca2, GRB.MAXIMIZE)
spca_gurobi2.params.NonConvex = 2 # Let gurobi know that constraint is non convex
spca_gurobi2.optimize()
weights_spca_gurobi2 = np.array([s_weights2[i].X for i in range(p)])
I am pretty new to Gurobi so I would really appreciate any advice here.
-
Official comment
-
Gurobi Staff
This post is more than three years old. Some information may not be up to date. For current information, please check the Gurobi Documentation or Knowledge Base. If you need more help, please create a new post in the community forum. Or why not try our AI Gurobot?. -
-
-
Hi Alvaro,
You need to write the L1 constraints like this:
for v, absv in zip(weights, abs_weights):
pca_gurobi.addConstr(abs_weights[absv] == gp.abs_(weights[v]))Please note that looping over the variable dicts will generate index objects as looping variables, so you need to use those to specify the actual decision variables in the respective dicts weights and abs_weights.
The error message here is quite misleading because the input data for the addConstr call doesn't contain a variable.
I hope that helps.
Cheers,
Matthias0 -
Post is closed for comments.
Comments
2 comments