In this example we carefully looked at official currencies of the United States of America and Republic of South Africa bought and sold on foreign exchange markets. More generally, we properly focused on the closing price of the currency pair. We intentionally used time series analysis to correctly identify consistent patterns and key trends, so as forecasting predictions of the insistent demand of the currency pair. We used closing price to assess changes in the currency pair. We compared closing price of USD/ZAR over a year (from 03 July 2019 to 03 July 2020) to measure market sentiments of the currency pair.
Time series analysis enables one to correctly identify and accurately predict the underlying nature of a phenomenon represented by a specific sequence of independent observations. This statistical model is carefully considered as a complex linear regression.
The ARIMA (p,d,q) model was used to compute a trend and predict future values of our time series. The term ARIMA is threefold.
- AR (Autoregressive) – linear combination of previous values’ influence.
- I (Integrative) – random walk.
- MA (Moving average) – linear combination of previous errors.
This model is flexible and easy to understand. Like any other statistical model, there are certain key assumptions that must be adequately fulfilled. Underneath are several assumptions associated with this model.
- There must be more than 50 independent observations.
- There must be no missing values.
- There must be no outliers in the series.
- The independent variable must not follow a normal distribution.
- The series must not be stationary.
- There must be no white noise.
- There must be a correlation between variables with themselves.
Importing Python libraries
Underneath, key Python libraries were imported. These libraries gallantly helped us create arrays and profile tables, compute graphs and make accurate predictions.
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
%matplotlib inline
import seaborn as sns
sns.set("talk","ticks",font_scale=1,font="sans-serif",color_codes=True)
from pylab import rcParams
plt.rcParams["figure.figsize"] = [18,14]
import warnings
warnings.filterwarnings("ignore")
import statsmodels.api as sm
import scipy.stats as stats
import statsmodels.api as sm
from pandas.plotting import lag_plot, autocorrelation_plot
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
from sklearn import metrics
from statsmodels.tsa.arima_model import ARIMA
from statsmodels.tsa.stattools import adfuller
from statsmodels.tsa.seasonal import seasonal_decomposeLoading data
df = pd.read_csv(r"C:\Users\Tshepo\Desktop\Advanced_Data_Science\USDZAR.csv",index_col=[0],parse_dates=[0])
df.head()The pd.read_csv() function was called to load obtained csv file into a pandas data frame. See underneath table to examine the data structure. Bear in mind that time was parsed by calling the parse_dates function and then turned into an index.
| Open | High | Low | Close | Adj Close | Volume | |
|---|---|---|---|---|---|---|
| Date | ||||||
| 2019-07-03 | 14.0770 | 14.1474 | 14.0571 | 14.0743 | 14.0743 | 0 |
| 2019-07-04 | 14.0394 | 14.0789 | 13.9533 | 14.0529 | 14.0529 | 0 |
| 2019-07-05 | 14.0503 | 14.2717 | 14.0082 | 14.0385 | 14.0385 | 0 |
| 2019-07-08 | 14.1925 | 14.2318 | 14.0964 | 14.1952 | 14.1952 | 0 |
| 2019-07-09 | 14.1798 | 14.2567 | 14.1219 | 14.1795 | 14.1795 | 0 |
Data preprocessing
del df["Open"]
del df["High"]
del df["Low"]
del df["Adj Close"]
del df["Volume"]| Close | |
| Date | |
| 2019-07-03 | 14.0743 |
| 2019-07-04 | 14.0529 |
| 2019-07-05 | 14.0385 |
| 2019-07-08 | 14.1952 |
| 2019-07-09 | 14.1795 |
Above we deleted features not for use and remain only with the closing price of USD/ZAR from
Detecting missing values
The visible presence of missing values results in measurement errors that negatively impact conclusions. To studiously avoid such as a dilemma, data was carefully examined to unanimously determine whether there were any missing values. See underneath figure for a visual representation of the possible presence of missing values.
sns.heatmap(df.isnull(),cmap="Blues")
plt.title("Detect Missing Values")
plt.show()
Descriptive statistics
In the section underneath we summarize the data.
fig, (ax1,ax2) = plt.subplots(1,2,figsize=(16,7))
sns.boxplot(df["Close"],ax=ax1)
ax1.set_title("USD/ZAR closing price box plot")
sns.distplot(df["Close"],ax=ax2)
ax2.set_title("USD/ZAR closing price histogram")
ax2.set_ylabel("Close value")
plt.show()
The histogram above shows that data was skewed to the left. Independent observations were not spread symmetrically around the true value. The knowledge of mean and explained variance were not sufficient to summarize data from this distribution. Consequently, we have fulfilled the fourth assumption.
Testing of stationary
A time series is considered to be stationary when statistical properties (such as variance and autocorrelation) do not vary across time.
The constructed hypothesis was as follows:
- Null hypothesis: The series is stationary.
- Alternative hypothesis: The series is not stationary.
The Augmented Dickey Fuller (ADF) test statistics was used to determine whether or not the series was stationary. The hypothesis was translated as follows:
- Null hypothesis: There is a unit root.
- Alternative hypothesis: There is no unit root.
adfullerreport = adfuller(df["Close"])
adfullerreportdata = pd.DataFrame(adfullerreport[0:4], columns = ["Values"], index = ["ADF F% statistics","P-value","No. of lags used", "No. of observations"])
adfullerreportdata| Values | |
| ADF F% statistics | -1.031626 |
| P-value | 0.741555 |
| No. of lags used | 0.000000 |
| No. of observations | 262.000000 |
The p-value is greater than 0.05. We fail to reject the null hypothesis. There is no unit root. 262 independent observations were used to conduct a ADF test statistics. The ADF statistics: F% equals to -1.031626 (closer to 1 as required by the test).
Testing of white noise
A time-series has white noise when there is no correlation between variables, meaning all autocorrelations are equal to zero. White noise is an example of stationary.
The hypothesis was as follows: –
- Null hypothesis: There is white noise.
- Alternative hypothesis: There is no white noise.
randval = np.random.randn(1000)
fig, (ax1,ax2) = plt.subplots(1,2,figsize=(16,7))
ax1.plot(randval)
ax1.set_title("Random white noise plot")
ax1.set_ylabel("Values")
ax2.hist(randval)
ax2.set_title("Random white noise histogram")
ax2.set_ylabel("Related random white noise plot Frequency")
ax2.set_xlabel("Random white noise")
plt.show()
Random white noise is centered around the mean.
autocorrelation_plot(randval)
plt.title("Random white noise lag and autocorrelation Plot")
plt.show()
There are spikes above the 95% and 99% confidence interval.
Test of correlation
The autocorrelation plots for USD/ZAR daily returns shows most of the spikes are not statistically significant. Meaning that returns are not highly correlated. A lag plot can also be used to detect outliers in a time series. As you can see, there were no outliers detected.
lag_plot(df["Close"])
plt.title("USD/ZAR closing price lag Plot")
plt.show()
The lag plot does not show randomness. We can state that there is a positive autocorrelation in the series.
autocorrelation_plot(df["Close"])
plt.title("USD/ZAR closing price lag and autocorrelation Plot")
plt.show()
ACF and a PACF plots were used to determine parameters of estimating the model. ACF gives one an idea of series correlation of data over time. The plot basically depicts the extent to which present values of a series relate to previous values. An ACF plot takes into account the trend, seasonality, cyclic and residual components.
plot_acf(df["Close"])
plt.xlabel("Lag")
plt.ylabel("ACF")
plt.show()
PACF show the partial correlation coefficient that is not explained at low-level-lags. The PACF shows a significant spike at lag 1. This means that all higher order autocorrelation can be explained by lag 1. As such, fitting an autoregressive model with order AR(1). We did not immediately jump into conclusion and fit the model. We also used the itertools to determine best parameters for an ARIMA through testing different models and finding which order has the lowest AIC.
plot_pacf(df["Close"])
plt.xlabel("Lag")
plt.ylabel("PACF")
plt.show()
Original time series
Underneath is the original time series of USD/ZAR closing price from 03 July 2019 to 03 July 2020.
df.plot(kind="line", color="navy")
plt.title("USD/ZAR closing price (original) time series")
plt.xlabel("Date")
plt.ylabel("Closing price")
plt.show()
Time Series Smoothening
Techniques underneath were used to smoothen the time series.
- Moving averages smoothening.
- Standard deviation smoothening.
- Exponential smoothening.
MA5 = df["Close"].rolling(window=5).mean()
MA15 = df["Close"].rolling(window=15).mean()
MA30 = df["Close"].rolling(window=30).mean()
MA60 = df["Close"].rolling(window=60).mean()
MA90 = df["Close"].rolling(window=90).mean()
df.plot(kind="line",color="navy")
MA5.plot(kind="line",color="green",label="Rolling MA 5 days")
MA15.plot(kind="line",color="gray",label="Rolling MA 15 days")
MA30.plot(kind="line",color="orange",label="Rolling MA 30 days")
MA60.plot(kind="line",color="red",label="Rolling MA 60 days")
MA90.plot(kind="line",color="black",label="Rolling MA 90 days")
plt.title("USD/ZAR closing price smoothened (moving averages) time series")
plt.xlabel("Date")
plt.xticks(rotation=45)
plt.ylabel("Closing price")
plt.legend(loc=2)
plt.show()The figure underneath depicts the smoothened series. We used 5 days, 15 days, 30 days and, 60 days and 90 days moving averages.
The figure underneath also shows the smoothened series. However, this time around we used 5 days, 15 days, 30 days and, 60 days and 90 days standard deviations.
Std5 = df["Close"].rolling(window=5).std()
Std15 = df["Close"].rolling(window=15).std()
Std30 = df["Close"].rolling(window=30).std()
Std60 = df["Close"].rolling(window=60).std()
Std90 = df["Close"].rolling(window=90).std()
Std5.plot(kind="line",color="green",label="Rolling Std 5 days")
Std15.plot(kind="line",color="gray",label="Rolling Std 15 days")
Std30.plot(kind="line",color="orange",label="Rolling Std 30 days")
Std60.plot(kind="line",color="red",label="Rolling MA Std days")
Std90.plot(kind="line",color="black",label="Rolling MA 90 days")
plt.title("USD/ZAR close price smoothened (standard deviation) time series")
plt.xlabel("Date")
plt.xticks(rotation=45)
plt.ylabel("Closing price")
plt.legend(loc=2)
plt.show()
Running min and max
The figure underneath was used to determine whether all values of the specified last period fall in a particular range. It shows the running minimum and maximum closing price of USD/ZAR.
df_expanding = df
df_expanding["Running_min"] = df_expanding["Close"].expanding().min()
df_expanding["Running_max"] = df_expanding["Close"].expanding().max()
df_expanding.plot()
plt.ylabel("Closing price")
plt.title("USD/ZAR close price running min and max")
plt.show()
USD/ZAR closing price dropped to a low of 13.8531 on the 24th of July 2019. Thereafter, a psychological level was established. In the last quarter of 2019, the South African rand moderately steadied against the US dollar. However, those gain quickly wiped out in the beginning of the year. The US dollar edged higher, reaching a new 19.2486 on the 6th of April 2020. The second was not an easy one for the US dollar.
Rate of turn
The figure underneath illustrates the rate of return of the portfolio from 3 July 2019 to 3 July 2020. We found that in the last quarter there was moderate rate of return. In the first quarter of 2020, the rate of return tripled. Realized nominal rates was about 35%.
pr = df.pct_change()
pr_plus_one = pr.add(1)
cummulative_return = pr_plus_one.cumprod().sub(1)
cummulative_return.mul(100).plot()
plt.title("USD/ZAR closing price rate of return")
plt.ylabel("Return")
plt.show()
Rolling 90-days rate of return
Underneath we constructed a function for rolling 90 days rate of return.
def get_rate_of_return(peiod_return):
return np.prod(peiod_return + 1) - 1rolling_rate_of_return = df["Close"].rolling(window=90).apply(get_rate_of_return)
rolling_rate_of_return.plot()
plt.title("USD/ZAR closing price rolling rate of return")
plt.ylabel("Return")
plt.show()
The figure above depicts an exponential growth of rate of return of the currency pair in the last 90 days.
Training and testing data
The data was split into training data and test data using the code underneath. The training data was used to fit the data. Meanwhile, the test data was conveniently used to validate the ARIMA model.
train = df.iloc[:len(df) - 12]
test = df.iloc[:len(df) - 12]Determining best parameters for Seasonal ARIMA model
We used itertools to determine the best parameters of seasonal ARIMA model using the code underneath.
import itertools
p = d = q = range(0, 2)
pdq = list(itertools.product(p, d, q))
seasonal_pdq = [(x[0], x[1], x[2], 12) for x in list(itertools.product(p, d, q))]
print('SARIMAX: {} x {}'.format(pdq[1], seasonal_pdq[1]))
print('SARIMAX: {} x {}'.format(pdq[1], seasonal_pdq[2]))
print('SARIMAX: {} x {}'.format(pdq[2], seasonal_pdq[3]))
print('SARIMAX: {} x {}'.format(pdq[2], seasonal_pdq[4]))SARIMAX: (0, 0, 1) x (0, 0, 1, 12)
SARIMAX: (0, 0, 1) x (0, 1, 0, 12)
SARIMAX: (0, 1, 0) x (0, 1, 1, 12)
SARIMAX: (0, 1, 0) x (1, 0, 0, 12)
for param in pdq:
for param_seasonal in seasonal_pdq:
try:
mod = sm.tsa.statespace.SARIMAX(train["Close"],
order=param,
seasonal_order=param_seasonal,
enforce_stationarity=False,
enforce_invertibility=False)
results = mod.fit()
print('ARIMA{}x{}12 - AIC:{}'.format(param, param_seasonal, results.aic))
except:
continueWe found that the best parameters for the seasonal ARIMA model are as follows: –
ARIMA(1, 1, 1)x(1, 1, 1, 12)12 – AIC:-91.8349778459599
The criterion of selecting the parameters as based on which order has the lowest AIC. The best order was(1, 1, 1) and best seasonal order was (1, 1, 1, 12). As such, configured the ARIMA model using those parameters underneath.
timeseriesmodel1 = sm.tsa.statespace.SARIMAX(df["Close"],
order=(1, 1, 1),
seasonal_order=(1, 1, 1, 12),
enforce_stationarity=False,
enforce_invertibility=False)timeseriesmodel_fit1 = timeseriesmodel1.fit()
timeseriesmodel_fit1.summary()| Dep. Variable: | Close | No. Observations: | 263 |
|---|---|---|---|
| Model: | SARIMAX(1, 1, 1)x(1, 1, 1, 12) | Log Likelihood | 48.742 |
| Date: | Fri, 03 Jul 2020 | AIC | -87.485 |
| Time: | 20:16:46 | BIC | -70.166 |
| Sample: | 07-03-2019 | HQIC | -80.503 |
| – 07-03-2020 | |||
| Covariance Type: | opg |
| coef | std err | z | P>|z| | [0.025 | 0.975] | |
|---|---|---|---|---|---|---|
| ar.L1 | 0.4933 | 28.976 | 0.017 | 0.986 | -56.299 | 57.285 |
| ma.L1 | -0.4945 | 28.954 | -0.017 | 0.986 | -57.244 | 56.255 |
| ar.S.L12 | -0.0232 | 0.068 | -0.343 | 0.731 | -0.156 | 0.109 |
| ma.S.L12 | -0.9363 | 0.063 | -14.881 | 0.000 | -1.060 | -0.813 |
| sigma2 | 0.0351 | 0.003 | 12.619 | 0.000 | 0.030 | 0.041 |
| Ljung-Box (Q): | 64.42 | Jarque-Bera (JB): | 87.67 |
|---|---|---|---|
| Prob(Q): | 0.01 | Prob(JB): | 0.00 |
| Heteroskedasticity (H): | 4.66 | Skew: | 0.52 |
| Prob(H) (two-sided): | 0.00 | Kurtosis: | 5.80 |
Warnings:
[1] Covariance matrix calculated using the outer product of gradients (complex-step).
The profile table above shows that AIC is -87.485. Which low compared to other orders. 263 observations were used to compute the table. ar.L1 , ma.L1 and ar.S.L12 are not statistically significant. P-values are less greater than 0.005. For the purpose of this example, we are not going to try other orders.
NB: It is advisable to try all possible orders and seasonal orders to find the most parameters that produce a model with more accurate predictions.
Now let us proceed and make predictions.
start = len(train)
end = len(train) + len(train) -1Predictions
The future is unknown and can only be understood by previous price behavior. Underneath we attempted to make predictions about the future of the USD/ZAR pair. This involved using ARIMA model fit on historical data to predict future values. We used a medium term prediction. It is important to note that it is far easier to make predictions using a short-time horizon than medium and long-term predictions. To get the best out of forecasting, one must make forecasts more frequently.
predictions1 = timeseriesmodel_fit1.predict(start,end,typ="levels").rename("Prediction")test.plot()
predictions1.plot()
plt.title("USD/ZAR closing price predictions")
plt.ylabel("Close price")
plt.show()
The ARIMA model forecast a slight drop in USD/ZAR at the end of July 2020 and an emerging upward bullish for the next 8 months.
Advanced predictions using Prophet
For this section, we use prophet to develop a model and make predictions.
df["ds"] = df["Time"]
df["y"] = df["Close"]df.set_index("Date")Close | ds | y | |
|---|---|---|---|
| Date | |||
| 2019-07-03 | 14.074300 | 2019-07-03 | 14.074300 |
| 2019-07-04 | 14.052900 | 2019-07-04 | 14.052900 |
| 2019-07-05 | 14.038500 | 2019-07-05 | 14.038500 |
| 2019-07-08 | 14.195200 | 2019-07-08 | 14.195200 |
| 2019-07-09 | 14.179500 | 2019-07-09 | 14.179500 |
| … | … | … | … |
| 2020-06-29 | 17.298901 | 2020-06-29 | 17.298901 |
| 2020-06-30 | 17.219200 | 2020-06-30 | 17.219200 |
| 2020-07-01 | 17.341900 | 2020-07-01 | 17.341900 |
| 2020-07-02 | 17.039301 | 2020-07-02 | 17.039301 |
| 2020-07-03 | 17.037100 | 2020-07-03 | 17.037100 |
263 rows × 3 columns
model = Prophet(interval_width=0.95, yearly_seasonality=True, weekly_seasonality=True, daily_seasonality=True)
model.fit(df)
future = model.make_future_dataframe(periods=365)
fb_predictions = model.predict(future)
model.plot(fb_predictions)
plt.title("USD/ZAR closing price predictions")
plt.xlabel("Date")
plt.ylabel("Close price")
plt.show()
Prophet agrees with the ARIMA model developed above. Overall, there is an upward trend on the USD/ZAR closing prices.
model.plot_components(fb_predictions)
plt.show()The components underneath show a clear daily seasonality. USD/ZAR closing prices drop in the early hours of the day and start to peak at about 20:00 GMT +002. This is quick reasonable given time difference between the USA and South Africa. Most returns are realized on the first day of the week. In March, the US dollar rallied strongly, however, there were some corrections at the end of the month.
