Tuning penalty strength in scikit-learn logistic regression

Question

From scikit-learn's user guide, the loss function for logistic regression is expressed in this generalized form:

$$ \min_{w,c}\frac{1-\rho}{2}w^{T}w+\rho\lVert w\rVert_{1}+C\sum_{i=1}^{n}\log\left(\exp\left(-y_{i}\left(x_{i}^{T}w+c\right)\right)+1\right). $$

This is all fine if you are working with a static dataset. What I don't get is, once you have tuned your C using some cross-validation procedure, and then you go out and collect more data, you might have to proportionally adjust the optimal C or even re-tune C altogether.

In the extreme case, assume iid distribution of all samples, if we flood the original dataset with 100x more data, and we repeat our CV procedure, the new optimal C will surely look very different from the original one. This seems very unecessary to me. If we have a stable structure of the model, whether we fit the model on the original small sample or the flooded big sample, we should have similar betas and C's.

Instead, why don't we express penalty strength in terms of mean per-sample loss:

$$ \min_{w,c}\frac{1-\rho}{2}w^{T}w+\rho\lVert w\rVert_{1}+\frac{c}{n}\sum_{i=1}^{n}\log\left(\exp\left(-y_{i}\left(x_{i}^{T}w+c\right)\right)+1\right). $$

That way I feel more comfortable leaving the optimal c alone and only re-fit my betas when more data comes in.

Please share if you've encountered some discussion on this point.

Answer 1

This is indeed a reasonable approach from a machine learning perspective, and I did something similar in my Weighted Least-Squares Support Vector Machine implementation (see this paper) so that the range of hyper-parameter values that you need to search is more compact and the optimal value less dependent on the number of training samples.

I suspect the reason it is not as commonly seen in more statistically based models (rather than models from a more machine learning oriented source) is that Bayesian model selection schemes would require the overall loss, rather than the per-pattern loss, as might AIC or BIC.

However, different samples have different sampling variations (noise) so you may want to retune $C$ for different training samples, even if they are of the same size (as there is no telling that the value for the first sample of data was not an "outlier" in the distribution of optimal values).

Answer 2

Instead, why don't we express penalty strength in terms of mean per-sample loss?

Without the penalty term $\lambda_2||w||_2^2 + \lambda_1 ||w||_1$, the normalization of the loss would not matter. With penalty, the optimal values of the penalty strengths, $\lambda_1$ and $\lambda_2$, depend on the (size of the) training set. Hence, if a larger training set becomes available, one would usually again search for the (new) optimal $\lambda$ anyway.

In the end, it's a matter of convention. The one without $\frac{1}{n}$ is the one used in Eq. (4.31) of [1] and Eq. (1.8) of [2].

Interestingly, for Ridge regression (OLS + L2 penalty) on a Gaussian target with an orthogonal design matrix $X \in \mathbb{R}^{n,p}$, i.e., $Y \sim \mathcal{N}(X\beta, \sigma^2)$ , it turns out that the optimal penalty strength in terms of MSE of the estimated coefficients $\hat{\beta}$ is $\lambda = \frac{p\sigma^2}{||\beta||_2^2}$, see Example 1.7 of [2]. This does not shrink with larger sample size $n$.

[1] Hastie, T. and Tibshirani, R. and Friedman, "The Elements of Statistical Learning"

[2] Wessel N. van Wieringen, "Lecture notes on ridge regression", https://arxiv.org/pdf/1509.09169.pdf