<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/using-regression-trees-for-forecasting-double-seasonal-time-series-with-trend-in-r/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://petolau.github.io/Regression-trees-for-forecasting-time-series-in-R/"> Peter Laurinec</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p>After blogging break caused by writing research papers, I managed to secure time to write something new about time series forecasting. This time I want to share with you my experiences with seasonal-trend time series forecasting using simple regression trees. Classification and <strong>regression tree</strong> (or <a rel="nofollow" href="https://en.wikipedia.org/wiki/Decision_tree_learning" target="_blank">decision tree</a>) is broadly used <strong>machine learning</strong> method for modeling. They are favorite because of these factors:</p>
<ul>
<li>simple to understand (white box)</li>
<li>from a tree we can extract interpretable results and make simple decisions</li>
<li>they are helpful for exploratory analysis as binary structure of tree is simple to visualize</li>
<li>very good prediction accuracy performance</li>
<li>very fast</li>
<li>they can be simply tuned by ensemble learning techniques</li>
</ul>
<p>But! There is always some “but”, they poorly adapt when new unexpected situations (values) appears. In other words, they can not detect and adapt to change or <a rel="nofollow" href="https://en.wikipedia.org/wiki/Concept_drift" target="_blank">concept drift</a> well (absolutely not). This is due to the fact that tree creates during learning just simple rules based on training data. Simple decision tree does not compute any regression coefficients like linear regression, so trend modeling is not possible. You would ask now, so why we are talking about time series forecasting with regression tree together, right? I will explain how to deal with it in more detail further in this post.</p>
<p>You will learn in this post how to:</p>
<ul>
<li>decompose double-seasonal time series</li>
<li>detrend time series</li>
<li>model and forecast <strong>double-seasonal time series</strong> with trend</li>
<li>use two types of simple <strong>regression trees</strong></li>
<li>set important hyperparameters related to regression tree</li>
</ul>
<h2>Exploring time series data of electricity consumption</h2>
<p>As in previous posts, I will use smart meter data of electricity consumption for demonstrating forecasting of seasonal time series. I created a new dataset of aggregated electricity load of consumers from an anonymous area. Time series data have the length of 17 weeks.</p>
<p>Firstly, let’s scan all of the needed packages for data analysis, modeling and visualizing.</p>
<figure>
<pre><code data-lang="r"><span>library</span><span>(</span><span>feather</span><span>)</span><span> </span><span># data import
</span><span>library</span><span>(</span><span>data.table</span><span>)</span><span> </span><span># data handle
</span><span>library</span><span>(</span><span>rpart</span><span>)</span><span> </span><span># decision tree method
</span><span>library</span><span>(</span><span>rpart.plot</span><span>)</span><span> </span><span># tree plot
</span><span>library</span><span>(</span><span>party</span><span>)</span><span> </span><span># decision tree method
</span><span>library</span><span>(</span><span>forecast</span><span>)</span><span> </span><span># forecasting methods
</span><span>library</span><span>(</span><span>ggplot2</span><span>)</span><span> </span><span># visualizations
</span><span>library</span><span>(</span><span>ggforce</span><span>)</span><span> </span><span># visualization tools
</span><span>library</span><span>(</span><span>plotly</span><span>)</span><span> </span><span># interactive visualizations
</span><span>library</span><span>(</span><span>grid</span><span>)</span><span> </span><span># visualizations
</span><span>library</span><span>(</span><span>animation</span><span>)</span><span> </span><span>#</span><span> </span><span>gif</span></code></pre>
</figure>
<p>Now read the mentioned time series data by <code>read_feather</code> to one <code>data.table</code>.</p>
<figure>
<pre><code data-lang="r"><span>DT</span><span> </span><span>&lt;-</span><span> </span><span>as.data.table</span><span>(</span><span>read_feather</span><span>(</span><span>&quot;DT_load_17weeks&quot;</span><span>))</span></code></pre>
</figure>
<p>And store information of the date and period of time series that is 48.</p>
<figure>
<pre><code data-lang="r"><span>n_date</span><span> </span><span>&lt;-</span><span> </span><span>unique</span><span>(</span><span>DT</span><span>[,</span><span> </span><span>date</span><span>])</span><span>
</span><span>period</span><span> </span><span>&lt;-</span><span> </span><span>48</span></code></pre>
</figure>
<p>For data visualization needs, store my favorite ggplot theme settings by function <code>theme</code>.</p>
<figure>
<pre><code data-lang="r"><span>theme_ts</span><span> </span><span>&lt;-</span><span> </span><span>theme</span><span>(</span><span>panel.border</span><span> </span><span>=</span><span> </span><span>element_rect</span><span>(</span><span>fill</span><span> </span><span>=</span><span> </span><span>NA</span><span>,</span><span> </span><span>colour</span><span> </span><span>=</span><span> </span><span>&quot;grey10&quot;</span><span>),</span><span>
</span><span>panel.background</span><span> </span><span>=</span><span> </span><span>element_blank</span><span>(),</span><span>
</span><span>panel.grid.minor</span><span> </span><span>=</span><span> </span><span>element_line</span><span>(</span><span>colour</span><span> </span><span>=</span><span> </span><span>&quot;grey85&quot;</span><span>),</span><span>
</span><span>panel.grid.major</span><span> </span><span>=</span><span> </span><span>element_line</span><span>(</span><span>colour</span><span> </span><span>=</span><span> </span><span>&quot;grey85&quot;</span><span>),</span><span>
</span><span>panel.grid.major.x</span><span> </span><span>=</span><span> </span><span>element_line</span><span>(</span><span>colour</span><span> </span><span>=</span><span> </span><span>&quot;grey85&quot;</span><span>),</span><span>
</span><span>axis.text</span><span> </span><span>=</span><span> </span><span>element_text</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>13</span><span>,</span><span> </span><span>face</span><span> </span><span>=</span><span> </span><span>&quot;bold&quot;</span><span>),</span><span>
</span><span>axis.title</span><span> </span><span>=</span><span> </span><span>element_text</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>15</span><span>,</span><span> </span><span>face</span><span> </span><span>=</span><span> </span><span>&quot;bold&quot;</span><span>),</span><span>
</span><span>plot.title</span><span> </span><span>=</span><span> </span><span>element_text</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>16</span><span>,</span><span> </span><span>face</span><span> </span><span>=</span><span> </span><span>&quot;bold&quot;</span><span>),</span><span>
</span><span>strip.text</span><span> </span><span>=</span><span> </span><span>element_text</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>16</span><span>,</span><span> </span><span>face</span><span> </span><span>=</span><span> </span><span>&quot;bold&quot;</span><span>),</span><span>
</span><span>strip.background</span><span> </span><span>=</span><span> </span><span>element_rect</span><span>(</span><span>colour</span><span> </span><span>=</span><span> </span><span>&quot;black&quot;</span><span>),</span><span>
</span><span>legend.text</span><span> </span><span>=</span><span> </span><span>element_text</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>15</span><span>),</span><span>
</span><span>legend.title</span><span> </span><span>=</span><span> </span><span>element_text</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>16</span><span>,</span><span> </span><span>face</span><span> </span><span>=</span><span> </span><span>&quot;bold&quot;</span><span>),</span><span>
</span><span>legend.background</span><span> </span><span>=</span><span> </span><span>element_rect</span><span>(</span><span>fill</span><span> </span><span>=</span><span> </span><span>&quot;white&quot;</span><span>),</span><span>
</span><span>legend.key</span><span> </span><span>=</span><span> </span><span>element_rect</span><span>(</span><span>fill</span><span> </span><span>=</span><span> </span><span>&quot;white&quot;</span><span>))</span></code></pre>
</figure>
<p>Now, pick some dates of the length 3 weeks from dataset to split data on the train and test part. Test set has the length of only one day because we will perform one day ahead forecast of electricity consumption.</p>
<figure>
<pre><code data-lang="r"><span>data_train</span><span> </span><span>&lt;-</span><span> </span><span>DT</span><span>[</span><span>date</span><span> </span><span>%in%</span><span> </span><span>n_date</span><span>[</span><span>43</span><span>:</span><span>63</span><span>]]</span><span>
</span><span>data_test</span><span> </span><span>&lt;-</span><span> </span><span>DT</span><span>[</span><span>date</span><span> </span><span>%in%</span><span> </span><span>n_date</span><span>[</span><span>64</span><span>]]</span></code></pre>
</figure>
<p>Let’s plot the train set and corresponding average weekly values of electricity load.</p>
<figure>
<pre><code data-lang="r"><span>averages</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>value</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>sapply</span><span>(</span><span>0</span><span>:</span><span>2</span><span>,</span><span> </span><span>function</span><span>(</span><span>i</span><span>)</span><span>
</span><span>mean</span><span>(</span><span>data_train</span><span>[((</span><span>i</span><span>*</span><span>period</span><span>*</span><span>7</span><span>)</span><span>+1</span><span>)</span><span>:</span><span>((</span><span>i</span><span>+1</span><span>)</span><span>*</span><span>period</span><span>*</span><span>7</span><span>),</span><span> </span><span>value</span><span>])),</span><span>
</span><span>each</span><span> </span><span>=</span><span> </span><span>period</span><span> </span><span>*</span><span> </span><span>7</span><span>),</span><span>
</span><span>date_time</span><span> </span><span>=</span><span> </span><span>data_train</span><span>$</span><span>date_time</span><span>)</span><span>
</span><span>ggplot</span><span>(</span><span>data_train</span><span>,</span><span> </span><span>aes</span><span>(</span><span>date_time</span><span>,</span><span> </span><span>value</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>()</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>data</span><span> </span><span>=</span><span> </span><span>averages</span><span>,</span><span> </span><span>aes</span><span>(</span><span>date_time</span><span>,</span><span> </span><span>value</span><span>),</span><span>
</span><span>linetype</span><span> </span><span>=</span><span> </span><span>5</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.75</span><span>,</span><span> </span><span>size</span><span> </span><span>=</span><span> </span><span>1.2</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>&quot;firebrick2&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>&quot;Date&quot;</span><span>,</span><span> </span><span>y</span><span> </span><span>=</span><span> </span><span>&quot;Load (kW)&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-7" src="https://i1.wp.com/petolau.github.io/images/post_5/unnamed-chunk-7-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>We can see some trend increasing over time, maybe air conditioning is more used when gets hotter in summer. The double-seasonal (daily and weekly) character of time series is obvious.</p>
<p>A very useful method for visualization and analysis of time series is <a rel="nofollow" href="https://www.otexts.org/fpp/6/5" target="_blank">STL decomposition</a>.<br>
STL decomposition is based on Loess regression, and it decomposes time series to three parts: seasonal, trend and remainder.<br>
We will use results from the STL decomposition to model our data as well.<br>
I am using <code>stl()</code> from <code>stats</code> package and before computation we must define weekly seasonality to our time series object. Let’s look on results:</p>
<figure>
<pre><code data-lang="r"><span>data_ts</span><span> </span><span>&lt;-</span><span> </span><span>ts</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>freq</span><span> </span><span>=</span><span> </span><span>period</span><span> </span><span>*</span><span> </span><span>7</span><span>)</span><span>
</span><span>decomp_ts</span><span> </span><span>&lt;-</span><span> </span><span>stl</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>s.window</span><span> </span><span>=</span><span> </span><span>&quot;periodic&quot;</span><span>,</span><span> </span><span>robust</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>$</span><span>time.series</span><span>
</span><span>decomp_stl</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>as.numeric</span><span>(</span><span>decomp_ts</span><span>)),</span><span>
</span><span>Date</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>data_train</span><span>[,</span><span>date_time</span><span>],</span><span> </span><span>ncol</span><span>(</span><span>decomp_ts</span><span>)</span><span>+1</span><span>),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>factor</span><span>(</span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;original data&quot;</span><span>,</span><span> </span><span>colnames</span><span>(</span><span>decomp_ts</span><span>)),</span><span>
</span><span>each</span><span> </span><span>=</span><span> </span><span>nrow</span><span>(</span><span>decomp_ts</span><span>)),</span><span>
</span><span>levels</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>&quot;original data&quot;</span><span>,</span><span> </span><span>colnames</span><span>(</span><span>decomp_ts</span><span>))))</span><span>
</span><span>ggplot</span><span>(</span><span>decomp_stl</span><span>,</span><span> </span><span>aes</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>Date</span><span>,</span><span> </span><span>y</span><span> </span><span>=</span><span> </span><span>Load</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>()</span><span> </span><span>+</span><span> </span><span>facet_grid</span><span>(</span><span>Type</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>scales</span><span> </span><span>=</span><span> </span><span>&quot;free_y&quot;</span><span>,</span><span> </span><span>switch</span><span> </span><span>=</span><span> </span><span>&quot;y&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>&quot;Date&quot;</span><span>,</span><span> </span><span>y</span><span> </span><span>=</span><span> </span><span>NULL</span><span>,</span><span>
</span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Time Series Decomposition by STL&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-8" src="https://i2.wp.com/petolau.github.io/images/post_5/unnamed-chunk-8-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>As was expected from the previous picture, we can see that there is “slight” trend increasing and decreasing (by around 100 kW so slightly large <img alt="😉" src="https://s.w.org/images/core/emoji/2/72x72/1f609.png" class="wp-smiley"> ).<br>
Remainder part (noise) is very fluctuate and not seems like classical <a rel="nofollow" href="https://en.wikipedia.org/wiki/White_noise" target="_blank">white noise</a> (we obviously missing additional information like weather and other unexpected situations).</p>
<h2>Constructing features to model</h2>
<p>In this section I will do <strong>feature engineering</strong> for modeling double-seasonal <strong>time series</strong> with trend best as possible by just available historical values.</p>
<p>Classical way to handle seasonality is to add seasonal features to a model as vectors of form \( (1, \dots, DailyPeriod, 1, …, DailyPeriod,…) \) for daily season or \( (1, \dots, 1, 2, \dots, 2, \dots , 7, 1, \dots) \) for weekly season. I used it this way in my previous <a rel="nofollow" href="https://petolau.github.io/Analyzing-double-seasonal-time-series-with-GAM-in-R/" target="_blank">post about GAM</a> and somehow similar also with <a rel="nofollow" href="https://petolau.github.io/Forecast-double-seasonal-time-series-with-multiple-linear-regression-in-R/" target="_blank">multiple linear regression</a>.</p>
<p>A better way to model seasonal variables (features) with <strong>nonlinear regression</strong> methods like <strong>tree</strong> is to transform it to <a rel="nofollow" href="https://en.wikipedia.org/wiki/Fourier_transform" target="_blank">Fourier terms</a> (sinus and cosine). It is more effective to tree models and also other nonlinear machine learning methods. I will explain why it is like that further of this post.</p>
<p>Fourier daily signals (terms) are defined as:</p>
<p></p>
<p>where \( ds \) is number of daily seasonality Fourier pairs and Fourier weekly terms are defines as:</p>
<p></p>
<p>where \( ws \) is a number of weekly seasonality Fourier pairs.</p>
<p>Another great feature (most of the times most powerful) is a lag of original time series. We can use lag by one day, one week, etc…<br>
The lag of time series can be preprocessed by removing noise or trend for example by STL decomposition method to ensure stability.</p>
<p>As was earlier mentioned, regression trees can’t predict trend because they logically make rules and predict future values only by rules made by training set.<br>
Therefore original time series that inputs to regression tree as dependent variable must be <a rel="nofollow" href="https://en.wiktionary.org/wiki/detrending" target="_blank">detrended</a> (removing the trend part of the time series). The acquired trend part then can be forecasted by for example <a rel="nofollow" href="https://en.wikipedia.org/wiki/Autoregressive_integrated_moving_average" target="_blank">ARIMA</a> model.</p>
<p>Let’s go to constructing mentioned features and trend forecasting.</p>
<p>Double-seasonal Fourier terms can be simply extracted by <code>fourier</code> function from <code>forecast</code> package.<br>
Firstly, we must create multiple seasonal object with function <code>msts</code>.</p>
<figure>
<pre><code data-lang="r"><span>data_msts</span><span> </span><span>&lt;-</span><span> </span><span>msts</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>seasonal.periods</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>period</span><span>,</span><span> </span><span>period</span><span>*</span><span>7</span><span>))</span></code></pre>
</figure>
<p>Now use <code>fourier</code> function using two conditions for a number of K terms.<br>
Set K for example just to 2.</p>
<figure>
<pre><code data-lang="r"><span>K</span><span> </span><span>&lt;-</span><span> </span><span>2</span><span>
</span><span>fuur</span><span> </span><span>&lt;-</span><span> </span><span>fourier</span><span>(</span><span>data_msts</span><span>,</span><span> </span><span>K</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>K</span><span>,</span><span> </span><span>K</span><span>))</span></code></pre>
</figure>
<p>It made 2 pairs (sine and cosine) of daily and weekly seasonal signals.<br>
If we compare it with approach described in previous posts, so simple periodic vectors, it looks like this:</p>
<figure>
<pre><code data-lang="r"><span>Daily</span><span> </span><span>&lt;-</span><span> </span><span>rep</span><span>(</span><span>1</span><span>:</span><span>period</span><span>,</span><span> </span><span>21</span><span>)</span><span> </span><span># simple daily vector
</span><span>Weekly</span><span> </span><span>&lt;-</span><span> </span><span>data_train</span><span>[,</span><span> </span><span>week_num</span><span>]</span><span> </span><span># simple weekly vector
</span><span> </span><span>data_fuur_simple</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>value</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>scale</span><span>(</span><span>Daily</span><span>),</span><span> </span><span>fuur</span><span>[,</span><span>2</span><span>],</span><span> </span><span>scale</span><span>(</span><span>Weekly</span><span>),</span><span> </span><span>fuur</span><span>[,</span><span>6</span><span>]),</span><span>
</span><span>date</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>data_train</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>4</span><span>),</span><span>
</span><span>method</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;simple-daily&quot;</span><span>,</span><span> </span><span>&quot;four-daily&quot;</span><span>,</span><span>
</span><span>&quot;simple-weekly&quot;</span><span>,</span><span> </span><span>&quot;four-weekly&quot;</span><span>),</span><span>
</span><span>each</span><span> </span><span>=</span><span> </span><span>nrow</span><span>(</span><span>fuur</span><span>)),</span><span>
</span><span>type</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;Daily season&quot;</span><span>,</span><span> </span><span>&quot;Weekly season&quot;</span><span>),</span><span>
</span><span>each</span><span> </span><span>=</span><span> </span><span>nrow</span><span>(</span><span>fuur</span><span>)</span><span>*</span><span>2</span><span>))</span><span>
</span><span>ggplot</span><span>(</span><span>data_fuur_simple</span><span>,</span><span> </span><span>aes</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>date</span><span>,</span><span> </span><span>y</span><span> </span><span>=</span><span> </span><span>value</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>method</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>1.2</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.7</span><span>)</span><span> </span><span>+</span><span> </span><span>facet_grid</span><span>(</span><span>type</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>scales</span><span> </span><span>=</span><span> </span><span>&quot;free_y&quot;</span><span>,</span><span> </span><span>switch</span><span> </span><span>=</span><span> </span><span>&quot;y&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>&quot;Date&quot;</span><span>,</span><span> </span><span>y</span><span> </span><span>=</span><span> </span><span>NULL</span><span>,</span><span>
</span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Features Comparison&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-11" src="https://i2.wp.com/petolau.github.io/images/post_5/unnamed-chunk-11-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>where <em>four-daily</em> is the Fourier term for daily season, <em>simple-daily</em> is the simple feature for daily season, <em>four-weekly</em> is the Fourier term for weekly season, and <em>simple-weekly</em> is the simple feature for weekly season. The <strong>advantage</strong> of <strong>Fourier terms</strong> is that there is much more closeness between ending and starting of a day or a week, which is more natural.</p>
<p>Now, let’s use data from STL decomposition to forecast trend part of time series. I will use <code>auto.arima</code> procedure from the <code>forecast</code> package to perform this.</p>
<figure>
<pre><code data-lang="r"><span>trend_part</span><span> </span><span>&lt;-</span><span> </span><span>ts</span><span>(</span><span>decomp_ts</span><span>[,</span><span>2</span><span>])</span><span>
</span><span>trend_fit</span><span> </span><span>&lt;-</span><span> </span><span>auto.arima</span><span>(</span><span>trend_part</span><span>)</span><span>
</span><span>trend_for</span><span> </span><span>&lt;-</span><span> </span><span>forecast</span><span>(</span><span>trend_fit</span><span>,</span><span> </span><span>period</span><span>)</span><span>$</span><span>mean</span></code></pre>
</figure>
<p>Let’s plot it:</p>
<figure>
<pre><code data-lang="r"><span>trend_data</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>decomp_ts</span><span>[,</span><span>2</span><span>],</span><span> </span><span>trend_for</span><span>),</span><span>
</span><span>Date</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>data_test</span><span>$</span><span>date_time</span><span>),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>rep</span><span>(</span><span>&quot;Real&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)),</span><span> </span><span>rep</span><span>(</span><span>&quot;Forecast&quot;</span><span>,</span><span>
</span><span>nrow</span><span>(</span><span>data_test</span><span>))))</span><span>
</span><span>ggplot</span><span>(</span><span>trend_data</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Date</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>1.2</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>title</span><span> </span><span>=</span><span> </span><span>paste</span><span>(</span><span>trend_fit</span><span>))</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-13" src="https://i2.wp.com/petolau.github.io/images/post_5/unnamed-chunk-13-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>Function <code>auto.arima</code> chose ARIMA(0,2,0) model as best for trend forecasting.</p>
<p>Next, make the final feature to the model (lag) and construct train matrix (model matrix).<br>
I am creating lag by one day and just taking seasonal part from STL decomposition (for having smooth lag time series feature).</p>
<figure>
<pre><code data-lang="r"><span>N</span><span> </span><span>&lt;-</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)</span><span>
</span><span>window</span><span> </span><span>&lt;-</span><span> </span><span>(</span><span>N</span><span> </span><span>/</span><span> </span><span>period</span><span>)</span><span> </span><span>-</span><span> </span><span>1</span><span> </span><span># number of days in train set minus lag
</span><span> </span><span>new_load</span><span> </span><span>&lt;-</span><span> </span><span>rowSums</span><span>(</span><span>decomp_ts</span><span>[,</span><span> </span><span>c</span><span>(</span><span>1</span><span>,</span><span>3</span><span>)])</span><span> </span><span># detrended load
</span><span>lag_seas</span><span> </span><span>&lt;-</span><span> </span><span>decomp_ts</span><span>[</span><span>1</span><span>:</span><span>(</span><span>period</span><span>*</span><span>window</span><span>),</span><span> </span><span>1</span><span>]</span><span> </span><span># seasonal part of time series as lag feature
</span><span> </span><span>matrix_train</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>tail</span><span>(</span><span>new_load</span><span>,</span><span> </span><span>window</span><span>*</span><span>period</span><span>),</span><span>
</span><span>fuur</span><span>[(</span><span>period</span><span> </span><span>+</span><span> </span><span>1</span><span>)</span><span>:</span><span>N</span><span>,],</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>lag_seas</span><span>)</span></code></pre>
</figure>
<p>The accuracy of forecast (or fitted values of a model) will be measured by MAPE, let’s defined it:</p>
<figure>
<pre><code data-lang="r"><span>mape</span><span> </span><span>&lt;-</span><span> </span><span>function</span><span>(</span><span>real</span><span>,</span><span> </span><span>pred</span><span>){</span><span>
</span><span>return</span><span>(</span><span>100</span><span> </span><span>*</span><span> </span><span>mean</span><span>(</span><span>abs</span><span>((</span><span>real</span><span> </span><span>-</span><span> </span><span>pred</span><span>)</span><span>/</span><span>real</span><span>)))</span><span> </span><span># MAPE - Mean Absolute Percentage Error
</span><span>}</span></code></pre>
</figure>
<h2>RPART (CART) tree</h2>
<p>In the next two sections, I will describe two <strong>regression tree</strong> methods. The first is <strong>RPART</strong>, or CART (Classification and Regression Trees), the second will be CTREE. RPART is recursive partitioning type of binary tree for classification or regression tasks. It performs a search over all possible splits by maximizing an information measure of node impurity, selecting the covariate showing the best split.</p>
<p>I’m using <code>rpart</code> implementation from the same named package. Let’s go forward to modeling and try default settings of <code>rpart</code> function:</p>
<figure>
<pre><code data-lang="r"><span>tree_1</span><span> </span><span>&lt;-</span><span> </span><span>rpart</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train</span><span>)</span></code></pre>
</figure>
<p>It makes many interesting outputs to check, for example we can see a table of nodes and corresponding errors by <code>printcp(tree_1)</code> or see a detailed summary of created nodes by <code>summary(tree_1)</code>. We will check variable importance and number of created splits:</p>
<figure>
<pre><code data-lang="r"><span>tree_1</span><span>$</span><span>variable.importance</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">##       Lag     C2-48    C1-336     S1-48     C1-48    S1-336    C2-336 
## 100504751  45918330  44310331  36245736  32359598  27831258  25385506 
##    S2-336     S2-48 
##  15156041   7595266</code></pre>
</figure>
<figure>
<pre><code data-lang="r"><span>paste</span><span>(</span><span>&quot;Number of splits: &quot;</span><span>,</span><span> </span><span>tree_1</span><span>$</span><span>cptable</span><span>[</span><span>dim</span><span>(</span><span>tree_1</span><span>$</span><span>cptable</span><span>)[</span><span>1</span><span>],</span><span> </span><span>&quot;nsplit&quot;</span><span>])</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] &quot;Number of splits:  10&quot;</code></pre>
</figure>
<p>We can see that most important variables are Lag and cosine forms of the daily and weekly season. The number of splits is 10, ehm, is it enough for time series of length 1008 values?</p>
<p>Let’s plot created rules with fancy <code>rpart.plot</code> function from the same named package:</p>
<figure>
<pre><code data-lang="r"><span>rpart.plot</span><span>(</span><span>tree_1</span><span>,</span><span> </span><span>digits</span><span> </span><span>=</span><span> </span><span>2</span><span>,</span><span> </span><span>box.palette</span><span> </span><span>=</span><span> </span><span>viridis</span><span>::</span><span>viridis</span><span>(</span><span>10</span><span>,</span><span> </span><span>option</span><span> </span><span>=</span><span> </span><span>&quot;D&quot;</span><span>,</span><span> </span><span>begin</span><span> </span><span>=</span><span> </span><span>0.85</span><span>,</span><span> </span><span>end</span><span> </span><span>=</span><span> </span><span>0</span><span>),</span><span> </span><span>shadow.col</span><span> </span><span>=</span><span> </span><span>&quot;grey65&quot;</span><span>,</span><span> </span><span>col</span><span> </span><span>=</span><span> </span><span>&quot;grey99&quot;</span><span>)</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-18" src="https://i1.wp.com/petolau.github.io/images/post_5/unnamed-chunk-18-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>We can see values, rules, and percentage of values split each time. Pretty simple and interpretable.</p>
<p>Now plot fitted values to see results of the <code>tree_1</code> model.</p>
<figure>
<pre><code data-lang="r"><span>datas</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span>
</span><span>predict</span><span>(</span><span>tree_1</span><span>)),</span><span>
</span><span>Time</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>1</span><span>:</span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>),</span><span> </span><span>2</span><span>),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;Real&quot;</span><span>,</span><span> </span><span>&quot;RPART&quot;</span><span>),</span><span> </span><span>each</span><span> </span><span>=</span><span> </span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>)))</span><span>
</span><span>ggplot</span><span>(</span><span>datas</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Time</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.75</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>y</span><span> </span><span>=</span><span> </span><span>&quot;Detrended load&quot;</span><span>,</span><span> </span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Fitted values from RPART tree&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-19" src="https://i2.wp.com/petolau.github.io/images/post_5/unnamed-chunk-19-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>And see the error of fitted values against real values.</p>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span> </span><span>predict</span><span>(</span><span>tree_1</span><span>))</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 180.6669</code></pre>
</figure>
<p>Whups. It’s a little bit simple (rectangular) and not really accurate, but it’s logical result from a simple tree model.<br>
The key to achieving better results and have more accurate fit is to set manually control hyperparameters of rpart.<br>
Check <code>?rpart.control</code> to get more information.<br>
The “hack” is to change cp (complexity) parameter to very low to produce more splits (nodes). The <code>cp</code> is a threshold deciding if each branch fulfills conditions for further processing, so only nodes with fitness larger than factor cp are processed. Other important parameters are the minimum number of observations in needed in a node to split (<code>minsplit</code>) and the maximal depth of a tree (<code>maxdepth</code>).<br>
Set the <code>minsplit</code> to 2 and set the <code>maxdepth</code> to its maximal value – 30.</p>
<figure>
<pre><code data-lang="r"><span>tree_2</span><span> </span><span>&lt;-</span><span> </span><span>rpart</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train</span><span>,</span><span>
</span><span>control</span><span> </span><span>=</span><span> </span><span>rpart.control</span><span>(</span><span>minsplit</span><span> </span><span>=</span><span> </span><span>2</span><span>,</span><span>
</span><span>maxdepth</span><span> </span><span>=</span><span> </span><span>30</span><span>,</span><span>
</span><span>cp</span><span> </span><span>=</span><span> </span><span>0.000001</span><span>))</span></code></pre>
</figure>
<p>Now make simple plot to see depth of the created tree…</p>
<figure>
<pre><code data-lang="r"><span>plot</span><span>(</span><span>tree_2</span><span>,</span><span> </span><span>compress</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-22" src="https://i1.wp.com/petolau.github.io/images/post_5/unnamed-chunk-22-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>That’s little bit impressive difference than previous one, isn’t it?<br>
Check also number of splits.</p>
<figure>
<pre><code data-lang="r"><span>tree_2</span><span>$</span><span>cptable</span><span>[</span><span>dim</span><span>(</span><span>tree_2</span><span>$</span><span>cptable</span><span>)[</span><span>1</span><span>],</span><span> </span><span>&quot;nsplit&quot;</span><span>]</span><span> </span><span>#</span><span> </span><span>Number</span><span> </span><span>of</span><span> </span><span>splits</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 600</code></pre>
</figure>
<p>600 is higher than 10 <img alt="🙂" src="https://s.w.org/images/core/emoji/2/72x72/1f642.png" class="wp-smiley"></p>
<p>Let’s plot fitted values from the model <code>tree_2</code>:</p>
<figure>
<pre><code data-lang="r"><span>datas</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span>
</span><span>predict</span><span>(</span><span>tree_2</span><span>)),</span><span>
</span><span>Time</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>1</span><span>:</span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>),</span><span> </span><span>2</span><span>),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;Real&quot;</span><span>,</span><span> </span><span>&quot;RPART&quot;</span><span>),</span><span> </span><span>each</span><span> </span><span>=</span><span> </span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>)))</span><span>
</span><span>ggplot</span><span>(</span><span>datas</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Time</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.75</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>y</span><span> </span><span>=</span><span> </span><span>&quot;Detrended load&quot;</span><span>,</span><span> </span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Fitted values from RPART&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-24" src="https://i1.wp.com/petolau.github.io/images/post_5/unnamed-chunk-24-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>And see the error of fitted values against real values.</p>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span> </span><span>predict</span><span>(</span><span>tree_2</span><span>))</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 16.0639</code></pre>
</figure>
<p>Much better, but obviously the model can be overfitted now.</p>
<p>Add together everything that we got till now, so forecast load one day ahead.<br>
Let’s create testing data matrix:</p>
<figure>
<pre><code data-lang="r"><span>test_lag</span><span> </span><span>&lt;-</span><span> </span><span>decomp_ts</span><span>[((</span><span>period</span><span>*</span><span>window</span><span>)</span><span>+1</span><span>)</span><span>:</span><span>N</span><span>,</span><span> </span><span>1</span><span>]</span><span>
</span><span>fuur_test</span><span> </span><span>&lt;-</span><span> </span><span>fourier</span><span>(</span><span>data_msts</span><span>,</span><span> </span><span>K</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>K</span><span>,</span><span> </span><span>K</span><span>),</span><span> </span><span>h</span><span> </span><span>=</span><span> </span><span>period</span><span>)</span><span>
</span><span>matrix_test</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>fuur_test</span><span>,</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>test_lag</span><span>)</span></code></pre>
</figure>
<p>Predict detrended time series part with <code>tree_2</code> model + add the trend part of time series forecasted by ARIMA model.</p>
<figure>
<pre><code data-lang="r"><span>for_rpart</span><span> </span><span>&lt;-</span><span> </span><span>predict</span><span>(</span><span>tree_2</span><span>,</span><span> </span><span>matrix_test</span><span>)</span><span> </span><span>+</span><span> </span><span>trend_for</span></code></pre>
</figure>
<p>Let’s plot the results and compare it with real values from <code>data_test</code>.</p>
<figure>
<pre><code data-lang="r"><span>data_for</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_rpart</span><span>),</span><span>
</span><span>Date</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>rep</span><span>(</span><span>data_test</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>2</span><span>)),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>rep</span><span>(</span><span>&quot;Train data&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;Test data&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;Forecast&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>))))</span><span>
</span><span>ggplot</span><span>(</span><span>data_for</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Date</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.75</span><span>)</span><span> </span><span>+</span><span>
</span><span>facet_zoom</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>Date</span><span> </span><span>%in%</span><span> </span><span>data_test</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>zoom.size</span><span> </span><span>=</span><span> </span><span>1.2</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Forecast from RPART&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-28" src="https://i0.wp.com/petolau.github.io/images/post_5/unnamed-chunk-28-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>Not bad. For clarity, compare forecasting results with model without separate trend forecasting and detrending.</p>
<figure>
<pre><code data-lang="r"><span>matrix_train_sim</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>tail</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>window</span><span>*</span><span>period</span><span>),</span><span>
</span><span>fuur</span><span>[(</span><span>period</span><span>+1</span><span>)</span><span>:</span><span>N</span><span>,],</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>lag_seas</span><span>)</span><span>
</span><span>tree_sim</span><span> </span><span>&lt;-</span><span> </span><span>rpart</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train_sim</span><span>,</span><span>
</span><span>control</span><span> </span><span>=</span><span> </span><span>rpart.control</span><span>(</span><span>minsplit</span><span> </span><span>=</span><span> </span><span>2</span><span>,</span><span>
</span><span>maxdepth</span><span> </span><span>=</span><span> </span><span>30</span><span>,</span><span>
</span><span>cp</span><span> </span><span>=</span><span> </span><span>0.000001</span><span>))</span><span>
</span><span>for_rpart_sim</span><span> </span><span>&lt;-</span><span> </span><span>predict</span><span>(</span><span>tree_sim</span><span>,</span><span> </span><span>matrix_test</span><span>)</span><span>
</span><span>data_for</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_rpart</span><span>,</span><span> </span><span>for_rpart_sim</span><span>),</span><span>
</span><span>Date</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>rep</span><span>(</span><span>data_test</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>3</span><span>)),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>rep</span><span>(</span><span>&quot;Train data&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;Test data&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;Forecast with trend&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;Forecast simple&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>))))</span><span>
</span><span>ggplot</span><span>(</span><span>data_for</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Date</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>,</span><span> </span><span>linetype</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.7</span><span>)</span><span> </span><span>+</span><span>
</span><span>facet_zoom</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>Date</span><span> </span><span>%in%</span><span> </span><span>data_test</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>zoom.size</span><span> </span><span>=</span><span> </span><span>1.2</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Forecasts from RPARTs with and without trend forecasting&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>scale_linetype_manual</span><span>(</span><span>values</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>5</span><span>,</span><span>6</span><span>,</span><span>1</span><span>,</span><span>1</span><span>))</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-29" src="https://i0.wp.com/petolau.github.io/images/post_5/unnamed-chunk-29-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>We can see that RPART model without trend manipulation has higher values of the forecast.<br>
Evaluate results with MAPE forecasting measure.</p>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_rpart</span><span>)</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 3.727473</code></pre>
</figure>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_rpart_sim</span><span>)</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 6.976259</code></pre>
</figure>
<p>We can see the large difference in MAPE. So detrending original time series and forecasting separately trend part really works, but not generalize the result now. You can read more about RPART method in its great <a rel="nofollow" href="https://cran.r-project.org/web/packages/rpart/vignettes/longintro.pdf" target="_blank">package vignette</a>.</p>
<h2>CTREE</h2>
<p>The second simple regression tree method that will be used is <strong>CTREE</strong>. Conditional inference trees (CTREE) is a statistical approach to recursive partitioning, which takes into account the distributional properties of the data. CTREE performs multiple test procedures that are applied to determine whether no significant association between any of the feature and the response (load in the our case) can be stated and the recursion needs to stop.<br>
In <strong>R</strong> CTREE is implemented in the package <code>party</code> in the function <code>ctree</code>.</p>
<p>Let’s try fit simple <code>ctree</code> with a default values.</p>
<figure>
<pre><code data-lang="r"><span>ctree_1</span><span> </span><span>&lt;-</span><span> </span><span>ctree</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train</span><span>)</span></code></pre>
</figure>
<p>Constructed tree can be again simply plotted by <code>plot</code> function, but it made many splits so it’s disarranged.</p>
<p>Let’s plot fitted values from <code>ctree_1</code> model.</p>
<figure>
<pre><code data-lang="r"><span>datas</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span>
</span><span>predict</span><span>(</span><span>ctree_1</span><span>)),</span><span>
</span><span>Time</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>1</span><span>:</span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>),</span><span> </span><span>2</span><span>),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;Real&quot;</span><span>,</span><span> </span><span>&quot;CTREE&quot;</span><span>),</span><span> </span><span>each</span><span> </span><span>=</span><span> </span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>)))</span><span>
</span><span>ggplot</span><span>(</span><span>datas</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Time</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.75</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>y</span><span> </span><span>=</span><span> </span><span>&quot;Detrended load&quot;</span><span>,</span><span> </span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Fitted values from CTREE&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-32" src="https://i1.wp.com/petolau.github.io/images/post_5/unnamed-chunk-32-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>And see the error of fitted values against real values.</p>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span> </span><span>predict</span><span>(</span><span>ctree_1</span><span>))</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 87.85983</code></pre>
</figure>
<p>Actually, this is pretty nice, but again, it can be tuned.</p>
<p>For available hyperparameters tuning check <code>?ctree_control</code>. I changed hyperparameters <code>minsplit</code> and <code>minbucket</code> that have similar meaning like the cp parameter in RPART. The <code>mincriterion</code> can be tuned also, and it is significance level (1 – p-value) that must be exceeded in order to implement a split. Let’s plot results.</p>
<figure>
<pre><code data-lang="r"><span>ctree_2</span><span> </span><span>&lt;-</span><span> </span><span>ctree</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train</span><span>,</span><span>
</span><span>controls</span><span> </span><span>=</span><span> </span><span>party</span><span>::</span><span>ctree_control</span><span>(</span><span>teststat</span><span> </span><span>=</span><span> </span><span>&quot;quad&quot;</span><span>,</span><span> </span><span>testtype</span><span> </span><span>=</span><span> </span><span>&quot;Teststatistic&quot;</span><span>,</span><span> </span><span>mincriterion</span><span> </span><span>=</span><span> </span><span>0.925</span><span>,</span><span>
</span><span>minsplit</span><span> </span><span>=</span><span> </span><span>1</span><span>,</span><span>
</span><span>minbucket</span><span> </span><span>=</span><span> </span><span>1</span><span>))</span><span>
</span><span>datas</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span>
</span><span>predict</span><span>(</span><span>ctree_2</span><span>)),</span><span>
</span><span>Time</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>1</span><span>:</span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>),</span><span> </span><span>2</span><span>),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>rep</span><span>(</span><span>c</span><span>(</span><span>&quot;Real&quot;</span><span>,</span><span> </span><span>&quot;CTREE&quot;</span><span>),</span><span> </span><span>each</span><span> </span><span>=</span><span> </span><span>length</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>)))</span><span>
</span><span>ggplot</span><span>(</span><span>datas</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Time</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.75</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>y</span><span> </span><span>=</span><span> </span><span>&quot;Detrended load&quot;</span><span>,</span><span> </span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Fitted values from CTREE&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-34" src="https://i1.wp.com/petolau.github.io/images/post_5/unnamed-chunk-34-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<p>And see the error of fitted values against real values.</p>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>matrix_train</span><span>$</span><span>Load</span><span>,</span><span> </span><span>predict</span><span>(</span><span>ctree_2</span><span>))</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 39.70532</code></pre>
</figure>
<p>It’s better. Now forecast values with <code>ctree_2</code> model.</p>
<figure>
<pre><code data-lang="r"><span>for_ctree</span><span> </span><span>&lt;-</span><span> </span><span>predict</span><span>(</span><span>ctree_2</span><span>,</span><span> </span><span>matrix_test</span><span>)</span><span> </span><span>+</span><span> </span><span>trend_for</span></code></pre>
</figure>
<p>And compare CTREE with RPART model.</p>
<figure>
<pre><code data-lang="r"><span>data_for</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_rpart</span><span>,</span><span> </span><span>for_ctree</span><span>),</span><span>
</span><span>Date</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>data_train</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>rep</span><span>(</span><span>data_test</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>3</span><span>)),</span><span>
</span><span>Type</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>rep</span><span>(</span><span>&quot;Train data&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;Test data&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;RPART&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>)),</span><span>
</span><span>rep</span><span>(</span><span>&quot;CTREE&quot;</span><span>,</span><span> </span><span>nrow</span><span>(</span><span>data_test</span><span>))))</span><span>
</span><span>ggplot</span><span>(</span><span>data_for</span><span>,</span><span> </span><span>aes</span><span>(</span><span>Date</span><span>,</span><span> </span><span>Load</span><span>,</span><span> </span><span>color</span><span> </span><span>=</span><span> </span><span>Type</span><span>,</span><span> </span><span>linetype</span><span> </span><span>=</span><span> </span><span>Type</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_line</span><span>(</span><span>size</span><span> </span><span>=</span><span> </span><span>0.8</span><span>,</span><span> </span><span>alpha</span><span> </span><span>=</span><span> </span><span>0.7</span><span>)</span><span> </span><span>+</span><span>
</span><span>facet_zoom</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>Date</span><span> </span><span>%in%</span><span> </span><span>data_test</span><span>$</span><span>date_time</span><span>,</span><span> </span><span>zoom.size</span><span> </span><span>=</span><span> </span><span>1.2</span><span>)</span><span> </span><span>+</span><span>
</span><span>labs</span><span>(</span><span>title</span><span> </span><span>=</span><span> </span><span>&quot;Forecasts from RPART and CTREE models&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>scale_linetype_manual</span><span>(</span><span>values</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>5</span><span>,</span><span>6</span><span>,</span><span>1</span><span>,</span><span>1</span><span>))</span><span> </span><span>+</span><span>
</span><span>theme_ts</span></code></pre>
</figure>
<p><img alt="plot of chunk unnamed-chunk-37" src="https://i0.wp.com/petolau.github.io/images/post_5/unnamed-chunk-37-1.png?w=456&ssl=1" data-recalc-dims="1"></p>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_rpart</span><span>)</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 3.727473</code></pre>
</figure>
<figure>
<pre><code data-lang="r"><span>mape</span><span>(</span><span>data_test</span><span>$</span><span>value</span><span>,</span><span> </span><span>for_ctree</span><span>)</span></code></pre>
</figure>
<figure>
<pre><code data-lang="text">## [1] 4.020834</code></pre>
</figure>
<p>Slightly better MAPE value with RPART, but again now it can not be anything to generalize. You can read more about CTREE method in its great <a rel="nofollow" href="https://cran.r-project.org/web/packages/party/vignettes/party.pdf" target="_blank">package vignette</a>.<br>
Try to forecast future values with all available electricity load data with sliding window approach (window of the length of three weeks) for a period of more than three months (98 days).</p>
<h2>Comparison</h2>
<p>Define functions that produce forecasts, so add up everything that we learned so far.</p>
<figure>
<pre><code data-lang="r"><span>RpartTrend</span><span> </span><span>&lt;-</span><span> </span><span>function</span><span>(</span><span>data</span><span>,</span><span> </span><span>set_of_date</span><span>,</span><span> </span><span>K</span><span>,</span><span> </span><span>period</span><span> </span><span>=</span><span> </span><span>48</span><span>){</span><span>
</span><span>data_train</span><span> </span><span>&lt;-</span><span> </span><span>data</span><span>[</span><span>date</span><span> </span><span>%in%</span><span> </span><span>set_of_date</span><span>]</span><span>
</span><span>N</span><span> </span><span>&lt;-</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)</span><span>
</span><span>window</span><span> </span><span>&lt;-</span><span> </span><span>(</span><span>N</span><span> </span><span>/</span><span> </span><span>period</span><span>)</span><span> </span><span>-</span><span> </span><span>1</span><span>
</span><span>data_ts</span><span> </span><span>&lt;-</span><span> </span><span>msts</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>seasonal.periods</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>period</span><span>,</span><span> </span><span>period</span><span>*</span><span>7</span><span>))</span><span>
</span><span>fuur</span><span> </span><span>&lt;-</span><span> </span><span>fourier</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>K</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>K</span><span>,</span><span> </span><span>K</span><span>))</span><span>
</span><span>fuur_test</span><span> </span><span>&lt;-</span><span> </span><span>as.data.frame</span><span>(</span><span>fourier</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>K</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>K</span><span>,</span><span> </span><span>K</span><span>),</span><span> </span><span>h</span><span> </span><span>=</span><span> </span><span>period</span><span>))</span><span>
</span><span>data_ts</span><span> </span><span>&lt;-</span><span> </span><span>ts</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>freq</span><span> </span><span>=</span><span> </span><span>period</span><span>*</span><span>7</span><span>)</span><span>
</span><span>decomp_ts</span><span> </span><span>&lt;-</span><span> </span><span>stl</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>s.window</span><span> </span><span>=</span><span> </span><span>&quot;periodic&quot;</span><span>,</span><span> </span><span>robust</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>
</span><span>new_load</span><span> </span><span>&lt;-</span><span> </span><span>rowSums</span><span>(</span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[,</span><span> </span><span>c</span><span>(</span><span>1</span><span>,</span><span>3</span><span>)])</span><span>
</span><span>trend_part</span><span> </span><span>&lt;-</span><span> </span><span>ts</span><span>(</span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[,</span><span>2</span><span>])</span><span>
</span><span>trend_fit</span><span> </span><span>&lt;-</span><span> </span><span>auto.arima</span><span>(</span><span>trend_part</span><span>)</span><span>
</span><span>trend_for</span><span> </span><span>&lt;-</span><span> </span><span>as.vector</span><span>(</span><span>forecast</span><span>(</span><span>trend_fit</span><span>,</span><span> </span><span>period</span><span>)</span><span>$</span><span>mean</span><span>)</span><span>
</span><span>lag_seas</span><span> </span><span>&lt;-</span><span> </span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[</span><span>1</span><span>:</span><span>(</span><span>period</span><span>*</span><span>window</span><span>),</span><span> </span><span>1</span><span>]</span><span>
</span><span>matrix_train</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>tail</span><span>(</span><span>new_load</span><span>,</span><span> </span><span>window</span><span>*</span><span>period</span><span>),</span><span>
</span><span>fuur</span><span>[(</span><span>period</span><span>+1</span><span>)</span><span>:</span><span>N</span><span>,],</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>lag_seas</span><span>)</span><span>
</span><span>tree_1</span><span> </span><span>&lt;-</span><span> </span><span>rpart</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train</span><span>,</span><span>
</span><span>control</span><span> </span><span>=</span><span> </span><span>rpart.control</span><span>(</span><span>minsplit</span><span> </span><span>=</span><span> </span><span>2</span><span>,</span><span>
</span><span>maxdepth</span><span> </span><span>=</span><span> </span><span>30</span><span>,</span><span>
</span><span>cp</span><span> </span><span>=</span><span> </span><span>0.000001</span><span>))</span><span>
</span><span>test_lag</span><span> </span><span>&lt;-</span><span> </span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[((</span><span>period</span><span>*</span><span>(</span><span>window</span><span>))</span><span>+1</span><span>)</span><span>:</span><span>N</span><span>,</span><span> </span><span>1</span><span>]</span><span>
</span><span>matrix_test</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>fuur_test</span><span>,</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>test_lag</span><span>)</span><span>
</span><span># prediction
</span><span> </span><span>pred_tree</span><span> </span><span>&lt;-</span><span> </span><span>predict</span><span>(</span><span>tree_1</span><span>,</span><span> </span><span>matrix_test</span><span>)</span><span> </span><span>+</span><span> </span><span>trend_for</span><span>
</span><span>return</span><span>(</span><span>as.vector</span><span>(</span><span>pred_tree</span><span>))</span><span>
</span><span>}</span><span>
</span><span>CtreeTrend</span><span> </span><span>&lt;-</span><span> </span><span>function</span><span>(</span><span>data</span><span>,</span><span> </span><span>set_of_date</span><span>,</span><span> </span><span>K</span><span>,</span><span> </span><span>period</span><span> </span><span>=</span><span> </span><span>48</span><span>){</span><span>
</span><span># subsetting the dataset by dates
</span><span> </span><span>data_train</span><span> </span><span>&lt;-</span><span> </span><span>data</span><span>[</span><span>date</span><span> </span><span>%in%</span><span> </span><span>set_of_date</span><span>]</span><span>
</span><span>N</span><span> </span><span>&lt;-</span><span> </span><span>nrow</span><span>(</span><span>data_train</span><span>)</span><span>
</span><span>window</span><span> </span><span>&lt;-</span><span> </span><span>(</span><span>N</span><span> </span><span>/</span><span> </span><span>period</span><span>)</span><span> </span><span>-</span><span> </span><span>1</span><span>
</span><span>data_ts</span><span> </span><span>&lt;-</span><span> </span><span>msts</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>seasonal.periods</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>period</span><span>,</span><span> </span><span>period</span><span>*</span><span>7</span><span>))</span><span>
</span><span>fuur</span><span> </span><span>&lt;-</span><span> </span><span>fourier</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>K</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>K</span><span>,</span><span> </span><span>K</span><span>))</span><span>
</span><span>fuur_test</span><span> </span><span>&lt;-</span><span> </span><span>as.data.frame</span><span>(</span><span>fourier</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>K</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>K</span><span>,</span><span> </span><span>K</span><span>),</span><span> </span><span>h</span><span> </span><span>=</span><span> </span><span>period</span><span>))</span><span>
</span><span>data_ts</span><span> </span><span>&lt;-</span><span> </span><span>ts</span><span>(</span><span>data_train</span><span>$</span><span>value</span><span>,</span><span> </span><span>freq</span><span> </span><span>=</span><span> </span><span>period</span><span>*</span><span>7</span><span>)</span><span>
</span><span>decomp_ts</span><span> </span><span>&lt;-</span><span> </span><span>stl</span><span>(</span><span>data_ts</span><span>,</span><span> </span><span>s.window</span><span> </span><span>=</span><span> </span><span>&quot;periodic&quot;</span><span>,</span><span> </span><span>robust</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>
</span><span>new_load</span><span> </span><span>&lt;-</span><span> </span><span>rowSums</span><span>(</span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[,</span><span> </span><span>c</span><span>(</span><span>1</span><span>,</span><span>3</span><span>)])</span><span>
</span><span>trend_part</span><span> </span><span>&lt;-</span><span> </span><span>ts</span><span>(</span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[,</span><span>2</span><span>])</span><span>
</span><span>trend_fit</span><span> </span><span>&lt;-</span><span> </span><span>auto.arima</span><span>(</span><span>trend_part</span><span>)</span><span>
</span><span>trend_for</span><span> </span><span>&lt;-</span><span> </span><span>as.vector</span><span>(</span><span>forecast</span><span>(</span><span>trend_fit</span><span>,</span><span> </span><span>period</span><span>)</span><span>$</span><span>mean</span><span>)</span><span>
</span><span>lag_seas</span><span> </span><span>&lt;-</span><span> </span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[</span><span>1</span><span>:</span><span>(</span><span>period</span><span>*</span><span>window</span><span>),</span><span> </span><span>1</span><span>]</span><span>
</span><span>matrix_train</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>Load</span><span> </span><span>=</span><span> </span><span>tail</span><span>(</span><span>new_load</span><span>,</span><span> </span><span>window</span><span>*</span><span>period</span><span>),</span><span>
</span><span>fuur</span><span>[(</span><span>period</span><span>+1</span><span>)</span><span>:</span><span>N</span><span>,],</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>lag_seas</span><span>)</span><span>
</span><span>tree_2</span><span> </span><span>&lt;-</span><span> </span><span>party</span><span>::</span><span>ctree</span><span>(</span><span>Load</span><span> </span><span>~</span><span> </span><span>.</span><span>,</span><span> </span><span>data</span><span> </span><span>=</span><span> </span><span>matrix_train</span><span>,</span><span>
</span><span>controls</span><span> </span><span>=</span><span> </span><span>party</span><span>::</span><span>ctree_control</span><span>(</span><span>teststat</span><span> </span><span>=</span><span> </span><span>&quot;quad&quot;</span><span>,</span><span>
</span><span>testtype</span><span> </span><span>=</span><span> </span><span>&quot;Teststatistic&quot;</span><span>,</span><span>
</span><span>mincriterion</span><span> </span><span>=</span><span> </span><span>0.925</span><span>,</span><span>
</span><span>minsplit</span><span> </span><span>=</span><span> </span><span>1</span><span>,</span><span>
</span><span>minbucket</span><span> </span><span>=</span><span> </span><span>1</span><span>))</span><span>
</span><span>test_lag</span><span> </span><span>&lt;-</span><span> </span><span>decomp_ts</span><span>$</span><span>time.series</span><span>[((</span><span>period</span><span>*</span><span>(</span><span>window</span><span>))</span><span>+1</span><span>)</span><span>:</span><span>N</span><span>,</span><span> </span><span>1</span><span>]</span><span>
</span><span>matrix_test</span><span> </span><span>&lt;-</span><span> </span><span>data.table</span><span>(</span><span>fuur_test</span><span>,</span><span>
</span><span>Lag</span><span> </span><span>=</span><span> </span><span>test_lag</span><span>)</span><span>
</span><span>pred_tree</span><span> </span><span>&lt;-</span><span> </span><span>predict</span><span>(</span><span>tree_2</span><span>,</span><span> </span><span>matrix_test</span><span>)</span><span> </span><span>+</span><span> </span><span>trend_for</span><span>
</span><span>return</span><span>(</span><span>as.vector</span><span>(</span><span>pred_tree</span><span>))</span><span>
</span><span>}</span></code></pre>
</figure>
<p>I created <code>plotly</code> boxplots graph of MAPE values from four models – CTREE simple, CTREE with detrending, RPART simple and RPART with detrending. Whole evaluation can be seen in the script that is stored in <a rel="nofollow" href="https://github.com/PetoLau/petolau.github.io/tree/master/_Rscripts" target="_blank">my GitHub repository</a>.</p>
<p><iframe width="450" src="https://plot.ly/~PetoLau2/56.embed"></iframe></p>
<p>We can see that <strong>detrending time series</strong> of electricity consumption improves the accuracy of the forecast with the combination of both <strong>regression tree</strong> methods – <strong>RPART</strong> and <strong>CTREE</strong>. My approach works as expected.</p>
<p>The habit of my posts is that animation must appear. So, I prepared for you two animations (animated dashboards) using <code>animation</code>, <code>grid</code>, <code>ggplot</code> and <code>ggforce</code> (for zooming) packages that visualize results of forecasting.</p>
<p><img alt="" src="https://i0.wp.com/petolau.github.io/images/post_5/trendRPART.gif?w=456&ssl=1" data-recalc-dims="1"></p>
<p><img alt="" src="https://i1.wp.com/petolau.github.io/images/post_5/trendCtree.gif?w=456&ssl=1" data-recalc-dims="1"></p>
<p>We can see that in many days it is almost perfect forecast, but on some days it has some potential for improving.</p>
<h2>Conclusion</h2>
<p>In this post, I showed you how to solve <strong>trend</strong> appearance in <strong>seasonal time series</strong> with using a <strong>regression tree</strong> model. <strong>Detrending time series</strong> for regression tree methods is a important (must) procedure due to the character of decision trees. The trend part of a time series was acquired by STL decomposition and separately forecasted by a simple ARIMA model. I evaluated this approach on the dataset from smart meters measurements of electricity consumption. The regression (decision) tree is a great technique for getting simple and interpretable results in very fast computational time.</p>
<p>In the future post, I will focus on enhancing the predictive performance of simple regression tree methods by ensemble learning methods like Bagging, Random Forest, and similar.</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://petolau.github.io/Regression-trees-for-forecasting-time-series-in-R/"> Peter Laurinec</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/n22byVOjaos" height="1"></a></div>