diff --git a/figures/cost_function_comparison.pdf b/figures/cost_function_comparison.pdf new file mode 100644 index 0000000..db9eea5 Binary files /dev/null and b/figures/cost_function_comparison.pdf differ diff --git a/figures/data_scatter.png b/figures/data_scatter.png index b85489b..8162739 100644 Binary files a/figures/data_scatter.png and b/figures/data_scatter.png differ diff --git a/figures/gradient_descent_convergence.pdf b/figures/gradient_descent_convergence.pdf index 0e19a08..6fa93b5 100644 Binary files a/figures/gradient_descent_convergence.pdf and b/figures/gradient_descent_convergence.pdf differ diff --git a/figures/ols_mse_r2.pdf b/figures/ols_mse_r2.pdf index bbaf21a..6df088f 100644 Binary files a/figures/ols_mse_r2.pdf and b/figures/ols_mse_r2.pdf differ diff --git a/figures/ols_parameter_plot.pdf b/figures/ols_parameter_plot.pdf index d7f21b0..96f13f6 100644 Binary files a/figures/ols_parameter_plot.pdf and b/figures/ols_parameter_plot.pdf differ diff --git a/figures/optimizer_comparison.pdf b/figures/optimizer_comparison.pdf new file mode 100644 index 0000000..c93b799 Binary files /dev/null and b/figures/optimizer_comparison.pdf differ diff --git a/figures/ridge_mse_r2.pdf b/figures/ridge_mse_r2.pdf index b5d8bed..942e495 100644 Binary files a/figures/ridge_mse_r2.pdf and b/figures/ridge_mse_r2.pdf differ diff --git a/figures/ridge_mse_r2_lambda.pdf b/figures/ridge_mse_r2_lambda.pdf index 411b9f1..b116957 100644 Binary files a/figures/ridge_mse_r2_lambda.pdf and b/figures/ridge_mse_r2_lambda.pdf differ diff --git a/figures/ridge_parameter_plot.pdf b/figures/ridge_parameter_plot.pdf index f63190f..0a36a8b 100644 Binary files a/figures/ridge_parameter_plot.pdf and b/figures/ridge_parameter_plot.pdf differ diff --git a/src/main.ipynb b/src/main.ipynb index eecd1e4..55562da 100644 --- a/src/main.ipynb +++ b/src/main.ipynb @@ -6,6 +6,17 @@ "id": "0", "metadata": {}, "outputs": [], + "source": [ + "%load_ext autoreload\n", + "%autoreload 2" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "1", + "metadata": {}, + "outputs": [], "source": [ "import optimizers\n", "import datamanip\n", @@ -20,7 +31,7 @@ { "cell_type": "code", "execution_count": null, - "id": "1", + "id": "2", "metadata": {}, "outputs": [], "source": [ @@ -33,7 +44,7 @@ { "cell_type": "code", "execution_count": null, - "id": "2", + "id": "3", "metadata": {}, "outputs": [], "source": [ @@ -53,7 +64,7 @@ { "cell_type": "code", "execution_count": null, - "id": "3", + "id": "4", "metadata": {}, "outputs": [], "source": [ @@ -88,7 +99,7 @@ { "cell_type": "code", "execution_count": null, - "id": "4", + "id": "5", "metadata": {}, "outputs": [], "source": [ @@ -101,7 +112,7 @@ { "cell_type": "code", "execution_count": null, - "id": "5", + "id": "6", "metadata": {}, "outputs": [], "source": [ @@ -112,7 +123,7 @@ { "cell_type": "code", "execution_count": null, - "id": "6", + "id": "7", "metadata": {}, "outputs": [], "source": [ @@ -148,7 +159,7 @@ { "cell_type": "code", "execution_count": null, - "id": "7", + "id": "8", "metadata": {}, "outputs": [], "source": [ @@ -165,7 +176,7 @@ { "cell_type": "code", "execution_count": null, - "id": "8", + "id": "9", "metadata": {}, "outputs": [], "source": [ @@ -202,7 +213,7 @@ { "cell_type": "code", "execution_count": null, - "id": "9", + "id": "10", "metadata": {}, "outputs": [], "source": [ @@ -222,7 +233,7 @@ { "cell_type": "code", "execution_count": null, - "id": "10", + "id": "11", "metadata": {}, "outputs": [], "source": [ @@ -249,7 +260,7 @@ { "cell_type": "code", "execution_count": null, - "id": "11", + "id": "12", "metadata": {}, "outputs": [], "source": [ @@ -291,21 +302,136 @@ { "cell_type": "code", "execution_count": null, - "id": "12", + "id": "13", "metadata": {}, "outputs": [], "source": [ - "RidgeMom = optimizers.RidgeMomentum(\n", - " learning_rate=0.1, num_iterations=1000, lam=0.1, delta=0.5\n", - ")\n", - "RidgeMom.fit(X_tr, y_tr)\n", - "plt.plot(RidgeMom.cost_history)" + "def get_cost_history(optimizer, X, y, **kwargs):\n", + " opt = optimizer(**kwargs)\n", + " opt.fit(X, y)\n", + " return opt.cost_history" ] }, { "cell_type": "code", "execution_count": null, - "id": "13", + "id": "14", + "metadata": {}, + "outputs": [], + "source": [ + "fig, (ax1, ax2) = plt.subplots(2, 1, figsize=plotting.get_figsize(0.8))\n", + "\n", + "X_train = datamanip.polynomial_features(x_train, 10, False)\n", + "X_test = datamanip.polynomial_features(x_test, 10, False)\n", + "X_tr, X_te = datamanip.scale_data(X_train, X_test)\n", + "y_tr, y_te = datamanip.scale_data(y_train, y_test)\n", + "\n", + "num_iters = 250\n", + "learning_rate_ols = 0.1\n", + "learning_rate_ridge = 0.01\n", + "lam = 0.1\n", + "\n", + "optimizers_ols = [\n", + " (optimizers.OLSGradientDescent, {\"learning_rate\": learning_rate_ols}),\n", + " (optimizers.OLSMomentum, {\"learning_rate\": learning_rate_ols, \"delta\": 0.9}),\n", + " (optimizers.OLSAdaGrad, {\"learning_rate\": learning_rate_ols}),\n", + " (optimizers.OLSRMSProp, {\"learning_rate\": learning_rate_ols, \"gamma\": 0.9}),\n", + " (\n", + " optimizers.OLSAdam,\n", + " {\"learning_rate\": learning_rate_ols, \"beta1\": 0.9, \"beta2\": 0.999},\n", + " ),\n", + "]\n", + "\n", + "optimizers_ridge = [\n", + " (\n", + " optimizers.RidgeGradientDescent,\n", + " {\"learning_rate\": learning_rate_ridge, \"lam\": lam},\n", + " ),\n", + " (\n", + " optimizers.RidgeMomentum,\n", + " {\"learning_rate\": learning_rate_ridge, \"delta\": 0.9, \"lam\": lam},\n", + " ),\n", + " (optimizers.RidgeAdaGrad, {\"learning_rate\": learning_rate_ridge, \"lam\": lam}),\n", + " (\n", + " optimizers.RidgeRMSProp,\n", + " {\"learning_rate\": learning_rate_ridge, \"gamma\": 0.9, \"lam\": lam},\n", + " ),\n", + " (\n", + " optimizers.RidgeAdam,\n", + " {\n", + " \"learning_rate\": learning_rate_ridge,\n", + " \"beta1\": 0.9,\n", + " \"beta2\": 0.999,\n", + " \"lam\": lam,\n", + " },\n", + " ),\n", + "]\n", + "\n", + "\n", + "def plot_optimizers(ax, optimizer_list, X, y, num_iters, ylabel=\"Cost\"):\n", + " for opt_class, params in optimizer_list:\n", + " history = get_cost_history(opt_class, X, y, num_iterations=num_iters, **params)\n", + " ax.plot(history, label=f\"{opt_class.__name__}\")\n", + " ax.set_xlabel(\"Iteration\")\n", + " ax.set_ylabel(ylabel)\n", + " ax.legend()\n", + "\n", + "\n", + "plot_optimizers(ax1, optimizers_ols, X_tr, y_tr, num_iters, ylabel=\"Cost (OLS)\")\n", + "ax1.set_ylim(bottom=0.465, top=0.505)\n", + "\n", + "plot_optimizers(ax2, optimizers_ridge, X_tr, y_tr, num_iters, ylabel=\"Cost (Ridge)\")\n", + "ax2.set_ylim(bottom=0.4775, top=0.505)\n", + "fig.tight_layout()\n", + "fig.savefig(os.path.join(FIG_DIR, \"optimizer_comparison.pdf\"))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "15", + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=plotting.get_figsize(0.5))\n", + "\n", + "ols_adam = optimizers.OLSAdam(learning_rate=0.01, beta1=0.9, beta2=0.999)\n", + "ridge_adam = optimizers.RidgeAdam(learning_rate=0.01, beta1=0.9, beta2=0.999, lam=0.1)\n", + "lasso_adam = optimizers.LASSOAdam(learning_rate=0.01, beta1=0.9, beta2=0.999, lam=0.1)\n", + "ols_gd = optimizers.OLSGradientDescent(learning_rate=0.1)\n", + "ridge_gd = optimizers.RidgeGradientDescent(learning_rate=0.01, lam=0.1)\n", + "lasso_gd = optimizers.LASSOGradientDescent(learning_rate=0.01, lam=0.1)\n", + "\n", + "for i, (optimizer, label) in enumerate(\n", + " zip(\n", + " [ols_adam, ridge_adam, lasso_adam, ols_gd, ridge_gd, lasso_gd],\n", + " [\"OLS Adam\", \"Ridge Adam\", \"LASSO Adam\", \"OLS GD\", \"Ridge GD\", \"LASSO GD\"],\n", + " )\n", + "):\n", + " optimizer.fit(X_tr, y_tr)\n", + " cost_history = optimizer.cost_history\n", + " ls = \"-\" if i < 3 else \"--\"\n", + " c = \"C\" + str(i % 3)\n", + " label = label.split()[0] if i < 3 else \"\"\n", + " ax.plot(cost_history, label=label, linestyle=ls, color=c)\n", + "ax.set_xlabel(\"Iteration\")\n", + "ax.set_ylabel(\"Cost\")\n", + "# Add linestyles manually to legend\n", + "handles, labels = ax.get_legend_handles_labels()\n", + "handles += [\n", + " plt.Line2D([0], [0], color=\"k\", linestyle=\"-\"),\n", + " plt.Line2D([0], [0], color=\"k\", linestyle=\"--\"),\n", + "]\n", + "labels += [\"Adam Opt.\", \"Grad. Desc.\"]\n", + "ax.legend(handles, labels)\n", + "fig.tight_layout()\n", + "fig.savefig(os.path.join(FIG_DIR, \"cost_function_comparison.pdf\"))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "16", "metadata": {}, "outputs": [], "source": [] diff --git a/src/optimizers.py b/src/optimizers.py index 011ee57..8685473 100644 --- a/src/optimizers.py +++ b/src/optimizers.py @@ -123,3 +123,124 @@ class OLSMomentum(OLSGradientDescent): class RidgeMomentum(OLSMomentum, RidgeGradientDescent): def __init__(self, *args, **kwargs): super().__init__(*args, **kwargs) + + +class OLSAdaGrad(OLSGradientDescent): + def __init__(self, *args, learning_rate=0.1, num_iterations=1000, **kwargs): + super().__init__( + *args, learning_rate=learning_rate, num_iterations=num_iterations, **kwargs + ) + + def _precomp(self): + self.quad_sum = np.zeros_like(self.theta) + return super()._precomp() + + def _update_theta(self): + grad = self._compute_grad() + self.quad_sum += np.square(grad) + self.theta -= ( + self.learning_rate * grad / (np.sqrt(self.quad_sum) + 1e-10) + ) # to avoid division by zero + + +class RidgeAdaGrad(OLSAdaGrad, RidgeGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + +class OLSRMSProp(OLSGradientDescent): + def __init__( + self, *args, learning_rate=0.1, num_iterations=1000, gamma=0.9, **kwargs + ): + super().__init__( + *args, learning_rate=learning_rate, num_iterations=num_iterations, **kwargs + ) + self.gamma = gamma + + def _precomp(self): + self.prev_v = np.zeros_like(self.theta) + return super()._precomp() + + def _update_theta(self): + grad = self._compute_grad() + v = self.gamma * self.prev_v + (1 - self.gamma) * np.square(grad) + self.theta -= ( + self.learning_rate * grad / (np.sqrt(v) + 1e-10) + ) # to avoid division by zero + self.prev_v = v + + +class RidgeRMSProp(OLSRMSProp, RidgeGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + +class OLSAdam(OLSGradientDescent): + def __init__( + self, + *args, + learning_rate=0.1, + num_iterations=1000, + beta1=0.9, + beta2=0.999, + **kwargs, + ): + super().__init__( + *args, learning_rate=learning_rate, num_iterations=num_iterations, **kwargs + ) + self.beta1 = beta1 + self.beta2 = beta2 + + def _precomp(self): + self.m = np.zeros_like(self.theta) + self.v = np.zeros_like(self.theta) + return super()._precomp() + + def _update_theta(self): + grad = self._compute_grad() + self.m = self.beta1 * self.m + (1 - self.beta1) * grad + self.v = self.beta2 * self.v + (1 - self.beta2) * np.square(grad) + + m_hat = self.m / (1 - self.beta1 ** (self.num_iterations)) + v_hat = self.v / (1 - self.beta2 ** (self.num_iterations)) + + self.theta -= ( + self.learning_rate * m_hat / (np.sqrt(v_hat) + 1e-10) + ) # to avoid division by zero + + +class RidgeAdam(OLSAdam, RidgeGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + +class LASSOGradientDescent(OLSGradientDescent): + def __init__(self, *args, lam: float = 0.1, **kwargs): + super().__init__(*args, **kwargs) + self.lam = lam + + def _compute_cost(self) -> float: + return super()._compute_cost() + self.lam * np.sum(np.abs(self.theta)) + + def _compute_grad(self) -> np.ndarray: + return super()._compute_grad() + self.lam * np.sign(self.theta) + + +class LASSOMomentum(OLSMomentum, LASSOGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + +class LASSOAdaGrad(OLSAdaGrad, LASSOGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + +class LASSORMSProp(OLSRMSProp, LASSOGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) + + +class LASSOAdam(OLSAdam, LASSOGradientDescent): + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs)