eBook for Financial Time Series Analysis with R

eBook - Financial Time Series Analysis with R

Data and R Code Files for Financial Time Series Analysis with R

Data and R Code

<p>Welcome to this course on financial time series analysis using R. In this course, we will learn about financial time series data analysis in R. You will learn about how to explore and build time series data, calculate its key statistics, and plot time series charts. You will also learn about how to use the important time series models such as White Noise, Random Walk, Autoregression and Moving Average, learn how to simulate these models, fit these models into financial time series data and use the models to predict the future. So, let's get started...</p>



<h3>What is Time Series?</h3>



<p>Time series refers to a series of data in a chronological order. A lot of data in this world is recorded sequentially, over time, in the form of time series. Some common examples include the weather in a city over time, the prices of a listed stock, the commodity prices and so on. While studying financial assets, the asset prices as well as asset returns are represented as time series. Investors generally prefer to use asset returns, over asset prices, in their analysis. This is primarily for two reasons: 1) the asset returns provide a complete and scale-free summary of asset returns and 2) the asset returns are easier to analyze compared to asset prices because of their statistical properties. In our lessons, we will use examples of both asset prices and asset returns.</p>



<h3>Financial Time Series Data</h3>



<p>There are plenty of financial time series data sources on the Internet. One popular source is Quandl, which contains thousands of datasets including financial and economic datasets.</p>



<p>For the purpose of this course, I would suggest you to signup for a free account on&nbsp;<a href="https://www.quandl.com/" target="_blank" rel="noreferrer noopener">Quandl.com</a>. Once you signup, you will get an API key that you can use to fetch data directly in R from Quandl.</p>



<h3>Install Quandl R Package</h3>



<p>In order to work with Quandl datasets in R, you need to install and load the Quandl package.</p>



<pre class="wp-block-preformatted">install.packages("Quandl")
library(Quandl)
</pre>



<h3>Load Sample Financial Data</h3>



<p>Now that you have your Quandl API key and have loaded Quandl package in your R session, you can download some financial data. In the following example, I am loading the daily closing prices of the Microsoft stock (Symbol: MSFT), for the year 2016 from the&nbsp;<code>WIKI/PRICES</code>&nbsp;dataset in Quandl.</p>



<p>In the following code, add your own api key and then run it in your R console. The code will create an R time series object containing Microsoft stock prices for the year 2016. If you don't understand what's happening in this code, don't worry about it. We will be studying in detail about how to create time series objects in the coming lessons.</p>



<pre class="wp-block-preformatted">Quandl.api_key("YOUR API KEY")
msft_data <- Quandl.datatable("WIKI/PRICES" ,
                           qopts.columns=c("date", "close"),
                           ticker=c("MSFT"),
                           date.gte=c("2016-01-01"),
                           date.lte=c("2016-12-31"))
ZOO <- zoo(msft_data$close, order.by=as.Date(as.character(msft_data$date), format='%Y-%m-%d'))
msft_ts <- ts(ZOO)
</pre>



<h3>View Time Series</h3>



<p>The easiest way to explore a time series is to view or print it as shown below:</p>



<pre class="wp-block-preformatted">< print(msft_ts)
Time Series:
Start = 1
End = 252
Frequency = 1
  [1] 54.800 55.050 54.050 52.170 52.330 52.300 52.780 51.640 53.110 50.990 50.560 50.790 50.480
 [14] 52.290 51.790 52.170 51.220 52.055 55.090 54.710 53.000 52.160 52.000 50.160 49.410 49.280
 [27] 49.710 49.690 50.500 51.090 52.420 52.190 51.820 52.650 51.180 51.360 52.100 51.300 50.880
 [40] 52.580 52.950 52.350 52.030 51.030 51.650 52.840 52.050 53.070 53.170 53.590 54.350 54.660
 [53] 53.490 53.860 54.070 53.970 54.210 53.540 54.710 55.050 55.230 55.570 55.430 54.560 55.120
 [66] 54.460 54.420 54.310 54.650 55.350 55.360 55.650 56.460 56.390 55.590 55.780 51.780 52.110
 ...
 ...
 ...
</pre>



<p>As you can see, the result explicitly tells us that it is a time series, which starts with 1 and ends at 252, i.e., there are a total of 252 observations with a frequency of 1.</p>



<p>In the next lesson, we will learn more about exploring the time series.</p>



<div class="wp-block-file"><a href="https://financetrain.com/wp-content/uploads/msft_data.csv">msft_data.csv</a><a href="https://financetrain.com/wp-content/uploads/msft_data.csv" class="wp-block-file__button" download>Download Data</a></div>



<pre class="wp-block-preformatted">install.packages("readr") 
library(readr) 
msft_data <- read_csv("msft_data.csv") 
View(msft_data)</pre>

Financial Time Series Data

<p>Let's look at a few commands that we will frequently use while exploring time series data.</p>



<h3>length()</h3>



<p>The&nbsp;<code>length()</code>&nbsp;function tells us the number of elements in out time series dataset.</p>



<pre class="wp-block-code"><code>> length(msft_ts)
[1] 252
></code></pre>



<h3>head()</h3>



<p>The&nbsp;<code>head()</code>&nbsp;function displays the top n elements of the dataset. This is useful while exploring large datasets.</p>



<pre class="wp-block-code"><code>> head(msft_ts,n=10)
 [1] 54.80 55.05 54.05 52.17 52.33 52.30 52.78 51.64 53.11 50.99
></code></pre>



<h3>tail()</h3>



<p>The&nbsp;<code>tail()</code>&nbsp;function displays the last n elements of the dataset. This is useful while exploring large datasets.</p>



<pre class="wp-block-code"><code>> tail(msft_ts,n=10)
 [1] 62.30 63.62 63.54 63.54 63.55 63.24 63.28 62.99 62.90 62.14
></code></pre>



<h3>Example 2 - Quarterly GDP Data</h3>



<p>Let's take one more example, this time quarterly GDP data. We will load this data from Quandl. The following command will load the quarterly GDP data from Quandl for the years 2014 to 2016</p>



<pre class="wp-block-code"><code>> GDP_data = Quandl("FRED/GDP", start_date="2014-01-01", end_date="2016-12-31",type="ts")</code></pre>



<p>Let's explore this data.</p>



<h3>print()</h3>



<p>We earlier saw that we can also use the&nbsp;<code>print()</code>&nbsp;function to display the time series.</p>



<pre class="wp-block-code"><code>> print(GDP_data)
        Qtr1    Qtr2    Qtr3    Qtr4
2014 17025.2 17285.6 17569.4 17692.2
2015 17783.6 17998.3 18141.9 18222.8
2016 18281.6 18450.1 18675.3 18869.4
></code></pre>



<p>As we can see, this data is presented as yearly data with 4 observations in each year.</p>



<h3>start() and end()</h3>



<pre class="wp-block-code"><code>> start(GDP_data)
[1] 2014    1
> end(GDP_data)
[1] 2016    4
></code></pre>



<h3>frequency()</h3>



<p>This tells us number of observations per unit of time.</p>



<pre class="wp-block-code"><code>> frequency(GDP_data)
[1] 4
></code></pre>



<p>The data is quarterly, so the lag between successive observations is 1 quarter. The dataset&nbsp;<code>GDP_data</code>&nbsp;has been set up so that the unit of time is 1 year (frequency=4).</p>



<h3>deltat()</h3>



<p>This is the fraction of the sampling period between successive observations; e.g., 1/12 for monthly data, and 1/4 for quarterly data. The function&nbsp;<code>deltat</code>&nbsp;uses this time unit to compute the lag by the formula ∆t= 1/frequency). Only one of&nbsp;<code>frequency</code>&nbsp;or&nbsp;<code>deltat</code>&nbsp;should be provided.</p>



<pre class="wp-block-code"><code>> deltat(GDP_data)
[1] 0.25
>
</code></pre>



<h3>time()</h3>



<p>The functions time and cycle create time series of the times at which the observations in a time series are taken and their "seasons".</p>



<pre class="wp-block-code"><code>> time(GDP_data)
        Qtr1    Qtr2    Qtr3    Qtr4
2014 2014.00 2014.25 2014.50 2014.75
2015 2015.00 2015.25 2015.50 2015.75
2016 2016.00 2016.25 2016.50 2016.75
></code></pre>



<pre class="wp-block-code"><code>> cycle(GDP_data)
     Qtr1 Qtr2 Qtr3 Qtr4
2014    1    2    3    4
2015    1    2    3    4
2016    1    2    3    4
></code></pre>

Exploring Time Series Data in R

<p>While we can explore time series data using commands such as&nbsp;<code>print()</code>,&nbsp;<code>head()</code>,&nbsp;<code>tail()</code>, etc in R, it can be very helpful to plot the time series data as a line chart and explore it visually.</p>



<p>In the following examples, we plot the Microsoft stock data and the quarterly GDP data in two different plots using the <code>plot()</code> function in R. You can alternatively also use <code>ts.plot()</code> function.</p>



<pre class="wp-block-code"><code>> plot(msft_ts,xlab="Day",ylab="Closing price", main="MSFT Closing Stock Prices")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ts-1.png" alt="" class="wp-image-26185"/></figure></div>



<pre class="wp-block-code"><code>> plot(GDP_data,xlab="Time",ylab="GDP", main="GDP Data",type="b")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ts-2.png" alt="" class="wp-image-26186"/></figure></div>



<p>Both the series show an upward trend over time. The trend for the MSFT stock is not very clear during the first half of the year, however during the send half, it is a clear uptrend. The GDP graph on hand, has an uptrend since the beginning.</p>

Plotting Time Series in R

<p>In the examples we saw earlier, we had good quality data with all values available for all time indexes. However, in real life, the data may contain missing values which will influence our analysis. Depending on the nature of data, we may choose to ignore missing values. However, in some cases it might be more suitable to estimate and fill the missing values. Data scientists use various techniques to estimate missing values. One common technique is to take the mean of the time series and replace NA with the mean value. Depending on the data this may or may not be suitable. For example, if the data is about loan borrowers and there are missing values in the loan interest rate, then the data scientist may decide to use the average interest rates for missing values, or if he clearly sees a pattern such as the interest rates being higher for self-employed individuals compared to salaries individuals, then the data scientist may decide to fill the missing values with means based on the categories based on their employment status.</p>



<h3>Handling Missing Values in R</h3>



<p>We will use our GDP data example to understand how we can estimate and fill missing values in R. Since we don't have access to a real dataset with missing values, we will make one. We will create a copy of our&nbsp;<code>GDP_data</code>&nbsp;dataset and then deliberately turn some values to NA. This is done by the following code:</p>



<pre class="wp-block-code"><code>> GDP_mod &lt;-GDP_data
> GDP_mod[7] &lt;- NA
> GDP_mod
        Qtr1    Qtr2    Qtr3    Qtr4
2014 17025.2 17285.6 17569.4 17692.2
2015 17783.6 17998.3      NA 18222.8
2016 18281.6 18450.1 18675.3 18869.4
></code></pre>



<p>As you can see, we have made the GDP value for 2015 Q3 blank.</p>



<h3>Plotting the Time Series</h3>



<p>We can plot the new time series using the&nbsp;<code>plot.ts()</code>&nbsp;function.</p>



<pre class="wp-block-code"><code>> plot.ts(GDP_mod)
></code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-missing-1.png" alt="" class="wp-image-26190"/></figure></div>



<p>As you can see the value of 2015 Q3 is missing.</p>



<h3>Calculate GDP Mean</h3>



<p>We will fill the missing value with the average of GDP over the three years. In R, we can calculate the mean using the&nbsp;<code>mean()</code>&nbsp;function. However, in this case, the&nbsp;<code>mean()</code>&nbsp;function will fail because of the presence of missing values. We can remove the missing values using the parameter&nbsp;<code>na.rm=TRUE</code>&nbsp;to calculate the mean by removing all missing values.</p>



<pre class="wp-block-code"><code>> #Calculating mean() will fail and return NA
> mean(GDP_mod)
[1] NA
> #Calculate mean by removing all missing values
> mean(GDP_mod,na.rm=TRUE)
[1] 17986.68
></code></pre>



<h3>Replace NA's with Mean</h3>



<p>Now that we have the mean of the series, we can replace the missing values with the mean value as shown below:</p>



<pre class="wp-block-code"><code>> #Replace missing values with mean
> GDP_mod[7] &lt;- mean(GDP_mod, na.rm = TRUE)
> print the series and notice that the missing value is now filled.
> print(GDP_mod)
         Qtr1     Qtr2     Qtr3     Qtr4
2014 17025.20 17285.60 17569.40 17692.20
2015 17783.60 17998.30 17986.68 18222.80
2016 18281.60 18450.10 18675.30 18869.40
></code></pre>



<h3>Plot Both Original and the Modified Series</h3>



<p>We will now plot both the original&nbsp;<code>GDP_data</code>&nbsp;and the modified&nbsp;<code>GDP_mod</code>&nbsp;time series to see how well the mean estimates the original value.</p>



<pre class="wp-block-code"><code>> plot(GDP_data)
> points(GDP_mod, type = "l", col = 2, lty = 3)
></code></pre>



<p>In the below chart, the black line represents the original data and the red dotted line represents the modified data. As we can see, the mean is not really a good indicator of the actual GDP in that quarter.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-missing-2.png" alt="" class="wp-image-26191"/></figure></div>

Handling Missing Values in Time Series

<p>In R, we can use the&nbsp;<code>ts()</code>&nbsp;function to create a time series object.</p>



<h3>Usage</h3>



<p>Below is a simplified format of the&nbsp;<code>ts</code>&nbsp;function. For complete details use&nbsp;<code>?ts</code>&nbsp;in your R console.</p>



<pre class="wp-block-code"><code>ts(data = NA, start = 1, end = numeric(), frequency = 1,  deltat = 1, names = )</code></pre>



<ul><li><strong>data:</strong>&nbsp;a vector or matrix of the observed time-series values.</li><li><strong>start:</strong>&nbsp;the time of the first observation.</li><li><strong>end:</strong>&nbsp;the time of the last observation, specified in the same way as start.</li><li><strong>frequency:</strong>&nbsp;the number of observations per unit of time.</li><li><strong>deltat:</strong>&nbsp;the fraction of the sampling period between successive observations; e.g., 1/12 for monthly data. Only one of frequency or deltat should be provided.</li><li><strong>names:</strong>&nbsp;a character vector of names for the series in a multiple series: defaults to the colnames of data, or Series 1, Series 2, ....</li></ul>



<p>While you can have data containing dates and corresponding values in an R object of any other class such as a dataframe, creating objects of ts class offers many benefits such as the time index information. Also, when you plot a ts object, it automatically creates a plot over time.</p>



<p>Let's look at how we can create a time series object using the&nbsp;<code>ts</code>&nbsp;class.</p>



<h3>Dataset</h3>



<p>For this example, we will use a dataset that consists of monthly prices of the S&amp;P Composite Index for the last three years starting from July 2014 to June 2017. We load the original data in a data vector called&nbsp;<code>sp_vector</code>.</p>



<pre class="wp-block-code"><code>sp_vector &lt;- c(1973.1,1961.53,1993.23,1937.27,2044.57,2054.27,2028.18,2082.2,2079.99,2094.86,2111.94,2099.29,2094.14,2039.87,1944.41,2024.81,2080.62,2054.08,1918.6,1904.42,2021.95,2075.54,2065.55,2083.89,2148.9,2170.95,2157.69,2143.02,2164.99,2246.63,2275.12,2329.91,2366.82,2359.31,2395.35,2433.99)</code></pre>



<h3>Print and Plot the Vector</h3>



<p>We can first try to directly print and plot the vector to observe the dataset.</p>



<pre class="wp-block-code"><code>> print(sp_vector)
 [1] 1973.10 1961.53 1993.23 1937.27 2044.57 2054.27 2028.18 2082.20 2079.99 2094.86 2111.94
[12] 2099.29 2094.14 2039.87 1944.41 2024.81 2080.62 2054.08 1918.60 1904.42 2021.95 2075.54
[23] 2065.55 2083.89 2148.90 2170.95 2157.69 2143.02 2164.99 2246.63 2275.12 2329.91 2366.82
[34] 2359.31 2395.35 2433.99
></code></pre>



<pre class="wp-block-code"><code>> plot(sp_vector)
></code></pre>



<p>As you can see, the plot does not contain any time index information.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/ts-ts-1.png" alt="" class="wp-image-26203"/></figure></div>



<h3>Creating Time Series Object</h3>



<p>We can use the&nbsp;<code>ts()</code>&nbsp;function to convert this vector into a time series object. Our data set contains monthly stock prices from July 2017 to June 2017. There are 12 observations per year starting from July 2014. So, we will set the frequency to 12 and the set start argument to c(2014,7) where the first value is the year and the second value is the month index (July).</p>



<pre class="wp-block-code"><code>> sp_ts &lt;- ts(sp_vector,start=c(2014,7),frequency=12)</code></pre>



<p>We can print to check that we have the correct time series.</p>



<pre class="wp-block-code"><code>> print(sp_ts)
         Jan     Feb     Mar     Apr     May     Jun     Jul     Aug     Sep     Oct     Nov     Dec
2014                                                 1973.10 1961.53 1993.23 1937.27 2044.57 2054.27
2015 2028.18 2082.20 2079.99 2094.86 2111.94 2099.29 2094.14 2039.87 1944.41 2024.81 2080.62 2054.08
2016 1918.60 1904.42 2021.95 2075.54 2065.55 2083.89 2148.90 2170.95 2157.69 2143.02 2164.99 2246.63
2017 2275.12 2329.91 2366.82 2359.31 2395.35 2433.99
>
</code></pre>



<h3>Plot the Time Series</h3>



<p>We can now use the&nbsp;<code>plot.ts()</code>&nbsp;function to plot the time series.</p>



<pre class="wp-block-code"><code>> plot.ts(sp_ts)
</code></pre>



<p>Note that it automatically recognized that it is a time series and plots time on the x-axis.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/ts-ts-2.png" alt="" class="wp-image-26204"/></figure></div>

Creating a Time Series Object in R

<p>In R, objects can be of different class such as vector, list, dataframe, ts, etc. When you load a dataset into R, it may not necessarily be a time series object. We can use the&nbsp;<code>is.ts()</code>&nbsp;function to test if the given object is a time series (ts) object or not.</p>



<p>Below we perform this test on all the objects we've created so far in this course.</p>



<pre class="wp-block-code"><code>> is.ts(msft_ts)
[1] TRUE
> is.ts(GDP_data)
[1] TRUE
> is.ts(sp_vector)
[1] FALSE
> is.ts(sp_ts)
[1] TRUE
></code></pre>



<p>As you can see, all objects return TRUE except&nbsp;<code>sp_vector</code>&nbsp;which is a vector.</p>

Check if an object is a time series object in R

<p>Let's take one more example of plotting financial time series data. This time we will use the&nbsp;<code>EuStockMarkets</code>&nbsp;data set that comes by default with R. It contains the daily closing prices of major European stock indices from 1991 to 1998.</p>



<h3>Check and Print the Data</h3>



<p>Let's first check if the data is a time series and then print a few values from the dataset to get a feel of the data.</p>



<pre class="wp-block-code"><code>> is.ts(EuStockMarkets)
[1] TRUE
> head(EuStockMarkets)
         DAX    SMI    CAC   FTSE
[1,] 1628.75 1678.1 1772.8 2443.6
[2,] 1613.63 1688.5 1750.5 2460.2
[3,] 1606.51 1678.6 1718.0 2448.2
[4,] 1621.04 1684.1 1708.1 2470.4
[5,] 1618.16 1686.6 1723.1 2484.7
[6,] 1610.61 1671.6 1714.3 2466.8
></code></pre>



<p>As you can see, the dataset is a time series object and contains prices for 4 indices, namely, DAX, SMI, CAC and FTSE.</p>



<h3>Plot the Data</h3>



<p>We can plot the data using the&nbsp;<code>plot()</code>&nbsp;function to create a simple plot which will print all 4 series as shown below:</p>



<pre class="wp-block-code"><code>> plot(EuStockMarkets)</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-tss-1.png" alt="" class="wp-image-26217"/></figure></div>



<p>Alternatively we can use the&nbsp;<code>ts.plot()</code>&nbsp;function to get the same plot but with a common y-axis. This time we are also adding a legend to the chart to make it more user-friendly.</p>



<pre class="wp-block-code"><code>> ts.plot(EuStockMarkets, col = 1:4, xlab = "Year", ylab = "Prices", main = "European Stock Indices")
> legend("topleft", colnames(EuStockMarkets), lty = 1, col = 1:4, bty = "n")
></code></pre>



<p>Notes: <code>col</code> sets the colors of lines, <code>lty</code> is the line type, and <code>bty</code> is for the type of box to be drawn around the legend. <code>n</code> represents no box.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-tss-2.png" alt="" class="wp-image-26218"/></figure></div>

Plotting Financial Time Series Data (Multiple Columns) in R

<p>Time series have several characteristics that make their analysis different from other types of data.</p>



<ul><li>The time series variable (for example, the stock price) may have a&nbsp;<strong>trend</strong>&nbsp;over time. This refers to the increasing or decreasing values in a given time series.</li><li>The variable may exhibit cyclicity or&nbsp;<strong>seasonality</strong>. This refers to the repeating cycle over a specific period (such as week, month, etc.) in the time series. A seasonal pattern exists when a series is influenced by seasonal factors (e.g., the quarter of the year, the month, or day of the week). Seasonality is always of a fixed and known period. As an example, if you study the quarterly earnings of a company such as John Deere, it will exhibit seasonality with the earnings peaking during the harvest season.</li><li>The data will have&nbsp;<strong>serial correlation</strong>&nbsp;between subsequent observations.</li><li>The data will almost always have an&nbsp;<strong>irregular component</strong>, which is referred to as the White Noise. This is the random variation not explained by any other factor. White Noise is a stationary process with a constant mean and variance.</li></ul>



<p>A cyclic pattern exists when data exhibit rises and falls that are&nbsp;<em>not of fixed period</em>. The duration of these fluctuations is usually of at least 2 years. Many people confuse cyclic behavior with seasonal behavior, but they are really quite different. If the fluctuations are not of fixed period then they are cyclic; if the period is unchanging and associated with some aspect of the calendar, then the pattern is seasonal. In general, the average length of cycles is longer than the length of a seasonal pattern, and the magnitude of cycles tends to be more variable than the magnitude of seasonal patterns.</p>



<p>Due to these factors, the time series data can exhibit a huge variety of patterns and it is therefore helpful to decompose the time series into its several components each representing an underlying pattern, namely, trend, cycle, seasonality, and random (stationary). Decomposing a time series into these components is an important part of time series analysis because this is the primary assumption to forecast the time series in the future. This is generally done to help understand the time series better, but it can also be used to improve forecasts.</p>



<h3>Examples</h3>



<p>The following two charts show examples of trend and seasonality. The first chart plots the quarterly GDP data which shows a clear uptrend over time. The second chart shows the quarterly earnings for John Deere which exhibits seasonality in earnings within a year.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-trends-1.png" alt="" class="wp-image-26221"/></figure></div>



<p>In the following lessons we will learn about how we can use time series analysis to develop mathematical models that provide plausible descriptions for sample data. Before we get to the models, we will learn about a few important concepts such as the stationary process, autocorrelation, and decomposition of time series in more detail.</p>



<h3>Exercise</h3>



<p>The following vector contains the quarterly John Deere earnings for the period Jan 2011 to July 2017.</p>



<pre class="wp-block-code"><code>c(1.87,1.70,0.50,0.36,0.95,1.46,0.71,0.74,1.47,1.57,0.83,1.58,2.19,2.40,1.51,1.89,2.16,2.74,1.41,1.88,2.32,2.55,1.25,1.44,1.68,2.06,1.02)</code></pre>



<p>The data is in reverse order, i.e., the first observation (1.87) is for July 2017, the second observation 1.70 is for Apr 2017 and so on.</p>



<ol><li>Load this vector into a variable in your R session.</li><li>Reverse the vector using the&nbsp;<code>rev()</code>&nbsp;function so that the resulting data is in increasing chronological order.</li><li>Convert this into time series data using&nbsp;<code>ts()</code>&nbsp;function.</li><li>Check that the resulting data is of the class ts using the&nbsp;<code>is.ts()</code>&nbsp;function.</li><li>Use the&nbsp;<code>plot.ts()</code>&nbsp;function to plot the data as time series.</li></ol>

Characteristics of Time Series

<p>A common assumption made in time series analysis is that one of the components of the pattern exhibited by a time series is the stationary series. This is the random or irregular component we discussed earlier. This random variation is not explained by any other factor.</p>



<p>There are three characteristics of a stationary series:</p>



<ol><li>It has a constant mean over time, i.e., the mean does not change with the passage of time.</li><li>It has a constant variance over time, i.e., the variance/volatility in the data does not change over time.</li><li>It's autocorrelation remains the same over time. Autocorrelation is a situation in which time series data is influenced by its own historical values. We will learn about this shortly.</li></ol>



<p>The following are all examples of stationary series. They all have constant statistical properties such as mean, variation, autocorrelation, etc.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/stationary.png" alt="" class="wp-image-26225"/></figure></div>



<h3>Why Make Series Stationary?</h3>



<p>While performing times series analysis, data scientists will generally make the series stationary. There are two compelling reasons to stationarize a series:</p>



<ol><li>Most statistical forecasting models are based on the assumption that the underlying time series can be transformed to a stationary series using mathematical transformations. It is relatively easy to make prediction on a stationary series - the idea being that you can assume that its statistical properties will remain the same in the future as in the past! Once the prediction has been made with the stationary series, we need to untransform the series, that is, we reverse the mathematical transformations we applied to it in the past. This gives us the prediction of the original series. This is generally taken care of by the R software. Therefore for making a good prediction, it is important to find the most suitable sequence of transformations needed to stationarize a time series. This provides important clues in our search for the right forecasting model. Making a series stationary is an important part of fitting an ARIMA model, which we will be learning.</li><li>Stationarizing a time series also helps us get important sample statistics such as means, variances, and correlations with other variables. As described above, these statistics are useful to predict the future behavior only if the series is stationary. For example, if the series is consistently increasing over time, the sample mean and variance will grow with the size of the sample, and they will always underestimate the mean and variance in future periods. And if the mean and variance of a series are not well-defined, then neither are its correlations with other variables. For this reason you should be cautious about trying to extrapolate regression models fitted to nonstationary data.</li></ol>



<p><strong>What is meant by series is consistently increasing over time?</strong></p>



<p>This was an example of a non-stationary series for which both the mean and covariance are increasing over time. You would observe two trends in such a time series. 1) There will be a clear upward trend, i.e., the mean is increasing over time. 2) the variability of the data is increasing over time. As an example, think of a stock whose price is consistently increasing over time and the volatility of the stock price is also going up with the passage of time.</p>

Stationary Process in Time Series

<p>Most financial and economic times series are not stationary. Even when you adjust them for seasonal variations, they will exhibit trends, cycles, random walk and other non-stationary behavior. We can use a variety of techniques to make a non-stationary series stationary depending on the kind of non-stationary behavior present in the series. The two techniques we will learn are called differencing and logarithmic transformations.</p>



<h3>Differencing</h3>



<p>Under this technique, we difference the data. That is, given the series Z(t), we create the new series:</p>



<p><strong>Y(i) = Z(i) - Z(i-1)</strong></p>



<p>The differenced data will contain one less point than the original data. Usually, one differencing is sufficient to stationarize the data. However, you can difference the data more than once, if needed. In R, differencing is done using the&nbsp;<code>diff()</code>&nbsp;function.</p>



<p>Differencing a time series can remove a linear trend from it. In finance or stock markets, a series of asset returns or stock returns is the differenced time series which is calculated by taking differences of prices on consecutive time intervals.</p>



<h3>Seasonal Differencing</h3>



<p>If a series has seasonality present in it, then we can use seasonal differencing to remove these periodic patterns. For monthly data, in which there are 12 periods in a season, the seasonal difference of Y at period t is&nbsp;<code>Y(t) - Y(t-12)</code>. for quarterly data, the difference will be based on a lag of 4 data points.</p>



<h3>Log Transformation</h3>



<p>Log transformation can be used to stabilize the variance of a series with non-constant variance. This is done using the&nbsp;<code>log()</code>&nbsp;function. One limitation of log transformation is that it can be applied only to positively valued time series. Taking a log shrinks the values towards 0. For values that are close to 1, the shrinking is less and for the values that are higher, the shrinking is more, thus reducing the variance.</p>



<p>For negative data, you can add a suitable constant to make all the data positive before applying the transformation. This constant can then be subtracted from the model to obtain predicted (i.e., the fitted) values and forecasts for future points.</p>



<p>This also makes the case for why usage of log returns is popular in finance, rather than the prices or raw returns.</p>



<p>In the next lesson, we will learn about how we can apply these transformations.</p>

Transforming a Series to Stationary

<p>We will now learn about how we can perform the mathematical transformations in R in order to make a non-stationary series stationary.</p>



<p>We have the daily stock prices of an imaginary stock exhibiting rapid growth stored in a csv file. You can download the csv file and follow along this tutorial to perform the transformations in your R console. There are two hundred observations starting from 1st Jan , 2016.</p>



<div class="wp-block-file"><a href="https://financetrain.com/wp-content/uploads/stock_daily.csv">stock_daily.csv</a><a href="https://financetrain.com/wp-content/uploads/stock_daily.csv" class="wp-block-file__button" download>Download</a></div>



<p>We have placed the csv file in our working directory as shown below:</p>



<pre class="wp-block-code"><code>> setwd("C:/Users/Manish/Dropbox/Finance Train/Courses/Data")
> getwd()
[1] "C:/Users/Manish/Dropbox/Finance Train/Courses/Data"
> list.files()
 [1] "stock_daily.csv"</code></pre>



<p>As you can see, the file&nbsp;<code>stock_daily.csv</code>&nbsp;in available in our working directory.</p>



<h3>Load Data in R</h3>



<p>We will read this csv file in R using&nbsp;<code>read.csv()</code>&nbsp;function and store the dataset in a variable called&nbsp;<code>sp</code>. The data will be stored as a dataframe.</p>



<pre class="wp-block-code"><code>> sp &lt;- read.csv("stock_daily.csv", header=FALSE)
>
> str(sp)
'data.frame':    200 obs. of  1 variable:
 $ V1: num  254 229 277 262 264 ...
></code></pre>



<h3>Convert dataframe to time series using&nbsp;<code>ts()</code>&nbsp;function</h3>



<p>Now that we have the data in a dataframe, we can convert it into a time series object using the&nbsp;<code>ts()</code>&nbsp;function as shown below:</p>



<pre class="wp-block-code"><code>> sp_ts &lt;- ts(sp, start=2016, frequency=365)</code></pre>



<h3>Plot the series</h3>



<p>Let's now plot the time series using the&nbsp;<code>plot.ts()</code>&nbsp;function.</p>



<pre class="wp-block-code"><code>> plot.ts(sp_ts,col=4, main="Daily Stock Prices", ylab="Price")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-trans-1.png" alt="" class="wp-image-26230"/></figure></div>



<p>We can observe two things here:</p>



<ol><li>There is a clearly visible upward trend in the data</li><li>The variability in the data is increasing over time.</li></ol>



<p>Let's now look at how we can remove both these patterns to stationarize the series.</p>

Time Series Transformation in R

<h3>Removing Variability Using Logarithmic Transformation</h3>



<p>Since the data shows changing variance over time, the first thing we will do is stabilize the variance by applying log transformation using the&nbsp;<code>log()</code>&nbsp;function. The resulting series will be a linear time series.</p>



<pre class="wp-block-code"><code>> sp_linear&lt;-log(sp_ts)
> plot.ts(sp_linear, main="Daily Stock Prices (log)", ylab="Price", col=4)</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-transs-1.png" alt="" class="wp-image-26234"/></figure></div>



<h3>Removing Linear Trend</h3>



<p>We will now perform the first difference transformation&nbsp;<code>[z(t) - z(t-1)]</code>&nbsp;to our series to remove the linear trend.</p>



<pre class="wp-block-code"><code>> sp_linear_diff &lt;- diff(sp_linear)
> plot.ts(sp_linear_diff, main="Daily Stock Prices (log)", ylab="Price", col=4)
></code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-transs-2.png" alt="" class="wp-image-26235"/></figure></div>



<h3>Removing Seasonal Differencing</h3>



<p>Let's take another example to understand how we can use the&nbsp;<code>diff()</code>&nbsp;function to remove seasonal differencing from data. We will use the John Deer's Quarterly earnings data we used earlier as it exhibits seasonality.</p>



<pre class="wp-block-code"><code>> par(mfrow = c(1,2))
> de_earnings_diff &lt;- diff(johndeere_earnings,lag=4)
> plot.ts(johndeere_earnings, main="Earnings (Quarterly)")
> plot.ts(de_earnings_diff, main="Earnings (Differenced, lag=4)")
</code></pre>



<p>The chart on the left shows the original earnings. The chart on the right shows the difference in earnings with a lag of 4.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-transs-3.png" alt="" class="wp-image-26236"/></figure></div>

Differencing and Log Transformation

<p>Autocorrelation is an important part of time series analysis. It helps us understand how each observation in a time series is related to its recent past observations. When autocorrelation is high in a time series, it becomes easy to predict their future observations.</p>



<p>Let us consider the Microsoft stock prices for the year 2016, which we have as a time series object&nbsp;<code>msft_ts</code>. Below we have plotted the stock prices in the left chart and a scatter plot of the stock prices with a lag of 1 on the right hand side. We can clearly see a strong positive correlation between the two.</p>



<pre class="wp-block-code"><code>> plot.ts(msft_ts,main="MSFT Stock Prices",ylab="Prices")
> plot(msft_ts[-252],msft_ts[-1],main="Scatterplot (lag=1)")
> abline(lm(msft_ts[-1] ~ msft_ts[-252]),col=4)
></code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-auto-1.png" alt="" class="wp-image-26241"/></figure></div>



<p>We can also calculate the correlation between the actual series and the lagged series using the&nbsp;<code>cor()</code>&nbsp;function in R.</p>



<pre class="wp-block-code"><code>> #Correlation of stock price today and 1 day earlier
> cor(msft_ts[-252],msft_ts[-1])
[1] 0.9797061
></code></pre>



<p>Similarly we can calculate the correlation with different lags, such as between the stock price today and the stock price two days earlier. In our example, the correlation between these two pairs is 0.9623, which is still quite high but less compared to the correlation with a time lag of 1.</p>



<pre class="wp-block-code"><code>> cor(msft_ts[-(251:252)],msft_ts[-(1:2)])
[1] 0.9623822
></code></pre>



<h3>ACF Function</h3>



<p>Instead of calculating the correlation with each time lag manually, we can use the&nbsp;<code>acf()</code>&nbsp;function in R. The function acf computes (and by default plots) estimates of the autocovariance or autocorrelation function for different time lags.</p>



<p>Below we get the autocorrelations for lag 1 to 10. Notice that the correlation keeps reducing as the lag increases.</p>



<pre class="wp-block-code"><code>> acf(msft_ts, lag.max=10, plot=FALSE)</code></pre>



<p>Autocorrelations of series&nbsp;<code>msft_ts</code>, by lag:</p>



<pre class="wp-block-preformatted">    0     1     2     3     4     5     6     7     8     9    10
1.000 0.973 0.948 0.922 0.897 0.872 0.849 0.826 0.805 0.787 0.772
></pre>



<p>We can also get the same information in an acf plot as shown below. This is also called a correlogram, also known as an autocorrelation plot.</p>



<pre class="wp-block-code"><code>> acf(msft_ts, lag.max=10)</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-acf-1.png" alt="" class="wp-image-26242"/></figure></div>



<p>The x-axis donates the time lag, while the y-axis displays the estimated autocorrelation. Looking at this data, we can say that each observation is positively related to its recent past observations. However, the correlation decreases as the lag increases.</p>



<h3>Exercise</h3>



<p>Provided below is a csv file that contains the daily stock prices of 5 US stocks for 251 days. The five stocks are AAPL, MSFT, GOOG, IBM, and INTC.</p>



<ol><li>Load the data in R in a variable called&nbsp;<code>stock_data</code>.</li><li>Extract only the AAPL stock prices in another variable called&nbsp;<code>aapl_prices</code>.</li><li>Convert the&nbsp;<code>aapl_prices</code>&nbsp;into a time series using the&nbsp;<code>ts()</code>&nbsp;function. Store this in a variable called&nbsp;<code>aapl_prices_ts</code>.</li><li>Calculate the autocorrelation in&nbsp;<code>aapl_prices_ts</code>&nbsp;with 1 and 2 lags using the&nbsp;<code>cor()</code>&nbsp;function.</li><li>Use the acf function to find the autocorrelations in the&nbsp;<code>aapl_prices_ts</code>&nbsp;with 1 to 10 time lags.</li></ol>



<div class="wp-block-file"><a href="https://financetrain.com/wp-content/uploads/stock_data.csv">stock_data.csv</a><a href="https://financetrain.com/wp-content/uploads/stock_data.csv" class="wp-block-file__button" download>Download</a></div>



<p><strong>Why there is difference in the values of correlation coefficients with that of values given by  acf. Lag 2 gives a value 0.962 while acf shows 0.948</strong></p>



<p>The two estimates differ slightly as they use slightly different scalings in their calculation of sample covariance, 1/(n-l) in case of cor() versus 1/n in case of acf(). Even though the acf() method provides a biased estimate, it is preferred in time series analysis. The autocorrelation estimates differ by a factor of (n-l)/n.</p>



<p></p>

Autocorrelation in R

<p>By now we have a strong foundational understanding of various concepts essential for time series analysis. The rest of the course will focus on the following:</p>



<ol><li>A theoretical understanding of the important time series models (White Noise, AutoRegressive (AR), Moving Average (MA), ARMA.</li><li>The ARIMA model and how various time series processes can be explained by ARIMA.</li><li>Simulating and estimating these time series models in R.</li><li>Box-Jenkins (B-J) methodology for time series forecasting.</li><li>A comprehensive case study for time series analysis in R (With code and plots).</li></ol>



<h3>Time Series Models</h3>



<p>The objective of the following text is to provide a theoretical understanding of the time series models. It can be a bit complex to grasp, however, as we move to practical implementation in the next lessons, things will start to make more sense.</p>



<h3>White Noise</h3>



<p>White Noise is the simplest example of a stationary process and is also the foundation of the other models we will discuss. A white noise process is a random process of random variables that are uncorrelated, have mean zero, and a finite variance. In other words, a series is called white noise if it is purely random in nature. Random noise is donated by ε<sub>t</sub></p>



<p>Plots of white noise series exhibit a very erratic, jumpy, unpredictable behavior. Since the ε<sub>t</sub>&nbsp;are uncorrelated, previous values do not help us to forecast future values. White noise series themselves are quite uninteresting from a forecasting standpoint (they are not linearly forecastable), but they form the building blocks for more general models.</p>



<h3>Simple Time Series Models</h3>



<p>A simple trend model can be expressed as follows:</p>



<p><strong>y<sub>t</sub>&nbsp;= b<sub>0</sub>&nbsp;+ b<sub>1</sub>t+ ε</strong><sub><strong>t</strong></sub></p>



<ul><li>b<sub>0</sub>&nbsp;= the y-intercept; where t = 0.</li><li>b<sub>1</sub>&nbsp;= the slope coefficient of the time trend.</li><li>t = the time period.</li><li>ŷ<sub>t</sub>&nbsp;= the estimated value for time t based on the model.</li><li>ε<sub>t</sub>&nbsp;= the random error of the time trend.</li></ul>



<p>The big validity pit-fall for simple trend models is serial correlation; if this problem is present, then you will see an artificially high R<sup>2</sup>&nbsp;and your slope coefficient may falsely appear to be significant.</p>



<p>There is a visual way to detect serial correlation (<em>correlogram</em>) or you can perform a Dubin-Watson test.</p>



<h3>Autoregressive Models</h3>



<p>This type of time series model utilizes a time period lagged observation as the independent variable to predict the dependent variable, which is the value in the next time period.</p>



<p><strong>x<sub>t</sub>&nbsp;= b<sub>0</sub>&nbsp;+ b<sub>1</sub>x<sub>t-1</sub>&nbsp;+ ε</strong><sub><strong>t</strong></sub></p>



<p>There can be more than one time period lag independent variables. The model assumes that Xt depends only on its own past values, such as p past values - AR(p). So, AR(1) is the first order autoregression meaning that the current value is based on the immediately preceding value. An AR(2) process has the current value based on the previous two values.</p>



<h3>Random Walk</h3>



<p>A random walk is an example of a non-stationary process. It has no specific mean or variance, however, it exhibits very strong dependence over time, with each observation closely related to its immediate past observations. However, the changes in value or increments follow a white noise process which is stationary.</p>



<p>Random walk is the case of an AR time series model where the predicted value is expected to equal the previous period plus a random error:</p>



<p><strong>x<sub>t</sub>&nbsp;= b<sub>0</sub>&nbsp;+ x<sub>t-1</sub>&nbsp;+ ε</strong><sub><strong>t</strong></sub></p>



<p>When b<sub>0</sub>&nbsp;is not equal to zero, the model is a random walk with a drift, but the key characteristic is a b<sub>1</sub>&nbsp;= 1.</p>



<p>ε<sub>t</sub> is the mean zero white noise. The expected value of the error is still zero.</p>



<p>Some important points about Random Walk Process (Advanced material. Not mandatory):</p>



<ul><li>The mean reverting level for a random walk is not covariance stationary and the technique of first differencing is frequently used to transform an AR model with one time lag variable (AR1) into a model that is covariance stationary.</li><li>If an AR time series is covariance stationary, then the serial correlations for the lag variables are insignificant or they rapidly drop to zero as the number of time period lags rises.</li><li>When the lag coefficient is not statistically different from 1, a unit root exists. Dickey-Fuller test is applied to AR1 model to test for a unit root. If a unit root is present, then the model is not covariance stationary; if this is the case, the independent variable must be transformed, so you can re-model.</li></ul>



<h3>When to Use Autoregressive Models?</h3>



<p>If serial correlation exists in a simple time series model, the analyst can create an auto-regressive time series with the sample data, where the independent variable is a lagged (prior period) value. The AR model is appropriate where the prior period value is the best predictor for the future period dependent variable value.</p>



<h3>Moving Average Models</h3>



<p>A moving average model is one when x<sub>t</sub>&nbsp;depends only on the random error terms which follow a white noise process. A simple moving average is a series x<sub>t</sub>&nbsp;generated from a white noise series ε<sub>t</sub>&nbsp;by the rule:</p>



<p><strong>x<sub>t</sub>&nbsp;= ε<sub>t</sub>&nbsp;+ βε</strong><sub><strong>t-1</strong></sub></p>



<p>Rather than use past values of the forecast variable in a regression, a moving average model uses past forecast errors in a regression-like model.</p>



<p>The above equation is an example of a first order moving average model which is written as MA(1), as we are considering lag only for 1 past value. The general form is MA(q), where xt depends on q past values.</p>



<p>Note that unless&nbsp;<strong>β = 0</strong>, we will have a nontrivial correlation structure.</p>

Time Series Models

<p>If we combine differencing with autoregression and a moving average model, we obtain a non-seasonal ARIMA model. ARIMA is an acronym for&nbsp;<strong><a href="https://financetrain.com/arma-models-and-arch-testing/">AutoRegressive Integrated Moving Average</a></strong>&nbsp;model. The term "integration" in this context is the reverse of differencing.</p>



<p>ARIMA model is represented as&nbsp;<strong>ARIMA(p,d,q)</strong></p>



<p>Where:</p>



<ul><li>p = order of the autoregressive part</li><li>d = degree of first differencing involved</li><li>q = order of the moving average part</li></ul>



<p>A series which is stationary after being differentiated once is said to be integrated of order 1 and is denoted by I(1). In general a series which is <a href="https://financetrain.com/stationary-process-in-time-series/">stationary</a> after being differentiated d times is said to be integrated of order d, denoted I(d). A series which is stationary without differencing is said to be I(0).</p>



<p>While using ARIMA modeling for time series forecasting, selecting appropriate values for p, d and q can be difficult. We make use of ACF, PACF and other functions to help us determine these values.</p>



<p>Many of the models we previously discussed can be easily explained by the ARIMA model as shown below:</p>



<ul><li>White noise: ARIMA(0,0,0)</li><li>Random walk: ARIMA(0,1,0) with no constant</li><li>Random walk with drift: ARIMA(0,1,0) with a constant</li><li>Autoregression: ARIMA(p,0,0)</li><li>Moving average: ARIMA(0,0,q)</li></ul>



<h3>ARIMA Modeling in R</h3>



<p>In the next few lessons, we will learn about how we can simulate and estimate these models in R using the ARIMA class of functions.</p>

ARIMA Modeling

<p>The function <code>arima.sim()</code> can be used to simulate data from a variety of <a href="https://financetrain.com/time-series-models/">time series models</a>. Based on the model we want to apply, we specify the appropriate values for p, d and q to the model <code>ARIMA(p,d,q)</code>.</p>



<p>The general format for the&nbsp;<code>arima.sim()</code>&nbsp;function is as follows:</p>



<pre class="wp-block-code"><code>arima.sim(model, n)
model is where we specify the AR, differencing order and MA components.
n is the length of output series.</code></pre>



<p>For the White Noise model, all p, d and q in arima model are 0. So,&nbsp;<code>ARIMA(0,0,0)</code>&nbsp;is simply the White Noise(WN) model.</p>



<h3>Simulate White Noise Model in R</h3>



<p>To simulate WN model in R, we will set all, p, d and q to 0. To generate 200 observation series, we will set the n argument to 200.</p>



<pre class="wp-block-code"><code>WN &lt;- arima.sim(model = list(order = c(0, 0, 0)), n = 200)</code></pre>



<p>This will create a time series object that follows White Noise model.</p>



<p>You can plot the newly generated time series using the&nbsp;<code>plot.ts()</code>&nbsp;function.</p>



<pre class="wp-block-code"><code>> plot.ts(WN,col=4, main="White Noise Series")
</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/white-noise-1.png" alt="" class="wp-image-26251"/></figure></div>



<p>You can calculate the mean and <a href="https://financetrain.com/variance-standard-deviation/">standard deviation</a> of this series and notice that the series will have a mean close to 0 and a standard deviation close to 1.</p>



<pre class="wp-block-code"><code>> mean(WN)
[1] -0.08357315
> var(WN)
[1] 0.91063
> sd(WN)
[1] 0.9542694</code></pre>



<p>While generating the White Noise series, we can also specify the mean and standard deviation. In the following example, we simulate White Noise model with mean=50 and sd=10.</p>



<pre class="wp-block-code"><code>> WN_2 &lt;- arima.sim(model = list(order = c(0, 0, 0)), n = 200, mean=50, sd=10)
> plot.ts(WN_2,col=4, main="White Noise Series (mean=50, sd=10)")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/white-noise-2.png" alt="" class="wp-image-26252"/></figure></div>



<p>Again, we can verify that the statistics for this time series are close to what we fed to it.</p>



<pre class="wp-block-code"><code>> mean(WN_2)
[1] 50.64454
> sd(WN_2)
[1] 10.0867
></code></pre>



<h3>Estimating the White Noise Model (Model Fitting)</h3>



<p>If we have a series (where we feel White Noise model is suitable), then we can fit the White Noise model using the <code><a href="https://financetrain.com/arima-modeling/">arima()</a></code> function. When applied, the <code>arima()</code> function returns the important information about the estimated model, such as its estimated mean and variance.</p>



<p>We will use the White Noise series, WN, that we created above to fit the model to it.</p>



<pre class="wp-block-code"><code>> arima(WN, order = c(0, 0, 0))
Call:
arima(x = WN, order = c(0, 0, 0))
Coefficients:
      intercept
        50.6445
s.e.     0.7115
sigma^2 estimated as 101.2:  log likelihood = -745.53,  aic = 1495.06
></code></pre>



<p>The intercept 50.6445 is the mean and sigma^2 101.2 is the variance, which gives a standard deviation, sigma of 10.05982. Compare this with the mean and standard deviation we calculated using the standard&nbsp;<code>mean()</code>&nbsp;and&nbsp;<code>sd()</code>&nbsp;functions.</p>

Simulate White Noise (WN) in R

<p>When a series follows a random walk model, it is said to be non-stationary. We can stationarize it by taking a first-order difference of the time series, which will produce a stationary series, that is, a Zero Mean White Noise series. For example, the stock prices of a stock follow a random walk model, and the series of returns (differencing of pricing series) will follow White Noise model.</p>



<p>Let's understand this phenomenon in more detail.</p>



<p>Since the differencing of a Random Walk series is a White Noise series, we can say that a Random Walk series is a cumulative sum (i.e., Integration) of a zero mean White Noise series. With this information, we can define Random Walk series in the form of ARIMA model as follows:</p>



<pre class="wp-block-code"><code>ARIMA(0,1,0)
Where:
- Augoregressive part, p = 0
- Integration, d = 1
- Moving average part, q = 0
</code></pre>



<h3>Simulate Random Walk Series</h3>



<p>We can now simulate a random walk series in R by supplying the appropriate parameters to the&nbsp;<code>arima.sim()</code>&nbsp;function as shown below:</p>



<pre class="wp-block-code"><code>RW &lt;- arima.sim(model= list(order = c(0, 1, 0)), n=200)
</code></pre>



<p>We can plot the newly generated series using the&nbsp;<code>plot.ts()</code>&nbsp;function.</p>



<pre class="wp-block-code"><code>> plot.ts(RW,main="Random Walk", col=4)</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/random-walk-1.png" alt="" class="wp-image-26256"/></figure></div>



<p>As we can clearly observe, this is a non-stationary series and its mean and standard deviation is not constant over time.</p>



<h3>First Difference Series</h3>



<p>In order to make the series stationary, we take the first order difference of the series.</p>



<pre class="wp-block-code"><code>> RW_diff &lt;- diff(RW)</code></pre>



<p>When plotted you will notice that the difference series resembles White Noise.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/random-walk-2.png" alt="" class="wp-image-26257"/></figure></div>



<p>The statistics for the&nbsp;<code>RW_diff</code>&nbsp;series are calculated below:</p>



<pre class="wp-block-code"><code>> mean(RW_diff)
[1] -0.00903282
> sd(RW_diff)
[1] 1.020447
></code></pre>



<h3>Random Walk with Drift</h3>



<p>The above Random Walk series that we simulated wanders up and down around the mean. However, we can have the Random Walk series follow an up or a down trend, called drift. To do so, we provide an additional argument mean/intercept to the&nbsp;<code>arima.sim()</code>&nbsp;function. This intercept is the slope for the model. We can also change the standard deviation of the simulated series. In the following code, we supply a mean of 1 and a standard deviation of 5.</p>



<pre class="wp-block-code"><code>> RW_drift &lt;- arima.sim(model= list(order = c(0, 1, 0)), n=200, mean=1,sd=5)
> plot.ts(RW_drift, main="Random Walk with Drift")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/random-walk-3.png" alt="" class="wp-image-26258"/></figure></div>



<h3>Estimating Random Walk Model</h3>



<p>To fit a random walk model with a drift to a time series, we will follow the following steps</p>



<ol><li>Take the first order difference of the data.</li><li>Fit the white noise model to the differenced data using&nbsp;<code>arima()</code>&nbsp;function with order of&nbsp;<code>c(0,0,0)</code>.</li><li>Plot the original time series plot.</li><li>Add the estimated trend using the&nbsp;<code>abline()</code>&nbsp;function by supplying the intercept we got by fitting the White Noise model as the slope.</li></ol>



<h4>1. First Order Difference</h4>



<p>In order to make this series stationary, we will take the difference of the series.</p>



<pre class="wp-block-code"><code>> RW_drift_diff &lt;-diff(RW_drift)
> plot.ts(RW_drift_diff, col=4, main="First Order Difference")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/random-walk-4.png" alt="" class="wp-image-26259"/></figure></div>



<h4>2. Fit White Noise model to differenced data</h4>



<p>We can now use the&nbsp;<code>arima()</code>&nbsp;function to fit the White Noise model to the differenced data.</p>



<pre class="wp-block-code"><code>> whitenoise_model &lt;- arima(RW_drift_diff, order = c(0, 0, 0))
> whitenoise_model
Call:
arima(x = RW_drift_diff, order = c(0, 0, 0))
Coefficients:
      intercept
         0.9189
s.e.     0.3631
sigma^2 estimated as 26.37:  log likelihood = -611.01,  aic = 1226.01
></code></pre>



<p>We can see that the fitted White Noise model has an intercept of 0.9189.</p>



<h4>3. Plot the original Random Walk data</h4>



<p>This can be done using the following command:</p>



<pre class="wp-block-code"><code>> plot.ts(RW_drift,col=4, main="Random Walk with Drift")
</code></pre>



<h4>4. Add the estimated trend</h4>



<p>Now on the same plot we want to add the estimated trend. At the beginning of this lesson we explained how a Random Walk series is a cumulative sum (i.e., Integration) of a zero mean White Noise series. So, the intercept, in effect, is actually the slope for our random walk series.</p>



<p>We can plot the trend line using the&nbsp;<code>abline(a,b)</code>&nbsp;function, where a is the intercept and b is the slope of the line. In our case, we will specify 'a=0' and 'b=intercept' of the White Noise model.</p>



<pre class="wp-block-code"><code>> abline(0, whitenoise_model$coef,col=2)</code></pre>



<p>The estimated trend line will be added to our plot.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/random-walk-5.png" alt="" class="wp-image-26260"/></figure></div>

Simulate Random Walk (RW) in R

<p>AutoRegressive (AR) model is one of the most popular time series model. In this model, each value is regressed to its previous observations. AR(1) is the first order autoregression meaning that the current value is based on the immediately preceding value.</p>



<p>We can use the <code>arima.sim()</code> function to simulate the AutoRegressive (AR) model.</p>



<p>Note that model argument is meant to be a list giving the ARMA order, not an actual arima model. So, for the AutoRegressive model, we will specify model as&nbsp;<code>list(ar = phi)</code>&nbsp;, in which phi is a slope parameter from the interval (-1, 1).</p>



<p>Below we create two sets of simulations with AR model, one with a slope of 0.5 and another with a slope of 0.8.</p>



<pre class="wp-block-code"><code># Simulate AutoRegressive model with 0.5 slope
AR_1 &lt;- arima.sim(model = list(ar = 0.5), n = 200)
# Simulate AutoRegressive model with 0.8 slope
#
AR_2 &lt;- arima.sim(model = list(ar = 0.9), n = 200)
plot.ts(cbind(AR_1 , AR_2 ), main="AR Model Simulated Data")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ar-1.png" alt="" class="wp-image-26265"/></figure></div>



<p>It is important to note that Random Walk is a special case of AutoRegressive models, where the slope parameter is equal to 1.</p>



<h3>ACF and PACF of Autoregressive Model</h3>



<p>We can calculate the Autocorrelation and Partial Autocorrelation Functions of the Autoregressive model using the&nbsp;<code>acf()</code>&nbsp;and the&nbsp;<code>pacf()</code>&nbsp;functions.</p>



<p>The following are the respective ACF and PACF plots for the&nbsp;<code>AR_1</code>&nbsp;series.</p>



<pre class="wp-block-code"><code>> acf(AR_1)
> pacf(AR_1)</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ar-2.png" alt="" class="wp-image-26266"/></figure></div>



<h3>Characteristics of AutoRegressive Model</h3>



<p><strong>Persistence:</strong>&nbsp;The slope in an AR model can range from -1 to 1. As the slope gets closer to 1, the model shows higher persistence, i.e., it shows higher correlation with previous values. Also, the higher the slope, the slower is the decay of ACF to 0.</p>



<p><strong>Oscillatory behavior:</strong>&nbsp;This refers to a large amount of variation between an observation and its lag. As the slope reduces, the AR model exhibits oscillatory behavior.</p>

AutoRegressive (AR) Model in R

<p>We will now see how we can fit an <a href="https://financetrain.com/autoregressive-ar-model-in-r/">AR model</a> to a given time series using the <code>arima()</code> function in R. Recall that AR model is an <code>ARIMA(1, 0, 0)</code> model.</p>



<p>We can use the&nbsp;<code>arima()</code>&nbsp;function in R to fit the AR model by specifying the&nbsp;<code>order = c(1, 0, 0)</code>.</p>



<p>We will perform the estimation using the&nbsp;<code>msft_ts</code>&nbsp;time series that we created earlier in the first lesson. If you don't have the&nbsp;<code>msft_ts</code>&nbsp;loaded in your R session, please follow the steps to create it as specified in the first lesson.</p>



<p>Let's start by creating a plot of the original data using the&nbsp;<code>plot.ts()</code>&nbsp;function.</p>



<pre class="wp-block-code"><code>> plot.ts(msft_ts, main="MSFT prices", ylab="Prices")</code></pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/ar-estimate-1.png" alt="" class="wp-image-26269"/></figure></div>



<p>As of now we are not worried about whether an AR model is best suited for this data or not. Our objective is to understand the process of fitting the AR model to this data.</p>



<p>We will fit the AR model to this data using the following command:</p>



<pre class="wp-block-code"><code>msft_ar &lt;- arima(msft_ts , order = c(1, 0, 0))</code></pre>



<p>The output contains many things including the estimated slope (ar1), mean (intercept), and innovation variance (sigma^2) as shown below:</p>



<pre class="wp-block-code"><code>> msft_ar
Call:
arima(x = msft_ts, order = c(1, 0, 0))
Coefficients:
         ar1  intercept
      0.9815    56.2260
s.e.  0.0113     2.2359
sigma^2 estimated as 0.5891:  log likelihood = -292.55,  aic = 591.1
></code></pre>



<p>The&nbsp;<code>msft_ar</code>&nbsp;object also contains the residuals (ε<sub>t</sub>&nbsp;). Using the&nbsp;<code>summary()</code>&nbsp;function, you can see that the object contains a time series of residuals.</p>



<pre class="wp-block-preformatted">&gt; summary(msft_ar)
          Length Class  Mode
coef        2    -none- numeric
sigma2      1    -none- numeric
var.coef    4    -none- numeric
mask        2    -none- logical
loglik      1    -none- numeric
aic         1    -none- numeric
arma        7    -none- numeric
residuals 252    ts     numeric
call        3    -none- call
series      1    -none- character
code        1    -none- numeric
n.cond      1    -none- numeric
nobs        1    -none- numeric
model      10    -none- list
&gt;
</pre>



<p>You can extract the residuals using the&nbsp;<code>residuals()</code>&nbsp;function in R.</p>



<pre class="wp-block-preformatted">residuals &lt;- residuals(msft_ar)</pre>



<p>Once you find the residuals ε<sub>t</sub>, the fitted values are just X̂<sub>t</sub>=X<sub>t</sub>−ε<sub>t</sub>. In R, we can do it as follows:</p>



<pre class="wp-block-preformatted">msft_fitted &lt;- msft_ts - residuals</pre>



<p>We can now plot both the original and the fitted time series to see how close the fit is.</p>



<pre class="wp-block-preformatted">ts.plot(msft_ts)
points(msft_fitted, type = "l", col = 2, lty = 2)</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ar-3.png" alt="" class="wp-image-26270"/></figure></div>



<h3>Calculating Residuals Manually from the Fitted Model (Advanced Content)</h3>



<p>The AR(1) model is expressed as follows:</p>



<p><strong>Xt = Ar * Xt-1 + εt</strong></p>



<p>So,</p>



<p><strong>εt = Xt - Ar * Xt-1</strong></p>



<p>Where Ar is the estimated autoregressive part in the fitted model.</p>



<pre class="wp-block-preformatted">n=252
e=rep(1,n)
e[1]=0
for (t in (2 : n)){
  e[t] = msft_ts [t]-coef(msft_ar )[1]*msft_ts [t-1]
}</pre>



<hr class="wp-block-separator"/>



<p>Suppose the model was ARIMA(1,0,1)</p>



<p><strong>εt = Xt - Ar * Xt-1 - Ma * Xt-1</strong></p>



<pre class="wp-block-preformatted">&gt; n=252
&gt; e=rep(1,n)
&gt; e[1]=0
for (t in (2 : n)){
  e[t] = msft_ts [t]-coef(msft_ar )[1]*msft_ts [t-1]-coef(msft_ar )[2]*e[t-1]
}</pre>

Estimating AutoRegressive (AR) Model in R

<p>Now that we know how to estimate the AR model using ARIMA, we can create a simple forecast based on the model.</p>



<h3>Step 1: Fit the model</h3>



<p>The first step is to fit the model as&nbsp;<code>ARIMA(1, 0, 0)</code>. We have already seen this in the previous lesson.</p>



<pre class="wp-block-preformatted">&gt; msft_ar&lt;-arima(msft_ts,c(1,0,0))</pre>



<h3>Step 2: Create Forecast</h3>



<p>We can now use the&nbsp;<code>predict()</code>&nbsp;function to create a forecast using the fitted AR model. It takes as its inputs, the model object that we created in step 1, and an additional parameter&nbsp;<code>n.ahead</code>&nbsp;which establishes the forecast horizon, that is, how many steps (periods) in the future we want to create the forecast. In our example, we will provide&nbsp;<code>n.ahead=20</code>, which will create forecast for next 20 steps which corresponds to 20 days for our daily data.</p>



<pre class="wp-block-preformatted">&gt; msft_forecast &lt;- predict(msft_ar, n.ahead = 20)</pre>



<p>The object generated by the&nbsp;<code>predict()</code>&nbsp;command contains two time series: 1)&nbsp;<code>$pred</code>&nbsp;which contains the forecasted values and 2)&nbsp;<code>$se</code>&nbsp;which contains the standard error for the forecast. We will use the&nbsp;<code>$pred</code>&nbsp;time series to plot the forecast and the&nbsp;<code>$se</code>&nbsp;time series to add confidence intervals to our plot.</p>



<pre class="wp-block-preformatted">&gt; msft_forecast
$pred
Time Series:
Start = 253
End = 272
Frequency = 1
 [1] 62.03064 61.92331 61.81796 61.71456 61.61307 61.51346 61.41569
 [8] 61.31973 61.22554 61.13309 61.04236 60.95330 60.86588 60.78009
[15] 60.69588 60.61323 60.53210 60.45248 60.37432 60.29762
$se
Time Series:
Start = 253
End = 272
Frequency = 1
 [1] 0.7675185 1.0754465 1.3051026 1.4933081 1.6544932 1.7961449
 [7] 1.9227626 2.0373125 2.1418796 2.2379995 2.3268456 2.4093400
[13] 2.4862246 2.5581077 2.6254964 2.6888192 2.7484426 2.8046829
[19] 2.8578161 2.9080845
</pre>



<p>For simplicity sake, let's also extract the two series in their own respective variables.</p>



<pre class="wp-block-preformatted">&gt; msft_forecast_values &lt;- msft_forecast$pred
&gt; msft_forecast_se &lt;- msft_forecast$se
</pre>



<h3>Step 3: Plot the Forecast</h3>



<p>We can now use the&nbsp;<code>plot.ts()</code>&nbsp;function to first plot the original data and then add points for the forecasted values using the&nbsp;<code>points()</code>&nbsp;function as shown below:</p>



<pre class="wp-block-preformatted">&gt; plot.ts(msft_ts, xlim = c(0, 300), ylim = c(40,80))
&gt; points(msft_forecast_values , type = "l", col = 2)</pre>



<p>Notice that while creating the initial plot, we've specified scale limits for x and y axis in order to provision for the forecast values.</p>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-forecast-1.png" alt="" class="wp-image-26291"/></figure></div>



<h3>Step 4: Add Confidence Intervals to Forecast</h3>



<p>We can add a 95% confidence interval to our forecast using the standard error values.</p>



<pre class="wp-block-preformatted">&gt; points(msft_forecast_values - 2*msft_forecast_se, type = "l", col = 4, lty = 2)
&gt; points(msft_forecast_values + 2*msft_forecast_se, type = "l", col = 4, lty = 2)
&gt;</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-forecast-2.png" alt="" class="wp-image-26290"/></figure></div>



<p>You have successfully created your first forecast.</p>

Forecasting with AutoRegressive (AR) Model in R

<p>A Moving Average is a process where each value is a function of the noise in the past observations. These are the random error terms which follow a white noise process. The general form is MA(q), where x<sub>t</sub>&nbsp;depends on q past values.</p>



<p>Just like AR models, this also has a regression like structure, however, instead of the actual values, we are regressing each value on the noise/error in the previous observations.</p>



<p>We can use the&nbsp;<code>arima.sim()</code>&nbsp;function to simulate the MA model. For the MA model, we will specify model as&nbsp;<code>list(ma = theta)</code>&nbsp;, where theta is the slope parameter from the interval (-1, 1).</p>



<p>Below we create two sets of simulations with MA model, one with a slope of 0.5 and another with a slope of 0.8.</p>



<pre class="wp-block-preformatted"># Simulate AutoRegressive model with 0.5 slope
MA_1 &lt;- arima.sim(model = list(ma = 0.5), n = 200)
# Simulate AutoRegressive model with 0.8 slope
MA_2 &lt;- arima.sim(model = list(ma = 0.9), n = 200)
# Simulate AutoRegressive model with -0.6 slope
MA_3 &lt;- arima.sim(model = list(ma = -0.6), n = 200)
plot.ts(cbind(MA_1 , MA_2, MA_3 ), main="MA Model Simulated Data")
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ma-1.png" alt="" class="wp-image-26295"/></figure></div>



<h3>ACF and PACF of Autoregressive Model</h3>



<p>We can calculate the Autocorrelation and Partial Autocorrelation functions of the Moving Average model using the&nbsp;<code>acf()</code>&nbsp;and the&nbsp;<code>pacf()</code>&nbsp;functions.</p>



<p>The following are the respective ACF and PACF plots for the&nbsp;<code>MA_1</code>&nbsp;series.</p>



<pre class="wp-block-preformatted">&gt; acf(MA_1)
&gt; pacf(MA_1)
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-ma-2.png" alt="" class="wp-image-26296"/></figure></div>

Moving Average (MA) Model in R

<p>We will now see how we can fit an MA model to a given time series using the&nbsp;<code>arima()</code>&nbsp;function in R. Recall that MA model is an&nbsp;<code>ARIMA(0, 0, 1)</code>&nbsp;model.</p>



<p>We can use the&nbsp;<code>arima()</code>&nbsp;function in R to fit the MA model by specifying&nbsp;<code>order = c(0, 0, 1)</code>.</p>



<p>We will perform the estimation using the&nbsp;<code>msft_ts</code>&nbsp;time series that we created earlier in the first lesson. If you don't have the&nbsp;<code>msft_ts</code>&nbsp;loaded in your R session, please follow the steps to create it as specified in the first lesson.</p>



<p>Let's start by creating a plot of the original data using the&nbsp;<code>plot.ts()</code>&nbsp;function.</p>



<pre class="wp-block-preformatted">&gt; plot.ts(msft_ts, main="MSFT prices", ylab="Prices")
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/ar-estimate-1-1.png" alt="" class="wp-image-26300"/></figure></div>



<p>We will fit the MA model to this data using the following command:</p>



<pre class="wp-block-preformatted">msft_ma &lt;- arima(msft_ts, order = c(0, 0, 1))</pre>



<p>The output contains many things including the estimated slope (ma1), mean (intercept), and innovation variance (sigma^2) as shown below:</p>



<pre class="wp-block-preformatted">> msft_ma
Call:
arima(x = msft_ts, order = c(0, 0, 1))
Coefficients:
         ma1  intercept
      0.8602    55.2735
s.e.  0.0259     0.2603
sigma^2 estimated as 4.952:  log likelihood = -559.83,  aic = 1125.66
></pre>



<p>The&nbsp;<code>msft_ma</code>&nbsp;object also contains the residuals (ε<sub>t</sub>&nbsp;). We can extract the residuals, using the&nbsp;<code>residuals()</code>&nbsp;function in R.</p>



<pre class="wp-block-preformatted">residuals &lt;- residuals(msft_ma)</pre>



<p>Once you find the residuals ε<sub>t</sub>, the fitted values are just X̂t=X<sub>t</sub>−ε<sub>t</sub>. In R, we can do it as follows:</p>



<pre class="wp-block-preformatted">msft_fitted &lt;- msft_ts - residuals</pre>



<p>We can now plot both the original and the fitted time series to see how close the fit is.</p>



<pre class="wp-block-preformatted">ts.plot(msft_ts)
points(msft_fitted, type = "l", col = 2, lty = 2)</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/estimate-ma-1.png" alt="" class="wp-image-26301"/></figure></div>



<p>As we can see, the model does not provide a good fit for the original series. One possible reason is that our original series was not stationary. We can fit the ARIMA model with first order differencing by passing the parameter I=1.</p>



<pre class="wp-block-preformatted">> msft_ma &lt;- arima(msft_ts, order = c(0, 1, 1))
> residuals &lt;- residuals(msft_ma)
> msft_fitted &lt;- msft_ts - residuals
> ts.plot(msft_ts)
> points(msft_fitted, type = "l", col = 2, lty = 2)</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/plot-estimate-2.png" alt="" class="wp-image-26302"/></figure></div>

Estimating Moving Average (MA) Model in R

<p>We now have a fair idea about how we can use <a href="https://financetrain.com/arima-modeling/">ARIMA modelling</a> in R to estimate and forecast a time series.</p>



<p>This is also called the Box–Jenkins method, named after the statisticians George Box and Gwilym Jenkins, that applies <em>autoregressive</em> moving average (ARMA) or <em>autoregressive</em> integrated moving average (<em>ARIMA</em>) models to find the best fit of a time-series model to past values of a time series.</p>



<p>Box-Jenkins models are quite flexible due to the inclusion of both <a href="https://financetrain.com/autoregressive-ar-model-in-r/">autoregressive</a> and <a href="https://financetrain.com/moving-average-ma-model-in-r/">moving average</a> terms. The Box-Jenkins model assumes that the time series is stationary. Box and Jenkins recommend differencing non-stationary series one or more times to achieve stationarity. Doing so produces an ARIMA model, with the "I" standing for "Integrated".</p>



<p>We will now go through all the steps involved in ARIMA modelling.</p>



<p>Understanding the abbreviations:</p>



<ul><li>When a model only involves autoregressive terms it may be referred to as an AR model. When a model only involves moving average terms, it may be referred to as an MA model.</li><li>When no differencing is involved, the abbreviation ARMA may be used.</li></ul>



<p>The examples we have considered are all for non-seasonal data. Box-Jenkins models can be extended to include seasonal autoregressive and seasonal moving average terms. Although this complicates the notation and mathematics of the model, the underlying concepts for seasonal autoregressive and seasonal moving average terms are similar to the non-seasonal autoregressive and moving average terms.</p>



<p><strong>Specifying the Elements of the Model</strong></p>



<p>In an ARIMA model, the elements are specified in the order (AR order, differencing, MA order).</p>



<ul><li>A model with one AR term would be specified as an ARIMA of order (1,0,0).</li><li>A model with two AR terms would be specified as an ARIMA of order (2,0,0)</li><li>A model with two MA terms (MA(2)) would be specified as an ARIMA of order (0,0,2).</li><li>A model with one AR term, a first difference, and one MA term would be specified as ARIMA of order (1,1,1).</li></ul>



<p>In the last model,&nbsp;<code>ARIMA (1,1,1)</code>, we are applying a model with one AR term and one MA term to the variable&nbsp;<em>z<sub>t</sub>=x<sub>t</sub>-x<sub>t-1</sub></em>. This is the first difference which is used to account for a linear trend in the data.</p>

ARIMA Modelling in R

<p>The first step is to identify a possible model for a given time series. To do so, we need three things: a time series plot of the data, ACF plot and the ACF plot. Analysis of these three plots can help us fairly identify the suitable model.</p>



<h3>Observing the Time Series Plot</h3>



<ul><li>The very first thing that we should do with our data is to create a time series plot and look for patterns such as trends, seasonality, outliers, and non-constant variance.</li><li>If we spot an obvious uptrend or down trend, we may need to consider first order differencing. If the time series exhibits a quadratic trend we may want to go for 2nd order differencing.</li><li>Generally, we don't go beyond two levels of differencing. We don't want to introduce overdifferencing.</li><li>If the series is still not stationary after two levels of differencing, then we may want to go for smoothing.</li><li>If the data exhibits a curved upward trend along with non-constant variance such as increasing variance, we can transform the data series by either a logarithm or a square root. Data series with non-constant variance may also call for more advanced models such as ARCH (Not covered in this course).</li></ul>



<h3>ACF and PACF Plots</h3>



<p>We should consider ACF and PACF plots together to identify the order (i.e., the p and q) of the autoregressive and moving average terms.</p>



<p><strong>Important Note:</strong>&nbsp;If the ACF and PACF do not tail off, but instead have values that stay close to 1 over many lags, the series is non-stationary and differencing will be needed. Try a first difference and then look at the ACF and PACF of the differenced data.</p>



<h3>Identification of an MA model is often best done with the ACF rather than the PACF.</h3>



<ul><li>For an MA model, the theoretical PACF does not shut off, but instead tapers toward 0 in some manner. A clearer pattern for an MA model is in the ACF. The ACF will have non-zero autocorrelations only at lags involved in the model.</li><li>A sample ACF with a significant autocorrelation only at lag 1 is an indicator of a possible MA(1) model.</li><li>A sample ACF with significant autocorrelations at lags 1 and 2, but non-significant autocorrelations for higher lags indicates a possible MA(2) model.</li><li>A property of MA(q) models in general is that there are non-zero autocorrelations for the first q lags and autocorrelations = 0 for all lags &gt; q.</li></ul>



<p>So, to identify the order of <a href="https://financetrain.com/moving-average-ma-model-in-r/">moving average</a> process, we examine the sample <a href="https://financetrain.com/autocorrelation-in-r/">autocorrelation</a> function to see where it essentially becomes zero.</p>



<h3>Identification of an AR model is often best done with the PACF.</h3>



<p>For an AR model, the theoretical PACF “shuts off” past the order of the model. The phrase “shuts off” means that in theory the partial autocorrelations are equal to 0 beyond that point. Put another way, the number of non-zero partial autocorrelations gives the order of the AR model. By the “order of the model” we mean the most extreme lag of x that is used as a predictor.</p>



<p><strong>Summary:</strong>&nbsp;For AR models, the ACF will dampen exponentially and the PACF plot will be used to identify the order (p) of the AR model. For MA models, the PACF will dampen exponentially and the ACF plot will be used to identify the order (q) of the MA model.</p>



<h3>ARMA Models</h3>



<p>ARMA models (including both AR and MA terms) have ACFs and PACFs that both tail off to 0. These are the trickiest because the order will not be particularly obvious. Basically you just have to guess that one or two terms of each type may be needed and then see what happens when you estimate the model.</p>



<table class="wp-block-table"><tbody><tr><th>SHAPE</th><th>INDICATED MODEL</th></tr><tr><td>Exponential, decaying to zero</td><td>Autoregressive model. Use the partial autocorrelation plot to identify the order of the autoregressive model.</td></tr><tr><td>Alternating positive and negative, decaying to zero</td><td>Autoregressive model. Use the partial autocorrelation plot to help identify the order.</td></tr><tr><td>One or more spikes, rest are essentially zero</td><td>Moving average model, order identified by where plot becomes zero.</td></tr><tr><td>Decay, starting after a few lags</td><td>Mixed autoregressive and moving average model.</td></tr><tr><td>All zero or close to zero</td><td>Data is essentially random.</td></tr><tr><td>High values at fixed intervals</td><td>Include seasonal autoregressive term.</td></tr><tr><td>No decay to zero</td><td>Series is not stationary.</td></tr></tbody></table>



<h3>Difficulty in Identifying Mixed Model</h3>



<p>In practice, the autocorrelation and partial autocorrelation functions may not give a very clear picture as suggested in theory. Sometimes both MA and AR model will seem okay. In such a case, we can conduct a few more tests to identify the best suited model.</p>



<ul><li>Analyze the standard errors of forecasted values. We should forecast with both models and analyze the standard errors of the forecasted values. Finally pickup the model which has lowest standard errors for predictions.</li><li>Compare AIC and BIC statistics: We should compare the models w.r.t. statistics such as AIC (Akaike Information Criterion) and BIC. The model with least values for these statistics should be selected.</li></ul>

ARIMA Modelling - Identify Model for a Time Series

<p>In this lesson, we will take a new dataset (stock prices) and use all that we have learned to create a forecast using the ARIMA Models.</p>



<p>We will take the closing prices of Facebook stock for this example.</p>



<h3>Step 1: Load the Data</h3>



<p>We will load Facebook daily closing prices for the past 3 years from the WIKI/PRICES dataset in Quandl.</p>



<p>The following code loads the data from Quandl and converts it into a time series object called&nbsp;<code>fb_ts</code>.</p>



<pre class="wp-block-preformatted">Quandl.api_key("YOUR API KEY")
fb_data &lt;- Quandl.datatable("WIKI/PRICES" ,
                              qopts.columns=c("date", "close"),
                              ticker=c("FB"),
                              date.gte=c("2014-01-01"),
                              date.lte=c("2016-12-31"))
ZOO &lt;- zoo(fb_data$close, order.by=as.Date(as.character(fb_data$date), format='%Y-%m-%d'))
fb_ts &lt;- ts(ZOO)
</pre>



<h3>Step 2: Plot the Time Series and Test for Stationarity</h3>



<p>We can now plot the time series and look for common patterns.</p>



<pre class="wp-block-preformatted">&gt; plot.ts(fb_ts, main="Facebook Stock Prices", col=4)
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/case-1.png" alt="" class="wp-image-26309"/></figure></div>



<p>The above chart exhibits a clear upward trend in the stock prices. This means that we will be differencing to make the series stationary.</p>



<p>Sometimes a visual inspection may not clearly identify such pattern. In such cases, we can use the Augmented Dickey-Fuller unit root test to check whether the series is stationary or not. The p-value resulting from the ADF test must be less than 0.05 or 5% for a time series to be stationary. If the p-value is greater than 0.05 or 5%, you conclude that the time series has a unit root which means that it is a non-stationary process.</p>



<p>In R, we can do this using&nbsp;<code>adf.test()</code>&nbsp;function available in the tseries package. The following code loads the required package and performs the test.</p>



<pre class="wp-block-preformatted">&gt; install.packages("tseries")
&gt; library(tseries)
&gt; adf.test(fb_ts)
    Augmented Dickey-Fuller Test
data:  fb_ts
Dickey-Fuller = -3.1368, Lag order = 9, p-value = 0.09883
alternative hypothesis: stationary
&gt;
</pre>



<p>The test results confirm our observation that series is non-stationary (p-value &gt;0.05) and will need differencing to make it stationary.</p>



<h3>Step 3: Identify the Model</h3>



<p>The next step is to identify the model, i.e., the appropriate order of Autoregressive (AR) and Moving Average (MA) processes p, and q. We will do so using the Autocorrelation function (ACF) and Partial Autocorrelation function (PACF).</p>



<p>Let's create the ACF and PACF plots.</p>



<pre class="wp-block-preformatted">&gt; acf(fb_ts)
&gt; pacf(fb_ts)
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/case-2.png" alt="" class="wp-image-26310"/></figure></div>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/case-3.png" alt="" class="wp-image-26311"/></figure></div>



<p>Recall our analysis of these two functions. The ACF plot shows slow decay of lag to 0 indicating an AR model. The PACF plot suggests AR model of the order 1 AR(1) as PACF number is close to 0 after lag 1.</p>



<p>Identifying the best fit model is a complex process and we may want to test multiple models to check what best fits our data. Both experience and knowledge of advanced topics can be helpful. However, based on our limited analysis, let's say we will go with p=1, d=1, and q=0.</p>



<p>Our suggested model is then ARIMA(1,1,0).</p>



<h3>Step 4: Estimate the Model</h3>



<p>We can now estimate the model for our data as shown below:</p>



<pre class="wp-block-preformatted">fb_fit &lt;- arima(fb_ts, order = c(1, 1, 0))
</pre>



<pre class="wp-block-preformatted">&gt; fb_fit
Call:
arima(x = fb_ts, order = c(1, 1, 0))
Coefficients:
         ar1
      0.0306
s.e.  0.0364
sigma^2 estimated as 2.642:  log likelihood = -1438.1,  aic = 2880.2
&gt;
</pre>



<h3>Step 5: Create Forecast</h3>



<p>The next step is to create a forecast using the fitted model. We are predicting the stock prices for the next 20 days.</p>



<pre class="wp-block-preformatted">&gt; fb_forecast &lt;- predict(fb_fit , n.ahead = 20)
</pre>



<pre class="wp-block-preformatted">&gt; fb_forecast
$pred
Time Series:
Start = 757
End = 776
Frequency = 1
 [1] 115.0102 115.0090 115.0089 115.0089 115.0089 115.0089 115.0089 115.0089 115.0089
[10] 115.0089 115.0089 115.0089 115.0089 115.0089 115.0089 115.0089 115.0089 115.0089
[19] 115.0089 115.0089
$se
Time Series:
Start = 757
End = 776
Frequency = 1
 [1] 1.625520 2.334309 2.874161 3.327575 3.726221 4.086157 4.416859 4.724469 5.013239
[10] 5.286259 5.545854 5.793829 6.031618 6.260382 6.481076 6.694498 6.901324 7.102128
[19] 7.297410 7.487600
&gt;
</pre>



<p>Notice that all predicted values for the next 20 periods look the same. This is because the forecast will be flat with no drift. There is a possibility that our model has a drift also along with a trend.</p>



<p>For simplicity sake, let's also extract the two series in their own respective variables.</p>



<pre class="wp-block-preformatted">&gt; fb_forecast_values &lt;- fb_forecast$pred
&gt; fb_forecast_se &lt;- fb_forecast$se
</pre>



<h3>Step 6: Plot the Forecast</h3>



<p>We can now plot the original series and the forecast as shown below:</p>



<pre class="wp-block-preformatted">&gt; plot.ts(fb_ts, xlim = c(0, 900), ylim = c(50,160))
&gt; points(fb_forecast_values , type = "l", col = 2)
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/case-4.png" alt="" class="wp-image-26312"/></figure></div>



<p>We can also plot the fitted data for the current time period on the chart to show how well the model fits the data.</p>



<pre class="wp-block-preformatted">&gt; fb_fitted_data &lt;- fb_ts - residuals(fb_fit)
&gt; points(fb_fitted_data, type = "l", col = 4, lty = 2)</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/case-5.png" alt="" class="wp-image-26313"/></figure></div>

Forecasting with ARIMA Modeling in R - Case Study

<h3>auto.arima() Function<br></h3>



<p>R also has a package called forecast, which contains many forecasting functions for time series and linear models. It also contains a very useful function called&nbsp;<code>auto.arima</code>, which returns the best ARIMA model according to either AIC, AICc or BIC value. The function conducts a search over possible models within the order constraints provided.</p>



<p>Let's try the&nbsp;<code>auto.arima</code>&nbsp;function on our time series.</p>



<pre class="wp-block-preformatted">> fb_fit_auto &lt;- auto.arima(fb_ts)
> fb_fit_auto
Series: fb_ts
ARIMA(1,1,4) with drift
Coefficients:
          ar1     ma1      ma2      ma3      ma4   drift
      -0.8405  0.8626  -0.0123  -0.1632  -0.1768  0.0800
s.e.   0.1080  0.1111   0.0488   0.0509   0.0381  0.0478
sigma^2 estimated as 2.581:  log likelihood=-1426.29
AIC=2866.58   AICc=2866.73   BIC=2898.97
></pre>



<p>The&nbsp;<code>auto.arima</code>&nbsp;function suggests the best fit model as ARIMA(1,1,4) with drift.</p>



<p>The following code creates the forecast for the FB stock prices using the suggested model:</p>



<pre class="wp-block-preformatted">fb_fit &lt;- arima(fb_ts, order = c(1, 1, 4))
fb_forecast &lt;- predict(fb_fit , n.ahead = 20)
fb_forecast_values &lt;- fb_forecast$pred
plot.ts(fb_ts, xlim = c(0, 900), ylim = c(50,160))
points(fb_forecast_values , type = "l", col = 2)
</pre>



<div class="wp-block-image"><figure class="aligncenter"><img src="https://financetrain.com/wp-content/uploads/case-6.png" alt="" class="wp-image-26316"/></figure></div>

Automatic Identification of Model Using auto.arima() Function in R

<p>This course provided an overview of the fundamentals of time series analysis and how we can perform time series analysis in R. We reviewed some of the most important concepts of time series analysis and looked at the process involved in modeling a time series using the ARIMA models.</p>



<p>Time series analysis is a complex subject and this course focused on building a strong foundation for someone interested in time series analysis and forecasting. The skills you learned here will go a long way in building time series models in real-life scenarios. If you are interested in further building your expertise in time series analysis, there is practically end-less material available for you to explore. It will help to take up lots of sample time series data sets and try to forecast them. Some of the topics you should research next include how to go about validating your model, smoothing and decomposition methods, modeling seasonal time series (SARIMA - Seasonal ARIMA), ARCH models and spectral analysis.</p>

Financial Time Series in R - Course Conclusion

Financial Time Series Analysis in R

This frontmatter description will appear below the title

Financial Time Series Analysis in R

Learn the fundamentals of analyzing a financial time series in R

What's Included

Course Resources

Lessons

Free Guides - Getting Started with R and Python

Accelerate your finance career with cutting-edge data skills.