<p>Feature selection is a critical component of many machine learning tasks. It helps prevent overfitting, improves inference efficiency, and makes a given model more interpretable. In this paper, we address the problem of optimal feature selection in support vector machines (SVMs) by proposing an NP-hard bilevel optimization problem that performs embedded cross-validation to identify the feature set that minimizes the average hinge loss across folds. We reformulate the bilevel problem as a mixed-binary quadratic program with quadratic constraints, making our model amenable to commercial solvers. We also derive variable bounds for the constraints in our model, which we show do not exclude any optimal solutions. Finally, we conduct a series of computational experiments using both synthetic datasets and a publicly available credit risk classification dataset to demonstrate that the proposed model outperforms other feature selection methodologies for linear SVMs. Notably, even without enforcing a hard cardinality constraint, our model effectively identifies sparse solutions comparable to those obtained by methods that impose such constraints.</p>