\subsection{Influence of the Number of Parameters on Model Performance}
\label{sec:results_num_params}
\subsubsection{Ordinary Least Squares}
\label{sec:results_ols}
Even the best optimization algorithm cannot turn a bad model into a good one. If a model is not suited to describe the trend in the data irrelevant of the parameters, it is logical to assume that no optimization algorithm can find a good fit. A common approach to counteract this problem is to increase the number of parameters in a model and use a very flexible model, e.g. a polynomial of high degree. In this section the influence of the number of parameters on the model performance is studied in the context of ordinary least squares (OLS) regression.

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/ols_mse_r2.pdf}
    \caption{Mean squared error (MSE) and $R^2$ score for polynomial fits of different degrees to noisy data from the Runge function using ordinary least squares regression. $\N=\num{80000}, \Nt = \num{20000}$.}
    \label{fig:ols_mse_r2}
\end{figure}

How the mean squared error (MSE) and the $R^2$ score depend on the polynomial degree of the model is shown in \cref{fig:ols_mse_r2}. As expected, the MSE decreases and the $R^2$ score increases with increasing polynomial degree. However, it is also visible that for very high polynomial degrees the performance does not improve significantly anymore. Using a dataset with \num{80000} samples, no negative side effects of overfitting are visible even for polynomial degrees as high as 20. To highlight the danger of negative impact due to high parameter count it is also interesting to study the paramter values of the fitted models. \Cref{fig:ols_params} shows the parameter values for polynomial fits of different degrees to the same dataset. It is visible that for high polynomial degrees the parameter values become very large. Especially it is observable that the parameter values flip sign for every second parameter, i.e. compensating negative effects of one parameter with the next. While such behavior is not necessarily a problem within the feature space of the training data, it can lead to very bad performance outside of this feature space. This is especially relevant for extrapolation tasks, where the model is used to predict values outside of the feature space of the training data. In such cases, the model can produce completely meaningless results.

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/ols_parameter_plot.pdf}
    \caption{Parameter values for polynomial fits of different degrees to noisy data from the Runge function using ordinary least squares regression. $\N=\num{80000}$.}
    \label{fig:ols_params}
\end{figure}

In the limit of small datasets this issue becomes even more relevant. As the number of parameters approaches the number of data points, the model can fit the data perfectly, leading to a MSE of zero. However, this is usually not a desirable outcome, as the model will not generalize well to new data, as it has simply \textquote{memorized} the noise in the training data. This is a classic example of overfitting, where the model performs well on the training data but poorly on unseen data. In such cases, it is crucial to use techniques to prevent overfitting, such as regularization or cross-validation.

\subsubsection{Regularization Techniques}
\label{sec:results_reg_techniques}
To counteract the negative effects of high parameter count, regularization techniques such as Ridge and Lasso regression can be used. These techniques add a penalty term to the cost function that discourages large parameter values. In Ridge regression, the penalty term is proportional to the square of the parameter values, while in Lasso regression, it is proportional to the absolute value of the parameter values. This leads to a trade-off between fitting the data well and keeping the parameter values small.
With an increasing regularization strength $\lambda$ the parameter values are pushed towards zero, as visible in \cref{fig:ridge_parameter_plot}. This effect is especially strong for $\lambda > 0$. For very high regularization strengths the parameter values become very small, leading to a model that is almost constant. This is also visible in \cref{fig:ridge_mse_heatmap}, where the MSE is shown for polynomial fits of different degrees to noisy data from the Runge function using Ridge regression with different regularization strengths $\lambda$. For very high regularization strengths the MSE increases significantly, as the model is not flexible enough to fit the data well anymore. The optimal regularization strength is promoted to a hyperparameter that needs to be tuned to the problem at hand. However, with a suitable choice of $\lambda$ it is possible to use models with a high number of parameters without suffering from the negative effects of overfitting.

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/ridge_parameter_plot.pdf}
    \caption{Parameter values for polynomial fits of different degrees to noisy data from the Runge function using Ridge regression with different regularization strengths $\lambda$. $\N=\num{80000}$.}
    \label{fig:ridge_parameter_plot}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/ridge_mse_heatmap.pdf}
    \caption{Mean squared error (MSE) for polynomial fits of different degrees to noisy data from the Runge function using Ridge regression with different regularization strengths $\lambda$. $\N=\num{80000}, \Nt = \num{20000}$.}
    \label{fig:ridge_mse_heatmap}
\end{figure}





\subsection{Stochastic Gradient Descent}
\label{sec:results_sgd}
\Cref{fig:optimization_performance} shows the large decrease in computational cost, i.e. the decrease in computation time for \num{1000} epochs when using SGD instead of GD. The computation time decreases from \qty{14.82}{\s} in the case of OLS to \qty{3.71}{\s}. This is a decrease by over a factor of 3. For Ridge regression the computation time decreases from \qty{15.04}{\s} to \qty{4.29}{\s}. It is furthermore visible that the convergence of SGD is significantly faster during the first epochs as the parameters are updated more frequently. This leads to a significant decrease in the cost function during the first epochs, as the parameters are updated more frequently. 


\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/optimization_performance.pdf}
    \caption{Convergence of GD and SGD for polynomial fits of degree 10 to noisy data from the Runge function using OLS and Ridge regression with $\lambda = \num{0.1}$. $\N=\num{800000}, B=\num{512}, N_\mathrm{batches} = \num{100}$ and $\sigma = \num{0.1}$.}
    \label{fig:optimization_performance}
\end{figure}



However a second aspect of SGD is visible in \cref{fig:optimization_performance}. While GD converges smoothly to the optimal solution, SGD shows significant oscillations around the optimal solution. This is a direct consequence of the increased variance in the gradient estimate, as the gradient is only computed over a small subset of the data. This can lead to slower convergence and even no convergence at all, if the batch size is too small or the noise in the data is too high. The intrinsic noise in the gradient estimate leads also to the issue, that SGD may be unable to reach an optimal solution with a large amount of noise in the training data. To counteract this problem, the noise in the dataset has been reduced by one order of magnitude compared to the previous sections. However, even with a reduced noise level, SGD is unable to reach the optimal solution. Further fine tuning of the new hyperparameters introduced by SGD, i.e. the batch size and the batches per epoch, could lead to better results. However, this is a non-trivial task and is highly dependent on the problem at hand. To arrive at a robust solution it might be necessary to increase the batch size until reaching the limit of full-batch GD, which would defeat the purpose of using SGD in the first place.

\subsection{Effect of Learning Rate on Convergence}
\label{sec:results_learning_rate}
A crucial hyperparameter in gradient-based optimization algorithms is the learning rate $\eta$. The learning rate determines the step size in the direction of the negative gradient and thus has a significant impact on the convergence behavior of the optimization algorithm. A learning rate that is too high can lead to divergence, as the optimization algorithm overshoots the optimal solution. On the other hand, a learning rate that is too low can lead to very slow convergence, as the optimization algorithm takes very small steps towards the optimal solution. In this section the effect of the learning rate on the convergence behavior of GD is studied.

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/gradient_descent_convergence.pdf}
    \caption{Convergence of GD for polynomial fits of degree 10 to noisy data from the Runge function using OLS and Ridge regression with $\lambda = \num{0.1}$. $\N=\num{80000}, N_\mathrm{iterations} = \num{1000}$ with different learning rates $\eta$.}
    \label{fig:learning_rates}
\end{figure}

In \cref{fig:learning_rates} it is easily observable, that the convergence towards the minimum of the cost function accelerates with increasing learning rate. However, it is also visible that for very high learning rates, i.e. $\eta \geq 1$, the optimization algorithm diverges, as the cost function increases instead of decreasing. This is a direct consequence of the optimization algorithm overshooting the optimal solution and thus moving away from it instead of towards it. For very low learning rates, i.e. $\eta \leq 0.001$, the convergence is very slow, as the optimization algorithm takes very small steps towards the optimal solution. In such cases, it can be beneficial to increase the learning rate to speed up convergence. However, this comes at the risk of divergence, if the learning rate is increased too much. Thus, it is crucial to find a suitable learning rate that balances convergence speed and stability.


\subsection{Alternative Optimization Algorithms}
\label{sec:results_alt_opt_algos}
A different approach to limit the computational cost of fitting models with a high number of parameters is to use alternative optimization algorithms that reduce the number of gradient evaluations before reaching convergence. The optimized algorithms as introduced in \Cref{sec:numerical_optimization} are evaluated in regards to their convergence behavior on the same dataset and model.

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/optimizer_comparison.pdf}
    \caption{Minimization of OLS and Ridge regression cost functions using GD, Momentum GD, Adagrad, RMSProp and Adam optimization algorithms. $\N=\num{80000}, N_\mathrm{iterations}=250, p=10$ and learning rates $\eta_\mathrm{OLS} = \num{0.1}, \eta_\mathrm{Ridge} = \num{0.01}$.}
    \label{fig:optimizer_comparison}
\end{figure}

\Cref{fig:optimizer_comparison} shows the convergence behavior of GD, Momentum GD, Adagrad, RMSProp and Adam optimization algorithms for OLS and Ridge regression. It is visible that all optimized algorithms converge faster than GD, with Adam showing the best performance. While GD is unable to reach perfect convergence within \num{250} iterations. Especially for the Ridge regression cost function, the optimized algorithms converge within a very short frame. On the contrary it is also visible that in multiple optimized algorithms there are stark oscillations in the convergence behavior. This is especially visible for the Momentum GD algorithm and the Adam algorithm. While these oscillations can be reduced by fine-tuning the hyperparameters of the algorithms, they are an intrinsic property of the algorithms and cannot be completely eliminated. However, as the overall convergence behavior is significantly improved compared to GD, these oscillations are usually not a problem in practice.

\subsection{Bias-Variance Tradeoff and Resampling Techniques}
\label{sec:results_bias_variance}
\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/bias_variance_tradeoff_combined_plot.pdf}
    \caption{Bias-variance decomposition of the mean squared error (MSE) for polynomial fits of different degrees to noisy data from the Runge function using OLS. The Bias-Variance decomposition is calculated from bootstraped samples. The MSE is also calculated using 5-fold cross-validation and bootstrap resampling techniques. $\N=\num{80}, \Nt = \num{20}$.}
    \label{fig:bias_variance_tradeoff}
\end{figure}

To evaluate the bias-variance tradeoff of the optimized models, the MSE is decomposed into its bias and variance components as introduced in \cref{sec:bias_variance_tradeoff}. The results are shown in \cref{fig:bias_variance_tradeoff}. To accurately depict the bias-variance tradeoff, there were two sets of test data used. To compute the bias, no noise was added to the test data, while for the rest of the computations the same noise level as in the training data was used. It is visible that for increasing polynomial degrees the variance increases, while the bias stays almost constant. This is a direct consequence of the increased flexibility of the model with increasing polynomial degree. To ensure numerical stability, Ridge regression with a very small regularization strength of $\lambda = \num{1e-10}$ has been used to compute the bias-variance decomposition in an approximation of OLS. This ensures that the matrix inversion in the computation of the model parameters is numerically stable, while the regularization term has almost no effect on the model performance.
The optimal model complexity is usually found at the point where the sum of bias and variance is minimal. Usually this coincides with the point of the minimum of the mean squared error. How the MSE behaves on average for $k$-fold cross-validation and bootstrap resampling techniques is also shown in \cref{fig:bias_variance_tradeoff}.
Both resampling techniques show a minimum in the MSE for for the polynomial degree of 2. The crossvalidation technique shows a much more constant MSE value over the entire range of polynomial degrees, while the bootstrap technique shows a more pronounced minimum. This indicates that the bootstrap technique is more sensitive to the choice of polynomial degree and thus provides a better estimate in the case of hyperparameter tuning at the cost of a much higher computational cost.

This trend of optimal model complexity being achieved somewhere between the high bias case of low parameter count and the high variance case of high parameter count is a general property of machine learning models and is not limited to polynomial regression with OLS cost functions. As seen in \cref{fig:kfold_mse_comp} this trend of a minimum in the MSE is also observable for Ridge and Lasso regression. 

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{../figures/kfold_mse_comparison_per_cost_function.pdf}
    \caption{Mean squared error (MSE) for polynomial fits of different degrees to noisy data from the Runge function using OLS, Ridge and Lasso regression cost functions with 5-fold cross-validation resampling technique. $\N=\num{240}, \Nt = \num{60}$ and $\lambda = \num{0.1}$ for Ridge and Lasso regression.}
    \label{fig:kfold_mse_comp}
\end{figure}

It is also observable, that the minimum of the MSE is shifted towards higher polynomial degrees for Ridge and Lasso regression with a slower increase in the mean squared error above the optimal polynomial degree compared to OLS regression. This is a direct consequence of the regularization term in the cost function, which decreases the effective degrees of freedom of the model. Thus, a model with a higher number of parameters is needed to achieve the same flexibility as a model with a lower number of parameters without regularization. This effect is especially strong for Lasso regression, where the regularization term can lead to some parameters being set to zero, effectively removing them from the model. This leads to a significant reduction in the effective degrees of freedom of the model and thus a shift of the optimal polynomial degree towards higher values. The reduction in effective degrees of freedom always depends on the data at hand and the chosen regularization strength $\lambda$ and thus cannot be quantified in a general way.
