diff --git a/doc/pub/Splines/html/._Splines-bs001.html b/doc/pub/Splines/html/._Splines-bs001.html index 49705e442..ad05fabe2 100644 --- a/doc/pub/Splines/html/._Splines-bs001.html +++ b/doc/pub/Splines/html/._Splines-bs001.html @@ -178,11 +178,11 @@ MathJax.Hub.Config({

Almost every problem in machine learning and data science starts with -a dataset \( X \), a model \( g(\theta) \), which is a function of the -parameters \( \theta \) and a cost function \( C(X, g(\theta)) \) that allows -us to judge how well the model \( g(\theta) \) explains the observations -\( X \). The model is fit by finding the values of \( \theta \) that minimize -the cost function. Ideally we would be able to solve for \( \theta \) +a dataset \( X \), a model \( g(\beta) \), which is a function of the +parameters \( \beta \) and a cost function \( C(X, g(\beta)) \) that allows +us to judge how well the model \( g(\beta) \) explains the observations +\( X \). The model is fit by finding the values of \( \beta \) that minimize +the cost function. Ideally we would be able to solve for \( \beta \) analytically, however this is not possible in general and we must use some approximative/numerical method to compute the minimum. diff --git a/doc/pub/Splines/html/._Splines-bs003.html b/doc/pub/Splines/html/._Splines-bs003.html index 3580baca9..90b01b5c8 100644 --- a/doc/pub/Splines/html/._Splines-bs003.html +++ b/doc/pub/Splines/html/._Splines-bs003.html @@ -180,7 +180,7 @@ MathJax.Hub.Config({ The previous observation is the basis of the method of steepest descent, which is also referred to as just gradient descent (GD). One starts with an initial guess \( \mathbf{x}_0 \) for a minimum of \( F \) and -compute new approximations according to +computes new approximations according to $$ \mathbf{x}_{k+1} = \mathbf{x}_k - \gamma_k \nabla F(\mathbf{x}_k), \ \ k \geq 0. @@ -188,7 +188,7 @@ $$

The parameter \( \gamma_k \) is often referred to as the step length or -the learning rate in the context of Machine Learning. +the learning rate within the context of Machine Learning.

diff --git a/doc/pub/Splines/html/._Splines-bs004.html b/doc/pub/Splines/html/._Splines-bs004.html index f85e6f20c..9484ecc0c 100644 --- a/doc/pub/Splines/html/._Splines-bs004.html +++ b/doc/pub/Splines/html/._Splines-bs004.html @@ -177,10 +177,21 @@ MathJax.Hub.Config({

The ideal

-Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global minimum of the function \( F \). In general we do not know if we are in a global or local minimum. In the special case when \( F \) is a convex function, all local minima are also global minima, so in this case gradient descent can converge to the global solution. The advantage of this scheme is that it is conceptually simple and straightforward to implement. However the method in this form has some severe limitations: +Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global +minimum of the function \( F \). In general we do not know if we are in a +global or local minimum. In the special case when \( F \) is a convex +function, all local minima are also global minima, so in this case +gradient descent can converge to the global solution. The advantage of +this scheme is that it is conceptually simple and straightforward to +implement. However the method in this form has some severe +limitations:

-In machine learing we are often faced with non-convex high dimensional cost functions with many local minimum. Since GD is deterministic we will get stuck in a local minimum, if the method converges, unless we have a very good intial guess. This also implies that the scheme is sensitive to the chosen initial condition. +In machine learing we are often faced with non-convex high dimensional +cost functions with many local minima. Since GD is deterministic we +will get stuck in a local minimum, if the method converges, unless we +have a very good intial guess. This also implies that the scheme is +sensitive to the chosen initial condition.

Note that the gradient is a function of \( \mathbf{x} = diff --git a/doc/pub/Splines/html/._Splines-bs005.html b/doc/pub/Splines/html/._Splines-bs005.html index a0db40746..f9a3208bd 100644 --- a/doc/pub/Splines/html/._Splines-bs005.html +++ b/doc/pub/Splines/html/._Splines-bs005.html @@ -180,14 +180,14 @@ MathJax.Hub.Config({ GD is sensitive to the choice of learning rate \( \gamma_k \). This is due to the fact that we are only guaranteed that \( F(\mathbf{x}_{k+1}) \leq F(\mathbf{x}_k) \) for sufficiently small \( \gamma_k \). The problem is to -determine an optimal learning rate. If the learning rate is chosen to -small the method will take a long to converge and if it is to large we -can experience erratic behavior. +determine an optimal learning rate. If the learning rate is chosen too +small the method will take a long time to converge and if it is too +large we can experience erratic behavior.

Many of these shortcomings can be alleviated by introducing randomness. One such method is that of Stochastic Gradient Descent -(SGD), see below +(SGD), see below.

diff --git a/doc/pub/Splines/html/._Splines-bs006.html b/doc/pub/Splines/html/._Splines-bs006.html index 9126bdff5..306b99cdb 100644 --- a/doc/pub/Splines/html/._Splines-bs006.html +++ b/doc/pub/Splines/html/._Splines-bs006.html @@ -175,6 +175,8 @@ MathJax.Hub.Config({

Gradient Descent Example

+ +

We revisit now our simple linear regression example with a linear polynomial.

@@ -192,30 +194,30 @@ x = 2*np y = 4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) -theta = np.random.randn(2,1) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) +beta = np.random.randn(2,1) eta = 0.1 Niterations = 1000 m = 100 for iter in range(Niterations): - gradients = 2.0/m*xb.T.dot(xb.dot(theta)-y) - theta -= eta*gradients + gradients = 2.0/m*xb.T.dot(xb.dot(beta)-y) + beta -= eta*gradients -print(theta) +print(beta) xnew = np.array([[0],[2]]) xbnew = np.c_[np.ones((2,1)), xnew] -ypredict = xbnew.dot(theta) -ypredict2 = xbnew.dot(theta_linreg) +ypredict = xbnew.dot(beta) +ypredict2 = xbnew.dot(beta_linreg) plt.plot(xnew, ypredict, "r-") plt.plot(xnew, ypredict2, "b-") plt.plot(x, y ,'ro') plt.axis([0,2.0,0, 15.0]) plt.xlabel(r'$x$') plt.ylabel(r'$y$') -plt.title(r'Random numbers ') +plt.title(r'Gradient descent example') plt.show()

diff --git a/doc/pub/Splines/html/._Splines-bs007.html b/doc/pub/Splines/html/._Splines-bs007.html index 686d4fb1a..6b2852c4c 100644 --- a/doc/pub/Splines/html/._Splines-bs007.html +++ b/doc/pub/Splines/html/._Splines-bs007.html @@ -189,8 +189,8 @@ x = 2*np y = 4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) sgdreg = SGDRegressor(n_iter = 50, penalty=None, eta0=0.1) sgdreg.fit(x,y.ravel()) print(sgdreg.intercept_, sgdreg.coef_) diff --git a/doc/pub/Splines/html/._Splines-bs008.html b/doc/pub/Splines/html/._Splines-bs008.html index b01eef450..86e9ac845 100644 --- a/doc/pub/Splines/html/._Splines-bs008.html +++ b/doc/pub/Splines/html/._Splines-bs008.html @@ -175,6 +175,8 @@ MathJax.Hub.Config({

Convex functions

+ +

Ideally we want our cost/loss function to be convex(concave).

diff --git a/doc/pub/Splines/html/._Splines-bs009.html b/doc/pub/Splines/html/._Splines-bs009.html index bf637803e..7ae79d754 100644 --- a/doc/pub/Splines/html/._Splines-bs009.html +++ b/doc/pub/Splines/html/._Splines-bs009.html @@ -177,7 +177,7 @@ MathJax.Hub.Config({

Convex function

-Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below. +Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below.

diff --git a/doc/pub/Splines/html/._Splines-bs010.html b/doc/pub/Splines/html/._Splines-bs010.html index 07406c768..c1bb98aa7 100644 --- a/doc/pub/Splines/html/._Splines-bs010.html +++ b/doc/pub/Splines/html/._Splines-bs010.html @@ -180,9 +180,7 @@ MathJax.Hub.Config({ In the following we state first and second-order conditions which ensures convexity of a function \( f \). We write \( D_f \) to denote the domain of \( f \), i.e the subset of \( R^n \) where \( f \) is defined. For more -details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex -Optimization. Cambridge University Press, http://stanford.edu/ -boyd/cvxbook/, 2004. +details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press.

diff --git a/doc/pub/Splines/html/._Splines-bs011.html b/doc/pub/Splines/html/._Splines-bs011.html index 36831b01e..8a0925577 100644 --- a/doc/pub/Splines/html/._Splines-bs011.html +++ b/doc/pub/Splines/html/._Splines-bs011.html @@ -180,9 +180,12 @@ MathJax.Hub.Config({ The next result is of great importance to us and the reason why we are going on about convex functions. In machine learning we frequently have to minimize a loss/cost function in order to find the best -parameters for the model we are considering. Ideally we want the -global minimum, however for high-dimensional models it is hard to know -if we have local or global minimum. However, if the cost/loss function +parameters for the model we are considering. + +

+Ideally we want the +global minimum (for high-dimensional models it is hard to know +if we have local or global minimum). However, if the cost/loss function is convex the following result provides invaluable information:

@@ -195,6 +198,8 @@ is minimal, where \( f \) is convex and differentiable. Then, any point

+ +

This result means that if we know that the cost/loss function is convex and we are able to find a minimum, we are guaranteed that it is a global minimum.

diff --git a/doc/pub/Splines/html/._Splines-bs012.html b/doc/pub/Splines/html/._Splines-bs012.html index 65fbf0e4a..29a039b11 100644 --- a/doc/pub/Splines/html/._Splines-bs012.html +++ b/doc/pub/Splines/html/._Splines-bs012.html @@ -177,27 +177,26 @@ MathJax.Hub.Config({

Some simple problems

    -
  1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.
  2. -
- -Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $$\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $$ - -
    +
  1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity. Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $
  2. Using the second order condition show that the following functions are convex on the specified domain.
  3. -
-\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \). -\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \). + -
  1. Let \( f(x) = x^2 \) and \( g(x) = e^x \). Show that \( f(g(x)) \) and \( g(f(x)) \) is convex for \( x \in \mathbb{R} \). Also show that if \( f(x) \) is any convex function than \( h(x) = e^{f(x)} \) is convex.
  2. A norm is any function that satisfy the following properties
  3. + + +
-\( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \). -\( f(x+y) \leq f(x) + f(y) \) -\( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \) -Using the definition of convexity, show that a function satisfying the properties above is convex (the third condition is not needed to show this). +Using the definition of convexity, try to show that a function satisfying the properties above is convex (the third condition is not needed to show this).

diff --git a/doc/pub/Splines/html/._Splines-bs013.html b/doc/pub/Splines/html/._Splines-bs013.html index 742f4f3de..2a7103a97 100644 --- a/doc/pub/Splines/html/._Splines-bs013.html +++ b/doc/pub/Splines/html/._Splines-bs013.html @@ -177,7 +177,10 @@ MathJax.Hub.Config({

Revisiting our first homework

-We will use linear regression as a case study for the gradient descent methods. Linear regression is a great test case for the gradient descent methods discussed in the lectures since it has several desirable properties such as: +We will use linear regression as a case study for the gradient descent +methods. Linear regression is a great test case for the gradient +descent methods discussed in the lectures since it has several +desirable properties such as:

  1. An analytical solution (recall homework set 1).
  2. @@ -193,12 +196,12 @@ $$ with \( x_i \in [0,1] \) chosen randomly with a uniform distribution. Additionally \( \xi_i \) represents stochastic noise chosen according to a normal distribution \( \cal {N}(0,1) \). The linear regression model is given by $$ -h_\theta(x) = \hat{y} = \theta_0 + \theta_1 x, +h_\beta(x) = \hat{y} = \beta_0 + \beta_1 x, $$ such that $$ -\hat{y}_i = \theta_0 + \theta_1 x_i. +\hat{y}_i = \beta_0 + \beta_1 x_i. $$

    diff --git a/doc/pub/Splines/html/._Splines-bs014.html b/doc/pub/Splines/html/._Splines-bs014.html index 18989a1d1..596f81cd4 100644 --- a/doc/pub/Splines/html/._Splines-bs014.html +++ b/doc/pub/Splines/html/._Splines-bs014.html @@ -177,10 +177,10 @@ MathJax.Hub.Config({

    Gradient descent example

    -Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \theta = (\theta_0, \theta_1)^T \) +Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \beta = (\beta_0, \beta_1)^T \)

    -t is convenient to write \( \mathbf{\hat{y}} = X\theta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by +t is convenient to write \( \mathbf{\hat{y}} = X\beta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by $$ \begin{equation} X \equiv \begin{bmatrix} @@ -194,10 +194,10 @@ $$ The loss function is given by $$ -C(\theta) = ||X\theta-\mathbf{y}||^2 = ||X\theta||^2 - 2 \mathbf{y}^T X\theta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\theta_0 + \theta_1 x_i)^2 - 2 y_i (\theta_0 + \theta_1 x_i) + y_i^2 +C(\beta) = ||X\beta-\mathbf{y}||^2 = ||X\beta||^2 - 2 \mathbf{y}^T X\beta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\beta_0 + \beta_1 x_i)^2 - 2 y_i (\beta_0 + \beta_1 x_i) + y_i^2 $$ -and we want to find \( \theta \) such that \( C(\theta) \) is minimized. +and we want to find \( \beta \) such that \( C(\beta) \) is minimized.

    diff --git a/doc/pub/Splines/html/._Splines-bs015.html b/doc/pub/Splines/html/._Splines-bs015.html index c0ed6b688..1778f7a66 100644 --- a/doc/pub/Splines/html/._Splines-bs015.html +++ b/doc/pub/Splines/html/._Splines-bs015.html @@ -177,11 +177,11 @@ MathJax.Hub.Config({

    The derivative of the cost/loss function

    -Computing \( \partial C(\theta) / \partial \theta_0 \) and \( \partial C(\theta) / \partial \theta_1 \) we can show that the gradient can be written as +Computing \( \partial C(\beta) / \partial \beta_0 \) and \( \partial C(\beta) / \partial \beta_1 \) we can show that the gradient can be written as $$ -\nabla_\theta C(\theta) = (\partial C(\theta) / \partial \theta_0, \partial C(\theta) / \partial \theta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} = 2X^T(X\theta - \mathbf{y}), +\nabla_\beta C(\beta) = (\partial C(\beta) / \partial \beta_0, \partial C(\beta) / \partial \beta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} = 2X^T(X\beta - \mathbf{y}), $$ where \( X \) is the design matrix defined above. diff --git a/doc/pub/Splines/html/._Splines-bs016.html b/doc/pub/Splines/html/._Splines-bs016.html index beaf9cc37..d4a0a6cb0 100644 --- a/doc/pub/Splines/html/._Splines-bs016.html +++ b/doc/pub/Splines/html/._Splines-bs016.html @@ -175,15 +175,15 @@ MathJax.Hub.Config({

    The Hessian matrix

    -The Hessian matrix of \( C(\theta) \) is given by +The Hessian matrix of \( C(\beta) \) is given by $$ \hat{H} \equiv \begin{bmatrix} -\frac{\partial^2 C(\theta)}{\partial \theta_0^2} & \frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} \\ -\frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} & \frac{\partial^2 C(\theta)}{\partial \theta_1^2} & \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0^2} & \frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} & \frac{\partial^2 C(\beta)}{\partial \beta_1^2} & \\ \end{bmatrix} = 2X^T X. $$ -This result implies that \( C(\theta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite. +This result implies that \( C(\beta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite.

    diff --git a/doc/pub/Splines/html/._Splines-bs017.html b/doc/pub/Splines/html/._Splines-bs017.html index 5eef83c31..d02aaae45 100644 --- a/doc/pub/Splines/html/._Splines-bs017.html +++ b/doc/pub/Splines/html/._Splines-bs017.html @@ -177,19 +177,19 @@ MathJax.Hub.Config({

    Simple program

    -We can now write a program that minimizes \( C(\theta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to +We can now write a program that minimizes \( C(\beta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to $$ -\theta_{k+1} = \theta_k - \gamma \nabla_\theta C(\theta_k), \ k=0,1,\cdots +\beta_{k+1} = \beta_k - \gamma \nabla_\beta C(\beta_k), \ k=0,1,\cdots $$

    We can use the expression we computed for the gradient and let use a -\( \theta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating -when \( ||\nabla_\theta C(\theta_k) || < \epsilon = 10^{-8} \). +\( \beta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating +when \( ||\nabla_\beta C(\beta_k) || < \epsilon = 10^{-8} \).

    -And finally we can compare our solution for \( \theta \) with the analytic result given by -\( \theta= (X^TX)^{-1} X^T \mathbf{y} \). +And finally we can compare our solution for \( \beta \) with the analytic result given by +\( \beta= (X^TX)^{-1} X^T \mathbf{y} \).

    @@ -206,11 +206,11 @@ x = np. y = 5*x**2 + 0.1*np.random.randn(N) X = np.c_[np.ones(N),x] #Construct design matrix -#Compute theta according to normal equations to compare with GD solution +#Compute beta according to normal equations to compare with GD solution Xt_X_inv = np.linalg.inv(np.dot(X.T,X)) Xt_y = np.dot(X.transpose(),y) -theta_NE = np.dot(Xt_X_inv,Xt_y) -print(theta_NE) +beta_NE = np.dot(Xt_X_inv,Xt_y) +print(beta_NE)

    diff --git a/doc/pub/Splines/html/._Splines-bs018.html b/doc/pub/Splines/html/._Splines-bs018.html index e35445845..878874987 100644 --- a/doc/pub/Splines/html/._Splines-bs018.html +++ b/doc/pub/Splines/html/._Splines-bs018.html @@ -177,27 +177,27 @@ MathJax.Hub.Config({

    Gradient descent and Ridge

    -We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \theta \), +We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \beta \), $$ -C_{\text{ridge}}(\theta) = ||X\theta -\mathbf{y}||^2 + \lambda ||\theta||^2, \ \lambda \geq 0. +C_{\text{ridge}}(\beta) = ||X\beta -\mathbf{y}||^2 + \lambda ||\beta||^2, \ \lambda \geq 0. $$

    -In order to minimize \( C_{\text{ridge}}(\theta) \) using GD we only have adjust the gradient as follows +In order to minimize \( C_{\text{ridge}}(\beta) \) using GD we only have adjust the gradient as follows $$ -\nabla_\theta C_{\text{ridge}}(\theta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} + 2\lambda\begin{bmatrix} \theta_0 \\ \theta_1\end{bmatrix} = 2 (X^T(X\theta - \mathbf{y})+\lambda \theta). +\nabla_\beta C_{\text{ridge}}(\beta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} + 2\lambda\begin{bmatrix} \beta_0 \\ \beta_1\end{bmatrix} = 2 (X^T(X\beta - \mathbf{y})+\lambda \beta). $$

    -We can now extend our program to minimize \( C_{\text{ridge}}(\theta) \) using gradient descent and compare with the analytical solution given by +We can now extend our program to minimize \( C_{\text{ridge}}(\beta) \) using gradient descent and compare with the analytical solution given by $$ -\theta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, +\beta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, $$ for \( \lambda = {0,1,10,50,100} \) (\( \lambda = 0 \) corresponds to ordinary least squares). -We can then compute \( ||\theta_{\text{ridge}}|| \) for each \( \lambda \). +We can then compute \( ||\beta_{\text{ridge}}|| \) for each \( \lambda \).

    @@ -215,7 +215,7 @@ x = np. y = 5*x**2 + 0.1*np.random.randn(N) -#Compute analytic theta for Ridge regression +#Compute analytic beta for Ridge regression X = np.c_[np.ones(N),x] XT_X = np.dot(X.T,X) @@ -223,10 +223,10 @@ l = 0.1 Id = np.eye(XT_X.shape[0]) Z = np.linalg.inv(XT_X+l*Id) -theta_ridge = np.dot(Z,np.dot(X.T,y)) +beta_ridge = np.dot(Z,np.dot(X.T,y)) -print(theta_ridge) -print(np.linalg.norm(theta_ridge)) #||theta|| +print(beta_ridge) +print(np.linalg.norm(beta_ridge)) #||beta||

    diff --git a/doc/pub/Splines/html/._Splines-bs019.html b/doc/pub/Splines/html/._Splines-bs019.html index 30a05cbaa..464a234f4 100644 --- a/doc/pub/Splines/html/._Splines-bs019.html +++ b/doc/pub/Splines/html/._Splines-bs019.html @@ -185,8 +185,8 @@ The underlying idea of SGD comes from the observation that the cost function, which we want to minimize, can almost always be written as a sum over \( n \) datapoints \( \{\mathbf{x}_i\}_{i=1}^n \), $$ -C(\mathbf{\theta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, -\mathbf{\theta}). +C(\mathbf{\beta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

    diff --git a/doc/pub/Splines/html/._Splines-bs020.html b/doc/pub/Splines/html/._Splines-bs020.html index c92093207..953411294 100644 --- a/doc/pub/Splines/html/._Splines-bs020.html +++ b/doc/pub/Splines/html/._Splines-bs020.html @@ -180,8 +180,8 @@ MathJax.Hub.Config({ This in turn means that the gradient can be computed as a sum over \( i \)-gradients $$ -\nabla_\theta C(\mathbf{\theta}) = \sum_i^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}). +\nabla_\beta C(\mathbf{\beta}) = \sum_i^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

    diff --git a/doc/pub/Splines/html/._Splines-bs021.html b/doc/pub/Splines/html/._Splines-bs021.html index 9014dbac7..062ea2901 100644 --- a/doc/pub/Splines/html/._Splines-bs021.html +++ b/doc/pub/Splines/html/._Splines-bs021.html @@ -189,10 +189,10 @@ The idea is now to approximate the gradient by replacing the sum over all datapoints with a sum over the datapoints in one the minibatches picked at random in each gradient descent step $$ -\nabla_\theta -C(\mathbf{\theta}) = \sum_{i=1}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) \rightarrow \sum_{i \in B_k}^n \nabla_\theta -c_i(\mathbf{x}_i, \mathbf{\theta}). +\nabla_\beta +C(\mathbf{\beta}) = \sum_{i=1}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) \rightarrow \sum_{i \in B_k}^n \nabla_\beta +c_i(\mathbf{x}_i, \mathbf{\beta}). $$

    diff --git a/doc/pub/Splines/html/._Splines-bs022.html b/doc/pub/Splines/html/._Splines-bs022.html index 141208bfc..6dc924616 100644 --- a/doc/pub/Splines/html/._Splines-bs022.html +++ b/doc/pub/Splines/html/._Splines-bs022.html @@ -179,8 +179,8 @@ MathJax.Hub.Config({

    Thus a gradient descent step now looks like $$ -\theta_{j+1} = \theta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) +\beta_{j+1} = \beta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) $$

    diff --git a/doc/pub/Splines/html/._Splines-bs023.html b/doc/pub/Splines/html/._Splines-bs023.html index 4e7d3649b..9c98ea73d 100644 --- a/doc/pub/Splines/html/._Splines-bs023.html +++ b/doc/pub/Splines/html/._Splines-bs023.html @@ -191,7 +191,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for j += 1

    diff --git a/doc/pub/Splines/html/._Splines-bs024.html b/doc/pub/Splines/html/._Splines-bs024.html index e2f98d13c..1653551d2 100644 --- a/doc/pub/Splines/html/._Splines-bs024.html +++ b/doc/pub/Splines/html/._Splines-bs024.html @@ -185,7 +185,7 @@ is zero is valid also for local minima, so this would only tell us that we are close to a local/global minimum. However, we could also evaluate the cost function at this point, store the result and continue the search. If the test kicks in at a later stage we can -compare the values of the cost function and keep the \( \theta \) that +compare the values of the cost function and keep the \( \beta \) that gave the lowest value.

    diff --git a/doc/pub/Splines/html/._Splines-bs025.html b/doc/pub/Splines/html/._Splines-bs025.html index 308a160c0..49c8b196a 100644 --- a/doc/pub/Splines/html/._Splines-bs025.html +++ b/doc/pub/Splines/html/._Splines-bs025.html @@ -185,10 +185,10 @@ reasonable time such that we do not move at all. As an example, let \( e = 0,1,2,3,\cdots \) denote the current epoch and let \( t_0, t_1 > 0 \) be two fixed numbers. Furthermore, let \( t = e \cdot m + i \) where \( m \) is the number of minibatches and \( i=0,\cdots,m-1 \). Then the function $$\gamma_j(t; t_0, t_1) = \frac{t_0}{t+t_1} $$ goes to zero as the number of epochs gets large. I.e. we start with a step length \( \gamma_j (0; t_0, t_1) = t_0/t_1 \) which decays in time \( t \).

    -In this way we can fix the number of epochs, compute \( \theta \) and +In this way we can fix the number of epochs, compute \( \beta \) and evaluate the cost function at the end. Repeating the computation will give a different result since the scheme is random by design. Then we -pick the final \( \theta \) that gives the lowest value of the cost +pick the final \( \beta \) that gives the lowest value of the cost function.

    @@ -212,7 +212,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for beta t = epoch*m+i gamma_j = step_length(t,t0,t1) j += 1 diff --git a/doc/pub/Splines/html/Splines-reveal.html b/doc/pub/Splines/html/Splines-reveal.html index 7834dc5b1..f78ef0d25 100644 --- a/doc/pub/Splines/html/Splines-reveal.html +++ b/doc/pub/Splines/html/Splines-reveal.html @@ -163,11 +163,11 @@ MathJax.Hub.Config({

    Almost every problem in machine learning and data science starts with -a dataset \( X \), a model \( g(\theta) \), which is a function of the -parameters \( \theta \) and a cost function \( C(X, g(\theta)) \) that allows -us to judge how well the model \( g(\theta) \) explains the observations -\( X \). The model is fit by finding the values of \( \theta \) that minimize -the cost function. Ideally we would be able to solve for \( \theta \) +a dataset \( X \), a model \( g(\beta) \), which is a function of the +parameters \( \beta \) and a cost function \( C(X, g(\beta)) \) that allows +us to judge how well the model \( g(\beta) \) explains the observations +\( X \). The model is fit by finding the values of \( \beta \) that minimize +the cost function. Ideally we would be able to solve for \( \beta \) analytically, however this is not possible in general and we must use some approximative/numerical method to compute the minimum. @@ -204,7 +204,7 @@ we are always moving towards smaller function values, i.e a minimum. The previous observation is the basis of the method of steepest descent, which is also referred to as just gradient descent (GD). One starts with an initial guess \( \mathbf{x}_0 \) for a minimum of \( F \) and -compute new approximations according to +computes new approximations according to

     
    $$ @@ -214,7 +214,7 @@ $$

    The parameter \( \gamma_k \) is often referred to as the step length or -the learning rate in the context of Machine Learning. +the learning rate within the context of Machine Learning. @@ -222,10 +222,21 @@ the learning rate in the context of Machine Learning.

    The ideal

    -Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global minimum of the function \( F \). In general we do not know if we are in a global or local minimum. In the special case when \( F \) is a convex function, all local minima are also global minima, so in this case gradient descent can converge to the global solution. The advantage of this scheme is that it is conceptually simple and straightforward to implement. However the method in this form has some severe limitations: +Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global +minimum of the function \( F \). In general we do not know if we are in a +global or local minimum. In the special case when \( F \) is a convex +function, all local minima are also global minima, so in this case +gradient descent can converge to the global solution. The advantage of +this scheme is that it is conceptually simple and straightforward to +implement. However the method in this form has some severe +limitations:

    -In machine learing we are often faced with non-convex high dimensional cost functions with many local minimum. Since GD is deterministic we will get stuck in a local minimum, if the method converges, unless we have a very good intial guess. This also implies that the scheme is sensitive to the chosen initial condition. +In machine learing we are often faced with non-convex high dimensional +cost functions with many local minima. Since GD is deterministic we +will get stuck in a local minimum, if the method converges, unless we +have a very good intial guess. This also implies that the scheme is +sensitive to the chosen initial condition.

    Note that the gradient is a function of \( \mathbf{x} = @@ -240,19 +251,21 @@ Note that the gradient is a function of \( \mathbf{x} = GD is sensitive to the choice of learning rate \( \gamma_k \). This is due to the fact that we are only guaranteed that \( F(\mathbf{x}_{k+1}) \leq F(\mathbf{x}_k) \) for sufficiently small \( \gamma_k \). The problem is to -determine an optimal learning rate. If the learning rate is chosen to -small the method will take a long to converge and if it is to large we -can experience erratic behavior. +determine an optimal learning rate. If the learning rate is chosen too +small the method will take a long time to converge and if it is too +large we can experience erratic behavior.

    Many of these shortcomings can be alleviated by introducing randomness. One such method is that of Stochastic Gradient Descent -(SGD), see below +(SGD), see below.

    Gradient Descent Example

    + +

    We revisit now our simple linear regression example with a linear polynomial.

    @@ -270,30 +283,30 @@ x = 2*np.random.rand(4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) -theta = np.random.randn(2,1) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) +beta = np.random.randn(2,1) eta = 0.1 Niterations = 1000 m = 100 for iter in range(Niterations): - gradients = 2.0/m*xb.T.dot(xb.dot(theta)-y) - theta -= eta*gradients + gradients = 2.0/m*xb.T.dot(xb.dot(beta)-y) + beta -= eta*gradients -print(theta) +print(beta) xnew = np.array([[0],[2]]) xbnew = np.c_[np.ones((2,1)), xnew] -ypredict = xbnew.dot(theta) -ypredict2 = xbnew.dot(theta_linreg) +ypredict = xbnew.dot(beta) +ypredict2 = xbnew.dot(beta_linreg) plt.plot(xnew, ypredict, "r-") plt.plot(xnew, ypredict2, "b-") plt.plot(x, y ,'ro') plt.axis([0,2.0,0, 15.0]) plt.xlabel(r'$x$') plt.ylabel(r'$y$') -plt.title(r'Random numbers ') +plt.title(r'Gradient descent example') plt.show()

    @@ -315,8 +328,8 @@ x = 2*np.random.rand(4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) sgdreg = SGDRegressor(n_iter = 50, penalty=None, eta0=0.1) sgdreg.fit(x,y.ravel()) print(sgdreg.intercept_, sgdreg.coef_) @@ -326,6 +339,8 @@ sgdreg.fit(x,y.ravel())

    Convex functions

    + +

    Ideally we want our cost/loss function to be convex(concave).

    @@ -346,7 +361,7 @@ regular polygons (triangles, rectangles, pentagons, etc...).

    Convex function

    -Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if

     
    +Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if

     
    $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$

     
    for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below.

    @@ -359,9 +374,7 @@ $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ In the following we state first and second-order conditions which ensures convexity of a function \( f \). We write \( D_f \) to denote the domain of \( f \), i.e the subset of \( R^n \) where \( f \) is defined. For more -details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex -Optimization. Cambridge University Press, http://stanford.edu/ -boyd/cvxbook/, 2004. +details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press.

    @@ -404,9 +417,12 @@ This condition is particularly useful since it gives us an procedure for determi The next result is of great importance to us and the reason why we are going on about convex functions. In machine learning we frequently have to minimize a loss/cost function in order to find the best -parameters for the model we are considering. Ideally we want the -global minimum, however for high-dimensional models it is hard to know -if we have local or global minimum. However, if the cost/loss function +parameters for the model we are considering. + +

    +Ideally we want the +global minimum (for high-dimensional models it is hard to know +if we have local or global minimum). However, if the cost/loss function is convex the following result provides invaluable information:

    @@ -418,6 +434,7 @@ is minimal, where \( f \) is convex and differentiable. Then, any point \( x^* \) that satisfies \( \nabla f(x^*) = 0 \) is a global minimum.

    +

    This result means that if we know that the cost/loss function is convex and we are able to find a minimum, we are guaranteed that it is a global minimum. @@ -426,32 +443,26 @@ This result means that if we know that the cost/loss function is convex and we a

    Some simple problems

      -

    1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.
    2. -
    -

    - -Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \)

     
    -$$\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $$ -

     
    - -

      +

    1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity. Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $
    2. Using the second order condition show that the following functions are convex on the specified domain.
    3. -
    -

    -\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \). -\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \). - -

      +

    1. Let \( f(x) = x^2 \) and \( g(x) = e^x \). Show that \( f(g(x)) \) and \( g(f(x)) \) is convex for \( x \in \mathbb{R} \). Also show that if \( f(x) \) is any convex function than \( h(x) = e^{f(x)} \) is convex.
    2. A norm is any function that satisfy the following properties
    3. + + +

    -\( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \). -\( f(x+y) \leq f(x) + f(y) \) -\( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \) -Using the definition of convexity, show that a function satisfying the properties above is convex (the third condition is not needed to show this). +Using the definition of convexity, try to show that a function satisfying the properties above is convex (the third condition is not needed to show this). @@ -459,7 +470,10 @@ Using the definition of convexity, show that a function satisfying the propertie

    Revisiting our first homework

    -We will use linear regression as a case study for the gradient descent methods. Linear regression is a great test case for the gradient descent methods discussed in the lectures since it has several desirable properties such as: +We will use linear regression as a case study for the gradient descent +methods. Linear regression is a great test case for the gradient +descent methods discussed in the lectures since it has several +desirable properties such as:

    1. An analytical solution (recall homework set 1).
    2. @@ -479,14 +493,14 @@ with \( x_i \in [0,1] \) chosen randomly with a uniform distribution. Additiona The linear regression model is given by

       
      $$ -h_\theta(x) = \hat{y} = \theta_0 + \theta_1 x, +h_\beta(x) = \hat{y} = \beta_0 + \beta_1 x, $$

       
      such that

       
      $$ -\hat{y}_i = \theta_0 + \theta_1 x_i. +\hat{y}_i = \beta_0 + \beta_1 x_i. $$

       
      @@ -496,10 +510,10 @@ $$

      Gradient descent example

      -Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \theta = (\theta_0, \theta_1)^T \) +Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \beta = (\beta_0, \beta_1)^T \)

      -t is convenient to write \( \mathbf{\hat{y}} = X\theta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by +t is convenient to write \( \mathbf{\hat{y}} = X\beta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by

       
      $$ \begin{equation} @@ -516,11 +530,11 @@ $$ The loss function is given by

       
      $$ -C(\theta) = ||X\theta-\mathbf{y}||^2 = ||X\theta||^2 - 2 \mathbf{y}^T X\theta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\theta_0 + \theta_1 x_i)^2 - 2 y_i (\theta_0 + \theta_1 x_i) + y_i^2 +C(\beta) = ||X\beta-\mathbf{y}||^2 = ||X\beta||^2 - 2 \mathbf{y}^T X\beta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\beta_0 + \beta_1 x_i)^2 - 2 y_i (\beta_0 + \beta_1 x_i) + y_i^2 $$

       
      -and we want to find \( \theta \) such that \( C(\theta) \) is minimized. +and we want to find \( \beta \) such that \( C(\beta) \) is minimized. @@ -528,12 +542,12 @@ and we want to find \( \theta \) such that \( C(\theta) \) is minimized.

      The derivative of the cost/loss function

      -Computing \( \partial C(\theta) / \partial \theta_0 \) and \( \partial C(\theta) / \partial \theta_1 \) we can show that the gradient can be written as +Computing \( \partial C(\beta) / \partial \beta_0 \) and \( \partial C(\beta) / \partial \beta_1 \) we can show that the gradient can be written as

       
      $$ -\nabla_\theta C(\theta) = (\partial C(\theta) / \partial \theta_0, \partial C(\theta) / \partial \theta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} = 2X^T(X\theta - \mathbf{y}), +\nabla_\beta C(\beta) = (\partial C(\beta) / \partial \beta_0, \partial C(\beta) / \partial \beta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} = 2X^T(X\beta - \mathbf{y}), $$

       
      @@ -543,17 +557,17 @@ where \( X \) is the design matrix defined above.

      The Hessian matrix

      -The Hessian matrix of \( C(\theta) \) is given by +The Hessian matrix of \( C(\beta) \) is given by

       
      $$ \hat{H} \equiv \begin{bmatrix} -\frac{\partial^2 C(\theta)}{\partial \theta_0^2} & \frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} \\ -\frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} & \frac{\partial^2 C(\theta)}{\partial \theta_1^2} & \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0^2} & \frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} & \frac{\partial^2 C(\beta)}{\partial \beta_1^2} & \\ \end{bmatrix} = 2X^T X. $$

       
      -This result implies that \( C(\theta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite. +This result implies that \( C(\beta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite.

      @@ -561,21 +575,21 @@ This result implies that \( C(\theta) \) is a convex function since the matrix \

      Simple program

      -We can now write a program that minimizes \( C(\theta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to +We can now write a program that minimizes \( C(\beta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to

       
      $$ -\theta_{k+1} = \theta_k - \gamma \nabla_\theta C(\theta_k), \ k=0,1,\cdots +\beta_{k+1} = \beta_k - \gamma \nabla_\beta C(\beta_k), \ k=0,1,\cdots $$

       

      We can use the expression we computed for the gradient and let use a -\( \theta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating -when \( ||\nabla_\theta C(\theta_k) || < \epsilon = 10^{-8} \). +\( \beta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating +when \( ||\nabla_\beta C(\beta_k) || < \epsilon = 10^{-8} \).

      -And finally we can compare our solution for \( \theta \) with the analytic result given by -\( \theta= (X^TX)^{-1} X^T \mathbf{y} \). +And finally we can compare our solution for \( \beta \) with the analytic result given by +\( \beta= (X^TX)^{-1} X^T \mathbf{y} \).

      @@ -592,11 +606,11 @@ x = np.random.rand(N) #Uniformly generated x-value y = 5*x**2 + 0.1*np.random.randn(N) X = np.c_[np.ones(N),x] #Construct design matrix -#Compute theta according to normal equations to compare with GD solution +#Compute beta according to normal equations to compare with GD solution Xt_X_inv = np.linalg.inv(np.dot(X.T,X)) Xt_y = np.dot(X.transpose(),y) -theta_NE = np.dot(Xt_X_inv,Xt_y) -print(theta_NE) +beta_NE = np.dot(Xt_X_inv,Xt_y) +print(beta_NE) @@ -605,33 +619,33 @@ theta_NE = np.dot(Xt_X_inv,Xt_y)

      Gradient descent and Ridge

      -We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \theta \), +We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \beta \),

       
      $$ -C_{\text{ridge}}(\theta) = ||X\theta -\mathbf{y}||^2 + \lambda ||\theta||^2, \ \lambda \geq 0. +C_{\text{ridge}}(\beta) = ||X\beta -\mathbf{y}||^2 + \lambda ||\beta||^2, \ \lambda \geq 0. $$

       

      -In order to minimize \( C_{\text{ridge}}(\theta) \) using GD we only have adjust the gradient as follows +In order to minimize \( C_{\text{ridge}}(\beta) \) using GD we only have adjust the gradient as follows

       
      $$ -\nabla_\theta C_{\text{ridge}}(\theta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} + 2\lambda\begin{bmatrix} \theta_0 \\ \theta_1\end{bmatrix} = 2 (X^T(X\theta - \mathbf{y})+\lambda \theta). +\nabla_\beta C_{\text{ridge}}(\beta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} + 2\lambda\begin{bmatrix} \beta_0 \\ \beta_1\end{bmatrix} = 2 (X^T(X\beta - \mathbf{y})+\lambda \beta). $$

       

      -We can now extend our program to minimize \( C_{\text{ridge}}(\theta) \) using gradient descent and compare with the analytical solution given by +We can now extend our program to minimize \( C_{\text{ridge}}(\beta) \) using gradient descent and compare with the analytical solution given by

       
      $$ -\theta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, +\beta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, $$

       
      for \( \lambda = {0,1,10,50,100} \) (\( \lambda = 0 \) corresponds to ordinary least squares). -We can then compute \( ||\theta_{\text{ridge}}|| \) for each \( \lambda \). +We can then compute \( ||\beta_{\text{ridge}}|| \) for each \( \lambda \).

      @@ -649,7 +663,7 @@ x = np.random.rand(N) y = 5*x**2 + 0.1*np.random.randn(N) -#Compute analytic theta for Ridge regression +#Compute analytic beta for Ridge regression X = np.c_[np.ones(N),x] XT_X = np.dot(X.T,X) @@ -657,10 +671,10 @@ l = 0.1 #Ridge Id = np.eye(XT_X.shape[0]) Z = np.linalg.inv(XT_X+l*Id) -theta_ridge = np.dot(Z,np.dot(X.T,y)) +beta_ridge = np.dot(Z,np.dot(X.T,y)) -print(theta_ridge) -print(np.linalg.norm(theta_ridge)) #||theta|| +print(beta_ridge) +print(np.linalg.norm(beta_ridge)) #||beta|| @@ -678,8 +692,8 @@ function, which we want to minimize, can almost always be written as a sum over \( n \) datapoints \( \{\mathbf{x}_i\}_{i=1}^n \),

       
      $$ -C(\mathbf{\theta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, -\mathbf{\theta}). +C(\mathbf{\beta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

       
      @@ -693,8 +707,8 @@ This in turn means that the gradient can be computed as a sum over \( i \)-gradients

       
      $$ -\nabla_\theta C(\mathbf{\theta}) = \sum_i^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}). +\nabla_\beta C(\mathbf{\beta}) = \sum_i^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

       
      @@ -724,10 +738,10 @@ all datapoints with a sum over the datapoints in one the minibatches picked at random in each gradient descent step

       
      $$ -\nabla_\theta -C(\mathbf{\theta}) = \sum_{i=1}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) \rightarrow \sum_{i \in B_k}^n \nabla_\theta -c_i(\mathbf{x}_i, \mathbf{\theta}). +\nabla_\beta +C(\mathbf{\beta}) = \sum_{i=1}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) \rightarrow \sum_{i \in B_k}^n \nabla_\beta +c_i(\mathbf{x}_i, \mathbf{\beta}). $$

       
      @@ -740,8 +754,8 @@ $$ Thus a gradient descent step now looks like

       
      $$ -\theta_{j+1} = \theta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) +\beta_{j+1} = \beta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) $$

       
      @@ -772,7 +786,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for j += 1

      @@ -798,7 +812,7 @@ is zero is valid also for local minima, so this would only tell us that we are close to a local/global minimum. However, we could also evaluate the cost function at this point, store the result and continue the search. If the test kicks in at a later stage we can -compare the values of the cost function and keep the \( \theta \) that +compare the values of the cost function and keep the \( \beta \) that gave the lowest value. @@ -817,10 +831,10 @@ $$\gamma_j(t; t_0, t_1) = \frac{t_0}{t+t_1} $$

       
      goes to zero as the number of epochs gets large. I.e. we start with a step length \( \gamma_j (0; t_0, t_1) = t_0/t_1 \) which decays in time \( t \).

      -In this way we can fix the number of epochs, compute \( \theta \) and +In this way we can fix the number of epochs, compute \( \beta \) and evaluate the cost function at the end. Repeating the computation will give a different result since the scheme is random by design. Then we -pick the final \( \theta \) that gives the lowest value of the cost +pick the final \( \beta \) that gives the lowest value of the cost function.

      @@ -844,7 +858,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for beta t = epoch*m+i gamma_j = step_length(t,t0,t1) j += 1 diff --git a/doc/pub/Splines/html/Splines-solarized.html b/doc/pub/Splines/html/Splines-solarized.html index c14baf537..fc0bfdc82 100644 --- a/doc/pub/Splines/html/Splines-solarized.html +++ b/doc/pub/Splines/html/Splines-solarized.html @@ -156,11 +156,11 @@ MathJax.Hub.Config({

      Almost every problem in machine learning and data science starts with -a dataset \( X \), a model \( g(\theta) \), which is a function of the -parameters \( \theta \) and a cost function \( C(X, g(\theta)) \) that allows -us to judge how well the model \( g(\theta) \) explains the observations -\( X \). The model is fit by finding the values of \( \theta \) that minimize -the cost function. Ideally we would be able to solve for \( \theta \) +a dataset \( X \), a model \( g(\beta) \), which is a function of the +parameters \( \beta \) and a cost function \( C(X, g(\beta)) \) that allows +us to judge how well the model \( g(\beta) \) explains the observations +\( X \). The model is fit by finding the values of \( \beta \) that minimize +the cost function. Ideally we would be able to solve for \( \beta \) analytically, however this is not possible in general and we must use some approximative/numerical method to compute the minimum. @@ -195,7 +195,7 @@ we are always moving towards smaller function values, i.e a minimum. The previous observation is the basis of the method of steepest descent, which is also referred to as just gradient descent (GD). One starts with an initial guess \( \mathbf{x}_0 \) for a minimum of \( F \) and -compute new approximations according to +computes new approximations according to $$ \mathbf{x}_{k+1} = \mathbf{x}_k - \gamma_k \nabla F(\mathbf{x}_k), \ \ k \geq 0. @@ -203,7 +203,7 @@ $$

      The parameter \( \gamma_k \) is often referred to as the step length or -the learning rate in the context of Machine Learning. +the learning rate within the context of Machine Learning.

      @@ -211,10 +211,21 @@ the learning rate in the context of Machine Learning.

      The ideal

      -Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global minimum of the function \( F \). In general we do not know if we are in a global or local minimum. In the special case when \( F \) is a convex function, all local minima are also global minima, so in this case gradient descent can converge to the global solution. The advantage of this scheme is that it is conceptually simple and straightforward to implement. However the method in this form has some severe limitations: +Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global +minimum of the function \( F \). In general we do not know if we are in a +global or local minimum. In the special case when \( F \) is a convex +function, all local minima are also global minima, so in this case +gradient descent can converge to the global solution. The advantage of +this scheme is that it is conceptually simple and straightforward to +implement. However the method in this form has some severe +limitations:

      -In machine learing we are often faced with non-convex high dimensional cost functions with many local minimum. Since GD is deterministic we will get stuck in a local minimum, if the method converges, unless we have a very good intial guess. This also implies that the scheme is sensitive to the chosen initial condition. +In machine learing we are often faced with non-convex high dimensional +cost functions with many local minima. Since GD is deterministic we +will get stuck in a local minimum, if the method converges, unless we +have a very good intial guess. This also implies that the scheme is +sensitive to the chosen initial condition.

      Note that the gradient is a function of \( \mathbf{x} = @@ -229,19 +240,21 @@ Note that the gradient is a function of \( \mathbf{x} = GD is sensitive to the choice of learning rate \( \gamma_k \). This is due to the fact that we are only guaranteed that \( F(\mathbf{x}_{k+1}) \leq F(\mathbf{x}_k) \) for sufficiently small \( \gamma_k \). The problem is to -determine an optimal learning rate. If the learning rate is chosen to -small the method will take a long to converge and if it is to large we -can experience erratic behavior. +determine an optimal learning rate. If the learning rate is chosen too +small the method will take a long time to converge and if it is too +large we can experience erratic behavior.

      Many of these shortcomings can be alleviated by introducing randomness. One such method is that of Stochastic Gradient Descent -(SGD), see below +(SGD), see below.











      Gradient Descent Example

      + +

      We revisit now our simple linear regression example with a linear polynomial.

      @@ -259,30 +272,30 @@ x = 2*np.random.rand(4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) -theta = np.random.randn(2,1) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) +beta = np.random.randn(2,1) eta = 0.1 Niterations = 1000 m = 100 for iter in range(Niterations): - gradients = 2.0/m*xb.T.dot(xb.dot(theta)-y) - theta -= eta*gradients + gradients = 2.0/m*xb.T.dot(xb.dot(beta)-y) + beta -= eta*gradients -print(theta) +print(beta) xnew = np.array([[0],[2]]) xbnew = np.c_[np.ones((2,1)), xnew] -ypredict = xbnew.dot(theta) -ypredict2 = xbnew.dot(theta_linreg) +ypredict = xbnew.dot(beta) +ypredict2 = xbnew.dot(beta_linreg) plt.plot(xnew, ypredict, "r-") plt.plot(xnew, ypredict2, "b-") plt.plot(x, y ,'ro') plt.axis([0,2.0,0, 15.0]) plt.xlabel(r'$x$') plt.ylabel(r'$y$') -plt.title(r'Random numbers ') +plt.title(r'Gradient descent example') plt.show()

      @@ -303,8 +316,8 @@ x = 2*np.random.rand(4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) sgdreg = SGDRegressor(n_iter = 50, penalty=None, eta0=0.1) sgdreg.fit(x,y.ravel()) print(sgdreg.intercept_, sgdreg.coef_) @@ -313,6 +326,8 @@ sgdreg.fit(x,y.ravel())

      Convex functions

      + +

      Ideally we want our cost/loss function to be convex(concave).

      @@ -333,7 +348,7 @@ regular polygons (triangles, rectangles, pentagons, etc...).

      Convex function

      -Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below. +Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below.











      @@ -344,9 +359,7 @@ Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that t In the following we state first and second-order conditions which ensures convexity of a function \( f \). We write \( D_f \) to denote the domain of \( f \), i.e the subset of \( R^n \) where \( f \) is defined. For more -details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex -Optimization. Cambridge University Press, http://stanford.edu/ -boyd/cvxbook/, 2004. +details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press.

      @@ -389,9 +402,12 @@ This condition is particularly useful since it gives us an procedure for determi The next result is of great importance to us and the reason why we are going on about convex functions. In machine learning we frequently have to minimize a loss/cost function in order to find the best -parameters for the model we are considering. Ideally we want the -global minimum, however for high-dimensional models it is hard to know -if we have local or global minimum. However, if the cost/loss function +parameters for the model we are considering. + +

      +Ideally we want the +global minimum (for high-dimensional models it is hard to know +if we have local or global minimum). However, if the cost/loss function is convex the following result provides invaluable information:

      @@ -403,6 +419,8 @@ is minimal, where \( f \) is convex and differentiable. Then, any point \( x^* \) that satisfies \( \nabla f(x^*) = 0 \) is a global minimum.

      + +

      This result means that if we know that the cost/loss function is convex and we are able to find a minimum, we are guaranteed that it is a global minimum.

      @@ -411,27 +429,26 @@ This result means that if we know that the cost/loss function is convex and we a

      Some simple problems

        -
      1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.
      2. -
      - -Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $$\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $$ - -
        +
      1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity. Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $
      2. Using the second order condition show that the following functions are convex on the specified domain.
      3. -
      -\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \). -\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \). + -
      1. Let \( f(x) = x^2 \) and \( g(x) = e^x \). Show that \( f(g(x)) \) and \( g(f(x)) \) is convex for \( x \in \mathbb{R} \). Also show that if \( f(x) \) is any convex function than \( h(x) = e^{f(x)} \) is convex.
      2. A norm is any function that satisfy the following properties
      3. + +
          +
        • \( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \).
        • +
        • \( f(x+y) \leq f(x) + f(y) \)
        • +
        • \( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \)
        • +
        +
      -\( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \). -\( f(x+y) \leq f(x) + f(y) \) -\( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \) -Using the definition of convexity, show that a function satisfying the properties above is convex (the third condition is not needed to show this). +Using the definition of convexity, try to show that a function satisfying the properties above is convex (the third condition is not needed to show this).

      @@ -439,7 +456,10 @@ Using the definition of convexity, show that a function satisfying the propertie

      Revisiting our first homework

      -We will use linear regression as a case study for the gradient descent methods. Linear regression is a great test case for the gradient descent methods discussed in the lectures since it has several desirable properties such as: +We will use linear regression as a case study for the gradient descent +methods. Linear regression is a great test case for the gradient +descent methods discussed in the lectures since it has several +desirable properties such as:

      1. An analytical solution (recall homework set 1).
      2. @@ -455,12 +475,12 @@ $$ with \( x_i \in [0,1] \) chosen randomly with a uniform distribution. Additionally \( \xi_i \) represents stochastic noise chosen according to a normal distribution \( \cal {N}(0,1) \). The linear regression model is given by $$ -h_\theta(x) = \hat{y} = \theta_0 + \theta_1 x, +h_\beta(x) = \hat{y} = \beta_0 + \beta_1 x, $$ such that $$ -\hat{y}_i = \theta_0 + \theta_1 x_i. +\hat{y}_i = \beta_0 + \beta_1 x_i. $$

        @@ -469,10 +489,10 @@ $$

        Gradient descent example

        -Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \theta = (\theta_0, \theta_1)^T \) +Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \beta = (\beta_0, \beta_1)^T \)

        -t is convenient to write \( \mathbf{\hat{y}} = X\theta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by +t is convenient to write \( \mathbf{\hat{y}} = X\beta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by $$ \begin{equation} X \equiv \begin{bmatrix} @@ -486,10 +506,10 @@ $$ The loss function is given by $$ -C(\theta) = ||X\theta-\mathbf{y}||^2 = ||X\theta||^2 - 2 \mathbf{y}^T X\theta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\theta_0 + \theta_1 x_i)^2 - 2 y_i (\theta_0 + \theta_1 x_i) + y_i^2 +C(\beta) = ||X\beta-\mathbf{y}||^2 = ||X\beta||^2 - 2 \mathbf{y}^T X\beta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\beta_0 + \beta_1 x_i)^2 - 2 y_i (\beta_0 + \beta_1 x_i) + y_i^2 $$ -and we want to find \( \theta \) such that \( C(\theta) \) is minimized. +and we want to find \( \beta \) such that \( C(\beta) \) is minimized.











        @@ -497,11 +517,11 @@ and we want to find \( \theta \) such that \( C(\theta) \) is minimized.

        The derivative of the cost/loss function

        -Computing \( \partial C(\theta) / \partial \theta_0 \) and \( \partial C(\theta) / \partial \theta_1 \) we can show that the gradient can be written as +Computing \( \partial C(\beta) / \partial \beta_0 \) and \( \partial C(\beta) / \partial \beta_1 \) we can show that the gradient can be written as $$ -\nabla_\theta C(\theta) = (\partial C(\theta) / \partial \theta_0, \partial C(\theta) / \partial \theta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} = 2X^T(X\theta - \mathbf{y}), +\nabla_\beta C(\beta) = (\partial C(\beta) / \partial \beta_0, \partial C(\beta) / \partial \beta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} = 2X^T(X\beta - \mathbf{y}), $$ where \( X \) is the design matrix defined above. @@ -510,15 +530,15 @@ where \( X \) is the design matrix defined above.









        The Hessian matrix

        -The Hessian matrix of \( C(\theta) \) is given by +The Hessian matrix of \( C(\beta) \) is given by $$ \hat{H} \equiv \begin{bmatrix} -\frac{\partial^2 C(\theta)}{\partial \theta_0^2} & \frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} \\ -\frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} & \frac{\partial^2 C(\theta)}{\partial \theta_1^2} & \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0^2} & \frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} & \frac{\partial^2 C(\beta)}{\partial \beta_1^2} & \\ \end{bmatrix} = 2X^T X. $$ -This result implies that \( C(\theta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite. +This result implies that \( C(\beta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite.











        @@ -526,19 +546,19 @@ This result implies that \( C(\theta) \) is a convex function since the matrix \

        Simple program

        -We can now write a program that minimizes \( C(\theta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to +We can now write a program that minimizes \( C(\beta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to $$ -\theta_{k+1} = \theta_k - \gamma \nabla_\theta C(\theta_k), \ k=0,1,\cdots +\beta_{k+1} = \beta_k - \gamma \nabla_\beta C(\beta_k), \ k=0,1,\cdots $$

        We can use the expression we computed for the gradient and let use a -\( \theta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating -when \( ||\nabla_\theta C(\theta_k) || < \epsilon = 10^{-8} \). +\( \beta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating +when \( ||\nabla_\beta C(\beta_k) || < \epsilon = 10^{-8} \).

        -And finally we can compare our solution for \( \theta \) with the analytic result given by -\( \theta= (X^TX)^{-1} X^T \mathbf{y} \). +And finally we can compare our solution for \( \beta \) with the analytic result given by +\( \beta= (X^TX)^{-1} X^T \mathbf{y} \).

        @@ -555,11 +575,11 @@ x = np.random.rand(N) #Uniformly generated x-value y = 5*x**2 + 0.1*np.random.randn(N) X = np.c_[np.ones(N),x] #Construct design matrix -#Compute theta according to normal equations to compare with GD solution +#Compute beta according to normal equations to compare with GD solution Xt_X_inv = np.linalg.inv(np.dot(X.T,X)) Xt_y = np.dot(X.transpose(),y) -theta_NE = np.dot(Xt_X_inv,Xt_y) -print(theta_NE) +beta_NE = np.dot(Xt_X_inv,Xt_y) +print(beta_NE)

        @@ -567,27 +587,27 @@ theta_NE = np.dot(Xt_X_inv,Xt_y)

        Gradient descent and Ridge

        -We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \theta \), +We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \beta \), $$ -C_{\text{ridge}}(\theta) = ||X\theta -\mathbf{y}||^2 + \lambda ||\theta||^2, \ \lambda \geq 0. +C_{\text{ridge}}(\beta) = ||X\beta -\mathbf{y}||^2 + \lambda ||\beta||^2, \ \lambda \geq 0. $$

        -In order to minimize \( C_{\text{ridge}}(\theta) \) using GD we only have adjust the gradient as follows +In order to minimize \( C_{\text{ridge}}(\beta) \) using GD we only have adjust the gradient as follows $$ -\nabla_\theta C_{\text{ridge}}(\theta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} + 2\lambda\begin{bmatrix} \theta_0 \\ \theta_1\end{bmatrix} = 2 (X^T(X\theta - \mathbf{y})+\lambda \theta). +\nabla_\beta C_{\text{ridge}}(\beta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} + 2\lambda\begin{bmatrix} \beta_0 \\ \beta_1\end{bmatrix} = 2 (X^T(X\beta - \mathbf{y})+\lambda \beta). $$

        -We can now extend our program to minimize \( C_{\text{ridge}}(\theta) \) using gradient descent and compare with the analytical solution given by +We can now extend our program to minimize \( C_{\text{ridge}}(\beta) \) using gradient descent and compare with the analytical solution given by $$ -\theta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, +\beta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, $$ for \( \lambda = {0,1,10,50,100} \) (\( \lambda = 0 \) corresponds to ordinary least squares). -We can then compute \( ||\theta_{\text{ridge}}|| \) for each \( \lambda \). +We can then compute \( ||\beta_{\text{ridge}}|| \) for each \( \lambda \).

        @@ -605,7 +625,7 @@ x = np.random.rand(N) y = 5*x**2 + 0.1*np.random.randn(N) -#Compute analytic theta for Ridge regression +#Compute analytic beta for Ridge regression X = np.c_[np.ones(N),x] XT_X = np.dot(X.T,X) @@ -613,10 +633,10 @@ l = 0.1 #Ridge Id = np.eye(XT_X.shape[0]) Z = np.linalg.inv(XT_X+l*Id) -theta_ridge = np.dot(Z,np.dot(X.T,y)) +beta_ridge = np.dot(Z,np.dot(X.T,y)) -print(theta_ridge) -print(np.linalg.norm(theta_ridge)) #||theta|| +print(beta_ridge) +print(np.linalg.norm(beta_ridge)) #||beta||











        @@ -632,8 +652,8 @@ The underlying idea of SGD comes from the observation that the cost function, which we want to minimize, can almost always be written as a sum over \( n \) datapoints \( \{\mathbf{x}_i\}_{i=1}^n \), $$ -C(\mathbf{\theta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, -\mathbf{\theta}). +C(\mathbf{\beta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

        @@ -645,8 +665,8 @@ $$ This in turn means that the gradient can be computed as a sum over \( i \)-gradients $$ -\nabla_\theta C(\mathbf{\theta}) = \sum_i^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}). +\nabla_\beta C(\mathbf{\beta}) = \sum_i^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

        @@ -674,10 +694,10 @@ The idea is now to approximate the gradient by replacing the sum over all datapoints with a sum over the datapoints in one the minibatches picked at random in each gradient descent step $$ -\nabla_\theta -C(\mathbf{\theta}) = \sum_{i=1}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) \rightarrow \sum_{i \in B_k}^n \nabla_\theta -c_i(\mathbf{x}_i, \mathbf{\theta}). +\nabla_\beta +C(\mathbf{\beta}) = \sum_{i=1}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) \rightarrow \sum_{i \in B_k}^n \nabla_\beta +c_i(\mathbf{x}_i, \mathbf{\beta}). $$

        @@ -688,8 +708,8 @@ $$

        Thus a gradient descent step now looks like $$ -\theta_{j+1} = \theta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) +\beta_{j+1} = \beta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) $$

        @@ -719,7 +739,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for j += 1

        @@ -745,7 +765,7 @@ is zero is valid also for local minima, so this would only tell us that we are close to a local/global minimum. However, we could also evaluate the cost function at this point, store the result and continue the search. If the test kicks in at a later stage we can -compare the values of the cost function and keep the \( \theta \) that +compare the values of the cost function and keep the \( \beta \) that gave the lowest value.

        @@ -762,10 +782,10 @@ reasonable time such that we do not move at all. As an example, let \( e = 0,1,2,3,\cdots \) denote the current epoch and let \( t_0, t_1 > 0 \) be two fixed numbers. Furthermore, let \( t = e \cdot m + i \) where \( m \) is the number of minibatches and \( i=0,\cdots,m-1 \). Then the function $$\gamma_j(t; t_0, t_1) = \frac{t_0}{t+t_1} $$ goes to zero as the number of epochs gets large. I.e. we start with a step length \( \gamma_j (0; t_0, t_1) = t_0/t_1 \) which decays in time \( t \).

        -In this way we can fix the number of epochs, compute \( \theta \) and +In this way we can fix the number of epochs, compute \( \beta \) and evaluate the cost function at the end. Repeating the computation will give a different result since the scheme is random by design. Then we -pick the final \( \theta \) that gives the lowest value of the cost +pick the final \( \beta \) that gives the lowest value of the cost function.

        @@ -789,7 +809,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for beta t = epoch*m+i gamma_j = step_length(t,t0,t1) j += 1 diff --git a/doc/pub/Splines/html/Splines.html b/doc/pub/Splines/html/Splines.html index f616f2f8c..b7f5a50f9 100644 --- a/doc/pub/Splines/html/Splines.html +++ b/doc/pub/Splines/html/Splines.html @@ -161,11 +161,11 @@ MathJax.Hub.Config({

        Almost every problem in machine learning and data science starts with -a dataset \( X \), a model \( g(\theta) \), which is a function of the -parameters \( \theta \) and a cost function \( C(X, g(\theta)) \) that allows -us to judge how well the model \( g(\theta) \) explains the observations -\( X \). The model is fit by finding the values of \( \theta \) that minimize -the cost function. Ideally we would be able to solve for \( \theta \) +a dataset \( X \), a model \( g(\beta) \), which is a function of the +parameters \( \beta \) and a cost function \( C(X, g(\beta)) \) that allows +us to judge how well the model \( g(\beta) \) explains the observations +\( X \). The model is fit by finding the values of \( \beta \) that minimize +the cost function. Ideally we would be able to solve for \( \beta \) analytically, however this is not possible in general and we must use some approximative/numerical method to compute the minimum. @@ -200,7 +200,7 @@ we are always moving towards smaller function values, i.e a minimum. The previous observation is the basis of the method of steepest descent, which is also referred to as just gradient descent (GD). One starts with an initial guess \( \mathbf{x}_0 \) for a minimum of \( F \) and -compute new approximations according to +computes new approximations according to $$ \mathbf{x}_{k+1} = \mathbf{x}_k - \gamma_k \nabla F(\mathbf{x}_k), \ \ k \geq 0. @@ -208,7 +208,7 @@ $$

        The parameter \( \gamma_k \) is often referred to as the step length or -the learning rate in the context of Machine Learning. +the learning rate within the context of Machine Learning.

        @@ -216,10 +216,21 @@ the learning rate in the context of Machine Learning.

        The ideal

        -Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global minimum of the function \( F \). In general we do not know if we are in a global or local minimum. In the special case when \( F \) is a convex function, all local minima are also global minima, so in this case gradient descent can converge to the global solution. The advantage of this scheme is that it is conceptually simple and straightforward to implement. However the method in this form has some severe limitations: +Ideally the sequence \( \{ \mathbf{x}_k \}_{k=0} \) converges to a global +minimum of the function \( F \). In general we do not know if we are in a +global or local minimum. In the special case when \( F \) is a convex +function, all local minima are also global minima, so in this case +gradient descent can converge to the global solution. The advantage of +this scheme is that it is conceptually simple and straightforward to +implement. However the method in this form has some severe +limitations:

        -In machine learing we are often faced with non-convex high dimensional cost functions with many local minimum. Since GD is deterministic we will get stuck in a local minimum, if the method converges, unless we have a very good intial guess. This also implies that the scheme is sensitive to the chosen initial condition. +In machine learing we are often faced with non-convex high dimensional +cost functions with many local minima. Since GD is deterministic we +will get stuck in a local minimum, if the method converges, unless we +have a very good intial guess. This also implies that the scheme is +sensitive to the chosen initial condition.

        Note that the gradient is a function of \( \mathbf{x} = @@ -234,19 +245,21 @@ Note that the gradient is a function of \( \mathbf{x} = GD is sensitive to the choice of learning rate \( \gamma_k \). This is due to the fact that we are only guaranteed that \( F(\mathbf{x}_{k+1}) \leq F(\mathbf{x}_k) \) for sufficiently small \( \gamma_k \). The problem is to -determine an optimal learning rate. If the learning rate is chosen to -small the method will take a long to converge and if it is to large we -can experience erratic behavior. +determine an optimal learning rate. If the learning rate is chosen too +small the method will take a long time to converge and if it is too +large we can experience erratic behavior.

        Many of these shortcomings can be alleviated by introducing randomness. One such method is that of Stochastic Gradient Descent -(SGD), see below +(SGD), see below.











        Gradient Descent Example

        + +

        We revisit now our simple linear regression example with a linear polynomial.

        @@ -264,30 +277,30 @@ x = 2*np y = 4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) -theta = np.random.randn(2,1) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) +beta = np.random.randn(2,1) eta = 0.1 Niterations = 1000 m = 100 for iter in range(Niterations): - gradients = 2.0/m*xb.T.dot(xb.dot(theta)-y) - theta -= eta*gradients + gradients = 2.0/m*xb.T.dot(xb.dot(beta)-y) + beta -= eta*gradients -print(theta) +print(beta) xnew = np.array([[0],[2]]) xbnew = np.c_[np.ones((2,1)), xnew] -ypredict = xbnew.dot(theta) -ypredict2 = xbnew.dot(theta_linreg) +ypredict = xbnew.dot(beta) +ypredict2 = xbnew.dot(beta_linreg) plt.plot(xnew, ypredict, "r-") plt.plot(xnew, ypredict2, "b-") plt.plot(x, y ,'ro') plt.axis([0,2.0,0, 15.0]) plt.xlabel(r'$x$') plt.ylabel(r'$y$') -plt.title(r'Random numbers ') +plt.title(r'Gradient descent example') plt.show()

        @@ -308,8 +321,8 @@ x = 2*np y = 4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) sgdreg = SGDRegressor(n_iter = 50, penalty=None, eta0=0.1) sgdreg.fit(x,y.ravel()) print(sgdreg.intercept_, sgdreg.coef_) @@ -318,6 +331,8 @@ sgdreg.fit(x,y.

        Convex functions

        + +

        Ideally we want our cost/loss function to be convex(concave).

        @@ -338,7 +353,7 @@ regular polygons (triangles, rectangles, pentagons, etc...).

        Convex function

        -Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below. +Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that the function \( f: X \rightarrow \mathbb{R} \) is continuous, then \( f \) is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all \( x_1, x_2 \in X \) and for all \( t \in [0,1] \). If \( \leq \) is replaced with a strict inequaltiy in the definition, we demand \( x_1 \neq x_2 \) and \( t\in(0,1) \) then \( f \) is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting \( f(x_1) \) and \( f(x_2) \), the value of the function on the interval \( [x_1,x_2] \) is always below the line as illustrated below.











        @@ -349,9 +364,7 @@ Convex function: Let \( X \subset \mathbb{R}^n \) be a convex set. Assume that t In the following we state first and second-order conditions which ensures convexity of a function \( f \). We write \( D_f \) to denote the domain of \( f \), i.e the subset of \( R^n \) where \( f \) is defined. For more -details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex -Optimization. Cambridge University Press, http://stanford.edu/ -boyd/cvxbook/, 2004. +details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press.

        @@ -394,9 +407,12 @@ This condition is particularly useful since it gives us an procedure for determi The next result is of great importance to us and the reason why we are going on about convex functions. In machine learning we frequently have to minimize a loss/cost function in order to find the best -parameters for the model we are considering. Ideally we want the -global minimum, however for high-dimensional models it is hard to know -if we have local or global minimum. However, if the cost/loss function +parameters for the model we are considering. + +

        +Ideally we want the +global minimum (for high-dimensional models it is hard to know +if we have local or global minimum). However, if the cost/loss function is convex the following result provides invaluable information:

        @@ -408,6 +424,8 @@ is minimal, where \( f \) is convex and differentiable. Then, any point \( x^* \) that satisfies \( \nabla f(x^*) = 0 \) is a global minimum.

        + +

        This result means that if we know that the cost/loss function is convex and we are able to find a minimum, we are guaranteed that it is a global minimum.

        @@ -416,27 +434,26 @@ This result means that if we know that the cost/loss function is convex and we a

        Some simple problems

          -
        1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.
        2. -
        - -Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $$\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $$ - -
          +
        1. Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity. Hint: If you re-write the definition, \( f \) is convex if the following holds for all \( x,y \in D_f \) and any \( \lambda \in [0,1] \) $\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $
        2. Using the second order condition show that the following functions are convex on the specified domain.
        3. -
        -\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \). -\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \). +
          +
        • \( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).
        • +
        • \( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).
        • +
        -
        1. Let \( f(x) = x^2 \) and \( g(x) = e^x \). Show that \( f(g(x)) \) and \( g(f(x)) \) is convex for \( x \in \mathbb{R} \). Also show that if \( f(x) \) is any convex function than \( h(x) = e^{f(x)} \) is convex.
        2. A norm is any function that satisfy the following properties
        3. + +
            +
          • \( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \).
          • +
          • \( f(x+y) \leq f(x) + f(y) \)
          • +
          • \( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \)
          • +
          +
        -\( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \). -\( f(x+y) \leq f(x) + f(y) \) -\( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \) -Using the definition of convexity, show that a function satisfying the properties above is convex (the third condition is not needed to show this). +Using the definition of convexity, try to show that a function satisfying the properties above is convex (the third condition is not needed to show this).

        @@ -444,7 +461,10 @@ Using the definition of convexity, show that a function satisfying the propertie

        Revisiting our first homework

        -We will use linear regression as a case study for the gradient descent methods. Linear regression is a great test case for the gradient descent methods discussed in the lectures since it has several desirable properties such as: +We will use linear regression as a case study for the gradient descent +methods. Linear regression is a great test case for the gradient +descent methods discussed in the lectures since it has several +desirable properties such as:

        1. An analytical solution (recall homework set 1).
        2. @@ -460,12 +480,12 @@ $$ with \( x_i \in [0,1] \) chosen randomly with a uniform distribution. Additionally \( \xi_i \) represents stochastic noise chosen according to a normal distribution \( \cal {N}(0,1) \). The linear regression model is given by $$ -h_\theta(x) = \hat{y} = \theta_0 + \theta_1 x, +h_\beta(x) = \hat{y} = \beta_0 + \beta_1 x, $$ such that $$ -\hat{y}_i = \theta_0 + \theta_1 x_i. +\hat{y}_i = \beta_0 + \beta_1 x_i. $$

          @@ -474,10 +494,10 @@ $$

          Gradient descent example

          -Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \theta = (\theta_0, \theta_1)^T \) +Let \( \mathbf{y} = (y_1,\cdots,y_n)^T \), \( \mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T \) and \( \beta = (\beta_0, \beta_1)^T \)

          -t is convenient to write \( \mathbf{\hat{y}} = X\theta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by +t is convenient to write \( \mathbf{\hat{y}} = X\beta \) where \( X \in \mathbb{R}^{100 \times 2} \) is the design matrix given by $$ \begin{equation} X \equiv \begin{bmatrix} @@ -491,10 +511,10 @@ $$ The loss function is given by $$ -C(\theta) = ||X\theta-\mathbf{y}||^2 = ||X\theta||^2 - 2 \mathbf{y}^T X\theta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\theta_0 + \theta_1 x_i)^2 - 2 y_i (\theta_0 + \theta_1 x_i) + y_i^2 +C(\beta) = ||X\beta-\mathbf{y}||^2 = ||X\beta||^2 - 2 \mathbf{y}^T X\beta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\beta_0 + \beta_1 x_i)^2 - 2 y_i (\beta_0 + \beta_1 x_i) + y_i^2 $$ -and we want to find \( \theta \) such that \( C(\theta) \) is minimized. +and we want to find \( \beta \) such that \( C(\beta) \) is minimized.











          @@ -502,11 +522,11 @@ and we want to find \( \theta \) such that \( C(\theta) \) is minimized.

          The derivative of the cost/loss function

          -Computing \( \partial C(\theta) / \partial \theta_0 \) and \( \partial C(\theta) / \partial \theta_1 \) we can show that the gradient can be written as +Computing \( \partial C(\beta) / \partial \beta_0 \) and \( \partial C(\beta) / \partial \beta_1 \) we can show that the gradient can be written as $$ -\nabla_\theta C(\theta) = (\partial C(\theta) / \partial \theta_0, \partial C(\theta) / \partial \theta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} = 2X^T(X\theta - \mathbf{y}), +\nabla_\beta C(\beta) = (\partial C(\beta) / \partial \beta_0, \partial C(\beta) / \partial \beta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} = 2X^T(X\beta - \mathbf{y}), $$ where \( X \) is the design matrix defined above. @@ -515,15 +535,15 @@ where \( X \) is the design matrix defined above.









          The Hessian matrix

          -The Hessian matrix of \( C(\theta) \) is given by +The Hessian matrix of \( C(\beta) \) is given by $$ \hat{H} \equiv \begin{bmatrix} -\frac{\partial^2 C(\theta)}{\partial \theta_0^2} & \frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} \\ -\frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} & \frac{\partial^2 C(\theta)}{\partial \theta_1^2} & \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0^2} & \frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} & \frac{\partial^2 C(\beta)}{\partial \beta_1^2} & \\ \end{bmatrix} = 2X^T X. $$ -This result implies that \( C(\theta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite. +This result implies that \( C(\beta) \) is a convex function since the matrix \( X^T X \) always is positive semi-definite.











          @@ -531,19 +551,19 @@ This result implies that \( C(\theta) \) is a convex function since the matrix \

          Simple program

          -We can now write a program that minimizes \( C(\theta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to +We can now write a program that minimizes \( C(\beta) \) using the gradient descent method with a constant learning rate \( \gamma \) according to $$ -\theta_{k+1} = \theta_k - \gamma \nabla_\theta C(\theta_k), \ k=0,1,\cdots +\beta_{k+1} = \beta_k - \gamma \nabla_\beta C(\beta_k), \ k=0,1,\cdots $$

          We can use the expression we computed for the gradient and let use a -\( \theta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating -when \( ||\nabla_\theta C(\theta_k) || < \epsilon = 10^{-8} \). +\( \beta_0 \) be chosen randomly and let \( \gamma = 0.001 \). Stop iterating +when \( ||\nabla_\beta C(\beta_k) || < \epsilon = 10^{-8} \).

          -And finally we can compare our solution for \( \theta \) with the analytic result given by -\( \theta= (X^TX)^{-1} X^T \mathbf{y} \). +And finally we can compare our solution for \( \beta \) with the analytic result given by +\( \beta= (X^TX)^{-1} X^T \mathbf{y} \).

          @@ -560,11 +580,11 @@ x = np. y = 5*x**2 + 0.1*np.random.randn(N) X = np.c_[np.ones(N),x] #Construct design matrix -#Compute theta according to normal equations to compare with GD solution +#Compute beta according to normal equations to compare with GD solution Xt_X_inv = np.linalg.inv(np.dot(X.T,X)) Xt_y = np.dot(X.transpose(),y) -theta_NE = np.dot(Xt_X_inv,Xt_y) -print(theta_NE) +beta_NE = np.dot(Xt_X_inv,Xt_y) +print(beta_NE)

          @@ -572,27 +592,27 @@ theta_NE = np.<

          Gradient descent and Ridge

          -We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \theta \), +We have also discussed Ridge regression where the loss function contains a regularized given by the \( L_2 \) norm of \( \beta \), $$ -C_{\text{ridge}}(\theta) = ||X\theta -\mathbf{y}||^2 + \lambda ||\theta||^2, \ \lambda \geq 0. +C_{\text{ridge}}(\beta) = ||X\beta -\mathbf{y}||^2 + \lambda ||\beta||^2, \ \lambda \geq 0. $$

          -In order to minimize \( C_{\text{ridge}}(\theta) \) using GD we only have adjust the gradient as follows +In order to minimize \( C_{\text{ridge}}(\beta) \) using GD we only have adjust the gradient as follows $$ -\nabla_\theta C_{\text{ridge}}(\theta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} + 2\lambda\begin{bmatrix} \theta_0 \\ \theta_1\end{bmatrix} = 2 (X^T(X\theta - \mathbf{y})+\lambda \theta). +\nabla_\beta C_{\text{ridge}}(\beta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} + 2\lambda\begin{bmatrix} \beta_0 \\ \beta_1\end{bmatrix} = 2 (X^T(X\beta - \mathbf{y})+\lambda \beta). $$

          -We can now extend our program to minimize \( C_{\text{ridge}}(\theta) \) using gradient descent and compare with the analytical solution given by +We can now extend our program to minimize \( C_{\text{ridge}}(\beta) \) using gradient descent and compare with the analytical solution given by $$ -\theta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, +\beta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, $$ for \( \lambda = {0,1,10,50,100} \) (\( \lambda = 0 \) corresponds to ordinary least squares). -We can then compute \( ||\theta_{\text{ridge}}|| \) for each \( \lambda \). +We can then compute \( ||\beta_{\text{ridge}}|| \) for each \( \lambda \).

          @@ -610,7 +630,7 @@ x = np. y = 5*x**2 + 0.1*np.random.randn(N) -#Compute analytic theta for Ridge regression +#Compute analytic beta for Ridge regression X = np.c_[np.ones(N),x] XT_X = np.dot(X.T,X) @@ -618,10 +638,10 @@ l = 0.1 Id = np.eye(XT_X.shape[0]) Z = np.linalg.inv(XT_X+l*Id) -theta_ridge = np.dot(Z,np.dot(X.T,y)) +beta_ridge = np.dot(Z,np.dot(X.T,y)) -print(theta_ridge) -print(np.linalg.norm(theta_ridge)) #||theta|| +print(beta_ridge) +print(np.linalg.norm(beta_ridge)) #||beta||











          @@ -637,8 +657,8 @@ The underlying idea of SGD comes from the observation that the cost function, which we want to minimize, can almost always be written as a sum over \( n \) datapoints \( \{\mathbf{x}_i\}_{i=1}^n \), $$ -C(\mathbf{\theta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, -\mathbf{\theta}). +C(\mathbf{\beta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

          @@ -650,8 +670,8 @@ $$ This in turn means that the gradient can be computed as a sum over \( i \)-gradients $$ -\nabla_\theta C(\mathbf{\theta}) = \sum_i^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}). +\nabla_\beta C(\mathbf{\beta}) = \sum_i^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}). $$

          @@ -679,10 +699,10 @@ The idea is now to approximate the gradient by replacing the sum over all datapoints with a sum over the datapoints in one the minibatches picked at random in each gradient descent step $$ -\nabla_\theta -C(\mathbf{\theta}) = \sum_{i=1}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) \rightarrow \sum_{i \in B_k}^n \nabla_\theta -c_i(\mathbf{x}_i, \mathbf{\theta}). +\nabla_\beta +C(\mathbf{\beta}) = \sum_{i=1}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) \rightarrow \sum_{i \in B_k}^n \nabla_\beta +c_i(\mathbf{x}_i, \mathbf{\beta}). $$

          @@ -693,8 +713,8 @@ $$

          Thus a gradient descent step now looks like $$ -\theta_{j+1} = \theta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) +\beta_{j+1} = \beta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) $$

          @@ -724,7 +744,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for j += 1

          @@ -750,7 +770,7 @@ is zero is valid also for local minima, so this would only tell us that we are close to a local/global minimum. However, we could also evaluate the cost function at this point, store the result and continue the search. If the test kicks in at a later stage we can -compare the values of the cost function and keep the \( \theta \) that +compare the values of the cost function and keep the \( \beta \) that gave the lowest value.

          @@ -767,10 +787,10 @@ reasonable time such that we do not move at all. As an example, let \( e = 0,1,2,3,\cdots \) denote the current epoch and let \( t_0, t_1 > 0 \) be two fixed numbers. Furthermore, let \( t = e \cdot m + i \) where \( m \) is the number of minibatches and \( i=0,\cdots,m-1 \). Then the function $$\gamma_j(t; t_0, t_1) = \frac{t_0}{t+t_1} $$ goes to zero as the number of epochs gets large. I.e. we start with a step length \( \gamma_j (0; t_0, t_1) = t_0/t_1 \) which decays in time \( t \).

          -In this way we can fix the number of epochs, compute \( \theta \) and +In this way we can fix the number of epochs, compute \( \beta \) and evaluate the cost function at the end. Repeating the computation will give a different result since the scheme is random by design. Then we -pick the final \( \theta \) that gives the lowest value of the cost +pick the final \( \beta \) that gives the lowest value of the cost function.

          @@ -794,7 +814,7 @@ j = 0 for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for beta t = epoch*m+i gamma_j = step_length(t,t0,t1) j += 1 diff --git a/doc/pub/Splines/ipynb/Splines.ipynb b/doc/pub/Splines/ipynb/Splines.ipynb index ca6d541d7..5328f5ba0 100644 --- a/doc/pub/Splines/ipynb/Splines.ipynb +++ b/doc/pub/Splines/ipynb/Splines.ipynb @@ -20,11 +20,11 @@ "## Optimization, the central part of any Machine Learning algortithm\n", "\n", "Almost every problem in machine learning and data science starts with\n", - "a dataset $X$, a model $g(\\theta)$, which is a function of the\n", - "parameters $\\theta$ and a cost function $C(X, g(\\theta))$ that allows\n", - "us to judge how well the model $g(\\theta)$ explains the observations\n", - "$X$. The model is fit by finding the values of $\\theta$ that minimize\n", - "the cost function. Ideally we would be able to solve for $\\theta$\n", + "a dataset $X$, a model $g(\\beta)$, which is a function of the\n", + "parameters $\\beta$ and a cost function $C(X, g(\\beta))$ that allows\n", + "us to judge how well the model $g(\\beta)$ explains the observations\n", + "$X$. The model is fit by finding the values of $\\beta$ that minimize\n", + "the cost function. Ideally we would be able to solve for $\\beta$\n", "analytically, however this is not possible in general and we must use\n", "some approximative/numerical method to compute the minimum.\n", "\n", @@ -61,7 +61,7 @@ "The previous observation is the basis of the method of steepest\n", "descent, which is also referred to as just gradient descent (GD). One\n", "starts with an initial guess $\\mathbf{x}_0$ for a minimum of $F$ and\n", - "compute new approximations according to" + "computes new approximations according to" ] }, { @@ -78,14 +78,25 @@ "metadata": {}, "source": [ "The parameter $\\gamma_k$ is often referred to as the step length or\n", - "the learning rate in the context of Machine Learning.\n", + "the learning rate within the context of Machine Learning.\n", "\n", "\n", "## The ideal\n", "\n", - "Ideally the sequence $\\{ \\mathbf{x}_k \\}_{k=0}$ converges to a global minimum of the function $F$. In general we do not know if we are in a global or local minimum. In the special case when $F$ is a convex function, all local minima are also global minima, so in this case gradient descent can converge to the global solution. The advantage of this scheme is that it is conceptually simple and straightforward to implement. However the method in this form has some severe limitations:\n", + "Ideally the sequence $\\{ \\mathbf{x}_k \\}_{k=0}$ converges to a global\n", + "minimum of the function $F$. In general we do not know if we are in a\n", + "global or local minimum. In the special case when $F$ is a convex\n", + "function, all local minima are also global minima, so in this case\n", + "gradient descent can converge to the global solution. The advantage of\n", + "this scheme is that it is conceptually simple and straightforward to\n", + "implement. However the method in this form has some severe\n", + "limitations:\n", "\n", - "In machine learing we are often faced with non-convex high dimensional cost functions with many local minimum. Since GD is deterministic we will get stuck in a local minimum, if the method converges, unless we have a very good intial guess. This also implies that the scheme is sensitive to the chosen initial condition.\n", + "In machine learing we are often faced with non-convex high dimensional\n", + "cost functions with many local minima. Since GD is deterministic we\n", + "will get stuck in a local minimum, if the method converges, unless we\n", + "have a very good intial guess. This also implies that the scheme is\n", + "sensitive to the chosen initial condition.\n", "\n", "Note that the gradient is a function of $\\mathbf{x} =\n", "(x_1,\\cdots,x_n)$ which makes it expensive to compute numerically.\n", @@ -97,15 +108,16 @@ "GD is sensitive to the choice of learning rate $\\gamma_k$. This is due\n", "to the fact that we are only guaranteed that $F(\\mathbf{x}_{k+1}) \\leq\n", "F(\\mathbf{x}_k)$ for sufficiently small $\\gamma_k$. The problem is to\n", - "determine an optimal learning rate. If the learning rate is chosen to\n", - "small the method will take a long to converge and if it is to large we\n", - "can experience erratic behavior.\n", + "determine an optimal learning rate. If the learning rate is chosen too\n", + "small the method will take a long time to converge and if it is too\n", + "large we can experience erratic behavior.\n", "\n", "Many of these shortcomings can be alleviated by introducing\n", "randomness. One such method is that of Stochastic Gradient Descent\n", - "(SGD), see below\n", + "(SGD), see below.\n", "\n", "## Gradient Descent Example\n", + "\n", "We revisit now our simple linear regression example with a linear polynomial." ] }, @@ -133,30 +145,30 @@ "y = 4+3*x+np.random.randn(100,1)\n", "\n", "xb = np.c_[np.ones((100,1)), x]\n", - "theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)\n", - "print(theta_linreg)\n", - "theta = np.random.randn(2,1)\n", + "beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)\n", + "print(beta_linreg)\n", + "beta = np.random.randn(2,1)\n", "\n", "eta = 0.1\n", "Niterations = 1000\n", "m = 100\n", "\n", "for iter in range(Niterations):\n", - " gradients = 2.0/m*xb.T.dot(xb.dot(theta)-y)\n", - " theta -= eta*gradients\n", + " gradients = 2.0/m*xb.T.dot(xb.dot(beta)-y)\n", + " beta -= eta*gradients\n", "\n", - "print(theta)\n", + "print(beta)\n", "xnew = np.array([[0],[2]])\n", "xbnew = np.c_[np.ones((2,1)), xnew]\n", - "ypredict = xbnew.dot(theta)\n", - "ypredict2 = xbnew.dot(theta_linreg)\n", + "ypredict = xbnew.dot(beta)\n", + "ypredict2 = xbnew.dot(beta_linreg)\n", "plt.plot(xnew, ypredict, \"r-\")\n", "plt.plot(xnew, ypredict2, \"b-\")\n", "plt.plot(x, y ,'ro')\n", "plt.axis([0,2.0,0, 15.0])\n", "plt.xlabel(r'$x$')\n", "plt.ylabel(r'$y$')\n", - "plt.title(r'Random numbers ')\n", + "plt.title(r'Gradient descent example')\n", "plt.show()" ] }, @@ -185,8 +197,8 @@ "y = 4+3*x+np.random.randn(100,1)\n", "\n", "xb = np.c_[np.ones((100,1)), x]\n", - "theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)\n", - "print(theta_linreg)\n", + "beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)\n", + "print(beta_linreg)\n", "sgdreg = SGDRegressor(n_iter = 50, penalty=None, eta0=0.1)\n", "sgdreg.fit(x,y.ravel())\n", "print(sgdreg.intercept_, sgdreg.coef_)" @@ -198,6 +210,7 @@ "source": [ "\n", "## Convex functions\n", + "\n", "Ideally we want our cost/loss function to be convex(concave).\n", "\n", "First we give the definition of a convex set: A set $C$ in\n", @@ -212,16 +225,14 @@ "\n", "## Convex function\n", "\n", - "Convex function: Let $X \\subset \\mathbb{R}^n$ be a convex set. Assume that the function $f: X \\rightarrow \\mathbb{R}$ is continuous, then $f$ is said to be convex if $$f(tx_1 + (1-t)x_2) \\leq tf(x_1) + (1-t)f(x_2) $$ for all $x_1, x_2 \\in X$ and for all $t \\in [0,1]$. If $\\leq$ is replaced with a strict inequaltiy in the definition, we demand $x_1 \\neq x_2$ and $t\\in(0,1)$ then $f$ is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting $f(x_1)$ and $f(x_2)$, the value of the function on the interval $[x_1,x_2]$ is always below the line as illustrated below.\n", + "**Convex function**: Let $X \\subset \\mathbb{R}^n$ be a convex set. Assume that the function $f: X \\rightarrow \\mathbb{R}$ is continuous, then $f$ is said to be convex if $$f(tx_1 + (1-t)x_2) \\leq tf(x_1) + (1-t)f(x_2) $$ for all $x_1, x_2 \\in X$ and for all $t \\in [0,1]$. If $\\leq$ is replaced with a strict inequaltiy in the definition, we demand $x_1 \\neq x_2$ and $t\\in(0,1)$ then $f$ is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting $f(x_1)$ and $f(x_2)$, the value of the function on the interval $[x_1,x_2]$ is always below the line as illustrated below.\n", "\n", "## Conditions on convex functions\n", "\n", "In the following we state first and second-order conditions which\n", "ensures convexity of a function $f$. We write $D_f$ to denote the\n", "domain of $f$, i.e the subset of $R^n$ where $f$ is defined. For more\n", - "details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex\n", - "Optimization. Cambridge University Press, \n", - "boyd/cvxbook/, 2004.\n", + "details and proofs we refer to: [S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press](http://stanford.edu/boyd/cvxbook/, 2004).\n", "\n", "**First order condition.**\n", "\n", @@ -255,9 +266,11 @@ "The next result is of great importance to us and the reason why we are\n", "going on about convex functions. In machine learning we frequently\n", "have to minimize a loss/cost function in order to find the best\n", - "parameters for the model we are considering. Ideally we want the\n", - "global minimum, however for high-dimensional models it is hard to know\n", - "if we have local or global minimum. However, if the cost/loss function\n", + "parameters for the model we are considering. \n", + "\n", + "Ideally we want the\n", + "global minimum (for high-dimensional models it is hard to know\n", + "if we have local or global minimum). However, if the cost/loss function\n", "is convex the following result provides invaluable information:\n", "\n", "**Any minimum is global for convex functions.**\n", @@ -267,30 +280,40 @@ "$x^*$ that satisfies $\\nabla f(x^*) = 0$ is a global minimum.\n", "\n", "\n", + "\n", "This result means that if we know that the cost/loss function is convex and we are able to find a minimum, we are guaranteed that it is a global minimum.\n", "\n", "## Some simple problems\n", - "1. Show that $f(x)=x^2$ is convex for $x \\in \\mathbb{R}$ using the definition of convexity.\n", "\n", - "Hint: If you re-write the definition, $f$ is convex if the following holds for all $x,y \\in D_f$ and any $\\lambda \\in [0,1] $ $$\\lambda f(x) + (1-\\lambda)f(y) - f(\\lambda x + (1-\\lambda) y ) \\geq 0. $$\n", + "1. Show that $f(x)=x^2$ is convex for $x \\in \\mathbb{R}$ using the definition of convexity. Hint: If you re-write the definition, $f$ is convex if the following holds for all $x,y \\in D_f$ and any $\\lambda \\in [0,1] $ $\\lambda f(x) + (1-\\lambda)f(y) - f(\\lambda x + (1-\\lambda) y ) \\geq 0. $\n", "\n", - "1. Using the second order condition show that the following functions are convex on the specified domain.\n", + "2. Using the second order condition show that the following functions are convex on the specified domain.\n", "\n", - "$f(x) = e^x$ is convex for $x \\in \\mathbb{R}$.\n", - "$g(x) = -\\ln(x)$ is convex for $x \\in (0,\\infty)$.\n", - "1. Let $f(x) = x^2$ and $g(x) = e^x$. Show that $f(g(x))$ and $g(f(x))$ is convex for $x \\in \\mathbb{R}$. Also show that if $f(x)$ is any convex function than $h(x) = e^{f(x)}$ is convex.\n", + " * $f(x) = e^x$ is convex for $x \\in \\mathbb{R}$.\n", "\n", - "2. A norm is any function that satisfy the following properties\n", + " * $g(x) = -\\ln(x)$ is convex for $x \\in (0,\\infty)$.\n", "\n", - "$f(\\alpha x) = |\\alpha| f(x)$ for all $\\alpha \\in \\mathbb{R}$.\n", - "$f(x+y) \\leq f(x) + f(y)$\n", - "$f(x) \\leq 0$ for all $x \\in \\mathbb{R}^n$ with equality if and only if $x = 0$\n", - "Using the definition of convexity, show that a function satisfying the properties above is convex (the third condition is not needed to show this).\n", + "\n", + "3. Let $f(x) = x^2$ and $g(x) = e^x$. Show that $f(g(x))$ and $g(f(x))$ is convex for $x \\in \\mathbb{R}$. Also show that if $f(x)$ is any convex function than $h(x) = e^{f(x)}$ is convex.\n", + "\n", + "4. A norm is any function that satisfy the following properties\n", + "\n", + " * $f(\\alpha x) = |\\alpha| f(x)$ for all $\\alpha \\in \\mathbb{R}$.\n", + "\n", + " * $f(x+y) \\leq f(x) + f(y)$\n", + "\n", + " * $f(x) \\leq 0$ for all $x \\in \\mathbb{R}^n$ with equality if and only if $x = 0$\n", + "\n", + "\n", + "Using the definition of convexity, try to show that a function satisfying the properties above is convex (the third condition is not needed to show this).\n", "\n", "\n", "## Revisiting our first homework\n", "\n", - "We will use linear regression as a case study for the gradient descent methods. Linear regression is a great test case for the gradient descent methods discussed in the lectures since it has several desirable properties such as:\n", + "We will use linear regression as a case study for the gradient descent\n", + "methods. Linear regression is a great test case for the gradient\n", + "descent methods discussed in the lectures since it has several\n", + "desirable properties such as:\n", "\n", "1. An analytical solution (recall homework set 1).\n", "\n", @@ -323,7 +346,7 @@ "metadata": {}, "source": [ "$$\n", - "h_\\theta(x) = \\hat{y} = \\theta_0 + \\theta_1 x,\n", + "h_\\beta(x) = \\hat{y} = \\beta_0 + \\beta_1 x,\n", "$$" ] }, @@ -339,7 +362,7 @@ "metadata": {}, "source": [ "$$\n", - "\\hat{y}_i = \\theta_0 + \\theta_1 x_i.\n", + "\\hat{y}_i = \\beta_0 + \\beta_1 x_i.\n", "$$" ] }, @@ -350,9 +373,9 @@ "\n", "## Gradient descent example\n", "\n", - "Let $\\mathbf{y} = (y_1,\\cdots,y_n)^T$, $\\mathbf{\\hat{y}} = (\\hat{y}_1,\\cdots,\\hat{y}_n)^T$ and $\\theta = (\\theta_0, \\theta_1)^T$\n", + "Let $\\mathbf{y} = (y_1,\\cdots,y_n)^T$, $\\mathbf{\\hat{y}} = (\\hat{y}_1,\\cdots,\\hat{y}_n)^T$ and $\\beta = (\\beta_0, \\beta_1)^T$\n", "\n", - "t is convenient to write $\\mathbf{\\hat{y}} = X\\theta$ where $X \\in \\mathbb{R}^{100 \\times 2} $ is the design matrix given by" + "t is convenient to write $\\mathbf{\\hat{y}} = X\\beta$ where $X \\in \\mathbb{R}^{100 \\times 2} $ is the design matrix given by" ] }, { @@ -386,7 +409,7 @@ "metadata": {}, "source": [ "$$\n", - "C(\\theta) = ||X\\theta-\\mathbf{y}||^2 = ||X\\theta||^2 - 2 \\mathbf{y}^T X\\theta + ||\\mathbf{y}||^2 = \\sum_{i=1}^{100} (\\theta_0 + \\theta_1 x_i)^2 - 2 y_i (\\theta_0 + \\theta_1 x_i) + y_i^2\n", + "C(\\beta) = ||X\\beta-\\mathbf{y}||^2 = ||X\\beta||^2 - 2 \\mathbf{y}^T X\\beta + ||\\mathbf{y}||^2 = \\sum_{i=1}^{100} (\\beta_0 + \\beta_1 x_i)^2 - 2 y_i (\\beta_0 + \\beta_1 x_i) + y_i^2\n", "$$" ] }, @@ -394,11 +417,11 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "and we want to find $\\theta$ such that $C(\\theta)$ is minimized.\n", + "and we want to find $\\beta$ such that $C(\\beta)$ is minimized.\n", "\n", "## The derivative of the cost/loss function\n", "\n", - "Computing $\\partial C(\\theta) / \\partial \\theta_0$ and $\\partial C(\\theta) / \\partial \\theta_1$ we can show that the gradient can be written as" + "Computing $\\partial C(\\beta) / \\partial \\beta_0$ and $\\partial C(\\beta) / \\partial \\beta_1$ we can show that the gradient can be written as" ] }, { @@ -406,9 +429,9 @@ "metadata": {}, "source": [ "$$\n", - "\\nabla_\\theta C(\\theta) = (\\partial C(\\theta) / \\partial \\theta_0, \\partial C(\\theta) / \\partial \\theta_1)^T = 2\\begin{bmatrix} \\sum_{i=1}^{100} \\left(\\theta_0+\\theta_1x_i-y_i\\right) \\\\\n", - "\\sum_{i=1}^{100}\\left( x_i (\\theta_0+\\theta_1x_i)-y_ix_i\\right) \\\\\n", - "\\end{bmatrix} = 2X^T(X\\theta - \\mathbf{y}),\n", + "\\nabla_\\beta C(\\beta) = (\\partial C(\\beta) / \\partial \\beta_0, \\partial C(\\beta) / \\partial \\beta_1)^T = 2\\begin{bmatrix} \\sum_{i=1}^{100} \\left(\\beta_0+\\beta_1x_i-y_i\\right) \\\\\n", + "\\sum_{i=1}^{100}\\left( x_i (\\beta_0+\\beta_1x_i)-y_ix_i\\right) \\\\\n", + "\\end{bmatrix} = 2X^T(X\\beta - \\mathbf{y}),\n", "$$" ] }, @@ -419,7 +442,7 @@ "where $X$ is the design matrix defined above.\n", "\n", "## The Hessian matrix\n", - "The Hessian matrix of $C(\\theta)$ is given by" + "The Hessian matrix of $C(\\beta)$ is given by" ] }, { @@ -428,8 +451,8 @@ "source": [ "$$\n", "\\hat{H} \\equiv \\begin{bmatrix}\n", - "\\frac{\\partial^2 C(\\theta)}{\\partial \\theta_0^2} & \\frac{\\partial^2 C(\\theta)}{\\partial \\theta_0 \\partial \\theta_1} \\\\\n", - "\\frac{\\partial^2 C(\\theta)}{\\partial \\theta_0 \\partial \\theta_1} & \\frac{\\partial^2 C(\\theta)}{\\partial \\theta_1^2} & \\\\\n", + "\\frac{\\partial^2 C(\\beta)}{\\partial \\beta_0^2} & \\frac{\\partial^2 C(\\beta)}{\\partial \\beta_0 \\partial \\beta_1} \\\\\n", + "\\frac{\\partial^2 C(\\beta)}{\\partial \\beta_0 \\partial \\beta_1} & \\frac{\\partial^2 C(\\beta)}{\\partial \\beta_1^2} & \\\\\n", "\\end{bmatrix} = 2X^T X.\n", "$$" ] @@ -438,11 +461,11 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "This result implies that $C(\\theta)$ is a convex function since the matrix $X^T X$ always is positive semi-definite.\n", + "This result implies that $C(\\beta)$ is a convex function since the matrix $X^T X$ always is positive semi-definite.\n", "\n", "## Simple program\n", "\n", - "We can now write a program that minimizes $C(\\theta)$ using the gradient descent method with a constant learning rate $\\gamma$ according to" + "We can now write a program that minimizes $C(\\beta)$ using the gradient descent method with a constant learning rate $\\gamma$ according to" ] }, { @@ -450,7 +473,7 @@ "metadata": {}, "source": [ "$$\n", - "\\theta_{k+1} = \\theta_k - \\gamma \\nabla_\\theta C(\\theta_k), \\ k=0,1,\\cdots\n", + "\\beta_{k+1} = \\beta_k - \\gamma \\nabla_\\beta C(\\beta_k), \\ k=0,1,\\cdots\n", "$$" ] }, @@ -459,11 +482,11 @@ "metadata": {}, "source": [ "We can use the expression we computed for the gradient and let use a\n", - "$\\theta_0$ be chosen randomly and let $\\gamma = 0.001$. Stop iterating\n", - "when $||\\nabla_\\theta C(\\theta_k) || < \\epsilon = 10^{-8}$. \n", + "$\\beta_0$ be chosen randomly and let $\\gamma = 0.001$. Stop iterating\n", + "when $||\\nabla_\\beta C(\\beta_k) || < \\epsilon = 10^{-8}$. \n", "\n", - "And finally we can compare our solution for $\\theta$ with the analytic result given by \n", - "$\\theta= (X^TX)^{-1} X^T \\mathbf{y}$." + "And finally we can compare our solution for $\\beta$ with the analytic result given by \n", + "$\\beta= (X^TX)^{-1} X^T \\mathbf{y}$." ] }, { @@ -487,11 +510,11 @@ "y = 5*x**2 + 0.1*np.random.randn(N)\n", "X = np.c_[np.ones(N),x] #Construct design matrix\n", "\n", - "#Compute theta according to normal equations to compare with GD solution\n", + "#Compute beta according to normal equations to compare with GD solution\n", "Xt_X_inv = np.linalg.inv(np.dot(X.T,X))\n", "Xt_y = np.dot(X.transpose(),y)\n", - "theta_NE = np.dot(Xt_X_inv,Xt_y)\n", - "print(theta_NE)" + "beta_NE = np.dot(Xt_X_inv,Xt_y)\n", + "print(beta_NE)" ] }, { @@ -501,7 +524,7 @@ "\n", "## Gradient descent and Ridge\n", "\n", - "We have also discussed Ridge regression where the loss function contains a regularized given by the $L_2$ norm of $\\theta$," + "We have also discussed Ridge regression where the loss function contains a regularized given by the $L_2$ norm of $\\beta$," ] }, { @@ -509,7 +532,7 @@ "metadata": {}, "source": [ "$$\n", - "C_{\\text{ridge}}(\\theta) = ||X\\theta -\\mathbf{y}||^2 + \\lambda ||\\theta||^2, \\ \\lambda \\geq 0.\n", + "C_{\\text{ridge}}(\\beta) = ||X\\beta -\\mathbf{y}||^2 + \\lambda ||\\beta||^2, \\ \\lambda \\geq 0.\n", "$$" ] }, @@ -517,7 +540,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "In order to minimize $C_{\\text{ridge}}(\\theta)$ using GD we only have adjust the gradient as follows" + "In order to minimize $C_{\\text{ridge}}(\\beta)$ using GD we only have adjust the gradient as follows" ] }, { @@ -525,9 +548,9 @@ "metadata": {}, "source": [ "$$\n", - "\\nabla_\\theta C_{\\text{ridge}}(\\theta) = 2\\begin{bmatrix} \\sum_{i=1}^{100} \\left(\\theta_0+\\theta_1x_i-y_i\\right) \\\\\n", - "\\sum_{i=1}^{100}\\left( x_i (\\theta_0+\\theta_1x_i)-y_ix_i\\right) \\\\\n", - "\\end{bmatrix} + 2\\lambda\\begin{bmatrix} \\theta_0 \\\\ \\theta_1\\end{bmatrix} = 2 (X^T(X\\theta - \\mathbf{y})+\\lambda \\theta).\n", + "\\nabla_\\beta C_{\\text{ridge}}(\\beta) = 2\\begin{bmatrix} \\sum_{i=1}^{100} \\left(\\beta_0+\\beta_1x_i-y_i\\right) \\\\\n", + "\\sum_{i=1}^{100}\\left( x_i (\\beta_0+\\beta_1x_i)-y_ix_i\\right) \\\\\n", + "\\end{bmatrix} + 2\\lambda\\begin{bmatrix} \\beta_0 \\\\ \\beta_1\\end{bmatrix} = 2 (X^T(X\\beta - \\mathbf{y})+\\lambda \\beta).\n", "$$" ] }, @@ -535,7 +558,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "We can now extend our program to minimize $C_{\\text{ridge}}(\\theta)$ using gradient descent and compare with the analytical solution given by" + "We can now extend our program to minimize $C_{\\text{ridge}}(\\beta)$ using gradient descent and compare with the analytical solution given by" ] }, { @@ -543,7 +566,7 @@ "metadata": {}, "source": [ "$$\n", - "\\theta_{\\text{ridge}} = \\left(X^T X + \\lambda I_{2 \\times 2} \\right)^{-1} X^T \\mathbf{y},\n", + "\\beta_{\\text{ridge}} = \\left(X^T X + \\lambda I_{2 \\times 2} \\right)^{-1} X^T \\mathbf{y},\n", "$$" ] }, @@ -552,7 +575,7 @@ "metadata": {}, "source": [ "for $\\lambda = {0,1,10,50,100}$ ($\\lambda = 0$ corresponds to ordinary least squares). \n", - "We can then compute $||\\theta_{\\text{ridge}}||$ for each $\\lambda$." + "We can then compute $||\\beta_{\\text{ridge}}||$ for each $\\lambda$." ] }, { @@ -576,7 +599,7 @@ "y = 5*x**2 + 0.1*np.random.randn(N)\n", "\n", "\n", - "#Compute analytic theta for Ridge regression \n", + "#Compute analytic beta for Ridge regression \n", "X = np.c_[np.ones(N),x]\n", "XT_X = np.dot(X.T,X)\n", "\n", @@ -584,10 +607,10 @@ "Id = np.eye(XT_X.shape[0])\n", "\n", "Z = np.linalg.inv(XT_X+l*Id)\n", - "theta_ridge = np.dot(Z,np.dot(X.T,y))\n", + "beta_ridge = np.dot(Z,np.dot(X.T,y))\n", "\n", - "print(theta_ridge)\n", - "print(np.linalg.norm(theta_ridge)) #||theta||" + "print(beta_ridge)\n", + "print(np.linalg.norm(beta_ridge)) #||beta||" ] }, { @@ -609,8 +632,8 @@ "metadata": {}, "source": [ "$$\n", - "C(\\mathbf{\\theta}) = \\sum_{i=1}^n c_i(\\mathbf{x}_i,\n", - "\\mathbf{\\theta}).\n", + "C(\\mathbf{\\beta}) = \\sum_{i=1}^n c_i(\\mathbf{x}_i,\n", + "\\mathbf{\\beta}).\n", "$$" ] }, @@ -629,8 +652,8 @@ "metadata": {}, "source": [ "$$\n", - "\\nabla_\\theta C(\\mathbf{\\theta}) = \\sum_i^n \\nabla_\\theta c_i(\\mathbf{x}_i,\n", - "\\mathbf{\\theta}).\n", + "\\nabla_\\beta C(\\mathbf{\\beta}) = \\sum_i^n \\nabla_\\beta c_i(\\mathbf{x}_i,\n", + "\\mathbf{\\beta}).\n", "$$" ] }, @@ -664,10 +687,10 @@ "metadata": {}, "source": [ "$$\n", - "\\nabla_\\theta\n", - "C(\\mathbf{\\theta}) = \\sum_{i=1}^n \\nabla_\\theta c_i(\\mathbf{x}_i,\n", - "\\mathbf{\\theta}) \\rightarrow \\sum_{i \\in B_k}^n \\nabla_\\theta\n", - "c_i(\\mathbf{x}_i, \\mathbf{\\theta}).\n", + "\\nabla_\\beta\n", + "C(\\mathbf{\\beta}) = \\sum_{i=1}^n \\nabla_\\beta c_i(\\mathbf{x}_i,\n", + "\\mathbf{\\beta}) \\rightarrow \\sum_{i \\in B_k}^n \\nabla_\\beta\n", + "c_i(\\mathbf{x}_i, \\mathbf{\\beta}).\n", "$$" ] }, @@ -685,8 +708,8 @@ "metadata": {}, "source": [ "$$\n", - "\\theta_{j+1} = \\theta_j - \\gamma_j \\sum_{i \\in B_k}^n \\nabla_\\theta c_i(\\mathbf{x}_i,\n", - "\\mathbf{\\theta})\n", + "\\beta_{j+1} = \\beta_j - \\gamma_j \\sum_{i \\in B_k}^n \\nabla_\\beta c_i(\\mathbf{x}_i,\n", + "\\mathbf{\\beta})\n", "$$" ] }, @@ -723,7 +746,7 @@ " for i in range(m):\n", " k = np.random.randint(m) #Pick the k-th minibatch at random\n", " #Compute the gradient using the data in minibatch Bk\n", - " #Compute new suggestion for theta\n", + " #Compute new suggestion for \n", " j += 1" ] }, @@ -749,7 +772,7 @@ "that we are close to a local/global minimum. However, we could also\n", "evaluate the cost function at this point, store the result and\n", "continue the search. If the test kicks in at a later stage we can\n", - "compare the values of the cost function and keep the $\\theta$ that\n", + "compare the values of the cost function and keep the $\\beta$ that\n", "gave the lowest value.\n", "\n", "## Slightly different approach\n", @@ -760,10 +783,10 @@ "\n", "As an example, let $e = 0,1,2,3,\\cdots$ denote the current epoch and let $t_0, t_1 > 0$ be two fixed numbers. Furthermore, let $t = e \\cdot m + i$ where $m$ is the number of minibatches and $i=0,\\cdots,m-1$. Then the function $$\\gamma_j(t; t_0, t_1) = \\frac{t_0}{t+t_1} $$ goes to zero as the number of epochs gets large. I.e. we start with a step length $\\gamma_j (0; t_0, t_1) = t_0/t_1$ which decays in *time* $t$.\n", "\n", - "In this way we can fix the number of epochs, compute $\\theta$ and\n", + "In this way we can fix the number of epochs, compute $\\beta$ and\n", "evaluate the cost function at the end. Repeating the computation will\n", "give a different result since the scheme is random by design. Then we\n", - "pick the final $\\theta$ that gives the lowest value of the cost\n", + "pick the final $\\beta$ that gives the lowest value of the cost\n", "function." ] }, @@ -793,7 +816,7 @@ " for i in range(m):\n", " k = np.random.randint(m) #Pick the k-th minibatch at random\n", " #Compute the gradient using the data in minibatch Bk\n", - " #Compute new suggestion for theta\n", + " #Compute new suggestion for beta\n", " t = epoch*m+i\n", " gamma_j = step_length(t,t0,t1)\n", " j += 1\n", diff --git a/doc/pub/Splines/ipynb/ipynb-Splines-src.tar.gz b/doc/pub/Splines/ipynb/ipynb-Splines-src.tar.gz index a32e88d05..a1b16d968 100644 Binary files a/doc/pub/Splines/ipynb/ipynb-Splines-src.tar.gz and b/doc/pub/Splines/ipynb/ipynb-Splines-src.tar.gz differ diff --git a/doc/pub/Splines/pdf/Splines-minted.pdf b/doc/pub/Splines/pdf/Splines-minted.pdf index cc3b8bc5f..c5ac6574a 100644 Binary files a/doc/pub/Splines/pdf/Splines-minted.pdf and b/doc/pub/Splines/pdf/Splines-minted.pdf differ diff --git a/doc/src/Splines/Splines.do.txt b/doc/src/Splines/Splines.do.txt index 9eb3dfa37..a1465d76d 100644 --- a/doc/src/Splines/Splines.do.txt +++ b/doc/src/Splines/Splines.do.txt @@ -7,11 +7,11 @@ DATE: today ===== Optimization, the central part of any Machine Learning algortithm ===== Almost every problem in machine learning and data science starts with -a dataset $X$, a model $g(\theta)$, which is a function of the -parameters $\theta$ and a cost function $C(X, g(\theta))$ that allows -us to judge how well the model $g(\theta)$ explains the observations -$X$. The model is fit by finding the values of $\theta$ that minimize -the cost function. Ideally we would be able to solve for $\theta$ +a dataset $X$, a model $g(\beta)$, which is a function of the +parameters $\beta$ and a cost function $C(X, g(\beta))$ that allows +us to judge how well the model $g(\beta)$ explains the observations +$X$. The model is fit by finding the values of $\beta$ that minimize +the cost function. Ideally we would be able to solve for $\beta$ analytically, however this is not possible in general and we must use some approximative/numerical method to compute the minimum. @@ -40,7 +40,7 @@ we are always moving towards smaller function values, i.e a minimum. The previous observation is the basis of the method of steepest descent, which is also referred to as just gradient descent (GD). One starts with an initial guess $\mathbf{x}_0$ for a minimum of $F$ and -compute new approximations according to +computes new approximations according to !bt \[ @@ -49,14 +49,25 @@ compute new approximations according to !et The parameter $\gamma_k$ is often referred to as the step length or -the learning rate in the context of Machine Learning. +the learning rate within the context of Machine Learning. !split ===== The ideal ===== -Ideally the sequence $\{ \mathbf{x}_k \}_{k=0}$ converges to a global minimum of the function $F$. In general we do not know if we are in a global or local minimum. In the special case when $F$ is a convex function, all local minima are also global minima, so in this case gradient descent can converge to the global solution. The advantage of this scheme is that it is conceptually simple and straightforward to implement. However the method in this form has some severe limitations: +Ideally the sequence $\{ \mathbf{x}_k \}_{k=0}$ converges to a global +minimum of the function $F$. In general we do not know if we are in a +global or local minimum. In the special case when $F$ is a convex +function, all local minima are also global minima, so in this case +gradient descent can converge to the global solution. The advantage of +this scheme is that it is conceptually simple and straightforward to +implement. However the method in this form has some severe +limitations: -In machine learing we are often faced with non-convex high dimensional cost functions with many local minimum. Since GD is deterministic we will get stuck in a local minimum, if the method converges, unless we have a very good intial guess. This also implies that the scheme is sensitive to the chosen initial condition. +In machine learing we are often faced with non-convex high dimensional +cost functions with many local minima. Since GD is deterministic we +will get stuck in a local minimum, if the method converges, unless we +have a very good intial guess. This also implies that the scheme is +sensitive to the chosen initial condition. Note that the gradient is a function of $\mathbf{x} = (x_1,\cdots,x_n)$ which makes it expensive to compute numerically. @@ -68,16 +79,17 @@ Note that the gradient is a function of $\mathbf{x} = GD is sensitive to the choice of learning rate $\gamma_k$. This is due to the fact that we are only guaranteed that $F(\mathbf{x}_{k+1}) \leq F(\mathbf{x}_k)$ for sufficiently small $\gamma_k$. The problem is to -determine an optimal learning rate. If the learning rate is chosen to -small the method will take a long to converge and if it is to large we -can experience erratic behavior. +determine an optimal learning rate. If the learning rate is chosen too +small the method will take a long time to converge and if it is too +large we can experience erratic behavior. Many of these shortcomings can be alleviated by introducing randomness. One such method is that of Stochastic Gradient Descent -(SGD), see below +(SGD), see below. !split ===== Gradient Descent Example ===== + We revisit now our simple linear regression example with a linear polynomial. !bc pycod @@ -94,30 +106,30 @@ x = 2*np.random.rand(100,1) y = 4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) -theta = np.random.randn(2,1) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) +beta = np.random.randn(2,1) eta = 0.1 Niterations = 1000 m = 100 for iter in range(Niterations): - gradients = 2.0/m*xb.T.dot(xb.dot(theta)-y) - theta -= eta*gradients + gradients = 2.0/m*xb.T.dot(xb.dot(beta)-y) + beta -= eta*gradients -print(theta) +print(beta) xnew = np.array([[0],[2]]) xbnew = np.c_[np.ones((2,1)), xnew] -ypredict = xbnew.dot(theta) -ypredict2 = xbnew.dot(theta_linreg) +ypredict = xbnew.dot(beta) +ypredict2 = xbnew.dot(beta_linreg) plt.plot(xnew, ypredict, "r-") plt.plot(xnew, ypredict2, "b-") plt.plot(x, y ,'ro') plt.axis([0,2.0,0, 15.0]) plt.xlabel(r'$x$') plt.ylabel(r'$y$') -plt.title(r'Random numbers ') +plt.title(r'Gradient descent example') plt.show() !ec @@ -136,8 +148,8 @@ x = 2*np.random.rand(100,1) y = 4+3*x+np.random.randn(100,1) xb = np.c_[np.ones((100,1)), x] -theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) -print(theta_linreg) +beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y) +print(beta_linreg) sgdreg = SGDRegressor(n_iter = 50, penalty=None, eta0=0.1) sgdreg.fit(x,y.ravel()) print(sgdreg.intercept_, sgdreg.coef_) @@ -146,6 +158,7 @@ print(sgdreg.intercept_, sgdreg.coef_) !split ===== Convex functions ===== + Ideally we want our cost/loss function to be convex(concave). First we give the definition of a convex set: A set $C$ in @@ -161,7 +174,7 @@ regular polygons (triangles, rectangles, pentagons, etc...). !split ===== Convex function ===== -Convex function: Let $X \subset \mathbb{R}^n$ be a convex set. Assume that the function $f: X \rightarrow \mathbb{R}$ is continuous, then $f$ is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all $x_1, x_2 \in X$ and for all $t \in [0,1]$. If $\leq$ is replaced with a strict inequaltiy in the definition, we demand $x_1 \neq x_2$ and $t\in(0,1)$ then $f$ is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting $f(x_1)$ and $f(x_2)$, the value of the function on the interval $[x_1,x_2]$ is always below the line as illustrated below. +_Convex function_: Let $X \subset \mathbb{R}^n$ be a convex set. Assume that the function $f: X \rightarrow \mathbb{R}$ is continuous, then $f$ is said to be convex if $$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$ for all $x_1, x_2 \in X$ and for all $t \in [0,1]$. If $\leq$ is replaced with a strict inequaltiy in the definition, we demand $x_1 \neq x_2$ and $t\in(0,1)$ then $f$ is said to be strictly convex. For a single variable function, convexity means that if you draw a straight line connecting $f(x_1)$ and $f(x_2)$, the value of the function on the interval $[x_1,x_2]$ is always below the line as illustrated below. !split ===== Conditions on convex functions ===== @@ -169,9 +182,7 @@ Convex function: Let $X \subset \mathbb{R}^n$ be a convex set. Assume that the f In the following we state first and second-order conditions which ensures convexity of a function $f$. We write $D_f$ to denote the domain of $f$, i.e the subset of $R^n$ where $f$ is defined. For more -details and proofs we refer to: S. Boyd and L. Vandenberghe. Convex -Optimization. Cambridge University Press, http://stanford.edu/ -boyd/cvxbook/, 2004. +details and proofs we refer to: "S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press":"http://stanford.edu/boyd/cvxbook/, 2004". !bblock First order condition Suppose $f$ is differentiable (i.e $\nabla f(x)$ is well defined for @@ -202,9 +213,11 @@ This condition is particularly useful since it gives us an procedure for determi The next result is of great importance to us and the reason why we are going on about convex functions. In machine learning we frequently have to minimize a loss/cost function in order to find the best -parameters for the model we are considering. Ideally we want the -global minimum, however for high-dimensional models it is hard to know -if we have local or global minimum. However, if the cost/loss function +parameters for the model we are considering. + +Ideally we want the +global minimum (for high-dimensional models it is hard to know +if we have local or global minimum). However, if the cost/loss function is convex the following result provides invaluable information: !bblock Any minimum is global for convex functions @@ -212,30 +225,33 @@ Consider the problem of finding $x \in \mathbb{R}^n$ such that $f(x)$ is minimal, where $f$ is convex and differentiable. Then, any point $x^*$ that satisfies $\nabla f(x^*) = 0$ is a global minimum. !eblock + This result means that if we know that the cost/loss function is convex and we are able to find a minimum, we are guaranteed that it is a global minimum. !split ===== Some simple problems ===== -o Show that $f(x)=x^2$ is convex for $x \in \mathbb{R}$ using the definition of convexity. -Hint: If you re-write the definition, $f$ is convex if the following holds for all $x,y \in D_f$ and any $\lambda \in [0,1] $ $$\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $$ + +o Show that $f(x)=x^2$ is convex for $x \in \mathbb{R}$ using the definition of convexity. Hint: If you re-write the definition, $f$ is convex if the following holds for all $x,y \in D_f$ and any $\lambda \in [0,1] $ $\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $ o Using the second order condition show that the following functions are convex on the specified domain. - -$f(x) = e^x$ is convex for $x \in \mathbb{R}$. -$g(x) = -\ln(x)$ is convex for $x \in (0,\infty)$. + * $f(x) = e^x$ is convex for $x \in \mathbb{R}$. + * $g(x) = -\ln(x)$ is convex for $x \in (0,\infty)$. o Let $f(x) = x^2$ and $g(x) = e^x$. Show that $f(g(x))$ and $g(f(x))$ is convex for $x \in \mathbb{R}$. Also show that if $f(x)$ is any convex function than $h(x) = e^{f(x)}$ is convex. o A norm is any function that satisfy the following properties + * $f(\alpha x) = |\alpha| f(x)$ for all $\alpha \in \mathbb{R}$. + * $f(x+y) \leq f(x) + f(y)$ + * $f(x) \leq 0$ for all $x \in \mathbb{R}^n$ with equality if and only if $x = 0$ -$f(\alpha x) = |\alpha| f(x)$ for all $\alpha \in \mathbb{R}$. -$f(x+y) \leq f(x) + f(y)$ -$f(x) \leq 0$ for all $x \in \mathbb{R}^n$ with equality if and only if $x = 0$ -Using the definition of convexity, show that a function satisfying the properties above is convex (the third condition is not needed to show this). +Using the definition of convexity, try to show that a function satisfying the properties above is convex (the third condition is not needed to show this). !split ===== Revisiting our first homework ===== -We will use linear regression as a case study for the gradient descent methods. Linear regression is a great test case for the gradient descent methods discussed in the lectures since it has several desirable properties such as: +We will use linear regression as a case study for the gradient descent +methods. Linear regression is a great test case for the gradient +descent methods discussed in the lectures since it has several +desirable properties such as: o An analytical solution (recall homework set 1). o The gradient can be computed analytically. @@ -251,22 +267,22 @@ with $x_i \in [0,1] $ chosen randomly with a uniform distribution. Additionally The linear regression model is given by !bt \[ -h_\theta(x) = \hat{y} = \theta_0 + \theta_1 x, +h_\beta(x) = \hat{y} = \beta_0 + \beta_1 x, \] !et such that !bt \[ -\hat{y}_i = \theta_0 + \theta_1 x_i. +\hat{y}_i = \beta_0 + \beta_1 x_i. \] !et !split ===== Gradient descent example ===== -Let $\mathbf{y} = (y_1,\cdots,y_n)^T$, $\mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T$ and $\theta = (\theta_0, \theta_1)^T$ +Let $\mathbf{y} = (y_1,\cdots,y_n)^T$, $\mathbf{\hat{y}} = (\hat{y}_1,\cdots,\hat{y}_n)^T$ and $\beta = (\beta_0, \beta_1)^T$ -t is convenient to write $\mathbf{\hat{y}} = X\theta$ where $X \in \mathbb{R}^{100 \times 2} $ is the design matrix given by +t is convenient to write $\mathbf{\hat{y}} = X\beta$ where $X \in \mathbb{R}^{100 \times 2} $ is the design matrix given by !bt \[ \begin{equation} @@ -281,53 +297,53 @@ X \equiv \begin{bmatrix} The loss function is given by !bt \[ -C(\theta) = ||X\theta-\mathbf{y}||^2 = ||X\theta||^2 - 2 \mathbf{y}^T X\theta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\theta_0 + \theta_1 x_i)^2 - 2 y_i (\theta_0 + \theta_1 x_i) + y_i^2 +C(\beta) = ||X\beta-\mathbf{y}||^2 = ||X\beta||^2 - 2 \mathbf{y}^T X\beta + ||\mathbf{y}||^2 = \sum_{i=1}^{100} (\beta_0 + \beta_1 x_i)^2 - 2 y_i (\beta_0 + \beta_1 x_i) + y_i^2 \] !et -and we want to find $\theta$ such that $C(\theta)$ is minimized. +and we want to find $\beta$ such that $C(\beta)$ is minimized. !split ===== The derivative of the cost/loss function ===== -Computing $\partial C(\theta) / \partial \theta_0$ and $\partial C(\theta) / \partial \theta_1$ we can show that the gradient can be written as +Computing $\partial C(\beta) / \partial \beta_0$ and $\partial C(\beta) / \partial \beta_1$ we can show that the gradient can be written as !bt \[ -\nabla_\theta C(\theta) = (\partial C(\theta) / \partial \theta_0, \partial C(\theta) / \partial \theta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} = 2X^T(X\theta - \mathbf{y}), +\nabla_\beta C(\beta) = (\partial C(\beta) / \partial \beta_0, \partial C(\beta) / \partial \beta_1)^T = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} = 2X^T(X\beta - \mathbf{y}), \] !et where $X$ is the design matrix defined above. !split ===== The Hessian matrix ===== -The Hessian matrix of $C(\theta)$ is given by +The Hessian matrix of $C(\beta)$ is given by !bt \[ \hat{H} \equiv \begin{bmatrix} -\frac{\partial^2 C(\theta)}{\partial \theta_0^2} & \frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} \\ -\frac{\partial^2 C(\theta)}{\partial \theta_0 \partial \theta_1} & \frac{\partial^2 C(\theta)}{\partial \theta_1^2} & \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0^2} & \frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} \\ +\frac{\partial^2 C(\beta)}{\partial \beta_0 \partial \beta_1} & \frac{\partial^2 C(\beta)}{\partial \beta_1^2} & \\ \end{bmatrix} = 2X^T X. \] !et -This result implies that $C(\theta)$ is a convex function since the matrix $X^T X$ always is positive semi-definite. +This result implies that $C(\beta)$ is a convex function since the matrix $X^T X$ always is positive semi-definite. !split ===== Simple program ===== -We can now write a program that minimizes $C(\theta)$ using the gradient descent method with a constant learning rate $\gamma$ according to +We can now write a program that minimizes $C(\beta)$ using the gradient descent method with a constant learning rate $\gamma$ according to !bt \[ -\theta_{k+1} = \theta_k - \gamma \nabla_\theta C(\theta_k), \ k=0,1,\cdots +\beta_{k+1} = \beta_k - \gamma \nabla_\beta C(\beta_k), \ k=0,1,\cdots \] !et We can use the expression we computed for the gradient and let use a -$\theta_0$ be chosen randomly and let $\gamma = 0.001$. Stop iterating -when $||\nabla_\theta C(\theta_k) || < \epsilon = 10^{-8}$. +$\beta_0$ be chosen randomly and let $\gamma = 0.001$. Stop iterating +when $||\nabla_\beta C(\beta_k) || < \epsilon = 10^{-8}$. -And finally we can compare our solution for $\theta$ with the analytic result given by -$\theta= (X^TX)^{-1} X^T \mathbf{y}$. +And finally we can compare our solution for $\beta$ with the analytic result given by +$\beta= (X^TX)^{-1} X^T \mathbf{y}$. !bc pycod import numpy as np @@ -342,40 +358,40 @@ x = np.random.rand(N) #Uniformly generated x-values in [0,1] y = 5*x**2 + 0.1*np.random.randn(N) X = np.c_[np.ones(N),x] #Construct design matrix -#Compute theta according to normal equations to compare with GD solution +#Compute beta according to normal equations to compare with GD solution Xt_X_inv = np.linalg.inv(np.dot(X.T,X)) Xt_y = np.dot(X.transpose(),y) -theta_NE = np.dot(Xt_X_inv,Xt_y) -print(theta_NE) +beta_NE = np.dot(Xt_X_inv,Xt_y) +print(beta_NE) !ec !split ===== Gradient descent and Ridge ===== -We have also discussed Ridge regression where the loss function contains a regularized given by the $L_2$ norm of $\theta$, +We have also discussed Ridge regression where the loss function contains a regularized given by the $L_2$ norm of $\beta$, !bt \[ -C_{\text{ridge}}(\theta) = ||X\theta -\mathbf{y}||^2 + \lambda ||\theta||^2, \ \lambda \geq 0. +C_{\text{ridge}}(\beta) = ||X\beta -\mathbf{y}||^2 + \lambda ||\beta||^2, \ \lambda \geq 0. \] !et -In order to minimize $C_{\text{ridge}}(\theta)$ using GD we only have adjust the gradient as follows +In order to minimize $C_{\text{ridge}}(\beta)$ using GD we only have adjust the gradient as follows !bt \[ -\nabla_\theta C_{\text{ridge}}(\theta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\theta_0+\theta_1x_i-y_i\right) \\ -\sum_{i=1}^{100}\left( x_i (\theta_0+\theta_1x_i)-y_ix_i\right) \\ -\end{bmatrix} + 2\lambda\begin{bmatrix} \theta_0 \\ \theta_1\end{bmatrix} = 2 (X^T(X\theta - \mathbf{y})+\lambda \theta). +\nabla_\beta C_{\text{ridge}}(\beta) = 2\begin{bmatrix} \sum_{i=1}^{100} \left(\beta_0+\beta_1x_i-y_i\right) \\ +\sum_{i=1}^{100}\left( x_i (\beta_0+\beta_1x_i)-y_ix_i\right) \\ +\end{bmatrix} + 2\lambda\begin{bmatrix} \beta_0 \\ \beta_1\end{bmatrix} = 2 (X^T(X\beta - \mathbf{y})+\lambda \beta). \] !et -We can now extend our program to minimize $C_{\text{ridge}}(\theta)$ using gradient descent and compare with the analytical solution given by +We can now extend our program to minimize $C_{\text{ridge}}(\beta)$ using gradient descent and compare with the analytical solution given by !bt \[ -\theta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, +\beta_{\text{ridge}} = \left(X^T X + \lambda I_{2 \times 2} \right)^{-1} X^T \mathbf{y}, \] !et for $\lambda = {0,1,10,50,100}$ ($\lambda = 0$ corresponds to ordinary least squares). -We can then compute $||\theta_{\text{ridge}}||$ for each $\lambda$. +We can then compute $||\beta_{\text{ridge}}||$ for each $\lambda$. !bc pycod import numpy as np @@ -391,7 +407,7 @@ x = np.random.rand(N) y = 5*x**2 + 0.1*np.random.randn(N) -#Compute analytic theta for Ridge regression +#Compute analytic beta for Ridge regression X = np.c_[np.ones(N),x] XT_X = np.dot(X.T,X) @@ -399,10 +415,10 @@ l = 0.1 #Ridge parameter lambda Id = np.eye(XT_X.shape[0]) Z = np.linalg.inv(XT_X+l*Id) -theta_ridge = np.dot(Z,np.dot(X.T,y)) +beta_ridge = np.dot(Z,np.dot(X.T,y)) -print(theta_ridge) -print(np.linalg.norm(theta_ridge)) #||theta|| +print(beta_ridge) +print(np.linalg.norm(beta_ridge)) #||beta|| !ec @@ -419,8 +435,8 @@ function, which we want to minimize, can almost always be written as a sum over $n$ datapoints $\{\mathbf{x}_i\}_{i=1}^n$, !bt \[ -C(\mathbf{\theta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, -\mathbf{\theta}). +C(\mathbf{\beta}) = \sum_{i=1}^n c_i(\mathbf{x}_i, +\mathbf{\beta}). \] !et @@ -431,8 +447,8 @@ This in turn means that the gradient can be computed as a sum over $i$-gradients !bt \[ -\nabla_\theta C(\mathbf{\theta}) = \sum_i^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}). +\nabla_\beta C(\mathbf{\beta}) = \sum_i^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}). \] !et @@ -458,10 +474,10 @@ all datapoints with a sum over the datapoints in one the minibatches picked at random in each gradient descent step !bt \[ -\nabla_\theta -C(\mathbf{\theta}) = \sum_{i=1}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) \rightarrow \sum_{i \in B_k}^n \nabla_\theta -c_i(\mathbf{x}_i, \mathbf{\theta}). +\nabla_\beta +C(\mathbf{\beta}) = \sum_{i=1}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) \rightarrow \sum_{i \in B_k}^n \nabla_\beta +c_i(\mathbf{x}_i, \mathbf{\beta}). \] !et @@ -471,8 +487,8 @@ c_i(\mathbf{x}_i, \mathbf{\theta}). Thus a gradient descent step now looks like !bt \[ -\theta_{j+1} = \theta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\theta c_i(\mathbf{x}_i, -\mathbf{\theta}) +\beta_{j+1} = \beta_j - \gamma_j \sum_{i \in B_k}^n \nabla_\beta c_i(\mathbf{x}_i, +\mathbf{\beta}) \] !et @@ -498,7 +514,7 @@ for epoch in range(1,n_epochs+1): for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for j += 1 !ec @@ -521,7 +537,7 @@ is zero is valid also for local minima, so this would only tell us that we are close to a local/global minimum. However, we could also evaluate the cost function at this point, store the result and continue the search. If the test kicks in at a later stage we can -compare the values of the cost function and keep the $\theta$ that +compare the values of the cost function and keep the $\beta$ that gave the lowest value. !split @@ -533,10 +549,10 @@ reasonable time such that we do not move at all. As an example, let $e = 0,1,2,3,\cdots$ denote the current epoch and let $t_0, t_1 > 0$ be two fixed numbers. Furthermore, let $t = e \cdot m + i$ where $m$ is the number of minibatches and $i=0,\cdots,m-1$. Then the function $$\gamma_j(t; t_0, t_1) = \frac{t_0}{t+t_1} $$ goes to zero as the number of epochs gets large. I.e. we start with a step length $\gamma_j (0; t_0, t_1) = t_0/t_1$ which decays in *time* $t$. -In this way we can fix the number of epochs, compute $\theta$ and +In this way we can fix the number of epochs, compute $\beta$ and evaluate the cost function at the end. Repeating the computation will give a different result since the scheme is random by design. Then we -pick the final $\theta$ that gives the lowest value of the cost +pick the final $\beta$ that gives the lowest value of the cost function. !bc pycod @@ -558,7 +574,7 @@ for epoch in range(1,n_epochs+1): for i in range(m): k = np.random.randint(m) #Pick the k-th minibatch at random #Compute the gradient using the data in minibatch Bk - #Compute new suggestion for theta + #Compute new suggestion for beta t = epoch*m+i gamma_j = step_length(t,t0,t1) j += 1