{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from pandas.plotting import register_matplotlib_converters\n",
    "from statsmodels.graphics.tsaplots import plot_acf, plot_pacf\n",
    "from statsmodels.tsa.stattools import acf, pacf\n",
    "from statsmodels.tsa.arima_model import ARIMA\n",
    "from datetime import datetime, timedelta\n",
    "register_matplotlib_converters()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Generate Some Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# $y_t = 50 + 0.4\\varepsilon_{t-1} + 0.3\\varepsilon_{t-2} + \\varepsilon_t$\n",
    "# $\\varepsilon_t \\sim N(0,1)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 126,
   "metadata": {},
   "outputs": [],
   "source": [
    "errors = np.random.normal(0, 1, 400)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 127,
   "metadata": {},
   "outputs": [],
   "source": [
    "date_index = pd.date_range(start='9/1/2019', end='1/1/2020')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 128,
   "metadata": {},
   "outputs": [],
   "source": [
    "mu = 50\n",
    "series = []\n",
    "for t in range(1,len(date_index)+1):\n",
    "    series.append(mu + 0.4*errors[t-1] + 0.3*errors[t-2] + errors[t])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 129,
   "metadata": {},
   "outputs": [],
   "source": [
    "series = pd.Series(series, date_index)\n",
    "series = series.asfreq(pd.infer_freq(series.index))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 130,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x1bbb09d94c8>"
      ]
     },
     "execution_count": 130,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 43280 characters]\n",
      "text/plain": [
       "<Figure size 720x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,4))\n",
    "plt.plot(series)\n",
    "plt.axhline(mu, linestyle='--', color='grey')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 131,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calc_corr(series, lag):\n",
    "    return pearsonr(series[:-lag], series[lag:])[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# ACF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 151,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<BarContainer object of 10 artists>"
      ]
     },
     "execution_count": 151,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 5168 characters]\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "acf_vals = acf(series)\n",
    "num_lags = 10\n",
    "plt.bar(range(num_lags), acf_vals[:num_lags])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PACF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 147,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<BarContainer object of 25 artists>"
      ]
     },
     "execution_count": 147,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 5376 characters]\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "pacf_vals = pacf(series)\n",
    "num_lags = 25\n",
    "plt.bar(range(num_lags), pacf_vals[:num_lags])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Get training and testing sets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 134,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_end = datetime(2019,12,30)\n",
    "test_end = datetime(2020,1,1)\n",
    "\n",
    "train_data = series[:train_end]\n",
    "test_data = series[train_end + timedelta(days=1):test_end]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fit ARIMA Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 135,
   "metadata": {},
   "outputs": [],
   "source": [
    "#create the model\n",
    "model = ARIMA(train_data, order=(0,0,2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 136,
   "metadata": {},
   "outputs": [],
   "source": [
    "#fit the model\n",
    "model_fit = model.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 137,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                              ARMA Model Results                              \n",
      "==============================================================================\n",
      "Dep. Variable:                      y   No. Observations:                  121\n",
      "Model:                     ARMA(0, 2)   Log Likelihood                -179.864\n",
      "Method:                       css-mle   S.D. of innovations              1.069\n",
      "Date:                Sun, 29 Mar 2020   AIC                            367.727\n",
      "Time:                        17:33:01   BIC                            378.910\n",
      "Sample:                    09-01-2019   HQIC                           372.269\n",
      "                         - 12-30-2019                                         \n",
      "==============================================================================\n",
      "                 coef    std err          z      P>|z|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "const         50.0478      0.157    318.883      0.000      49.740      50.355\n",
      "ma.L1.y        0.3686      0.088      4.176      0.000       0.196       0.542\n",
      "ma.L2.y        0.2538      0.086      2.949      0.004       0.085       0.422\n",
      "                                    Roots                                    \n",
      "=============================================================================\n",
      "                  Real          Imaginary           Modulus         Frequency\n",
      "-----------------------------------------------------------------------------\n",
      "MA.1           -0.7261           -1.8474j            1.9850           -0.3096\n",
      "MA.2           -0.7261           +1.8474j            1.9850            0.3096\n",
      "-----------------------------------------------------------------------------\n"
     ]
    }
   ],
   "source": [
    "#summary of the model\n",
    "print(model_fit.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Predicted Model:\n",
    "# $\\hat{y}_t = 50 + 0.37\\varepsilon_{t-1} + 0.25\\varepsilon_{t-2}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 138,
   "metadata": {},
   "outputs": [],
   "source": [
    "#get prediction start and end dates\n",
    "pred_start_date = test_data.index[0]\n",
    "pred_end_date = test_data.index[-1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 139,
   "metadata": {},
   "outputs": [],
   "source": [
    "#get the predictions and residuals\n",
    "predictions = model_fit.predict(start=pred_start_date, end=pred_end_date)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 140,
   "metadata": {},
   "outputs": [],
   "source": [
    "residuals = test_data - predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 141,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x1bbb0b50bc8>"
      ]
     },
     "execution_count": 141,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 28860 characters]\n",
      "text/plain": [
       "<Figure size 720x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,4))\n",
    "\n",
    "plt.plot(series[-14:])\n",
    "plt.plot(predictions)\n",
    "\n",
    "plt.legend(('Data', 'Predictions'), fontsize=16)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 142,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Mean Absolute Percent Error: 0.0053\n"
     ]
    }
   ],
   "source": [
    "print('Mean Absolute Percent Error:', round(np.mean(abs(residuals/test_data)),4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 109,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root Mean Squared Error: 1.8254981729090671\n"
     ]
    }
   ],
   "source": [
    "print('Root Mean Squared Error:', np.sqrt(np.mean(residuals**2)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
