<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/classification-from-scratch-boosting-11-8/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://freakonometrics.hypotheses.org/52782"> R-english – Freakonometrics</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p>Eleventh post of our series on <a rel="nofollow" href="https://freakonometrics.hypotheses.org/52731" target="_blank">classification from scratch</a>. Today, that should be the last one… unless I forgot something important. So today, we discuss <a rel="nofollow" href="https://en.wikipedia.org/wiki/Boosting_(machine_learning)" target="_blank">boosting</a>.</p>
<h2>An econometrician perspective</h2>
<p>I might start with a non-conventional introduction. But that’s actually how I understood what boosting was about. And I am quite sure it has to do with my background in econometrics.</p>
<p>The goal here is to solve something which looks like<img alt="m^\star=\underset{m\in\mathcal{M}}{\text{argmin}}\left\lbrace\sum_{i=1}^n \ell(y_i,m(\mathbf{x}_i))\right\rbrace" title="m^\star=\underset{m\in\mathcal{M}}{\text{argmin}}\left\lbrace\sum_{i=1}^n \ell(y_i,m(\mathbf{x}_i))\right\rbrace" src="https://s0.wp.com/latex.php?latex=m%5E%5Cstar%3D%5Cunderset%7Bm%5Cin%5Cmathcal%7BM%7D%7D%7B%5Ctext%7Bargmin%7D%7D%5Cleft%5Clbrace%5Csum_%7Bi%3D1%7D%5En+%5Cell%28y_i%2Cm%28%5Cmathbf%7Bx%7D_i%29%29%5Cright%5Crbrace&bg=ffffff&fg=000&s=0" class="latex">for some loss function <img alt="\ell" title="\ell" src="https://s0.wp.com/latex.php?latex=%5Cell&bg=ffffff&fg=000&s=0" class="latex">, and for some set of predictors <img alt="\mathcal{M}" title="\mathcal{M}" src="https://s0.wp.com/latex.php?latex=%5Cmathcal%7BM%7D&bg=ffffff&fg=000&s=0" class="latex">. This is an optimization problem. Well, optimization is here in a function space, but still, that’s simply an optimization problem. And from a numerical perspective, optimization is solve using gradient descent (this is why this technique is also called <a rel="nofollow" href="https://en.wikipedia.org/wiki/Gradient_boosting" target="_blank">gradient boosting</a>). And the gradient descent can be visualized like below</p>
<p><img src="https://i1.wp.com/freakonometrics.hypotheses.org/files/2018/06/gradient-desc.png?resize=274%2C243&ssl=1" alt="" width="274" data-recalc-dims="1" height="243"></p>
<p>Again, the optimum is not some some real value <img alt="x^\star" title="x^\star" src="https://s0.wp.com/latex.php?latex=x%5E%5Cstar&bg=ffffff&fg=000&s=0" class="latex">, but some function <img alt="m^\star" title="m^\star" src="https://s0.wp.com/latex.php?latex=m%5E%5Cstar&bg=ffffff&fg=000&s=0" class="latex">. Thus, here we will have something like<img alt="m^{(k)}=m^{(k-1)}+\underset{h\in\mathcal{H}}{\text{argmin}}\left\lbrace \sum_{i=1}^n \ell(y_i,m^{(k-1)}(\mathbf{x}_i)+h(\mathbf{x}_i))\right\rbrace" title="m^{(k)}=m^{(k-1)}+\underset{h\in\mathcal{H}}{\text{argmin}}\left\lbrace \sum_{i=1}^n \ell(y_i,m^{(k-1)}(\mathbf{x}_i)+h(\mathbf{x}_i))\right\rbrace" src="https://s0.wp.com/latex.php?latex=m%5E%7B%28k%29%7D%3Dm%5E%7B%28k-1%29%7D%2B%5Cunderset%7Bh%5Cin%5Cmathcal%7BH%7D%7D%7B%5Ctext%7Bargmin%7D%7D%5Cleft%5Clbrace+%5Csum_%7Bi%3D1%7D%5En+%5Cell%28y_i%2Cm%5E%7B%28k-1%29%7D%28%5Cmathbf%7Bx%7D_i%29%2Bh%28%5Cmathbf%7Bx%7D_i%29%29%5Cright%5Crbrace&bg=ffffff&fg=000&s=0" class="latex">(as they write it is serious articles) where the term on the right can also be written<img alt="m^{(k)}=m^{(k-1)}+\underset{h\in\mathcal{H}}{\text{argmin}}\left\lbrace \sum_{i=1}^n \ell(\underbrace{y_i-m^{(k-1)}(\mathbf{x}_i)}_{\varepsilon_{k,i}},h(\mathbf{x}_i))\right\rbrace" title="m^{(k)}=m^{(k-1)}+\underset{h\in\mathcal{H}}{\text{argmin}}\left\lbrace \sum_{i=1}^n \ell(\underbrace{y_i-m^{(k-1)}(\mathbf{x}_i)}_{\varepsilon_{k,i}},h(\mathbf{x}_i))\right\rbrace" src="https://s0.wp.com/latex.php?latex=m%5E%7B%28k%29%7D%3Dm%5E%7B%28k-1%29%7D%2B%5Cunderset%7Bh%5Cin%5Cmathcal%7BH%7D%7D%7B%5Ctext%7Bargmin%7D%7D%5Cleft%5Clbrace+%5Csum_%7Bi%3D1%7D%5En+%5Cell%28%5Cunderbrace%7By_i-m%5E%7B%28k-1%29%7D%28%5Cmathbf%7Bx%7D_i%29%7D_%7B%5Cvarepsilon_%7Bk%2Ci%7D%7D%2Ch%28%5Cmathbf%7Bx%7D_i%29%29%5Cright%5Crbrace&bg=ffffff&fg=000&s=0" class="latex">I prefer the later, because we see clearly that <img alt="f" title="f" src="https://s0.wp.com/latex.php?latex=f&bg=ffffff&fg=000&s=0" class="latex"> is some model we fit on the remaining residuals.</p>
<p>We can rewrite it like that: define<img alt="r_{i,k}=-\left.\frac{\partial \ell(y_i,m(\mathbf{x}_i))}{\partial m(\mathbf{x}_i)}\right\vert_{m(\mathbf{x}_i)=m^{(k-1)}(\mathbf{x}_i)}" title="r_{i,k}=-\left.\frac{\partial \ell(y_i,m(\mathbf{x}_i))}{\partial m(\mathbf{x}_i)}\right\vert_{m(\mathbf{x}_i)=m^{(k-1)}(\mathbf{x}_i)}" src="https://s0.wp.com/latex.php?latex=r_%7Bi%2Ck%7D%3D-%5Cleft.%5Cfrac%7B%5Cpartial+%5Cell%28y_i%2Cm%28%5Cmathbf%7Bx%7D_i%29%29%7D%7B%5Cpartial+m%28%5Cmathbf%7Bx%7D_i%29%7D%5Cright%5Cvert_%7Bm%28%5Cmathbf%7Bx%7D_i%29%3Dm%5E%7B%28k-1%29%7D%28%5Cmathbf%7Bx%7D_i%29%7D&bg=ffffff&fg=000&s=0" class="latex">for all <img alt="i=1,\cdots,n" title="i=1,\cdots,n" src="https://s0.wp.com/latex.php?latex=i%3D1%2C%5Ccdots%2Cn&bg=ffffff&fg=000&s=0" class="latex">. The goal is to fit a model so that <img alt="r_{i,k}=h^\star(\mathbf{x}_i)" title="r_{i,k}=h^\star(\mathbf{x}_i)" src="https://s0.wp.com/latex.php?latex=r_%7Bi%2Ck%7D%3Dh%5E%5Cstar%28%5Cmathbf%7Bx%7D_i%29&bg=ffffff&fg=000&s=0" class="latex">, and when we have that optimal function, set <img alt="m_k(\mathbf{x})=m_{k-1}(\mathbf{x})+\gamma_k h^\star(\mathbf{x})" title="m_k(\mathbf{x})=m_{k-1}(\mathbf{x})+\gamma_k h^\star(\mathbf{x})" src="https://s0.wp.com/latex.php?latex=m_k%28%5Cmathbf%7Bx%7D%29%3Dm_%7Bk-1%7D%28%5Cmathbf%7Bx%7D%29%2B%5Cgamma_k+h%5E%5Cstar%28%5Cmathbf%7Bx%7D%29&bg=ffffff&fg=000&s=0" class="latex"> (yes, we can include some shrinkage here). </p>
<p>Two important comments here. First of all, the idea should be weird to any econometrician. First, we fit a model to explain <img alt="y" title="y" src="https://s0.wp.com/latex.php?latex=y&bg=ffffff&fg=000&s=0" class="latex"> by some covariates <img alt="\mathbf{x}" title="\mathbf{x}" src="https://s0.wp.com/latex.php?latex=%5Cmathbf%7Bx%7D&bg=ffffff&fg=000&s=0" class="latex">. Then consider the residuals <img alt="\widehat{\varepsilon}" title="\widehat{\varepsilon}" src="https://s0.wp.com/latex.php?latex=%5Cwidehat%7B%5Cvarepsilon%7D&bg=ffffff&fg=000&s=0" class="latex">, and to explain them with the same covariate <img alt="\mathbf{x}" title="\mathbf{x}" src="https://s0.wp.com/latex.php?latex=%5Cmathbf%7Bx%7D&bg=ffffff&fg=000&s=0" class="latex">. If you try that with a linear regression, you’d done at the end of step 1, since residuals <img alt="\widehat{\varepsilon}" title="\widehat{\varepsilon}" src="https://s0.wp.com/latex.php?latex=%5Cwidehat%7B%5Cvarepsilon%7D&bg=ffffff&fg=000&s=0" class="latex"> are orthogonal to covariates <img alt="\mathbf{x}" title="\mathbf{x}" src="https://s0.wp.com/latex.php?latex=%5Cmathbf%7Bx%7D&bg=ffffff&fg=000&s=0" class="latex">: no way that we can learn from them. Here it works because we consider simple non linear model. And actually, something that can be used is to add a shrinkage parameter. Do not consider <img alt="\widehat{\varepsilon}=y-\widehat{m}(\mathbf{x})" title="\widehat{\varepsilon}=y-\widehat{m}(\mathbf{x})" src="https://s0.wp.com/latex.php?latex=%5Cwidehat%7B%5Cvarepsilon%7D%3Dy-%5Cwidehat%7Bm%7D%28%5Cmathbf%7Bx%7D%29&bg=ffffff&fg=000&s=0" class="latex"> but <img alt="\widehat{\varepsilon}=y-\gamma\widehat{m}(\mathbf{x})" title="\widehat{\varepsilon}=y-\gamma\widehat{m}(\mathbf{x})" src="https://s0.wp.com/latex.php?latex=%5Cwidehat%7B%5Cvarepsilon%7D%3Dy-%5Cgamma%5Cwidehat%7Bm%7D%28%5Cmathbf%7Bx%7D%29&bg=ffffff&fg=000&s=0" class="latex">. The idea of <em>weak</em> learners is extremely important here. The more we shrink, the longer it will take, but that’s not (too) important. </p>
<p>I should also mention that it’s nice to keep learning from our mistakes. But somehow, we should stop, someday. I said that I will not mention this part in this series of posts, maybe later on. But heuristically, we should stop when we start to overfit. And this can be observed either using a split training/validation of the initial dataset or to use cross validation. I will get back on that issue later one in this post, but again, those ideas should probably be dedicated to another series of posts.</p>
<h2>Learning with splines</h2>
<p>Just to make sure we get it, let’s try to learn with splines. Because <a rel="nofollow" href="https://freakonometrics.hypotheses.org/52771" target="_blank">standard splines</a> have fixed knots, actually, we do not really “learn” here (and after a few iterations we get to what we would have with a standard spline regression). So here, we will (somehow) optimize knots locations. There is a package to do so. And just to illustrate, use a Gaussian regression here, not a classification (we will do that later on). Consider the following dataset (with only one covariate)</p>
<div><table><tr><td><pre>1
2
3
4
5
</pre></td><td><pre>n<span>=</span><span>300</span>
 <span>set.<span>seed</span></span><span>(</span><span>1</span><span>)</span>
 u<span>=</span><span>sort</span><span>(</span><span>runif</span><span>(</span>n<span>)</span><span>*</span><span>2</span><span>*</span><span>pi</span><span>)</span>
 y<span>=</span><span>sin</span><span>(</span>u<span>)</span><span>+</span><span>rnorm</span><span>(</span>n<span>)</span><span>/</span><span>4</span>
 <span>df</span><span>=</span><span>data.<span>frame</span></span><span>(</span>x<span>=</span>u,y<span>=</span>y<span>)</span></pre></td></tr></table></div>
<p>For an optimal choice of knot locations, we can use</p>
<div><table><tr><td><pre>1
2
</pre></td><td><pre><span>library</span><span>(</span>freeknotsplines<span>)</span>
xy.<span>freekt</span><span>=</span>freelsgen<span>(</span><span>df</span>$x, <span>df</span>$y, degree <span>=</span> <span>1</span>, numknot <span>=</span> <span>2</span>, <span>555</span><span>)</span></pre></td></tr></table></div>
<p>With 5% shrinkage, the code it simply the following</p>
<div><table><tr><td><pre>1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
</pre></td><td><pre>v<span>=</span>.05
 <span>library</span><span>(</span>splines<span>)</span>
 xy.<span>freekt</span><span>=</span>freelsgen<span>(</span><span>df</span>$x, <span>df</span>$y, degree <span>=</span> <span>1</span>, numknot <span>=</span> <span>2</span>, <span>555</span><span>)</span>
 fit<span>=</span><span>lm</span><span>(</span>y~bs<span>(</span>x,degree<span>=</span><span>1</span>,<span>knots</span><span>=</span>xy.<span>freekt</span>@optknot<span>)</span>,<span>data</span><span>=</span><span>df</span><span>)</span>
 yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span><span>df</span><span>)</span>
 <span>df</span>$yr<span>=</span><span>df</span>$y <span>-</span> v<span>*</span>yp
 YP<span>=</span>v<span>*</span>yp
 <span>for</span><span>(</span><span>t</span> <span>in</span> <span>1</span><span>:</span><span>200</span><span>)</span><span>{</span>
   xy.<span>freekt</span><span>=</span>freelsgen<span>(</span><span>df</span>$x, <span>df</span>$yr, degree <span>=</span> <span>1</span>, numknot <span>=</span> <span>2</span>, <span>555</span><span>)</span>
   fit<span>=</span><span>lm</span><span>(</span>yr~bs<span>(</span>x,degree<span>=</span><span>1</span>,<span>knots</span><span>=</span>xy.<span>freekt</span>@optknot<span>)</span>,<span>data</span><span>=</span><span>df</span><span>)</span>
   yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span><span>df</span><span>)</span>
   <span>df</span>$yr<span>=</span><span>df</span>$yr <span>-</span> v<span>*</span>yp
   YP<span>=</span><span>cbind</span><span>(</span>YP,v<span>*</span>yp<span>)</span><span>}</span>
 nd<span>=</span><span>data.<span>frame</span></span><span>(</span>x<span>=</span><span>seq</span><span>(</span><span>0</span>,<span>2</span><span>*</span><span>pi</span>,<span>by</span><span>=</span>.01<span>)</span><span>)</span>
 viz<span>=</span><span>function</span><span>(</span>M<span>)</span><span>{</span>
    <span>if</span><span>(</span>M<span>==</span><span>1</span><span>)</span>  y<span>=</span>YP<span>[</span>,<span>1</span><span>]</span>
    <span>if</span><span>(</span>M<span>&amp;</span>gt<span>;</span><span>1</span><span>)</span>   y<span>=</span><span>apply</span><span>(</span>YP<span>[</span>,<span>1</span><span>:</span>M<span>]</span>,<span>1</span>,<span>sum</span><span>)</span>
    <span>plot</span><span>(</span><span>df</span>$x,<span>df</span>$y,ylab<span>=</span><span>&quot;&quot;</span>,xlab<span>=</span><span>&quot;&quot;</span><span>)</span>
    <span>lines</span><span>(</span><span>df</span>$x,y,type<span>=</span><span>&quot;l&quot;</span>,<span>col</span><span>=</span><span>&quot;red&quot;</span>,lwd<span>=</span><span>3</span><span>)</span>
    fit<span>=</span><span>lm</span><span>(</span>y~bs<span>(</span>x,degree<span>=</span><span>1</span>,<span>df</span><span>=</span><span>3</span><span>)</span>,<span>data</span><span>=</span><span>df</span><span>)</span>
    yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span>nd<span>)</span>
    <span>lines</span><span>(</span>nd$x,yp,type<span>=</span><span>&quot;l&quot;</span>,<span>col</span><span>=</span><span>&quot;blue&quot;</span>,lwd<span>=</span><span>3</span><span>)</span>
    <span>lines</span><span>(</span>nd$x,<span>sin</span><span>(</span>nd$x<span>)</span>,lty<span>=</span><span>2</span><span>)</span><span>}</span></pre></td></tr></table></div>
<p>To visualize the ouput after 100 iterations, use</p>
<div><table><tr><td><pre>1
</pre></td><td><pre>viz<span>(</span><span>100</span><span>)</span></pre></td></tr></table></div>
<p><img alt="" src="https://i2.wp.com/freakonometrics.hypotheses.org/files/2015/07/boosting-algo-0_v2.gif?w=456&ssl=1" data-recalc-dims="1"><br>
Clearly, we see that we learn from the data here… Cool, isn’t it?</p>
<h2>Learning with stumps (and trees)</h2>
<p>Let us try something else. What if we consider at each step a regression tree, instead of a linear-by-parts regression (that was considered with linear splines).</p>
<div><table><tr><td><pre>1
2
3
4
5
6
7
8
9
10
11
</pre></td><td><pre><span>library</span><span>(</span>rpart<span>)</span>
v<span>=</span>.1 
fit<span>=</span>rpart<span>(</span>y~x,<span>data</span><span>=</span><span>df</span><span>)</span>
yp<span>=</span><span>predict</span><span>(</span>fit<span>)</span>
<span>df</span>$yr<span>=</span><span>df</span>$y <span>-</span> v<span>*</span>yp
YP<span>=</span>v<span>*</span>yp
<span>for</span><span>(</span><span>t</span> <span>in</span> <span>1</span><span>:</span><span>100</span><span>)</span><span>{</span>
  fit<span>=</span>rpart<span>(</span>yr~x,<span>data</span><span>=</span><span>df</span><span>)</span>
  yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span><span>df</span><span>)</span>
  <span>df</span>$yr<span>=</span><span>df</span>$yr <span>-</span> v<span>*</span>yp
  YP<span>=</span><span>cbind</span><span>(</span>YP,v<span>*</span>yp<span>)</span><span>}</span></pre></td></tr></table></div>
<p>Again, to visualise the learning process, use</p>
<div><table><tr><td><pre>1
2
3
4
5
6
7
8
</pre></td><td><pre>viz<span>=</span><span>function</span><span>(</span>M<span>)</span><span>{</span>
y<span>=</span><span>apply</span><span>(</span>YP<span>[</span>,<span>1</span><span>:</span>M<span>]</span>,<span>1</span>,<span>sum</span><span>)</span>
<span>plot</span><span>(</span><span>df</span>$x,<span>df</span>$y,ylab<span>=</span><span>&quot;&quot;</span>,xlab<span>=</span><span>&quot;&quot;</span><span>)</span>
<span>lines</span><span>(</span><span>df</span>$x,y,type<span>=</span><span>&quot;s&quot;</span>,<span>col</span><span>=</span><span>&quot;red&quot;</span>,lwd<span>=</span><span>3</span><span>)</span>
fit<span>=</span>rpart<span>(</span>y~x,<span>data</span><span>=</span><span>df</span><span>)</span>
yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span>nd<span>)</span>
<span>lines</span><span>(</span>nd$x,yp,type<span>=</span><span>&quot;s&quot;</span>,<span>col</span><span>=</span><span>&quot;blue&quot;</span>,lwd<span>=</span><span>3</span><span>)</span>
<span>lines</span><span>(</span>nd$x,<span>sin</span><span>(</span>nd$x<span>)</span>,lty<span>=</span><span>2</span><span>)</span><span>}</span></pre></td></tr></table></div>
<p><img alt="" src="https://i0.wp.com/freakonometrics.hypotheses.org/files/2015/06/boosting-algo-3.gif?w=456&ssl=1" data-recalc-dims="1"><br>
This time, with those trees, it looks like not only we have a good model, but also a different model from the one we can get using a single regression tree.</p>
<p>What if we change the shrinkage parameter?</p>
<div><table><tr><td><pre>1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
</pre></td><td><pre>viz<span>=</span><span>function</span><span>(</span>v<span>=</span><span>0.05</span><span>)</span><span>{</span>
  fit<span>=</span>rpart<span>(</span>y~x,<span>data</span><span>=</span><span>df</span><span>)</span>
  yp<span>=</span><span>predict</span><span>(</span>fit<span>)</span>
  <span>df</span>$yr<span>=</span><span>df</span>$y <span>-</span> v<span>*</span>yp
  YP<span>=</span>v<span>*</span>yp
  <span>for</span><span>(</span><span>t</span> <span>in</span> <span>1</span><span>:</span><span>100</span><span>)</span><span>{</span>
    fit<span>=</span>rpart<span>(</span>yr~x,<span>data</span><span>=</span><span>df</span><span>)</span>
    yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span><span>df</span><span>)</span>
    <span>df</span>$yr<span>=</span><span>df</span>$yr <span>-</span> v<span>*</span>yp
    YP<span>=</span><span>cbind</span><span>(</span>YP,v<span>*</span>yp<span>)</span><span>}</span>
  y<span>=</span><span>apply</span><span>(</span>YP,<span>1</span>,<span>sum</span><span>)</span>
    <span>plot</span><span>(</span><span>df</span>$x,<span>df</span>$y,xlab<span>=</span><span>&quot;&quot;</span>,ylab<span>=</span><span>&quot;&quot;</span><span>)</span>
    <span>lines</span><span>(</span><span>df</span>$x,y,type<span>=</span><span>&quot;s&quot;</span>,<span>col</span><span>=</span><span>&quot;red&quot;</span>,lwd<span>=</span><span>3</span><span>)</span>
    fit<span>=</span>rpart<span>(</span>y~x,<span>data</span><span>=</span><span>df</span><span>)</span>
    yp<span>=</span><span>predict</span><span>(</span>fit,newdata<span>=</span>nd<span>)</span>
    <span>lines</span><span>(</span>nd$x,yp,type<span>=</span><span>&quot;s&quot;</span>,<span>col</span><span>=</span><span>&quot;blue&quot;</span>,lwd<span>=</span><span>3</span><span>)</span>
    <span>lines</span><span>(</span>nd$x,<span>sin</span><span>(</span>nd$x<span>)</span>,lty<span>=</span><span>2</span><span>)</span><span>}</span></pre></td></tr></table></div>
<p><img alt="" src="https://i1.wp.com/freakonometrics.hypotheses.org/files/2015/06/boosting-algo-4.gif?w=456&ssl=1" data-recalc-dims="1"><br>
There is clearly an impact of that shrinkage parameter. It has to be small to get a good model. This is the idea of using <em>weak learners</em> to get a good prediction.</p>
<h2>Classification and Adaboost</h2>
<p>Now that we understand how bootsting works, let’s try to adapt it to classification. It will be more complicated because residuals are usually not very informative in a classification. And it will be hard to shrink. So let’s try something slightly different, to introduce the <a rel="nofollow" href="https://en.wikipedia.org/wiki/AdaBoost" target="_blank">adaboost</a> algorithm.</p>
<p>In our initial discussion, the goal was to minimize a convex loss function. Here, if we express classes as <img alt="\{-1,+1\}" title="\{-1,+1\}" src="https://s0.wp.com/latex.php?latex=%5C%7B-1%2C%2B1%5C%7D&bg=ffffff&fg=000&s=0" class="latex">, the loss function we consider is <img alt="e^{-y\cdot m(\mathbf{x})}" title="e^{-y\cdot m(\mathbf{x})}" src="https://s0.wp.com/latex.php?latex=e%5E%7B-y%5Ccdot+m%28%5Cmathbf%7Bx%7D%29%7D&bg=ffffff&fg=000&s=0" class="latex"> (this product <img alt="y\cdot m(\mathbf{x})" title="y\cdot m(\mathbf{x})" src="https://s0.wp.com/latex.php?latex=y%5Ccdot+m%28%5Cmathbf%7Bx%7D%29&bg=ffffff&fg=000&s=0" class="latex">) was already discussed when we’ve seen the SVM algorithm. Note that the loss function related to the logistic model would be <img alt="\log(1+e^{-y\cdot m(\mathbf{x})})" title="\log(1+e^{-y\cdot m(\mathbf{x})})" src="https://s0.wp.com/latex.php?latex=%5Clog%281%2Be%5E%7B-y%5Ccdot+m%28%5Cmathbf%7Bx%7D%29%7D%29&bg=ffffff&fg=000&s=0" class="latex">.</p>
<p>What we do here is related to gradient descent (or Newton algorithm). Previously, we were learning from our errors. At each iteration, the residuals are computed and a (weak) model is fitted to these residuals. The the contribution of this weak model is used in a gradient descent optimization process. Here things will be different, because (from my understanding) it is more difficult to play with residuals, because null residuals never exist in classifications. So we will add weights. Initially, all the observations will have the same weights. But iteratively, we ill change them. We will increase the weights of the wrongly predicted individuals and decrease the ones of the correctly predicted individuals. Somehow, we want to focus more on the difficult predictions. That’s the trick. And I guess that’s why it performs so well. This algorithm is well described in <a rel="nofollow" href="https://en.wikipedia.org/wiki/AdaBoost#Example_algorithm_(Discrete_AdaBoost)" target="_blank">wikipedia</a>, so we will use it.</p>
<p>We start with <img alt="\mathbf{\omega}_0=\mathbf{1}/n" title="\mathbf{\omega}_0=\mathbf{1}/n" src="https://s0.wp.com/latex.php?latex=%5Cmathbf%7B%5Comega%7D_0%3D%5Cmathbf%7B1%7D%2Fn&bg=ffffff&fg=000&s=0" class="latex">, then at each step fit a model (a classification tree) with weights <img alt="\mathbf{\omega}_k" title="\mathbf{\omega}_k" src="https://s0.wp.com/latex.php?latex=%5Cmathbf%7B%5Comega%7D_k&bg=ffffff&fg=000&s=0" class="latex">(we did not discuss weights in the algorithms of trees, but it is straigtforward in the formula actually). Let <img alt="\widehat{h}_{\mathbf{\omega}_k}" title="\widehat{h}_{\mathbf{\omega}_k}" src="https://s0.wp.com/latex.php?latex=%5Cwidehat%7Bh%7D_%7B%5Cmathbf%7B%5Comega%7D_k%7D&bg=ffffff&fg=000&s=0" class="latex"> denote that model (i.e. the probability in each leaves). Then consider the classifier <img alt="2~\mathbf{1}[\widehat{h}_{\mathbf{\omega}_k}(\cdot)&gt;0.5]-1" title="2~\mathbf{1}[\widehat{h}_{\mathbf{\omega}_k}(\cdot)&gt;0.5]-1" src="https://s0.wp.com/latex.php?latex=2%7E%5Cmathbf%7B1%7D%5B%5Cwidehat%7Bh%7D_%7B%5Cmathbf%7B%5Comega%7D_k%7D%28%5Ccdot%29%3E0.5%5D-1&bg=ffffff&fg=000&s=0" class="latex"> which returns a value in <img alt="\{-1,+1\}" title="\{-1,+1\}" src="https://s0.wp.com/latex.php?latex=%5C%7B-1%2C%2B1%5C%7D&bg=ffffff&fg=000&s=0" class="latex">. Then set <img alt="\varepsilon_k=\sum_{i\in\mathcal{I}_k}\omega_i " title="\varepsilon_k=\sum_{i\in\mathcal{I}_k}\omega_i " src="https://s0.wp.com/latex.php?latex=%5Cvarepsilon_k%3D%5Csum_%7Bi%5Cin%5Cmathcal%7BI%7D_k%7D%5Comega_i+&bg=ffffff&fg=000&s=0" class="latex">where <img alt="\mathcal{I}_k" title="\mathcal{I}_k" src="https://s0.wp.com/latex.php?latex=%5Cmathcal%7BI%7D_k&bg=ffffff&fg=000&s=0" class="latex"> is the set of misclassified individuals,<img alt="\mathcal{I}_k=\big\lbrace i:2~\mathbf{1}[\widehat{h}_{\mathbf{\omega}_k}(\mathbf{x}_i)&gt;0.5]-1\neq y_i\big\rbrace " title="\mathcal{I}_k=\big\lbrace i:2~\mathbf{1}[\widehat{h}_{\mathbf{\omega}_k}(\mathbf{x}_i)&gt;0.5]-1\neq y_i\big\rbrace " src="https://s0.wp.com/latex.php?latex=%5Cmathcal%7BI%7D_k%3D%5Cbig%5Clbrace+i%3A2%7E%5Cmathbf%7B1%7D%5B%5Cwidehat%7Bh%7D_%7B%5Cmathbf%7B%5Comega%7D_k%7D%28%5Cmathbf%7Bx%7D_i%29%3E0.5%5D-1%5Cneq+y_i%5Cbig%5Crbrace+&bg=ffffff&fg=000&s=0" class="latex">Then set <img alt=" \alpha_k = \frac{1}{2} \ln \left(\frac{1-\epsilon_k}{\epsilon_k}\right)" title=" \alpha_k = \frac{1}{2} \ln \left(\frac{1-\epsilon_k}{\epsilon_k}\right)" src="https://s0.wp.com/latex.php?latex=+%5Calpha_k+%3D+%5Cfrac%7B1%7D%7B2%7D+%5Cln+%5Cleft%28%5Cfrac%7B1-%5Cepsilon_k%7D%7B%5Cepsilon_k%7D%5Cright%29&bg=ffffff&fg=000&s=0" class="latex">and update finally the model using<img alt="m_{k=1}=m_k+\alpha_k\widehat{h}_{\mathbf{\omega}_k}" title="m_{k=1}=m_k+\alpha_k\widehat{h}_{\mathbf{\omega}_k}" src="https://s0.wp.com/latex.php?latex=m_%7Bk%3D1%7D%3Dm_k%2B%5Calpha_k%5Cwidehat%7Bh%7D_%7B%5Cmathbf%7B%5Comega%7D_k%7D&bg=ffffff&fg=000&s=0" class="latex">as well as the weights<img alt="\mathbf{\omega}_{k+1}=\mathbf{\omega}_k e^{-\mathbf{y} \alpha_k \widehat{h}_{\mathbf{\omega}_k}(\mathbf{x}_i)}" title="\mathbf{\omega}_{k+1}=\mathbf{\omega}_k e^{-\mathbf{y} \alpha_k \widehat{h}_{\mathbf{\omega}_k}(\mathbf{x}_i)}" src="https://s0.wp.com/latex.php?latex=%5Cmathbf%7B%5Comega%7D_%7Bk%2B1%7D%3D%5Cmathbf%7B%5Comega%7D_k+e%5E%7B-%5Cmathbf%7By%7D+%5Calpha_k+%5Cwidehat%7Bh%7D_%7B%5Cmathbf%7B%5Comega%7D_k%7D%28%5Cmathbf%7Bx%7D_i%29%7D&bg=ffffff&fg=000&s=0" class="latex">(of course, devide by the sum to insure that the total sum is then 1). And as previously, one can include some shrinkage. To visualize the convergence of the process, we will plot the total error on our dataset.</p>
<div><table><tr><td><pre>1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
</pre></td><td><pre>n_iter <span>=</span> <span>100</span>
y <span>=</span> <span>(</span>myocarde<span>[</span>,<span>&quot;PRONO&quot;</span><span>]</span><span>==</span><span>1</span><span>)</span><span>*</span><span>2</span><span>-</span><span>1</span>
x <span>=</span> myocarde<span>[</span>,<span>1</span><span>:</span><span>7</span><span>]</span>
error <span>=</span> <span>rep</span><span>(</span><span>0</span>,n_iter<span>)</span> 
f <span>=</span> <span>rep</span><span>(</span><span>0</span>,<span>length</span><span>(</span>y<span>)</span><span>)</span> 
w <span>=</span> <span>rep</span><span>(</span><span>1</span>,<span>length</span><span>(</span>y<span>)</span><span>)</span> <span>#</span>
alpha <span>=</span> <span>1</span>
<span>library</span><span>(</span>rpart<span>)</span>
<span>for</span><span>(</span>i <span>in</span> <span>1</span><span>:</span>n_iter<span>)</span><span>{</span>
  w <span>=</span> <span>exp</span><span>(</span><span>-</span>alpha<span>*</span>y<span>*</span>f<span>)</span> <span>*</span>w 
  w <span>=</span> w<span>/</span><span>sum</span><span>(</span>w<span>)</span>
  rfit <span>=</span> rpart<span>(</span>y~., x, w, method<span>=</span><span>&quot;class&quot;</span><span>)</span>
  g <span>=</span> <span>-</span><span>1</span> <span>+</span> <span>2</span><span>*</span><span>(</span><span>predict</span><span>(</span>rfit,x<span>)</span><span>[</span>,<span>2</span><span>]</span><span>&amp;</span>gt<span>;</span>.5<span>)</span> 
  e <span>=</span> <span>sum</span><span>(</span>w<span>*</span><span>(</span>y<span>*</span>g<span>&amp;</span>lt<span>;</span><span>0</span><span>)</span><span>)</span>
  alpha <span>=</span> .5<span>*</span><span>log</span> <span>(</span> <span>(</span><span>1</span><span>-</span>e<span>)</span> <span>/</span> e <span>)</span>
  alpha <span>=</span> <span>0.1</span><span>*</span>alpha 
  f <span>=</span> f <span>+</span> alpha<span>*</span>g
  error<span>[</span>i<span>]</span> <span>=</span> <span>mean</span><span>(</span><span>1</span><span>*</span>f<span>*</span>y<span>&amp;</span>lt<span>;</span><span>0</span><span>)</span>
<span>}</span>
<span>plot</span><span>(</span><span>seq</span><span>(</span><span>1</span>,n_iter<span>)</span>,error,type<span>=&amp;</span>quot<span>;</span>l<span>&amp;</span>quot<span>;</span>,
     ylim<span>=</span><span>c</span><span>(</span><span>0</span>,.25<span>)</span>,<span>col</span><span>=&amp;</span>quot<span>;</span>blue<span>&amp;</span>quot<span>;</span>,
     ylab<span>=&amp;</span>quot<span>;</span>Error Rate<span>&amp;</span>quot<span>;</span>,xlab<span>=&amp;</span>quot<span>;</span>Iterations<span>&amp;</span>quot<span>;</span>,lwd<span>=</span><span>2</span><span>)</span></pre></td></tr></table></div>
<p><img alt="" src="https://i0.wp.com/freakonometrics.hypotheses.org/files/2018/06/Capture-d%E2%80%99e%CC%81cran-2018-06-09-a%CC%80-05.44.17.png?w=456&ssl=1" data-recalc-dims="1"><br>
Here we face a classical problem in machine learning: we have a perfect model. With zero error. That is nice, but not interesting. It is also possible in econometrics, with polynomial fits: with 10 observations, and a polynomial of degree 9, we have a perfect fit. But a poor model. Here it is the same. So the trick is to split our dataset in two, a training dataset, and a validation one</p>
<div><table><tr><td><pre>1
2
3
4
</pre></td><td><pre><span>set.<span>seed</span></span><span>(</span><span>123</span><span>)</span>
id_train <span>=</span> <span>sample</span><span>(</span><span>1</span><span>:</span><span>nrow</span><span>(</span>myocarde<span>)</span>, size<span>=</span><span>45</span>, <span>replace</span><span>=</span>FALSE<span>)</span>
train_myocarde <span>=</span> myocarde<span>[</span>id_train,<span>]</span>
test_myocarde <span>=</span> myocarde<span>[</span><span>-</span>id_train,<span>]</span></pre></td></tr></table></div>
<p>We construct the model on the first one, and we check on the second one that it’s not that bad…</p>
<div><table><tr><td><pre>1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
</pre></td><td><pre>y_train <span>=</span> <span>(</span>train_myocarde<span>[</span>,<span>&quot;PRONO&quot;</span><span>]</span><span>==</span><span>1</span><span>)</span><span>*</span><span>2</span><span>-</span><span>1</span>
x_train <span>=</span>  train_myocarde<span>[</span>,<span>1</span><span>:</span><span>7</span><span>]</span>
y_test <span>=</span> <span>(</span>test_myocarde<span>[</span>,<span>&quot;PRONO&quot;</span><span>]</span><span>==</span><span>1</span><span>)</span><span>*</span><span>2</span><span>-</span><span>1</span>
x_test <span>=</span> test_myocarde<span>[</span>,<span>1</span><span>:</span><span>7</span><span>]</span>
train_error <span>=</span> <span>rep</span><span>(</span><span>0</span>,n_iter<span>)</span> 
test_error <span>=</span> <span>rep</span><span>(</span><span>0</span>,n_iter<span>)</span>
f_train <span>=</span> <span>rep</span><span>(</span><span>0</span>,<span>length</span><span>(</span>y_train<span>)</span><span>)</span>
f_test <span>=</span> <span>rep</span><span>(</span><span>0</span>,<span>length</span><span>(</span>y_test<span>)</span><span>)</span> 
w_train <span>=</span> <span>rep</span><span>(</span><span>1</span>,<span>length</span><span>(</span>y_train<span>)</span><span>)</span> 
alpha <span>=</span> <span>1</span>
<span>for</span><span>(</span>i <span>in</span> <span>1</span><span>:</span>n_iter<span>)</span><span>{</span>
  w_train <span>=</span> w_train<span>*</span><span>exp</span><span>(</span><span>-</span>alpha<span>*</span>y_train<span>*</span>f_train<span>)</span> 
  w_train <span>=</span> w_train<span>/</span><span>sum</span><span>(</span>w_train<span>)</span>
  rfit <span>=</span> rpart<span>(</span>y_train~., x_train, w_train, method<span>=</span><span>&quot;class&quot;</span><span>)</span>
  g_train <span>=</span> <span>-</span><span>1</span> <span>+</span> <span>2</span><span>*</span><span>(</span><span>predict</span><span>(</span>rfit,x_train<span>)</span><span>[</span>,<span>2</span><span>]</span><span>&amp;</span>gt<span>;</span>.5<span>)</span>
  g_test <span>=</span> <span>-</span><span>1</span> <span>+</span> <span>2</span><span>*</span><span>(</span><span>predict</span><span>(</span>rfit,x_test<span>)</span><span>[</span>,<span>2</span><span>]</span><span>&amp;</span>gt<span>;</span>.5<span>)</span>
  e_train <span>=</span> <span>sum</span><span>(</span>w_train<span>*</span><span>(</span>y_train<span>*</span>g_train<span>&amp;</span>lt<span>;</span><span>0</span><span>)</span><span>)</span>
  alpha <span>=</span> .5<span>*</span><span>log</span> <span>(</span> <span>(</span><span>1</span><span>-</span>e_train<span>)</span> <span>/</span> e_train <span>)</span>
  alpha <span>=</span> <span>0.1</span><span>*</span>alpha 
  f_train <span>=</span> f_train <span>+</span> alpha<span>*</span>g_train
  f_test <span>=</span> f_test <span>+</span> alpha<span>*</span>g_test
  train_error<span>[</span>i<span>]</span> <span>=</span> <span>mean</span><span>(</span><span>1</span><span>*</span>f_train<span>*</span>y_train<span>&amp;</span>lt<span>;</span><span>0</span><span>)</span>
  test_error<span>[</span>i<span>]</span> <span>=</span> <span>mean</span><span>(</span><span>1</span><span>*</span>f_test<span>*</span>y_test<span>&amp;</span>lt<span>;</span><span>0</span><span>)</span><span>}</span>
<span>plot</span><span>(</span><span>seq</span><span>(</span><span>1</span>,n_iter<span>)</span>,test_error,<span>col</span><span>=</span><span>'red'</span><span>)</span>
<span>lines</span><span>(</span>train_error,lwd<span>=</span><span>2</span>,<span>col</span><span>=</span><span>'blue'</span><span>)</span></pre></td></tr></table></div>
<p><img alt="" src="https://i1.wp.com/freakonometrics.hypotheses.org/files/2018/06/Capture-d%E2%80%99e%CC%81cran-2018-06-09-a%CC%80-05.44.52.png?w=456&ssl=1" data-recalc-dims="1"><br>
Here, as previously, after 80 iterations, we have a perfect model on the training dataset, but it behaves badly on the validation dataset. But with 20 iterations, it seems to be ok…</p>
<h2>R function</h2>
<p>Of course, it’s possible to use R functions,</p>
<div><table><tr><td><pre>1
2
3
4
</pre></td><td><pre><span>library</span><span>(</span>gbm<span>)</span>
gbmWithCrossValidation <span>=</span> gbm<span>(</span>PRONO ~ .,distribution <span>=</span> <span>&quot;bernoulli&quot;</span>,
<span>data</span> <span>=</span> myocarde,n.<span>trees</span> <span>=</span> <span>2000</span>,shrinkage <span>=</span> .01,cv.<span>folds</span> <span>=</span> <span>5</span>,n.<span>cores</span> <span>=</span> <span>1</span><span>)</span>
bestTreeForPrediction <span>=</span> gbm.<span>perf</span><span>(</span>gbmWithCrossValidation<span>)</span></pre></td></tr></table></div>
<p>Here cross-validation is considered, and not training/validation, as well as forests instead of single trees, but overall, the idea is the same… Off course, the output is much nicer (here the shrinkage is a very small parameter, and learning is extremely slow)<br>
<img alt="" src="https://i2.wp.com/freakonometrics.hypotheses.org/files/2018/06/Capture-d%E2%80%99e%CC%81cran-2018-06-09-a%CC%80-05.55.20.png?w=456&ssl=1" data-recalc-dims="1"></p>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="https://freakonometrics.hypotheses.org/52782"> R-english – Freakonometrics</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img width="1" alt="" src="http://feeds.feedburner.com/~r/RBloggers/~4/CIl9njyD0t4" height="1"></a></div>