{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Tomorrow's Electricity Demand\n",
    "\n",
    "power_demand.py\n",
    "\n",
    "Capstone project: forecasting tomorrow's peak electricity demand for a\n",
    "city, so the distribution company can buy and schedule the right amount\n",
    "of power - no blackouts, and no expensive power bought and wasted.\n",
    "\n",
    "CBSE Class XII, Artificial Intelligence (843), Part C - Capstone Project.\n",
    "SDG 7 (Affordable and Clean Energy) and SDG 13 (Climate Action).\n",
    "\n",
    "Data Science Methodology, 10 steps. REGRESSION on a time series. Three\n",
    "models are compared - Linear Regression, a Random Forest and a neural\n",
    "network (a multi-layer perceptron) - and, because this is a forecast,\n",
    "they are trained on the past and tested on a later year, then checked\n",
    "with time-series cross validation.\n",
    "\n",
    "Libraries: pandas, numpy, matplotlib, scikit-learn, joblib."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 1 - Business Understanding"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"=\" * 64)\n",
    "print(\"TOMORROW'S ELECTRICITY DEMAND\")\n",
    "print(\"=\" * 64)\n",
    "print(\"Problem : power for tomorrow is bought today; a bad guess means\")\n",
    "print(\"          load-shedding or paying for power nobody uses\")\n",
    "print(\"Goal    : forecast tomorrow's peak demand (MW) from weather and calendar\")\n",
    "print()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 2 - Analytic Approach"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Analytic approach : predictive - regression on a time series\")\n",
    "print()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 3 - Data Requirements"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.linear_model import LinearRegression\n",
    "from sklearn.ensemble import RandomForestRegressor\n",
    "from sklearn.neural_network import MLPRegressor\n",
    "from sklearn.pipeline import make_pipeline\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.model_selection import TimeSeriesSplit\n",
    "from sklearn.metrics import mean_absolute_error, mean_squared_error, r2_score\n",
    "import joblib\n",
    "\n",
    "pd.set_option(\"display.width\", 130)\n",
    "pd.set_option(\"display.max_columns\", 20)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 4 - Data Collection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv(\"power_demand.csv\", parse_dates=[\"date\"])\n",
    "print(\"Rows collected :\", len(data))\n",
    "print(\"Dates          :\", data[\"date\"].min().date(), \"to\", data[\"date\"].max().date())\n",
    "print(data.head(6))\n",
    "print()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 5 - Data Understanding"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"--- Step 5: Data Understanding ---\")\n",
    "print(\"Missing values:\")\n",
    "print(data.isnull().sum().to_string())\n",
    "print(\"Days with demand recorded as 0 :\", (data[\"demand_mw\"] == 0).sum())\n",
    "print()\n",
    "print(data.describe().round(1))\n",
    "print()\n",
    "\n",
    "ok = data[data[\"demand_mw\"] > 0]\n",
    "\n",
    "# Chart 1: demand every day for three years.\n",
    "plt.figure(figsize=(10, 4))\n",
    "plt.plot(ok[\"date\"], ok[\"demand_mw\"], color=\"#f28e2b\", linewidth=0.8)\n",
    "plt.title(\"Daily peak demand, 2022-2024\")\n",
    "plt.xlabel(\"Date\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.grid(True, linestyle=\"--\", alpha=0.3)\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart1_demand.png\")\n",
    "plt.close()\n",
    "\n",
    "# Chart 2: temperature against demand.\n",
    "plt.figure(figsize=(7, 5))\n",
    "plt.scatter(ok[\"max_temp_c\"], ok[\"demand_mw\"], s=8, alpha=0.5, color=\"#e15759\")\n",
    "plt.title(\"Maximum temperature and demand\")\n",
    "plt.xlabel(\"Maximum temperature (deg C)\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.grid(True, linestyle=\"--\", alpha=0.3)\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart2_temperature.png\")\n",
    "plt.close()\n",
    "\n",
    "# Chart 3: average demand on each day of the week.\n",
    "weekday = ok.groupby(ok[\"date\"].dt.dayofweek)[\"demand_mw\"].mean()\n",
    "names = [\"Mon\", \"Tue\", \"Wed\", \"Thu\", \"Fri\", \"Sat\", \"Sun\"]\n",
    "print(\"Average demand by day of the week (MW):\")\n",
    "for i in range(7):\n",
    "    print(\"  \" + names[i], round(weekday[i]))\n",
    "print()\n",
    "plt.figure(figsize=(7, 4))\n",
    "plt.bar(names, weekday.values, color=[\"#4e79a7\"] * 5 + [\"#9c9c9c\", \"#9c9c9c\"])\n",
    "plt.ylim(weekday.min() - 300, weekday.max() + 100)\n",
    "plt.title(\"Average demand by day of the week\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart3_weekday.png\")\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 6 - Data Preparation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"--- Step 6: Data Preparation ---\")\n",
    "before = len(data)\n",
    "data = data.drop_duplicates(subset=\"date\").sort_values(\"date\").reset_index(drop=True)\n",
    "print(\"Duplicate dates removed        :\", before - len(data))\n",
    "\n",
    "# A demand of 0 is a grid failure, not a real demand: it becomes missing.\n",
    "# Missing days are then filled by interpolating between the days either\n",
    "# side, which suits a smooth day-to-day series.\n",
    "data[\"demand_mw\"] = pd.to_numeric(data[\"demand_mw\"]).replace(0, np.nan)\n",
    "gaps = data[\"demand_mw\"].isnull().sum()\n",
    "data[\"demand_mw\"] = data[\"demand_mw\"].interpolate()\n",
    "print(\"Missing or zero days filled    :\", gaps, \"(interpolated)\")\n",
    "\n",
    "# Feature engineering - all from information known the day before.\n",
    "data[\"day_of_week\"] = data[\"date\"].dt.dayofweek\n",
    "data[\"is_sunday\"] = (data[\"day_of_week\"] == 6).astype(int)\n",
    "data[\"is_saturday\"] = (data[\"day_of_week\"] == 5).astype(int)\n",
    "data[\"month\"] = data[\"date\"].dt.month\n",
    "data[\"cooling_degrees\"] = np.maximum(0, data[\"max_temp_c\"] - 27)       # heat above 27 drives ACs\n",
    "data[\"demand_yesterday\"] = data[\"demand_mw\"].shift(1)\n",
    "data[\"demand_last_week\"] = data[\"demand_mw\"].shift(7)\n",
    "data = data.dropna().reset_index(drop=True)\n",
    "features = [\"max_temp_c\", \"cooling_degrees\", \"humidity_pct\", \"holiday\", \"is_saturday\", \"is_sunday\",\n",
    "            \"month\", \"demand_yesterday\", \"demand_last_week\"]\n",
    "print(\"Rows after making lag features :\", len(data), \"(first 7 days have no 'last week')\")\n",
    "print(\"Features used                  :\", features)\n",
    "print()\n",
    "\n",
    "# Chart 4: demand in each month, as box plots.\n",
    "plt.figure(figsize=(9, 4.3))\n",
    "plt.boxplot([data[data[\"month\"] == m][\"demand_mw\"] for m in range(1, 13)],\n",
    "            labels=[\"J\", \"F\", \"M\", \"A\", \"M\", \"J\", \"J\", \"A\", \"S\", \"O\", \"N\", \"D\"])\n",
    "plt.title(\"Demand in each month (2022-2024)\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.grid(True, axis=\"y\", linestyle=\"--\", alpha=0.3)\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart4_months.png\")\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 7 - Modelling"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"--- Step 7: Modelling ---\")\n",
    "# A forecast must be tested on a future the model has not seen, so the\n",
    "# rows are NOT shuffled: train on 2022-2023, test on 2024.\n",
    "train = data[data[\"date\"].dt.year <= 2023]\n",
    "test = data[data[\"date\"].dt.year == 2024]\n",
    "X_train, y_train = train[features], train[\"demand_mw\"]\n",
    "X_test, y_test = test[features], test[\"demand_mw\"]\n",
    "print(\"Training days (2022-2023) :\", len(train))\n",
    "print(\"Testing days (2024)       :\", len(test))\n",
    "print()\n",
    "\n",
    "models = {\n",
    "    \"Linear Regression\": LinearRegression(),\n",
    "    \"Random Forest\": RandomForestRegressor(n_estimators=300, min_samples_leaf=2, random_state=0),\n",
    "    \"Neural Network\": make_pipeline(StandardScaler(),\n",
    "                                    MLPRegressor(hidden_layer_sizes=(32, 16), max_iter=3000,\n",
    "                                                 learning_rate_init=0.01, random_state=0)),\n",
    "}\n",
    "for name in models:\n",
    "    models[name].fit(X_train, y_train)\n",
    "print(\"Neural network layers : input 9 -> 32 -> 16 -> output 1\")\n",
    "print()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 8 - Evaluation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"--- Step 8: Evaluation ---\")\n",
    "average = y_test.mean()\n",
    "predictions = {}\n",
    "for name in models:\n",
    "    pred = models[name].predict(X_test)\n",
    "    predictions[name] = pred\n",
    "    rmse = np.sqrt(mean_squared_error(y_test, pred))\n",
    "    print(name)\n",
    "    print(\"  \" + name + \" MAE  :\", round(mean_absolute_error(y_test, pred), 1), \"MW\")\n",
    "    print(\"  \" + name + \" RMSE :\", round(rmse, 1), \"MW\")\n",
    "    print(\"  \" + name + \" R2   :\", round(r2_score(y_test, pred), 3))\n",
    "    print(\"  \" + name + \" RMSE as % of average demand :\", round(rmse / average * 100, 2), \"%\")\n",
    "naive = test[\"demand_yesterday\"]\n",
    "print(\"Baseline 'same as yesterday' RMSE :\", round(np.sqrt(mean_squared_error(y_test, naive)), 1), \"MW\")\n",
    "print()\n",
    "\n",
    "# Time-series cross validation: each fold trains on everything before it\n",
    "# and tests on the block that follows - never on the past.\n",
    "tscv = TimeSeriesSplit(n_splits=5)\n",
    "cv_rmse = {name: [] for name in models}\n",
    "X_all, y_all = data[features], data[\"demand_mw\"]\n",
    "for train_idx, test_idx in tscv.split(X_all):\n",
    "    for name in models:\n",
    "        m = {\"Linear Regression\": LinearRegression(),\n",
    "             \"Random Forest\": RandomForestRegressor(n_estimators=300, min_samples_leaf=2, random_state=0),\n",
    "             \"Neural Network\": make_pipeline(StandardScaler(),\n",
    "                                             MLPRegressor(hidden_layer_sizes=(32, 16), max_iter=3000,\n",
    "                                                          learning_rate_init=0.01, random_state=0))}[name]\n",
    "        m.fit(X_all.iloc[train_idx], y_all.iloc[train_idx])\n",
    "        p = m.predict(X_all.iloc[test_idx])\n",
    "        cv_rmse[name].append(np.sqrt(mean_squared_error(y_all.iloc[test_idx], p)))\n",
    "print(\"Time-series cross validation, RMSE (MW) in each fold:\")\n",
    "for name in models:\n",
    "    print(\"  \" + name.ljust(18), np.round(cv_rmse[name], 1), \" mean\", round(np.mean(cv_rmse[name]), 1))\n",
    "print()\n",
    "\n",
    "# Chart 5: actual and predicted demand through 2024.\n",
    "plt.figure(figsize=(10, 4.3))\n",
    "plt.plot(test[\"date\"], y_test, color=\"#999999\", linewidth=1, label=\"Actual\")\n",
    "plt.plot(test[\"date\"], predictions[\"Linear Regression\"], color=\"#4e79a7\", linewidth=1, label=\"Linear Regression\")\n",
    "plt.title(\"2024: actual demand and the Linear Regression forecast\")\n",
    "plt.xlabel(\"Date\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.legend()\n",
    "plt.grid(True, linestyle=\"--\", alpha=0.3)\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart5_forecast_2024.png\")\n",
    "plt.close()\n",
    "\n",
    "# Chart 6: error in each time-series fold for the three models.\n",
    "x = np.arange(5)\n",
    "plt.figure(figsize=(8, 4.3))\n",
    "plt.bar(x - 0.27, cv_rmse[\"Linear Regression\"], width=0.27, label=\"Linear Regression\", color=\"#4e79a7\")\n",
    "plt.bar(x, cv_rmse[\"Random Forest\"], width=0.27, label=\"Random Forest\", color=\"#59a14f\")\n",
    "plt.bar(x + 0.27, cv_rmse[\"Neural Network\"], width=0.27, label=\"Neural Network\", color=\"#b07aa1\")\n",
    "plt.xticks(x, [\"Fold 1\", \"Fold 2\", \"Fold 3\", \"Fold 4\", \"Fold 5\"])\n",
    "plt.ylabel(\"RMSE (MW)\")\n",
    "plt.title(\"Time-series cross validation\")\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart6_tscv.png\")\n",
    "plt.close()\n",
    "\n",
    "# Chart 7: what the linear model's coefficients mean, per unit.\n",
    "coef = pd.Series(models[\"Linear Regression\"].coef_, index=features)\n",
    "print(\"Linear Regression coefficients:\")\n",
    "print(coef.round(2).to_string())\n",
    "print()\n",
    "summer = test[(test[\"date\"].dt.month >= 4) & (test[\"date\"].dt.month <= 6)]\n",
    "plt.figure(figsize=(10, 4.3))\n",
    "plt.plot(summer[\"date\"], summer[\"demand_mw\"], color=\"#999999\", label=\"Actual\")\n",
    "for name, col in [[\"Linear Regression\", \"#4e79a7\"], [\"Random Forest\", \"#59a14f\"], [\"Neural Network\", \"#b07aa1\"]]:\n",
    "    plt.plot(summer[\"date\"], models[name].predict(summer[features]), color=col, linewidth=1, label=name)\n",
    "plt.title(\"Summer 2024 close up: April to June\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.legend()\n",
    "plt.grid(True, linestyle=\"--\", alpha=0.3)\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart7_summer.png\")\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 9 - Deployment"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"--- Step 9: Deployment ---\")\n",
    "best = models[\"Linear Regression\"]\n",
    "joblib.dump(best, \"demand_model.pkl\")\n",
    "model = joblib.load(\"demand_model.pkl\")\n",
    "last = data.iloc[-1]\n",
    "last_week = data.iloc[-7]\n",
    "tomorrow = pd.DataFrame({\"max_temp_c\": [24.0], \"cooling_degrees\": [0.0], \"humidity_pct\": [58], \"holiday\": [0],\n",
    "                         \"is_saturday\": [0], \"is_sunday\": [0], \"month\": [1],\n",
    "                         \"demand_yesterday\": [last[\"demand_mw\"]], \"demand_last_week\": [last_week[\"demand_mw\"]]})\n",
    "print(\"Forecast for 1 January 2025 (24 deg C forecast) :\", round(model.predict(tomorrow)[0]), \"MW\")\n",
    "print()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 10 - Feedback"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"--- Step 10: Feedback ---\")\n",
    "# Each day of January 2025 is forecast the evening before, then compared\n",
    "# with what the meters actually recorded.\n",
    "jan = pd.read_csv(\"january_2025.csv\", parse_dates=[\"date\"])\n",
    "history = pd.concat([data[[\"date\", \"demand_mw\"]], jan[[\"date\", \"demand_mw\"]]]).reset_index(drop=True)\n",
    "jan[\"day_of_week\"] = jan[\"date\"].dt.dayofweek\n",
    "jan[\"is_sunday\"] = (jan[\"day_of_week\"] == 6).astype(int)\n",
    "jan[\"is_saturday\"] = (jan[\"day_of_week\"] == 5).astype(int)\n",
    "jan[\"month\"] = 1\n",
    "jan[\"cooling_degrees\"] = np.maximum(0, jan[\"max_temp_c\"] - 27)\n",
    "jan[\"demand_yesterday\"] = history[\"demand_mw\"].shift(1).iloc[-len(jan):].values\n",
    "jan[\"demand_last_week\"] = history[\"demand_mw\"].shift(7).iloc[-len(jan):].values\n",
    "jan_pred = model.predict(jan[features])\n",
    "jan_rmse = np.sqrt(mean_squared_error(jan[\"demand_mw\"], jan_pred))\n",
    "test_rmse = np.sqrt(mean_squared_error(y_test, predictions[\"Linear Regression\"]))\n",
    "jan_2024 = test[test[\"date\"].dt.month == 1]\n",
    "jan_2024_rmse = np.sqrt(mean_squared_error(jan_2024[\"demand_mw\"], model.predict(jan_2024[features])))\n",
    "print(\"Days forecast in January 2025     :\", len(jan))\n",
    "print(\"RMSE, January 2025                :\", round(jan_rmse, 1), \"MW\")\n",
    "print(\"RMSE, January 2024 (test year)    :\", round(jan_2024_rmse, 1), \"MW\")\n",
    "print(\"RMSE, all of 2024                 :\", round(test_rmse, 1), \"MW\")\n",
    "print(\"Average error in January 2025     :\", round((jan[\"demand_mw\"] - jan_pred).mean(), 1),\n",
    "      \"MW (actual minus forecast)\")\n",
    "if jan_rmse > 1.25 * jan_2024_rmse:\n",
    "    print(\"Feedback: January 2025 errors are well above last January's - retrain.\")\n",
    "else:\n",
    "    print(\"Feedback: January 2025 is in line with last January - keep the model, retrain each quarter.\")\n",
    "\n",
    "# Chart 8: January 2025, forecast against actual.\n",
    "plt.figure(figsize=(9, 4))\n",
    "plt.plot(jan[\"date\"], jan[\"demand_mw\"], marker=\"o\", markersize=3, color=\"#999999\", label=\"Actual\")\n",
    "plt.plot(jan[\"date\"], jan_pred, marker=\"s\", markersize=3, color=\"#4e79a7\", label=\"Forecast the day before\")\n",
    "plt.title(\"Feedback: January 2025, day-ahead forecasts against the meters\")\n",
    "plt.ylabel(\"Demand (MW)\")\n",
    "plt.legend()\n",
    "plt.grid(True, linestyle=\"--\", alpha=0.3)\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"chart8_january_2025.png\")\n",
    "plt.close()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}