Regression algorithms predict continuous numerical outputs from input features. MLPP provides efficient implementations with automatic solver selection and built-in regularization.
Linear regression
Ordinary least squares (OLS) regression solves:
min_w (1/2n) ||Xw - y||² + (λ/2) ||w||²
LinearRegression class automatically selects the optimal solver based on problem geometry.
Basic usage
#include <mlpp/regression/linear_regression.hpp>
using namespace mlpp::regression;
// Create model with intercept, no regularization
LinearRegression<double> model(true, 0.0);
// Fit to training data
LinearRegression<double>::Matrix X_train(100, 5);
LinearRegression<double>::Vector y_train(100);
// ... populate X_train and y_train ...
model.fit(X_train, y_train);
// Make predictions
auto y_pred = model.predict(X_test);
// Evaluate R² score
double r2 = model.score(X_test, y_test);
Constructor parameters
Whether to fit a bias term. Set to false if data is already centered.
L2 penalty λ ≥ 0. Use λ > 0 for ridge regression to prevent overfitting.
method
SolveMethod
default:"Auto"
Linear solver strategy:
Auto: Chosen automatically (recommended)Cholesky: Normal equations, fast for n >> dSVD: Stable for any problem shapeJacobiSVD: Maximally stable for ill-conditioned problems
Methods
void fit(const Matrix& X, const Vector& y)
X: Feature matrix, shape (n_samples, n_features)y: Target vector, length n_samples
Vector predict(const Matrix& X) const
Scalar score(const Matrix& X, const Vector& y) const
Vector residuals(const Matrix& X, const Vector& y) const
Vector gradient(const Matrix& X, const Vector& y) const
∇L = (1/n) Xᵀ(Xw - y) + λw
const Vector& coefficients() const
Intercept (bias) term. Returns 0 if fit_intercept == false. Scalar condition_number() const
Ridge regression
Ridge regression adds L2 regularization to prevent overfitting:
// Ridge regression with λ = 0.5
LinearRegression<double> ridge(true, 0.5);
ridge.fit(X_train, y_train);
min_w (1/2n) ||Xw - y||² + (λ/2) ||w||²
The intercept is never regularized. MLPP centers the data before solving so the bias term is unaffected by λ.
Polynomial regression
Fit non-linear relationships using polynomial feature expansion:
#include <mlpp/regression/polynomial_regression.hpp>
using namespace mlpp::regression;
// Degree-3 polynomial with interactions
PolynomialRegression<double> poly(
3, // degree
true, // include_interactions
true, // fit_intercept
0.01 // regularization
);
poly.fit(X_train, y_train);
auto y_pred = poly.predict(X_test);
Feature expansion modes
Pure powers
Full interactions
With include_interactions = false, only pure power terms are added:[1, x₁, x₁², x₁³, x₂, x₂², x₂³, ...]
(include_bias ? 1 : 0) + d·DPolynomialRegression<double> model(
3, // degree
false, // no interactions
true, // include bias
0.0
);
With include_interactions = true, all cross-terms are included:[1, x₁, x₂, x₁², x₁x₂, x₂², x₁³, x₁²x₂, x₁x₂², x₂³, ...]
PolynomialRegression<double> model(
3, // degree
true, // include interactions
true,
0.01 // regularization recommended
);
Constructor parameters
Maximum polynomial degree D ≥ 1. Higher degrees can fit more complex patterns but risk overfitting.
Whether to include cross-feature monomials (e.g., x₁x₂, x₁²x₂). Dramatically increases feature count.
Forwarded to the underlying LinearRegression model.
L2 penalty forwarded to LinearRegression. Strongly recommended for high-degree polynomials.
method
SolveMethod
default:"Auto"
Solve strategy forwarded to LinearRegression.
Methods
PolynomialRegression provides the same interface as LinearRegression:
fit(X, y) - Expand features and fit modelpredict(X) - Expand features and predictscore(X, y) - Compute R² scoreresiduals(X, y) - Compute residualscoefficients() - Coefficients in expanded spaceintercept() - Intercept term
Polynomial features
You can also use PolynomialFeatures independently:
PolynomialFeatures<double> poly(3, true, true);
// Transform features
auto X_expanded = poly.transform(X_train);
// Check output dimension
size_t out_dim = poly.output_dim(X_train.cols());
Example: Complete workflow
#include <mlpp/regression/linear_regression.hpp>
#include <mlpp/regression/polynomial_regression.hpp>
#include <Eigen/Dense>
#include <iostream>
using namespace mlpp::regression;
using Matrix = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>;
using Vector = Eigen::Matrix<double, Eigen::Dynamic, 1>;
int main() {
// Generate synthetic data
Matrix X_train(100, 3);
Vector y_train(100);
// ... populate data ...
// Try linear regression
LinearRegression<double> linear(true, 0.0);
linear.fit(X_train, y_train);
std::cout << "Linear R²: " << linear.score(X_test, y_test) << std::endl;
// Try ridge regression
LinearRegression<double> ridge(true, 0.1);
ridge.fit(X_train, y_train);
std::cout << "Ridge R²: " << ridge.score(X_test, y_test) << std::endl;
// Try polynomial regression
PolynomialRegression<double> poly(2, false, true, 0.01);
poly.fit(X_train, y_train);
std::cout << "Poly R²: " << poly.score(X_test, y_test) << std::endl;
// Inspect coefficients
std::cout << "Coefficients: " << poly.coefficients().transpose() << std::endl;
std::cout << "Intercept: " << poly.intercept() << std::endl;
return 0;
}
Solver selection
Let the Auto solver choose for you unless you have specific requirements:
- Use
Cholesky when n >> d and the problem is well-conditioned
- Use
SVD for rank-deficient or ill-conditioned problems
- Use
JacobiSVD as a last resort for maximum stability
Regularization
Always use regularization for polynomial regression with degree ≥ 3 or when including interactions. Start with λ = 0.01 and tune via cross-validation.
Feature scaling
MLPP handles standardization automatically - you don’t need to scale features manually. The model:
- Standardizes features internally (zero mean, unit variance)
- Solves the problem in standardized space
- Transforms coefficients back to original scale
This ensures numerical stability while maintaining interpretability.