Improved policy executer

Result-Reports/Policies.md

@@ -141,7 +141,7 @@ Test data: 01-01-2023 until 08-10–2023
 
 # TODO
 - [ ] Baseline 
-- [ ] Profit penalty parameter als over charge cycles voor een dag -> parameter bepalen op training data (convex probleem) (< 400 charge cycles per jaar) (over een dag kijken hoeveel charge cycles -> profit - penalty * charge cycles erover, (misschien belonen als eronder charge cycles))
+- [x] Profit penalty parameter als over charge cycles voor een dag -> parameter bepalen op training data (convex probleem) (< 400 charge cycles per jaar) (over een dag kijken hoeveel charge cycles -> profit - penalty * charge cycles erover, (misschien belonen als eronder charge cycles))
 
 - [ ] Meer verschil bekijken tussen GRU en diffusion
 - [ ] Andere lagen voor diffusion model (GRU, kijken naar TSDiff)

src/notebooks/diffusion-training.ipynb

@@ -2,12 +2,13 @@
  "cells": [
   {
    "cell_type": "code",
-   "execution_count": 1,
+   "execution_count": 3,
    "metadata": {},
    "outputs": [],
    "source": [
     "import sys\n",
-    "sys.path.append('../..')"
+    "sys.path.append('../..')\n",
+    "import torch"
    ]
   },
   {
@@ -68,57 +69,30 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 4,
+   "execution_count": 5,
    "metadata": {},
    "outputs": [
     {
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
+      "ClearML Task: created new task id=b71216825809432682ea3c7841c07612\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
+      "ClearML results page: http://192.168.1.182:8080/projects/2e46d4af6f1e4c399cf9f5aa30bc8795/experiments/b71216825809432682ea3c7841c07612/output/log\n"
-      "      dtype='object')\n"
      ]
     },
     {
      "name": "stderr",
      "output_type": "stream",
      "text": [
-      "InsecureRequestWarning: Certificate verification is disabled! Adding certificate verification is strongly advised. See: https://urllib3.readthedocs.io/en/latest/advanced-usage.html#ssl-warnings\n"
+      "500 model found when searching for `file:///workspaces/Thesis/src/notebooks/checkpoint.pt`\n",
-     ]
+      "Selected model `Autoregressive Non Linear Quantile Regression + Quarter + DoW + Net` (id=bc0cb0d7fc614e2e8b0edf5b85348646)\n"
-    },
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "ClearML Task: created new task id=348145474a1140a6bdf8c81553a358b2\n",
-      "ClearML results page: http://192.168.1.182:8080/projects/2e46d4af6f1e4c399cf9f5aa30bc8795/experiments/348145474a1140a6bdf8c81553a358b2/output/log\n",
-      "2023-12-28 20:57:29,259 - clearml.Task - INFO - Storing jupyter notebook directly as code\n",
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
-      "      dtype='object')\n",
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
-      "      dtype='object')\n",
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
-      "      dtype='object')\n",
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
-      "      dtype='object')\n",
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
-      "      dtype='object')\n",
-      "Index(['datetime', 'nrv', 'load_forecast', 'total_load', 'wind_forecast',\n",
-      "       'wind_history', 'nominal_net_position', 'quarter', 'day_of_week'],\n",
-      "      dtype='object')\n"
      ]
     }
    ],
    "source": [
     "inputDim = data_processor.get_input_size()\n",
     "learningRate = 0.0001\n",
-    "epochs=5000\n",
+    "epochs=150\n",
     "\n",
     "#### Model ####\n",
     "model = SimpleDiffusionModel(96, [512, 512, 512], other_inputs_dim=inputDim[1], time_dim=64)\n",
@@ -138,6 +112,116 @@
    "outputs": [],
    "source": []
   },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "<class 'src.models.diffusion_model.SimpleDiffusionModel'>\n"
+     ]
+    }
+   ],
+   "source": [
+    "new_model = torch.load(\"checkpoint.pt\")\n",
+    "print(type(new_model))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 8,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Determine threshold based on predictions\n",
+    "from src.models.diffusion_model import DiffusionModel\n",
+    "\n",
+    "\n",
+    "def get_predicted_NRV(date):\n",
+    "    idx = test_loader.dataset.get_idx_for_date(date.date())\n",
+    "    initial, _, samples, target = auto_regressive(test_loader.dataset, [idx]*500, 96)\n",
+    "    samples = samples.cpu().numpy()\n",
+    "    target = target.cpu().numpy()\n",
+    "\n",
+    "    # inverse using data_processor\n",
+    "    samples = data_processor.inverse_transform(samples)\n",
+    "    target = data_processor.inverse_transform(target)\n",
+    "\n",
+    "    return initial.cpu().numpy()[0][-1], samples, target\n",
+    "\n",
+    "device = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n",
+    "\n",
+    "def sample_diffusion(model: DiffusionModel, n: int, inputs: torch.tensor):\n",
+    "    noise_steps = 1000\n",
+    "    beta_start = 1e-4\n",
+    "    beta_end = 0.02\n",
+    "    ts_length = 96\n",
+    "    \n",
+    "    beta = torch.linspace(beta_start, beta_end, noise_steps).to(device)\n",
+    "    alpha = 1. - beta\n",
+    "    alpha_hat = torch.cumprod(alpha, dim=0)\n",
+    "\n",
+    "    inputs = inputs.repeat(n, 1).to(device)\n",
+    "    model.eval()\n",
+    "    with torch.no_grad():\n",
+    "        x = torch.randn(inputs.shape[0], ts_length).to(device)\n",
+    "        for i in reversed(range(1, noise_steps)):\n",
+    "            t = (torch.ones(inputs.shape[0]) * i).long().to(device)\n",
+    "            predicted_noise = model(x, t, inputs)\n",
+    "            _alpha = alpha[t][:, None]\n",
+    "            _alpha_hat = alpha_hat[t][:, None]\n",
+    "            _beta = beta[t][:, None]\n",
+    "\n",
+    "            if i > 1:\n",
+    "                noise = torch.randn_like(x)\n",
+    "            else:\n",
+    "                noise = torch.zeros_like(x)\n",
+    "\n",
+    "            x = 1/torch.sqrt(_alpha) * (x-((1-_alpha) / (torch.sqrt(1 - _alpha_hat))) * predicted_noise) + torch.sqrt(_beta) * noise\n",
+    "    return x\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 27,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "tensor([[-178.8835,  -47.2518, -103.9158,   -9.8302,   15.9751,  138.9138,\n",
+       "          -56.8392, -128.0629, -128.3637,  -83.1066,   56.6656, -200.4618,\n",
+       "           10.8563, -146.4262,  120.4816,  -60.1130,  -18.7972, -214.0427,\n",
+       "          148.1229,  136.0194,   33.7580,   85.7884, -164.5678,   53.8879,\n",
+       "          187.6217,  -77.5978,  153.7462, -129.1419, -149.8551,  118.4640,\n",
+       "          -29.4688,  -37.3348, -104.4318,  -16.1735,  -29.9716,   -1.4205,\n",
+       "         -130.6785,   23.8387,   75.6755,  113.8617,  -61.4832,  -81.3838,\n",
+       "          -15.3194,  -63.5703,  215.4112,    8.0719,   26.4597,   72.4347,\n",
+       "          -23.1216,   44.8453,  -12.2994,   94.7612, -162.2193,   18.0694,\n",
+       "           31.2402,   78.6964,   35.1892, -105.0744,   38.7805,  -27.5867,\n",
+       "           39.5985,  136.5500, -179.8039,  231.9039,  116.1411, -226.0043,\n",
+       "         -149.2595,  -14.5097,  123.5570,  162.4510,  -62.9467,  -82.3552,\n",
+       "          187.5180,   12.3145, -189.3492, -159.3642, -144.8646,  130.9768,\n",
+       "          -79.4541,   53.5424,   35.7119,  134.5416,  -87.5582,   70.4020,\n",
+       "          -44.0516,  111.3181,   17.0087,  -14.9322, -187.4202,  -41.7765,\n",
+       "           11.2264,  221.0164, -106.3083, -123.9814,  -12.2132, -121.7845]],\n",
+       "       device='cuda:0')"
+      ]
+     },
+     "execution_count": 27,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "inputs = torch.randn(1, 672).to(device)\n",
+    "sample_diffusion(new_model, 1, inputs)"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,

src/policies/policy_executer.py

@@ -9,6 +9,7 @@ import datetime
 from tqdm import tqdm
 from src.utils.imbalance_price_calculator import ImbalancePriceCalculator
 import time
+import plotly.express as px
 
 ### import functions ###
 from src.trainers.quantile_trainer import auto_regressive as quantile_auto_regressive
@@ -19,10 +20,12 @@ from src.utils.clearml import ClearMLHelper
 parser = argparse.ArgumentParser()
 parser.add_argument('--task_id', type=str, default=None)
 parser.add_argument('--model_type', type=str, default=None)
+parser.add_argument('--model_name', type=str, default=None)
 args = parser.parse_args()
 
 assert args.task_id is not None, "Please specify task id"
 assert args.model_type is not None, "Please specify model type"
+assert args.model_name is not None, "Please specify model name"
 
 battery = Battery(2, 1)
 baseline_policy = BaselinePolicy(battery, data_path="")
@@ -43,22 +46,28 @@ def load_model(task_id: str):
     """
     task = Task.get_task(task_id=task_id)
 
+    lstm = task.get_parameter("data_processor/lstm")
+    full_day_skip = task.get_parameter("data_processor/full_day_skip")
+    output_size = int(task.get_parameter("data_processor/output_size"))
+
+    print(f"lstm: {lstm}")
+    print(f"full_day_skip: {full_day_skip}")
+    print(f"output_size: {output_size}")
+
     configuration = task.get_parameters_as_dict()
     data_features = configuration['data_features']
 
     ### Data Config ###
     data_config = DataConfig()
     for key, value in data_features.items():
-        setattr(data_config, key, bool(value))
+        setattr(data_config, key, value == "True")
-    data_config.PV_FORECAST = False
+    print(data_config.__dict__)
-    data_config.PV_HISTORY = False
-    data_config.QUARTER = False
-    data_config.DAY_OF_WEEK = False
 
     ### Data Processor ###
-    data_processor = DataProcessor(data_config, path="", lstm=False)
+    data_processor = DataProcessor(data_config, path="", lstm=lstm=="True")
     data_processor.set_batch_size(8192)
-    data_processor.set_full_day_skip(True)
+    data_processor.set_full_day_skip(full_day_skip == "True")
+    data_processor.set_output_size(int(output_size))
 
     ### Model ###
     output_model_id = task.output_models_id["checkpoint"]
@@ -72,7 +81,7 @@ def load_model(task_id: str):
     model.eval()
 
     _, test_loader = data_processor.get_dataloaders(
-        predict_sequence_length=96
+        predict_sequence_length=output_size
     )
 
     return configuration, model, data_processor, test_loader
@@ -80,7 +89,7 @@ def load_model(task_id: str):
 
 def quantile_auto_regressive_predicted_NRV(model, date, data_processor, test_loader):
     idx = test_loader.dataset.get_idx_for_date(date.date())
-    initial, _, samples, target = quantile_auto_regressive(test_loader.dataset, model, [idx]*500, 96)
+    initial, _, samples, target = quantile_auto_regressive(test_loader.dataset, model, model.quantiles, [idx]*500, 96)
     samples = samples.cpu().numpy()
     target = target.cpu().numpy()
 
@@ -147,7 +156,7 @@ def get_next_day_profits_for_date(model, data_processor, test_loader, date, ipc,
     return predicted_nrv_profits_cycles, baseline_profits_cycles
 
 def next_day_test_set(model, data_processor, test_loader, ipc, predict_NRV: callable):
-    penalties = [0, 10, 50, 150, 250, 350, 500]
+    penalties = [0, 10, 50, 150, 300, 500, 600, 800, 1000, 1500, 2000, 2500]
     predicted_nrv_profits_cycles = {i: [0, 0] for i in penalties}
     baseline_profits_cycles = {i: [0, 0] for i in penalties}
 
@@ -169,7 +178,6 @@ def next_day_test_set(model, data_processor, test_loader, ipc, predict_NRV: call
                 baseline_profits_cycles[penalty][1] += new_baseline_profits_cycles[penalty][1]
 
         except Exception as e:
-            # raise e
             # print(f"Error for date {date}")
             continue
 
@@ -179,14 +187,16 @@ def main():
     clearml_helper = ClearMLHelper(project_name="Thesis/NrvForecast")
     task = clearml_helper.get_task(task_name="Policy Test")
 
-    task.connect(args, name="Arguments")
     task.execute_remotely(queue_name="default", exit_process=True)
 
     configuration, model, data_processor, test_loader = load_model(args.task_id)
 
-    if args.model_type == "quantile":
+    if args.model_type == "autoregressive_quantile":
+        quantiles = configuration["general"]["quantiles"]
+        quantiles = list(map(float, quantiles.strip('[]').split(',')))
+        model.quantiles = quantiles
         predict_NRV = quantile_auto_regressive_predicted_NRV
-        task.add_tags(["quantile"])
+        task.add_tags(["autoregressive_quantile"])
     elif args.model_type == "diffusion":
         predict_NRV = diffusion_predicted_NRV
         task.add_tags(["diffusion"])
@@ -203,7 +213,7 @@ def main():
     # use concat
     for penalty in predicted_nrv_profits_cycles.keys():
         new_rows = pd.DataFrame({
-            "name": [args.model_type, "baseline"],
+            "name": [f"{args.model_type} ({args.model_name})", "baseline"],
             "penalty": [penalty, penalty],
             "profit": [predicted_nrv_profits_cycles[penalty][0], baseline_profits_cycles[penalty][0]],
             "cycles": [predicted_nrv_profits_cycles[penalty][1], baseline_profits_cycles[penalty][1]]
@@ -214,7 +224,7 @@ def main():
     df = df.sort_values(by=["name", "penalty"])
 
     # calculate profit per cycle
-    df["profit_per_cycle"] = df["profit"] / df["cycles"]
+    df["profit_per_cycle"] = df.apply(lambda row: row["profit"] / row["cycles"] if row["cycles"] != 0 else 0, axis=1)
 
     task.get_logger().report_table(
         "Policy Results", 
@@ -223,6 +233,28 @@ def main():
         table_plot=df
     )
 
+    # plotly to show profit on y axis and cycles on x axis (show 2 lines, one for each model)
+    fig = px.line(
+        df, 
+        x="cycles", 
+        y="profit", 
+        color="name", 
+        title="Profit vs. Cycles for Each Model",
+        labels={"cycles": "Cycles", "profit": "Profit"},
+        markers=True,  # Adds markers to the lines
+        hover_data=["penalty"]  # Adds additional hover information
+    )
+
+    fig.update_xaxes(autorange="reversed")
+
+
+    task.get_logger().report_plotly(
+        "Policy Results", 
+        "Profit vs. Cycles for Each Model", 
+        iteration=0, 
+        figure=fig
+    )
+
     # close task
     task.close()
 

src/trainers/quantile_trainer.py

@@ -9,56 +9,38 @@ from src.losses import PinballLoss, NonAutoRegressivePinballLoss, CRPSLoss
 import plotly.graph_objects as go
 import numpy as np
 import matplotlib.pyplot as plt
+from scipy.interpolate import CubicSpline
 
 
-def sample_from_dist(quantiles, output_values):
+def sample_from_dist(quantiles, preds):
-    # check if tensor:
+    if isinstance(preds, torch.Tensor):
-    if isinstance(quantiles, torch.Tensor):
+        preds = preds.detach().cpu()
-        quantiles = quantiles.cpu().numpy()
 
-    if isinstance(output_values, torch.Tensor):
+    # if preds more than 2 dimensions, flatten to 2
-        output_values = output_values.cpu().numpy()
+    if len(preds.shape) > 2:
+        preds = preds.reshape(-1, preds.shape[-1])
+        # target will be reshaped from (1024, 96, 15) to (1024*96, 15)
+        # our target (1024, 96) also needs to be reshaped to (1024*96, 1)
+        target = target.reshape(-1, 1)
 
-    if isinstance(quantiles, list):
+    # preds and target as numpy
-        quantiles = np.array(quantiles)
+    preds = preds.numpy()
 
-    reshaped_values = output_values.reshape(-1, len(quantiles))
+    # random probabilities of (1000, 1)
+    import random
+    probs = np.array([random.random() for _ in range(1000)])
 
-    uniform_random_numbers = np.random.uniform(0, 1, (reshaped_values.shape[0], 1000))
+    spline = CubicSpline(quantiles, preds, axis=1)
     
-    idx_below = np.searchsorted(quantiles, uniform_random_numbers, side="right") - 1
+    samples = spline(probs)
-    idx_above = np.clip(idx_below + 1, 0, len(quantiles) - 1)
 
-    # handle edge case where idx_below is -1
+    # get the diagonal
-    idx_below = np.clip(idx_below, 0, len(quantiles) - 1)
+    samples = np.diag(samples)
 
-    y_below = reshaped_values[np.arange(reshaped_values.shape[0])[:, None], idx_below]
+    return samples
-    y_above = reshaped_values[np.arange(reshaped_values.shape[0])[:, None], idx_above]
 
-    # Calculate the slopes for interpolation
+def auto_regressive(dataset, model, quantiles, idx_batch, sequence_length: int = 96):
-    x_below = quantiles[idx_below]
+    device = next(model.parameters()).device
-    x_above = quantiles[idx_above]
-
-    # Interpolate
-    # Ensure all variables are NumPy arrays
-    x_below_np = x_below.cpu().numpy() if isinstance(x_below, torch.Tensor) else x_below
-    x_above_np = x_above.cpu().numpy() if isinstance(x_above, torch.Tensor) else x_above
-    y_below_np = y_below.cpu().numpy() if isinstance(y_below, torch.Tensor) else y_below
-    y_above_np = y_above.cpu().numpy() if isinstance(y_above, torch.Tensor) else y_above
-
-    # Compute slopes for interpolation
-    slopes_np = (y_above_np - y_below_np) / (
-        np.clip(x_above_np - x_below_np, 1e-6, np.inf)
-    )
-
-    # Perform the interpolation
-    new_samples = y_below_np + slopes_np * (uniform_random_numbers - x_below_np)
-
-    # Return the mean of the samples
-    return np.mean(new_samples, axis=1)
-
-def auto_regressive(dataset, model, idx_batch, sequence_length: int = 96):
-    device = model.device
     prev_features, targets = dataset.get_batch(idx_batch)
     prev_features = prev_features.to(device)
     targets = targets.to(device)
@@ -72,7 +54,7 @@ def auto_regressive(dataset, model, idx_batch, sequence_length: int = 96):
     with torch.no_grad():
         new_predictions_full = model(prev_features)  # (batch_size, quantiles)
     samples = (
-        torch.tensor(sample_from_dist( new_predictions_full))
+        torch.tensor(sample_from_dist(quantiles, new_predictions_full))
         .unsqueeze(1)
         .to(device)
     )  # (batch_size, 1)
@@ -125,7 +107,7 @@ def auto_regressive(dataset, model, idx_batch, sequence_length: int = 96):
         )  # (batch_size, sequence_length, quantiles)
 
         samples = (
-            torch.tensor(sample_from_dist(new_predictions_full))
+            torch.tensor(sample_from_dist(quantiles, new_predictions_full))
             .unsqueeze(-1)
             .to(device)
         )  # (batch_size, 1)
@@ -353,7 +335,7 @@ class AutoRegressiveQuantileTrainer(AutoRegressiveTrainer):
         return fig
 
     def auto_regressive(self, dataset, idx_batch, sequence_length: int = 96):
-        return auto_regressive(dataset, self.model, idx_batch, sequence_length)
+        return auto_regressive(dataset, self.model, self.quantiles, idx_batch, sequence_length)
 
     def plot_quantile_percentages(
         self, task, data_loader, train: bool = True, iteration: int = None, full_day: bool = False

src/trainers/trainer.py

@@ -249,7 +249,7 @@ class Trainer:
 
     def finish_training(self, task):
         if self.best_score is not None:
-            self.model.load_state_dict(torch.load("checkpoint.pt"))
+            self.model = torch.load("checkpoint.pt")
         self.model.eval()
 
         

src/training_scripts/autoregressive_quantiles.py

@@ -13,7 +13,7 @@ from src.models.time_embedding_layer import TimeEmbedding
 
 #### ClearML ####
 clearml_helper = ClearMLHelper(project_name="Thesis/NrvForecast")
-task = clearml_helper.get_task(task_name="Autoregressive Quantile Regression: GRU + Quarter + DoW + Load + Wind + Net")
+task = clearml_helper.get_task(task_name="Autoregressive Quantile Regression: Non Linear + Quarter + DoW + Load + Wind + Net")
 
 
 #### Data Processor ####
@@ -35,7 +35,7 @@ data_config.NOMINAL_NET_POSITION = True
 
 data_config = task.connect(data_config, name="data_features")
 
-data_processor = DataProcessor(data_config, path="", lstm=True)
+data_processor = DataProcessor(data_config, path="", lstm=False)
 data_processor.set_batch_size(512)
 data_processor.set_full_day_skip(False)
 
@@ -67,9 +67,10 @@ model_parameters = {
 model_parameters = task.connect(model_parameters, name="model_parameters")
 
 time_embedding = TimeEmbedding(data_processor.get_time_feature_size(), model_parameters["time_feature_embedding"])
-lstm_model = GRUModel(time_embedding.output_dim(inputDim), len(quantiles), hidden_size=model_parameters["hidden_size"], num_layers=model_parameters["num_layers"], dropout=model_parameters["dropout"])
+# lstm_model = GRUModel(time_embedding.output_dim(inputDim), len(quantiles), hidden_size=model_parameters["hidden_size"], num_layers=model_parameters["num_layers"], dropout=model_parameters["dropout"])
+non_linear_model = NonLinearRegression(time_embedding.output_dim(inputDim), len(quantiles), hiddenSize=model_parameters["hidden_size"], numLayers=model_parameters["num_layers"], dropout=model_parameters["dropout"])
 
-model = nn.Sequential(time_embedding, lstm_model)
+model = nn.Sequential(time_embedding, non_linear_model)
 optimizer = torch.optim.Adam(model.parameters(), lr=model_parameters["learning_rate"])
 
 #### Trainer ####