Update on splines
This commit is contained in:
mhjensen
@@ -178,11 +178,11 @@ MathJax.Hub.Config({
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
|
||||
@@ -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 @@ $$
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -177,10 +177,21 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec3" class="anchor">The ideal </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
Note that the gradient is a function of \( \mathbf{x} =
|
||||
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -175,6 +175,8 @@ MathJax.Hub.Config({
|
||||
<!-- !split -->
|
||||
|
||||
<h2 id="___sec5" class="anchor">Gradient Descent Example </h2>
|
||||
|
||||
<p>
|
||||
We revisit now our simple linear regression example with a linear polynomial.
|
||||
<p>
|
||||
|
||||
@@ -192,30 +194,30 @@ x <span style="color: #666666">=</span> <span style="color: #666666">2*</span>np
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">4+3*</span>x<span style="color: #666666">+</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)
|
||||
|
||||
xb <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones((<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)), x]
|
||||
theta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_linreg)
|
||||
theta <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">2</span>,<span style="color: #666666">1</span>)
|
||||
beta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_linreg)
|
||||
beta <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">2</span>,<span style="color: #666666">1</span>)
|
||||
|
||||
eta <span style="color: #666666">=</span> <span style="color: #666666">0.1</span>
|
||||
Niterations <span style="color: #666666">=</span> <span style="color: #666666">1000</span>
|
||||
m <span style="color: #666666">=</span> <span style="color: #666666">100</span>
|
||||
|
||||
<span style="color: #008000; font-weight: bold">for</span> <span style="color: #008000">iter</span> <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(Niterations):
|
||||
gradients <span style="color: #666666">=</span> <span style="color: #666666">2.0/</span>m<span style="color: #666666">*</span>xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>dot(theta)<span style="color: #666666">-</span>y)
|
||||
theta <span style="color: #666666">-=</span> eta<span style="color: #666666">*</span>gradients
|
||||
gradients <span style="color: #666666">=</span> <span style="color: #666666">2.0/</span>m<span style="color: #666666">*</span>xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>dot(beta)<span style="color: #666666">-</span>y)
|
||||
beta <span style="color: #666666">-=</span> eta<span style="color: #666666">*</span>gradients
|
||||
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta)
|
||||
xnew <span style="color: #666666">=</span> np<span style="color: #666666">.</span>array([[<span style="color: #666666">0</span>],[<span style="color: #666666">2</span>]])
|
||||
xbnew <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones((<span style="color: #666666">2</span>,<span style="color: #666666">1</span>)), xnew]
|
||||
ypredict <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(theta)
|
||||
ypredict2 <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(theta_linreg)
|
||||
ypredict <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(beta)
|
||||
ypredict2 <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(beta_linreg)
|
||||
plt<span style="color: #666666">.</span>plot(xnew, ypredict, <span style="color: #BA2121">"r-"</span>)
|
||||
plt<span style="color: #666666">.</span>plot(xnew, ypredict2, <span style="color: #BA2121">"b-"</span>)
|
||||
plt<span style="color: #666666">.</span>plot(x, y ,<span style="color: #BA2121">'ro'</span>)
|
||||
plt<span style="color: #666666">.</span>axis([<span style="color: #666666">0</span>,<span style="color: #666666">2.0</span>,<span style="color: #666666">0</span>, <span style="color: #666666">15.0</span>])
|
||||
plt<span style="color: #666666">.</span>xlabel(<span style="color: #BA2121">r'$x$'</span>)
|
||||
plt<span style="color: #666666">.</span>ylabel(<span style="color: #BA2121">r'$y$'</span>)
|
||||
plt<span style="color: #666666">.</span>title(<span style="color: #BA2121">r'Random numbers '</span>)
|
||||
plt<span style="color: #666666">.</span>title(<span style="color: #BA2121">r'Gradient descent example'</span>)
|
||||
plt<span style="color: #666666">.</span>show()
|
||||
</pre></div>
|
||||
<p>
|
||||
|
||||
@@ -189,8 +189,8 @@ x <span style="color: #666666">=</span> <span style="color: #666666">2*</span>np
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">4+3*</span>x<span style="color: #666666">+</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)
|
||||
|
||||
xb <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones((<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)), x]
|
||||
theta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_linreg)
|
||||
beta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_linreg)
|
||||
sgdreg <span style="color: #666666">=</span> SGDRegressor(n_iter <span style="color: #666666">=</span> <span style="color: #666666">50</span>, penalty<span style="color: #666666">=</span><span style="color: #008000">None</span>, eta0<span style="color: #666666">=0.1</span>)
|
||||
sgdreg<span style="color: #666666">.</span>fit(x,y<span style="color: #666666">.</span>ravel())
|
||||
<span style="color: #008000; font-weight: bold">print</span>(sgdreg<span style="color: #666666">.</span>intercept_, sgdreg<span style="color: #666666">.</span>coef_)
|
||||
|
||||
@@ -175,6 +175,8 @@ MathJax.Hub.Config({
|
||||
<!-- !split -->
|
||||
|
||||
<h2 id="___sec7" class="anchor">Convex functions </h2>
|
||||
|
||||
<p>
|
||||
Ideally we want our cost/loss function to be convex(concave).
|
||||
|
||||
<p>
|
||||
|
||||
@@ -177,7 +177,7 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec8" class="anchor">Convex function </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
<b>Convex function</b>: 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.
|
||||
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -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, <a href="http://stanford.edu/" target="_self"><tt>http://stanford.edu/</tt></a>
|
||||
boyd/cvxbook/, 2004.
|
||||
details and proofs we refer to: <a href="http://stanford.edu/boyd/cvxbook/, 2004" target="_self">S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press</a>.
|
||||
|
||||
<p>
|
||||
<div class="panel panel-default">
|
||||
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
@@ -195,6 +198,8 @@ is minimal, where \( f \) is convex and differentiable. Then, any point
|
||||
</div>
|
||||
</div>
|
||||
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
|
||||
@@ -177,27 +177,26 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec11" class="anchor">Some simple problems </h2>
|
||||
|
||||
<ol>
|
||||
<li> Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.</li>
|
||||
</ol>
|
||||
|
||||
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. $$
|
||||
|
||||
<ol>
|
||||
<li> 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. $</li>
|
||||
<li> Using the second order condition show that the following functions are convex on the specified domain.</li>
|
||||
</ol>
|
||||
|
||||
\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).
|
||||
\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).
|
||||
<ul>
|
||||
<li> \( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).</li>
|
||||
<li> \( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).</li>
|
||||
</ul>
|
||||
|
||||
<ol>
|
||||
<li> 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.</li>
|
||||
<li> A norm is any function that satisfy the following properties</li>
|
||||
|
||||
<ul>
|
||||
<li> \( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \).</li>
|
||||
<li> \( f(x+y) \leq f(x) + f(y) \)</li>
|
||||
<li> \( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \)</li>
|
||||
</ul>
|
||||
|
||||
</ol>
|
||||
|
||||
\( 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).
|
||||
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -177,7 +177,10 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec12" class="anchor">Revisiting our first homework </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<ol>
|
||||
<li> An analytical solution (recall homework set 1).</li>
|
||||
@@ -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.
|
||||
$$
|
||||
|
||||
<p>
|
||||
|
||||
@@ -177,10 +177,10 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec13" class="anchor">Gradient descent example </h2>
|
||||
|
||||
<p>
|
||||
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 \)
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -177,11 +177,11 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec14" class="anchor">The derivative of the cost/loss function </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
@@ -175,15 +175,15 @@ MathJax.Hub.Config({
|
||||
<!-- !split -->
|
||||
|
||||
<h2 id="___sec15" class="anchor">The Hessian matrix </h2>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -177,19 +177,19 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec16" class="anchor">Simple program </h2>
|
||||
|
||||
<p>
|
||||
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
|
||||
$$
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
<p>
|
||||
|
||||
<!-- code=python (!bc pycod) typeset with pygments style "default" -->
|
||||
@@ -206,11 +206,11 @@ x <span style="color: #666666">=</span> np<span style="color: #666666">.</span>
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">5*</span>x<span style="color: #666666">**2</span> <span style="color: #666666">+</span> <span style="color: #666666">0.1*</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(N)
|
||||
X <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones(N),x] <span style="color: #408080; font-style: italic">#Construct design matrix</span>
|
||||
|
||||
<span style="color: #408080; font-style: italic">#Compute theta according to normal equations to compare with GD solution</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute beta according to normal equations to compare with GD solution</span>
|
||||
Xt_X_inv <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,X))
|
||||
Xt_y <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>transpose(),y)
|
||||
theta_NE <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Xt_X_inv,Xt_y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_NE)
|
||||
beta_NE <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Xt_X_inv,Xt_y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_NE)
|
||||
</pre></div>
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -177,27 +177,27 @@ MathJax.Hub.Config({
|
||||
<h2 id="___sec17" class="anchor">Gradient descent and Ridge </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
$$
|
||||
|
||||
<p>
|
||||
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).
|
||||
$$
|
||||
|
||||
<p>
|
||||
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 \).
|
||||
|
||||
<p>
|
||||
|
||||
@@ -215,7 +215,7 @@ x <span style="color: #666666">=</span> np<span style="color: #666666">.</span>
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">5*</span>x<span style="color: #666666">**2</span> <span style="color: #666666">+</span> <span style="color: #666666">0.1*</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(N)
|
||||
|
||||
|
||||
<span style="color: #408080; font-style: italic">#Compute analytic theta for Ridge regression </span>
|
||||
<span style="color: #408080; font-style: italic">#Compute analytic beta for Ridge regression </span>
|
||||
X <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones(N),x]
|
||||
XT_X <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,X)
|
||||
|
||||
@@ -223,10 +223,10 @@ l <span style="color: #666666">=</span> <span style="color: #666666">0.1</span>
|
||||
Id <span style="color: #666666">=</span> np<span style="color: #666666">.</span>eye(XT_X<span style="color: #666666">.</span>shape[<span style="color: #666666">0</span>])
|
||||
|
||||
Z <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(XT_X<span style="color: #666666">+</span>l<span style="color: #666666">*</span>Id)
|
||||
theta_ridge <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Z,np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,y))
|
||||
beta_ridge <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Z,np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,y))
|
||||
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_ridge)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>norm(theta_ridge)) <span style="color: #408080; font-style: italic">#||theta||</span>
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_ridge)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>norm(beta_ridge)) <span style="color: #408080; font-style: italic">#||beta||</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
<p>
|
||||
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
|
||||
@@ -179,8 +179,8 @@ MathJax.Hub.Config({
|
||||
<p>
|
||||
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})
|
||||
$$
|
||||
|
||||
<p>
|
||||
|
||||
@@ -191,7 +191,7 @@ j <span style="color: #666666">=</span> <span style="color: #666666">0</span>
|
||||
<span style="color: #008000; font-weight: bold">for</span> i <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(m):
|
||||
k <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randint(m) <span style="color: #408080; font-style: italic">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for theta</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for </span>
|
||||
j <span style="color: #666666">+=</span> <span style="color: #666666">1</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
|
||||
@@ -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 <em>time</em> \( t \).
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
@@ -212,7 +212,7 @@ j <span style="color: #666666">=</span> <span style="color: #666666">0</span>
|
||||
<span style="color: #008000; font-weight: bold">for</span> i <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(m):
|
||||
k <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randint(m) <span style="color: #408080; font-style: italic">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for theta</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for beta</span>
|
||||
t <span style="color: #666666">=</span> epoch<span style="color: #666666">*</span>m<span style="color: #666666">+</span>i
|
||||
gamma_j <span style="color: #666666">=</span> step_length(t,t0,t1)
|
||||
j <span style="color: #666666">+=</span> <span style="color: #666666">1</span>
|
||||
|
||||
@@ -163,11 +163,11 @@ MathJax.Hub.Config({
|
||||
|
||||
<p>
|
||||
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.
|
||||
</section>
|
||||
@@ -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
|
||||
|
||||
<p> <br>
|
||||
$$
|
||||
@@ -214,7 +214,7 @@ $$
|
||||
|
||||
<p>
|
||||
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.
|
||||
</section>
|
||||
|
||||
|
||||
@@ -222,10 +222,21 @@ the learning rate in the context of Machine Learning.
|
||||
<h2 id="___sec3">The ideal </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
</section>
|
||||
|
||||
|
||||
<section>
|
||||
<h2 id="___sec5">Gradient Descent Example </h2>
|
||||
|
||||
<p>
|
||||
We revisit now our simple linear regression example with a linear polynomial.
|
||||
<p>
|
||||
|
||||
@@ -270,30 +283,30 @@ x = <span style="color: #B452CD">2</span>*np.random.rand(<span style="color: #B4
|
||||
y = <span style="color: #B452CD">4</span>+<span style="color: #B452CD">3</span>*x+np.random.randn(<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)
|
||||
|
||||
xb = np.c_[np.ones((<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)), x]
|
||||
theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_linreg)
|
||||
theta = np.random.randn(<span style="color: #B452CD">2</span>,<span style="color: #B452CD">1</span>)
|
||||
beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_linreg)
|
||||
beta = np.random.randn(<span style="color: #B452CD">2</span>,<span style="color: #B452CD">1</span>)
|
||||
|
||||
eta = <span style="color: #B452CD">0.1</span>
|
||||
Niterations = <span style="color: #B452CD">1000</span>
|
||||
m = <span style="color: #B452CD">100</span>
|
||||
|
||||
<span style="color: #8B008B; font-weight: bold">for</span> <span style="color: #658b00">iter</span> <span style="color: #8B008B">in</span> <span style="color: #658b00">range</span>(Niterations):
|
||||
gradients = <span style="color: #B452CD">2.0</span>/m*xb.T.dot(xb.dot(theta)-y)
|
||||
theta -= eta*gradients
|
||||
gradients = <span style="color: #B452CD">2.0</span>/m*xb.T.dot(xb.dot(beta)-y)
|
||||
beta -= eta*gradients
|
||||
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta)
|
||||
xnew = np.array([[<span style="color: #B452CD">0</span>],[<span style="color: #B452CD">2</span>]])
|
||||
xbnew = np.c_[np.ones((<span style="color: #B452CD">2</span>,<span style="color: #B452CD">1</span>)), xnew]
|
||||
ypredict = xbnew.dot(theta)
|
||||
ypredict2 = xbnew.dot(theta_linreg)
|
||||
ypredict = xbnew.dot(beta)
|
||||
ypredict2 = xbnew.dot(beta_linreg)
|
||||
plt.plot(xnew, ypredict, <span style="color: #CD5555">"r-"</span>)
|
||||
plt.plot(xnew, ypredict2, <span style="color: #CD5555">"b-"</span>)
|
||||
plt.plot(x, y ,<span style="color: #CD5555">'ro'</span>)
|
||||
plt.axis([<span style="color: #B452CD">0</span>,<span style="color: #B452CD">2.0</span>,<span style="color: #B452CD">0</span>, <span style="color: #B452CD">15.0</span>])
|
||||
plt.xlabel(<span style="color: #CD5555">r'$x$'</span>)
|
||||
plt.ylabel(<span style="color: #CD5555">r'$y$'</span>)
|
||||
plt.title(<span style="color: #CD5555">r'Random numbers '</span>)
|
||||
plt.title(<span style="color: #CD5555">r'Gradient descent example'</span>)
|
||||
plt.show()
|
||||
</pre></div>
|
||||
</section>
|
||||
@@ -315,8 +328,8 @@ x = <span style="color: #B452CD">2</span>*np.random.rand(<span style="color: #B4
|
||||
y = <span style="color: #B452CD">4</span>+<span style="color: #B452CD">3</span>*x+np.random.randn(<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)
|
||||
|
||||
xb = np.c_[np.ones((<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)), x]
|
||||
theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_linreg)
|
||||
beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_linreg)
|
||||
sgdreg = SGDRegressor(n_iter = <span style="color: #B452CD">50</span>, penalty=<span style="color: #658b00">None</span>, eta0=<span style="color: #B452CD">0.1</span>)
|
||||
sgdreg.fit(x,y.ravel())
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(sgdreg.intercept_, sgdreg.coef_)
|
||||
@@ -326,6 +339,8 @@ sgdreg.fit(x,y.ravel())
|
||||
|
||||
<section>
|
||||
<h2 id="___sec7">Convex functions </h2>
|
||||
|
||||
<p>
|
||||
Ideally we want our cost/loss function to be convex(concave).
|
||||
|
||||
<p>
|
||||
@@ -346,7 +361,7 @@ regular polygons (triangles, rectangles, pentagons, etc...).
|
||||
<h2 id="___sec8">Convex function </h2>
|
||||
|
||||
<p>
|
||||
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 <p> <br>
|
||||
<b>Convex function</b>: 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 <p> <br>
|
||||
$$f(tx_1 + (1-t)x_2) \leq tf(x_1) + (1-t)f(x_2) $$
|
||||
<p> <br> 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.
|
||||
</section>
|
||||
@@ -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, <a href="http://stanford.edu/" target="_blank"><tt>http://stanford.edu/</tt></a>
|
||||
boyd/cvxbook/, 2004.
|
||||
details and proofs we refer to: <a href="http://stanford.edu/boyd/cvxbook/, 2004" target="_blank">S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press</a>.
|
||||
|
||||
<p>
|
||||
<div class="alert alert-block alert-block alert-text-normal">
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
@@ -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.
|
||||
</div>
|
||||
|
||||
<p>
|
||||
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.
|
||||
</section>
|
||||
|
||||
@@ -426,32 +443,26 @@ This result means that if we know that the cost/loss function is convex and we a
|
||||
<h2 id="___sec11">Some simple problems </h2>
|
||||
|
||||
<ol>
|
||||
<p><li> Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.</li>
|
||||
</ol>
|
||||
<p>
|
||||
|
||||
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] \) <p> <br>
|
||||
$$\lambda f(x) + (1-\lambda)f(y) - f(\lambda x + (1-\lambda) y ) \geq 0. $$
|
||||
<p> <br>
|
||||
|
||||
<ol>
|
||||
<p><li> 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. $</li>
|
||||
<p><li> Using the second order condition show that the following functions are convex on the specified domain.</li>
|
||||
</ol>
|
||||
<p>
|
||||
|
||||
\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).
|
||||
\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).
|
||||
|
||||
<ol>
|
||||
<ul>
|
||||
<p><li> \( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).</li>
|
||||
<p><li> \( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).</li>
|
||||
</ul>
|
||||
<p><li> 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.</li>
|
||||
<p><li> A norm is any function that satisfy the following properties</li>
|
||||
|
||||
<ul>
|
||||
<p><li> \( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \).</li>
|
||||
<p><li> \( f(x+y) \leq f(x) + f(y) \)</li>
|
||||
<p><li> \( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \)</li>
|
||||
</ul>
|
||||
<p>
|
||||
</ol>
|
||||
<p>
|
||||
|
||||
\( 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).
|
||||
</section>
|
||||
|
||||
|
||||
@@ -459,7 +470,10 @@ Using the definition of convexity, show that a function satisfying the propertie
|
||||
<h2 id="___sec12">Revisiting our first homework </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<ol>
|
||||
<p><li> An analytical solution (recall homework set 1).</li>
|
||||
@@ -479,14 +493,14 @@ with \( x_i \in [0,1] \) chosen randomly with a uniform distribution. Additiona
|
||||
The linear regression model is given by
|
||||
<p> <br>
|
||||
$$
|
||||
h_\theta(x) = \hat{y} = \theta_0 + \theta_1 x,
|
||||
h_\beta(x) = \hat{y} = \beta_0 + \beta_1 x,
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
such that
|
||||
<p> <br>
|
||||
$$
|
||||
\hat{y}_i = \theta_0 + \theta_1 x_i.
|
||||
\hat{y}_i = \beta_0 + \beta_1 x_i.
|
||||
$$
|
||||
<p> <br>
|
||||
</section>
|
||||
@@ -496,10 +510,10 @@ $$
|
||||
<h2 id="___sec13">Gradient descent example </h2>
|
||||
|
||||
<p>
|
||||
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 \)
|
||||
|
||||
<p>
|
||||
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
|
||||
<p> <br>
|
||||
$$
|
||||
\begin{equation}
|
||||
@@ -516,11 +530,11 @@ $$
|
||||
The loss function is given by
|
||||
<p> <br>
|
||||
$$
|
||||
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
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
and we want to find \( \theta \) such that \( C(\theta) \) is minimized.
|
||||
and we want to find \( \beta \) such that \( C(\beta) \) is minimized.
|
||||
</section>
|
||||
|
||||
|
||||
@@ -528,12 +542,12 @@ and we want to find \( \theta \) such that \( C(\theta) \) is minimized.
|
||||
<h2 id="___sec14">The derivative of the cost/loss function </h2>
|
||||
|
||||
<p>
|
||||
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
|
||||
<p> <br>
|
||||
$$
|
||||
\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}),
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
@@ -543,17 +557,17 @@ where \( X \) is the design matrix defined above.
|
||||
|
||||
<section>
|
||||
<h2 id="___sec15">The Hessian matrix </h2>
|
||||
The Hessian matrix of \( C(\theta) \) is given by
|
||||
The Hessian matrix of \( C(\beta) \) is given by
|
||||
<p> <br>
|
||||
$$
|
||||
\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.
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
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.
|
||||
</section>
|
||||
|
||||
|
||||
@@ -561,21 +575,21 @@ This result implies that \( C(\theta) \) is a convex function since the matrix \
|
||||
<h2 id="___sec16">Simple program </h2>
|
||||
|
||||
<p>
|
||||
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
|
||||
<p> <br>
|
||||
$$
|
||||
\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
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
<p>
|
||||
|
||||
<!-- code=python (!bc pycod) typeset with pygments style "perldoc" -->
|
||||
@@ -592,11 +606,11 @@ x = np.random.rand(N) <span style="color: #228B22">#Uniformly generated x-value
|
||||
y = <span style="color: #B452CD">5</span>*x**<span style="color: #B452CD">2</span> + <span style="color: #B452CD">0.1</span>*np.random.randn(N)
|
||||
X = np.c_[np.ones(N),x] <span style="color: #228B22">#Construct design matrix</span>
|
||||
|
||||
<span style="color: #228B22">#Compute theta according to normal equations to compare with GD solution</span>
|
||||
<span style="color: #228B22">#Compute beta according to normal equations to compare with GD solution</span>
|
||||
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)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_NE)
|
||||
beta_NE = np.dot(Xt_X_inv,Xt_y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_NE)
|
||||
</pre></div>
|
||||
</section>
|
||||
|
||||
@@ -605,33 +619,33 @@ theta_NE = np.dot(Xt_X_inv,Xt_y)
|
||||
<h2 id="___sec17">Gradient descent and Ridge </h2>
|
||||
|
||||
<p>
|
||||
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 \),
|
||||
<p> <br>
|
||||
$$
|
||||
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.
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
<p>
|
||||
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
|
||||
<p> <br>
|
||||
$$
|
||||
\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).
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
<p>
|
||||
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
|
||||
<p> <br>
|
||||
$$
|
||||
\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},
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
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 \).
|
||||
|
||||
<p>
|
||||
|
||||
@@ -649,7 +663,7 @@ x = np.random.rand(N)
|
||||
y = <span style="color: #B452CD">5</span>*x**<span style="color: #B452CD">2</span> + <span style="color: #B452CD">0.1</span>*np.random.randn(N)
|
||||
|
||||
|
||||
<span style="color: #228B22">#Compute analytic theta for Ridge regression </span>
|
||||
<span style="color: #228B22">#Compute analytic beta for Ridge regression </span>
|
||||
X = np.c_[np.ones(N),x]
|
||||
XT_X = np.dot(X.T,X)
|
||||
|
||||
@@ -657,10 +671,10 @@ l = <span style="color: #B452CD">0.1</span> <span style="color: #228B22">#Ridge
|
||||
Id = np.eye(XT_X.shape[<span style="color: #B452CD">0</span>])
|
||||
|
||||
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))
|
||||
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_ridge)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(np.linalg.norm(theta_ridge)) <span style="color: #228B22">#||theta||</span>
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_ridge)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(np.linalg.norm(beta_ridge)) <span style="color: #228B22">#||beta||</span>
|
||||
</pre></div>
|
||||
</section>
|
||||
|
||||
@@ -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 \),
|
||||
<p> <br>
|
||||
$$
|
||||
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}).
|
||||
$$
|
||||
<p> <br>
|
||||
</section>
|
||||
@@ -693,8 +707,8 @@ This in turn means that the gradient can be
|
||||
computed as a sum over \( i \)-gradients
|
||||
<p> <br>
|
||||
$$
|
||||
\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}).
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
@@ -724,10 +738,10 @@ all datapoints with a sum over the datapoints in one the minibatches
|
||||
picked at random in each gradient descent step
|
||||
<p> <br>
|
||||
$$
|
||||
\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}).
|
||||
$$
|
||||
<p> <br>
|
||||
</section>
|
||||
@@ -740,8 +754,8 @@ $$
|
||||
Thus a gradient descent step now looks like
|
||||
<p> <br>
|
||||
$$
|
||||
\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})
|
||||
$$
|
||||
<p> <br>
|
||||
|
||||
@@ -772,7 +786,7 @@ j = <span style="color: #B452CD">0</span>
|
||||
<span style="color: #8B008B; font-weight: bold">for</span> i <span style="color: #8B008B">in</span> <span style="color: #658b00">range</span>(m):
|
||||
k = np.random.randint(m) <span style="color: #228B22">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #228B22">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for theta</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for </span>
|
||||
j += <span style="color: #B452CD">1</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
@@ -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.
|
||||
</section>
|
||||
|
||||
@@ -817,10 +831,10 @@ $$\gamma_j(t; t_0, t_1) = \frac{t_0}{t+t_1} $$
|
||||
<p> <br> 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 <em>time</em> \( t \).
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
@@ -844,7 +858,7 @@ j = <span style="color: #B452CD">0</span>
|
||||
<span style="color: #8B008B; font-weight: bold">for</span> i <span style="color: #8B008B">in</span> <span style="color: #658b00">range</span>(m):
|
||||
k = np.random.randint(m) <span style="color: #228B22">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #228B22">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for theta</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for beta</span>
|
||||
t = epoch*m+i
|
||||
gamma_j = step_length(t,t0,t1)
|
||||
j += <span style="color: #B452CD">1</span>
|
||||
|
||||
@@ -156,11 +156,11 @@ MathJax.Hub.Config({
|
||||
|
||||
<p>
|
||||
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 @@ $$
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split -->
|
||||
@@ -211,10 +211,21 @@ the learning rate in the context of Machine Learning.
|
||||
<h2 id="___sec3">The ideal </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
|
||||
<h2 id="___sec5">Gradient Descent Example </h2>
|
||||
|
||||
<p>
|
||||
We revisit now our simple linear regression example with a linear polynomial.
|
||||
<p>
|
||||
|
||||
@@ -259,30 +272,30 @@ x = <span style="color: #B452CD">2</span>*np.random.rand(<span style="color: #B4
|
||||
y = <span style="color: #B452CD">4</span>+<span style="color: #B452CD">3</span>*x+np.random.randn(<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)
|
||||
|
||||
xb = np.c_[np.ones((<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)), x]
|
||||
theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_linreg)
|
||||
theta = np.random.randn(<span style="color: #B452CD">2</span>,<span style="color: #B452CD">1</span>)
|
||||
beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_linreg)
|
||||
beta = np.random.randn(<span style="color: #B452CD">2</span>,<span style="color: #B452CD">1</span>)
|
||||
|
||||
eta = <span style="color: #B452CD">0.1</span>
|
||||
Niterations = <span style="color: #B452CD">1000</span>
|
||||
m = <span style="color: #B452CD">100</span>
|
||||
|
||||
<span style="color: #8B008B; font-weight: bold">for</span> <span style="color: #658b00">iter</span> <span style="color: #8B008B">in</span> <span style="color: #658b00">range</span>(Niterations):
|
||||
gradients = <span style="color: #B452CD">2.0</span>/m*xb.T.dot(xb.dot(theta)-y)
|
||||
theta -= eta*gradients
|
||||
gradients = <span style="color: #B452CD">2.0</span>/m*xb.T.dot(xb.dot(beta)-y)
|
||||
beta -= eta*gradients
|
||||
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta)
|
||||
xnew = np.array([[<span style="color: #B452CD">0</span>],[<span style="color: #B452CD">2</span>]])
|
||||
xbnew = np.c_[np.ones((<span style="color: #B452CD">2</span>,<span style="color: #B452CD">1</span>)), xnew]
|
||||
ypredict = xbnew.dot(theta)
|
||||
ypredict2 = xbnew.dot(theta_linreg)
|
||||
ypredict = xbnew.dot(beta)
|
||||
ypredict2 = xbnew.dot(beta_linreg)
|
||||
plt.plot(xnew, ypredict, <span style="color: #CD5555">"r-"</span>)
|
||||
plt.plot(xnew, ypredict2, <span style="color: #CD5555">"b-"</span>)
|
||||
plt.plot(x, y ,<span style="color: #CD5555">'ro'</span>)
|
||||
plt.axis([<span style="color: #B452CD">0</span>,<span style="color: #B452CD">2.0</span>,<span style="color: #B452CD">0</span>, <span style="color: #B452CD">15.0</span>])
|
||||
plt.xlabel(<span style="color: #CD5555">r'$x$'</span>)
|
||||
plt.ylabel(<span style="color: #CD5555">r'$y$'</span>)
|
||||
plt.title(<span style="color: #CD5555">r'Random numbers '</span>)
|
||||
plt.title(<span style="color: #CD5555">r'Gradient descent example'</span>)
|
||||
plt.show()
|
||||
</pre></div>
|
||||
<p>
|
||||
@@ -303,8 +316,8 @@ x = <span style="color: #B452CD">2</span>*np.random.rand(<span style="color: #B4
|
||||
y = <span style="color: #B452CD">4</span>+<span style="color: #B452CD">3</span>*x+np.random.randn(<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)
|
||||
|
||||
xb = np.c_[np.ones((<span style="color: #B452CD">100</span>,<span style="color: #B452CD">1</span>)), x]
|
||||
theta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_linreg)
|
||||
beta_linreg = np.linalg.inv(xb.T.dot(xb)).dot(xb.T).dot(y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_linreg)
|
||||
sgdreg = SGDRegressor(n_iter = <span style="color: #B452CD">50</span>, penalty=<span style="color: #658b00">None</span>, eta0=<span style="color: #B452CD">0.1</span>)
|
||||
sgdreg.fit(x,y.ravel())
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(sgdreg.intercept_, sgdreg.coef_)
|
||||
@@ -313,6 +326,8 @@ sgdreg.fit(x,y.ravel())
|
||||
<!-- !split -->
|
||||
|
||||
<h2 id="___sec7">Convex functions </h2>
|
||||
|
||||
<p>
|
||||
Ideally we want our cost/loss function to be convex(concave).
|
||||
|
||||
<p>
|
||||
@@ -333,7 +348,7 @@ regular polygons (triangles, rectangles, pentagons, etc...).
|
||||
<h2 id="___sec8">Convex function </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
<b>Convex function</b>: 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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -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, <a href="http://stanford.edu/" target="_blank"><tt>http://stanford.edu/</tt></a>
|
||||
boyd/cvxbook/, 2004.
|
||||
details and proofs we refer to: <a href="http://stanford.edu/boyd/cvxbook/, 2004" target="_blank">S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press</a>.
|
||||
|
||||
<p>
|
||||
<div class="alert alert-block alert-block alert-text-normal">
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
@@ -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.
|
||||
</div>
|
||||
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
@@ -411,27 +429,26 @@ This result means that if we know that the cost/loss function is convex and we a
|
||||
<h2 id="___sec11">Some simple problems </h2>
|
||||
|
||||
<ol>
|
||||
<li> Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.</li>
|
||||
</ol>
|
||||
|
||||
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. $$
|
||||
|
||||
<ol>
|
||||
<li> 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. $</li>
|
||||
<li> Using the second order condition show that the following functions are convex on the specified domain.</li>
|
||||
</ol>
|
||||
|
||||
\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).
|
||||
\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).
|
||||
<ul>
|
||||
<li> \( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).</li>
|
||||
<li> \( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).</li>
|
||||
</ul>
|
||||
|
||||
<ol>
|
||||
<li> 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.</li>
|
||||
<li> A norm is any function that satisfy the following properties</li>
|
||||
|
||||
<ul>
|
||||
<li> \( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \).</li>
|
||||
<li> \( f(x+y) \leq f(x) + f(y) \)</li>
|
||||
<li> \( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \)</li>
|
||||
</ul>
|
||||
|
||||
</ol>
|
||||
|
||||
\( 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).
|
||||
|
||||
<p>
|
||||
<!-- !split -->
|
||||
@@ -439,7 +456,10 @@ Using the definition of convexity, show that a function satisfying the propertie
|
||||
<h2 id="___sec12">Revisiting our first homework </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<ol>
|
||||
<li> An analytical solution (recall homework set 1).</li>
|
||||
@@ -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.
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -469,10 +489,10 @@ $$
|
||||
<h2 id="___sec13">Gradient descent example </h2>
|
||||
|
||||
<p>
|
||||
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 \)
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -497,11 +517,11 @@ and we want to find \( \theta \) such that \( C(\theta) \) is minimized.
|
||||
<h2 id="___sec14">The derivative of the cost/loss function </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
|
||||
<h2 id="___sec15">The Hessian matrix </h2>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -526,19 +546,19 @@ This result implies that \( C(\theta) \) is a convex function since the matrix \
|
||||
<h2 id="___sec16">Simple program </h2>
|
||||
|
||||
<p>
|
||||
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
|
||||
$$
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
<p>
|
||||
|
||||
<!-- code=python (!bc pycod) typeset with pygments style "perldoc" -->
|
||||
@@ -555,11 +575,11 @@ x = np.random.rand(N) <span style="color: #228B22">#Uniformly generated x-value
|
||||
y = <span style="color: #B452CD">5</span>*x**<span style="color: #B452CD">2</span> + <span style="color: #B452CD">0.1</span>*np.random.randn(N)
|
||||
X = np.c_[np.ones(N),x] <span style="color: #228B22">#Construct design matrix</span>
|
||||
|
||||
<span style="color: #228B22">#Compute theta according to normal equations to compare with GD solution</span>
|
||||
<span style="color: #228B22">#Compute beta according to normal equations to compare with GD solution</span>
|
||||
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)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_NE)
|
||||
beta_NE = np.dot(Xt_X_inv,Xt_y)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_NE)
|
||||
</pre></div>
|
||||
<p>
|
||||
<!-- !split -->
|
||||
@@ -567,27 +587,27 @@ theta_NE = np.dot(Xt_X_inv,Xt_y)
|
||||
<h2 id="___sec17">Gradient descent and Ridge </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
$$
|
||||
|
||||
<p>
|
||||
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).
|
||||
$$
|
||||
|
||||
<p>
|
||||
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 \).
|
||||
|
||||
<p>
|
||||
|
||||
@@ -605,7 +625,7 @@ x = np.random.rand(N)
|
||||
y = <span style="color: #B452CD">5</span>*x**<span style="color: #B452CD">2</span> + <span style="color: #B452CD">0.1</span>*np.random.randn(N)
|
||||
|
||||
|
||||
<span style="color: #228B22">#Compute analytic theta for Ridge regression </span>
|
||||
<span style="color: #228B22">#Compute analytic beta for Ridge regression </span>
|
||||
X = np.c_[np.ones(N),x]
|
||||
XT_X = np.dot(X.T,X)
|
||||
|
||||
@@ -613,10 +633,10 @@ l = <span style="color: #B452CD">0.1</span> <span style="color: #228B22">#Ridge
|
||||
Id = np.eye(XT_X.shape[<span style="color: #B452CD">0</span>])
|
||||
|
||||
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))
|
||||
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(theta_ridge)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(np.linalg.norm(theta_ridge)) <span style="color: #228B22">#||theta||</span>
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(beta_ridge)
|
||||
<span style="color: #8B008B; font-weight: bold">print</span>(np.linalg.norm(beta_ridge)) <span style="color: #228B22">#||beta||</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -688,8 +708,8 @@ $$
|
||||
<p>
|
||||
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})
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -719,7 +739,7 @@ j = <span style="color: #B452CD">0</span>
|
||||
<span style="color: #8B008B; font-weight: bold">for</span> i <span style="color: #8B008B">in</span> <span style="color: #658b00">range</span>(m):
|
||||
k = np.random.randint(m) <span style="color: #228B22">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #228B22">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for theta</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for </span>
|
||||
j += <span style="color: #B452CD">1</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
@@ -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 <em>time</em> \( t \).
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
@@ -789,7 +809,7 @@ j = <span style="color: #B452CD">0</span>
|
||||
<span style="color: #8B008B; font-weight: bold">for</span> i <span style="color: #8B008B">in</span> <span style="color: #658b00">range</span>(m):
|
||||
k = np.random.randint(m) <span style="color: #228B22">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #228B22">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for theta</span>
|
||||
<span style="color: #228B22">#Compute new suggestion for beta</span>
|
||||
t = epoch*m+i
|
||||
gamma_j = step_length(t,t0,t1)
|
||||
j += <span style="color: #B452CD">1</span>
|
||||
|
||||
@@ -161,11 +161,11 @@ MathJax.Hub.Config({
|
||||
|
||||
<p>
|
||||
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 @@ $$
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split -->
|
||||
@@ -216,10 +216,21 @@ the learning rate in the context of Machine Learning.
|
||||
<h2 id="___sec3">The ideal </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
|
||||
<h2 id="___sec5">Gradient Descent Example </h2>
|
||||
|
||||
<p>
|
||||
We revisit now our simple linear regression example with a linear polynomial.
|
||||
<p>
|
||||
|
||||
@@ -264,30 +277,30 @@ x <span style="color: #666666">=</span> <span style="color: #666666">2*</span>np
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">4+3*</span>x<span style="color: #666666">+</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)
|
||||
|
||||
xb <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones((<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)), x]
|
||||
theta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_linreg)
|
||||
theta <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">2</span>,<span style="color: #666666">1</span>)
|
||||
beta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_linreg)
|
||||
beta <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">2</span>,<span style="color: #666666">1</span>)
|
||||
|
||||
eta <span style="color: #666666">=</span> <span style="color: #666666">0.1</span>
|
||||
Niterations <span style="color: #666666">=</span> <span style="color: #666666">1000</span>
|
||||
m <span style="color: #666666">=</span> <span style="color: #666666">100</span>
|
||||
|
||||
<span style="color: #008000; font-weight: bold">for</span> <span style="color: #008000">iter</span> <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(Niterations):
|
||||
gradients <span style="color: #666666">=</span> <span style="color: #666666">2.0/</span>m<span style="color: #666666">*</span>xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>dot(theta)<span style="color: #666666">-</span>y)
|
||||
theta <span style="color: #666666">-=</span> eta<span style="color: #666666">*</span>gradients
|
||||
gradients <span style="color: #666666">=</span> <span style="color: #666666">2.0/</span>m<span style="color: #666666">*</span>xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>dot(beta)<span style="color: #666666">-</span>y)
|
||||
beta <span style="color: #666666">-=</span> eta<span style="color: #666666">*</span>gradients
|
||||
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta)
|
||||
xnew <span style="color: #666666">=</span> np<span style="color: #666666">.</span>array([[<span style="color: #666666">0</span>],[<span style="color: #666666">2</span>]])
|
||||
xbnew <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones((<span style="color: #666666">2</span>,<span style="color: #666666">1</span>)), xnew]
|
||||
ypredict <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(theta)
|
||||
ypredict2 <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(theta_linreg)
|
||||
ypredict <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(beta)
|
||||
ypredict2 <span style="color: #666666">=</span> xbnew<span style="color: #666666">.</span>dot(beta_linreg)
|
||||
plt<span style="color: #666666">.</span>plot(xnew, ypredict, <span style="color: #BA2121">"r-"</span>)
|
||||
plt<span style="color: #666666">.</span>plot(xnew, ypredict2, <span style="color: #BA2121">"b-"</span>)
|
||||
plt<span style="color: #666666">.</span>plot(x, y ,<span style="color: #BA2121">'ro'</span>)
|
||||
plt<span style="color: #666666">.</span>axis([<span style="color: #666666">0</span>,<span style="color: #666666">2.0</span>,<span style="color: #666666">0</span>, <span style="color: #666666">15.0</span>])
|
||||
plt<span style="color: #666666">.</span>xlabel(<span style="color: #BA2121">r'$x$'</span>)
|
||||
plt<span style="color: #666666">.</span>ylabel(<span style="color: #BA2121">r'$y$'</span>)
|
||||
plt<span style="color: #666666">.</span>title(<span style="color: #BA2121">r'Random numbers '</span>)
|
||||
plt<span style="color: #666666">.</span>title(<span style="color: #BA2121">r'Gradient descent example'</span>)
|
||||
plt<span style="color: #666666">.</span>show()
|
||||
</pre></div>
|
||||
<p>
|
||||
@@ -308,8 +321,8 @@ x <span style="color: #666666">=</span> <span style="color: #666666">2*</span>np
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">4+3*</span>x<span style="color: #666666">+</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)
|
||||
|
||||
xb <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones((<span style="color: #666666">100</span>,<span style="color: #666666">1</span>)), x]
|
||||
theta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_linreg)
|
||||
beta_linreg <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(xb<span style="color: #666666">.</span>T<span style="color: #666666">.</span>dot(xb))<span style="color: #666666">.</span>dot(xb<span style="color: #666666">.</span>T)<span style="color: #666666">.</span>dot(y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_linreg)
|
||||
sgdreg <span style="color: #666666">=</span> SGDRegressor(n_iter <span style="color: #666666">=</span> <span style="color: #666666">50</span>, penalty<span style="color: #666666">=</span><span style="color: #008000">None</span>, eta0<span style="color: #666666">=0.1</span>)
|
||||
sgdreg<span style="color: #666666">.</span>fit(x,y<span style="color: #666666">.</span>ravel())
|
||||
<span style="color: #008000; font-weight: bold">print</span>(sgdreg<span style="color: #666666">.</span>intercept_, sgdreg<span style="color: #666666">.</span>coef_)
|
||||
@@ -318,6 +331,8 @@ sgdreg<span style="color: #666666">.</span>fit(x,y<span style="color: #666666">.
|
||||
<!-- !split -->
|
||||
|
||||
<h2 id="___sec7">Convex functions </h2>
|
||||
|
||||
<p>
|
||||
Ideally we want our cost/loss function to be convex(concave).
|
||||
|
||||
<p>
|
||||
@@ -338,7 +353,7 @@ regular polygons (triangles, rectangles, pentagons, etc...).
|
||||
<h2 id="___sec8">Convex function </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
<b>Convex function</b>: 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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -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, <a href="http://stanford.edu/" target="_blank"><tt>http://stanford.edu/</tt></a>
|
||||
boyd/cvxbook/, 2004.
|
||||
details and proofs we refer to: <a href="http://stanford.edu/boyd/cvxbook/, 2004" target="_blank">S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press</a>.
|
||||
|
||||
<p>
|
||||
<div class="alert alert-block alert-block alert-text-normal">
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<p>
|
||||
@@ -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.
|
||||
</div>
|
||||
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
@@ -416,27 +434,26 @@ This result means that if we know that the cost/loss function is convex and we a
|
||||
<h2 id="___sec11">Some simple problems </h2>
|
||||
|
||||
<ol>
|
||||
<li> Show that \( f(x)=x^2 \) is convex for \( x \in \mathbb{R} \) using the definition of convexity.</li>
|
||||
</ol>
|
||||
|
||||
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. $$
|
||||
|
||||
<ol>
|
||||
<li> 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. $</li>
|
||||
<li> Using the second order condition show that the following functions are convex on the specified domain.</li>
|
||||
</ol>
|
||||
|
||||
\( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).
|
||||
\( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).
|
||||
<ul>
|
||||
<li> \( f(x) = e^x \) is convex for \( x \in \mathbb{R} \).</li>
|
||||
<li> \( g(x) = -\ln(x) \) is convex for \( x \in (0,\infty) \).</li>
|
||||
</ul>
|
||||
|
||||
<ol>
|
||||
<li> 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.</li>
|
||||
<li> A norm is any function that satisfy the following properties</li>
|
||||
|
||||
<ul>
|
||||
<li> \( f(\alpha x) = |\alpha| f(x) \) for all \( \alpha \in \mathbb{R} \).</li>
|
||||
<li> \( f(x+y) \leq f(x) + f(y) \)</li>
|
||||
<li> \( f(x) \leq 0 \) for all \( x \in \mathbb{R}^n \) with equality if and only if \( x = 0 \)</li>
|
||||
</ul>
|
||||
|
||||
</ol>
|
||||
|
||||
\( 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).
|
||||
|
||||
<p>
|
||||
<!-- !split -->
|
||||
@@ -444,7 +461,10 @@ Using the definition of convexity, show that a function satisfying the propertie
|
||||
<h2 id="___sec12">Revisiting our first homework </h2>
|
||||
|
||||
<p>
|
||||
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:
|
||||
|
||||
<ol>
|
||||
<li> An analytical solution (recall homework set 1).</li>
|
||||
@@ -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.
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -474,10 +494,10 @@ $$
|
||||
<h2 id="___sec13">Gradient descent example </h2>
|
||||
|
||||
<p>
|
||||
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 \)
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -502,11 +522,11 @@ and we want to find \( \theta \) such that \( C(\theta) \) is minimized.
|
||||
<h2 id="___sec14">The derivative of the cost/loss function </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
|
||||
<h2 id="___sec15">The Hessian matrix </h2>
|
||||
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.
|
||||
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -531,19 +551,19 @@ This result implies that \( C(\theta) \) is a convex function since the matrix \
|
||||
<h2 id="___sec16">Simple program </h2>
|
||||
|
||||
<p>
|
||||
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
|
||||
$$
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
|
||||
<p>
|
||||
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} \).
|
||||
<p>
|
||||
|
||||
<!-- code=python (!bc pycod) typeset with pygments style "default" -->
|
||||
@@ -560,11 +580,11 @@ x <span style="color: #666666">=</span> np<span style="color: #666666">.</span>
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">5*</span>x<span style="color: #666666">**2</span> <span style="color: #666666">+</span> <span style="color: #666666">0.1*</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(N)
|
||||
X <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones(N),x] <span style="color: #408080; font-style: italic">#Construct design matrix</span>
|
||||
|
||||
<span style="color: #408080; font-style: italic">#Compute theta according to normal equations to compare with GD solution</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute beta according to normal equations to compare with GD solution</span>
|
||||
Xt_X_inv <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,X))
|
||||
Xt_y <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>transpose(),y)
|
||||
theta_NE <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Xt_X_inv,Xt_y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_NE)
|
||||
beta_NE <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Xt_X_inv,Xt_y)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_NE)
|
||||
</pre></div>
|
||||
<p>
|
||||
<!-- !split -->
|
||||
@@ -572,27 +592,27 @@ theta_NE <span style="color: #666666">=</span> np<span style="color: #666666">.<
|
||||
<h2 id="___sec17">Gradient descent and Ridge </h2>
|
||||
|
||||
<p>
|
||||
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.
|
||||
$$
|
||||
|
||||
<p>
|
||||
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).
|
||||
$$
|
||||
|
||||
<p>
|
||||
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 \).
|
||||
|
||||
<p>
|
||||
|
||||
@@ -610,7 +630,7 @@ x <span style="color: #666666">=</span> np<span style="color: #666666">.</span>
|
||||
y <span style="color: #666666">=</span> <span style="color: #666666">5*</span>x<span style="color: #666666">**2</span> <span style="color: #666666">+</span> <span style="color: #666666">0.1*</span>np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randn(N)
|
||||
|
||||
|
||||
<span style="color: #408080; font-style: italic">#Compute analytic theta for Ridge regression </span>
|
||||
<span style="color: #408080; font-style: italic">#Compute analytic beta for Ridge regression </span>
|
||||
X <span style="color: #666666">=</span> np<span style="color: #666666">.</span>c_[np<span style="color: #666666">.</span>ones(N),x]
|
||||
XT_X <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,X)
|
||||
|
||||
@@ -618,10 +638,10 @@ l <span style="color: #666666">=</span> <span style="color: #666666">0.1</span>
|
||||
Id <span style="color: #666666">=</span> np<span style="color: #666666">.</span>eye(XT_X<span style="color: #666666">.</span>shape[<span style="color: #666666">0</span>])
|
||||
|
||||
Z <span style="color: #666666">=</span> np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>inv(XT_X<span style="color: #666666">+</span>l<span style="color: #666666">*</span>Id)
|
||||
theta_ridge <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Z,np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,y))
|
||||
beta_ridge <span style="color: #666666">=</span> np<span style="color: #666666">.</span>dot(Z,np<span style="color: #666666">.</span>dot(X<span style="color: #666666">.</span>T,y))
|
||||
|
||||
<span style="color: #008000; font-weight: bold">print</span>(theta_ridge)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>norm(theta_ridge)) <span style="color: #408080; font-style: italic">#||theta||</span>
|
||||
<span style="color: #008000; font-weight: bold">print</span>(beta_ridge)
|
||||
<span style="color: #008000; font-weight: bold">print</span>(np<span style="color: #666666">.</span>linalg<span style="color: #666666">.</span>norm(beta_ridge)) <span style="color: #408080; font-style: italic">#||beta||</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
<!-- !split --><br><br><br><br><br><br><br><br><br><br>
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -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}).
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -693,8 +713,8 @@ $$
|
||||
<p>
|
||||
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})
|
||||
$$
|
||||
|
||||
<p>
|
||||
@@ -724,7 +744,7 @@ j <span style="color: #666666">=</span> <span style="color: #666666">0</span>
|
||||
<span style="color: #008000; font-weight: bold">for</span> i <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(m):
|
||||
k <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randint(m) <span style="color: #408080; font-style: italic">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for theta</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for </span>
|
||||
j <span style="color: #666666">+=</span> <span style="color: #666666">1</span>
|
||||
</pre></div>
|
||||
<p>
|
||||
@@ -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.
|
||||
|
||||
<p>
|
||||
@@ -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 <em>time</em> \( t \).
|
||||
|
||||
<p>
|
||||
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.
|
||||
|
||||
<p>
|
||||
@@ -794,7 +814,7 @@ j <span style="color: #666666">=</span> <span style="color: #666666">0</span>
|
||||
<span style="color: #008000; font-weight: bold">for</span> i <span style="color: #AA22FF; font-weight: bold">in</span> <span style="color: #008000">range</span>(m):
|
||||
k <span style="color: #666666">=</span> np<span style="color: #666666">.</span>random<span style="color: #666666">.</span>randint(m) <span style="color: #408080; font-style: italic">#Pick the k-th minibatch at random</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute the gradient using the data in minibatch Bk</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for theta</span>
|
||||
<span style="color: #408080; font-style: italic">#Compute new suggestion for beta</span>
|
||||
t <span style="color: #666666">=</span> epoch<span style="color: #666666">*</span>m<span style="color: #666666">+</span>i
|
||||
gamma_j <span style="color: #666666">=</span> step_length(t,t0,t1)
|
||||
j <span style="color: #666666">+=</span> <span style="color: #666666">1</span>
|
||||
|
||||
@@ -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",
|
||||
"<!-- !split -->\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": [
|
||||
"<!-- !split -->\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, <http://stanford.edu/>\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",
|
||||
"<!-- !split -->\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 @@
|
||||
"<!-- !split -->\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 @@
|
||||
"<!-- !split -->\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",
|
||||
|
||||
Binary file not shown.
Binary file not shown.
+108
-92
@@ -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
|
||||
|
||||
Reference in New Issue
Block a user