{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import seaborn as sns\n",
    "from matplotlib import pyplot as plt\n",
    "from statsmodels.tsa.stattools import adfuller\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from statsmodels.tsa.stattools import adfuller"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "series=pd.read_csv(\"daily-total-female-births.csv\", index_col=0, parse_dates=True, squeeze=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "series.head(5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "series.tail()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "series.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One Year data recorded on daily basis  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "series.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "series.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. Check Stationarity "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "result = adfuller(series)\n",
    "print('ADF Statistic: %f' % result[0])\n",
    "print('p-value: %f' % result[1])\n",
    "print('Critical Values:')\n",
    "for key, value in result[4].items():\n",
    "    print('\\t%s: %.3f' % (key, value))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Data is stationary"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2. Persistence Model "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def persistence(x):\n",
    "    return x\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dataframe=pd.DataFrame(series.values)\n",
    "dataframe=pd.concat([series.shift(1),series],axis=1)\n",
    "dataframe.columns=['t+1','t']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dataframe.head(5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# split into train and test sets\n",
    "X=dataframe.values\n",
    "train_size=int(len(X)*0.66)\n",
    "train, test=X[1:train_size], X[train_size:]\n",
    "train_X, train_y=train[:,0],train[:,1]\n",
    "test_X,test_y=test[:,0],test[:,1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predictions=[x for x in test_X]\n",
    "rmse=np.sqrt(mean_squared_error(test_y,predictions))\n",
    "print('Test RMSE: %.3f' %rmse)\n",
    "residuals=[test_y[i]-predictions[i]for i in range(len(predictions))]\n",
    "residuals=pd.DataFrame(residuals)\n",
    "print(residuals.head())\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "RMSE on test data is 9.151"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "3. Auto Regression ( AR)  of Residual Errors "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from statsmodels.tsa.ar_model import AR"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Training a persistence model on Training data , calcuting residuals and training an AR model on residuals (MA Model)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# persistence model on training set\n",
    "train_pred = [x for x in train_X]\n",
    "# calculate residuals\n",
    "train_resid = [train_y[i]-train_pred[i] for i in range(len(train_pred))]\n",
    "# model the training set residuals\n",
    "model = AR(train_resid)\n",
    "model_fit = model.fit()\n",
    "window = model_fit.k_ar\n",
    "coef = model_fit.params\n",
    "print('Lag=%d, Coef=%s' % (window, coef))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "window"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# walk forward over time steps in test\n",
    "history = train_resid[len(train_resid)-window:]\n",
    "history = [history[i] for i in range(len(history))]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "history"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "4. Apply Persistence Model on Test Data "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# walk forward over time steps in test\n",
    "history = train_resid[len(train_resid)-window:]\n",
    "history = [history[i] for i in range(len(history))]\n",
    "predictions = list()\n",
    "expected_error = list()\n",
    "for t in range(len(test_y)):\n",
    "# persistence\n",
    "    yhat = test_X[t]\n",
    "    error = test_y[t] - yhat\n",
    "    expected_error.append(error)\n",
    "# predict error\n",
    "    length = len(history)\n",
    "    lag = [history[i] for i in range(length-window,length)]\n",
    "    pred_error = coef[0]\n",
    "    for d in range(window):\n",
    "        pred_error += coef[d+1] * lag[window-d-1]\n",
    "    predictions.append(pred_error)\n",
    "    history.append(error)\n",
    "    print('predicted error=%f, expected error=%f' % (pred_error, error))\n",
    "# plot predicted error\n",
    "plt.plot(expected_error)\n",
    "plt.plot(predictions, color='red')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "5. Correct the Predictions using Residual Forecast ( MA Model) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# walk forward over time steps in test\n",
    "history = train_resid[len(train_resid)-window:]\n",
    "history = [history[i] for i in range(len(history))]\n",
    "predictions = list()\n",
    "expected_error = list()\n",
    "for t in range(len(test_y)):\n",
    "# persistence\n",
    "    yhat = test_X[t]\n",
    "    error = test_y[t] - yhat\n",
    "    expected_error.append(error)\n",
    "# predict error\n",
    "    length = len(history)\n",
    "    lag = [history[i] for i in range(length-window,length)]\n",
    "    pred_error = coef[0]\n",
    "    for d in range(window):\n",
    "        pred_error += coef[d+1] * lag[window-d-1]\n",
    "    # correct the prediction\n",
    "    yhat = yhat + pred_error\n",
    "    predictions.append(yhat)\n",
    "    history.append(error)\n",
    "    print('predicted error=%f, expected error=%f' % (pred_error, error))\n",
    "rmse = sqrt(mean_squared_error(test_y, predictions))\n",
    "print('Test RMSE: %.3f' % rmse)\n",
    "# plot predicted error\n",
    "plt.plot(test_y)\n",
    "plt.plot(predictions, color='red')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Part 2: Random Walk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from random import random\n",
    "random()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create and plot a random walk\n",
    "from random import seed\n",
    "from random import random\n",
    "from matplotlib import pyplot\n",
    "seed(1)\n",
    "random_walk = list()\n",
    "random_walk.append(-1 if random() < 0.5 else 1)\n",
    "for i in range(1, 1000):\n",
    "    movement = -1 if random() < 0.5 else 1\n",
    "    value = random_walk[i-1] + movement\n",
    "    random_walk.append(value)\n",
    "pyplot.figure(figsize=(10,6))\n",
    "pyplot.plot(random_walk)\n",
    "pyplot.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pandas.plotting import autocorrelation_plot\n",
    "pyplot.figure(figsize=(10,6))\n",
    "autocorrelation_plot(random_walk)\n",
    "pyplot.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# statistical test\n",
    "result = adfuller(random_walk)\n",
    "print('ADF Statistic: %f' % result[0])\n",
    "print('p-value: %f' % result[1])\n",
    "print('Critical Values:')\n",
    "for key, value in result[4].items():\n",
    "    print('\\t%s: %.3f' % (key, value))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The null hypothesis of the test is that the time series is non-stationary. Running the example,\n",
    "we can see that the test statistic value was 0.341605. This is larger than all of the critical values\n",
    "at the 1%, 5%, and 10% con\f",
    "dence levels. Therefore, we can say that the time series does\n",
    "appear to be non-stationary with a low likelihood of the result being a statistical \n",
    "fluke."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# take difference\n",
    "diff = list()\n",
    "for i in range(1, len(random_walk)):\n",
    "    value = random_walk[i] - random_walk[i - 1]\n",
    "    diff.append(value)\n",
    "# line plot\n",
    "pyplot.figure(figsize=(10,6))\n",
    "pyplot.plot(diff)\n",
    "pyplot.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pyplot.figure(figsize=(10,6))\n",
    "autocorrelation_plot(diff)\n",
    "pyplot.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Takin the difference on Random walk leaves no meaningful structure in the time series ( A white noise). It means that his time sereis cannot be predicted."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Predicting a Random Walk : Persistence Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# persistence forecasts for a random walk\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from math import sqrt\n",
    "\n",
    "train_size = int(len(random_walk) * 0.66)\n",
    "train, test = random_walk[0:train_size], random_walk[train_size:]\n",
    "# persistence\n",
    "predictions = list()\n",
    "history = train[-1]\n",
    "for i in range(len(test)):\n",
    "    yhat = history\n",
    "    predictions.append(yhat)\n",
    "    history = test[i]\n",
    "rmse = sqrt(mean_squared_error(test, predictions))\n",
    "print('Persistence RMSE: %.3f' % rmse)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Adding Variance as was addded in Random Walk generation process"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predictions = list()\n",
    "history = train[-1]\n",
    "for i in range(len(test)):\n",
    "    yhat = history + (-1 if random() < 0.5 else 1)\n",
    "    predictions.append(yhat)\n",
    "    history = test[i]\n",
    "rmse = sqrt(mean_squared_error(test, predictions))\n",
    "print(rmse)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Conclusion : Persistence or Naive Forecast is the best performance that can be achieved on a Random Walk"
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
