Looking at pre-dispatch demand forecast errors in 2021#
In this example, we will take a look at (30-minute) pre-dispatch (PREDISPATCH) demand forecast “error” (the difference between “actual” - the demand used in dispatch - and forecasted demand) for 2021. Unlike 5PMD, pre-dispatch extends out to 39 hours ahead, so it’s a good dataset to use to look at day-ahead forecast errors.
Key imports#
# standard libraries
import logging
from datetime import datetime, timedelta
from pathlib import Path
# NEM data libraries
# NEMOSIS for actual demand data
# NEMSEER for forecast demand data
import nemosis
from nemseer import compile_data, download_raw_data, generate_runtimes
# data wrangling libraries
import numpy as np
import pandas as pd
# interactive plotting
import plotly.express as px
import plotly.io as pio
import plotly.graph_objects as go
# progress bar for error computation
from tqdm.autonotebook import tqdm
# silence NEMSEER and NEMOSIS logging
logging.getLogger("nemseer").setLevel(logging.ERROR)
logging.getLogger("nemosis").setLevel(logging.WARNING)
Plot styling#
nemseer_template = dict(
layout=go.Layout(
font_family="Source Sans 3",
title_font_size=24,
title_x=0.05,
plot_bgcolor="#f0f0f0",
colorway=px.colors.qualitative.Bold,
)
)
Defining our analysis start and end dates#
analysis_start = "2021/01/01 00:00:00"
analysis_end = "2022/01/01 00:00:00"
Obtaining actual demand data from NEMOSIS
#
We will download DISPATCHREGIONSUM
to access the TOTALDEMAND
field.
We’ll first download the data we need and cache it so that it’s ready for computation.
nemosis_cache = Path("nemosis_cache/")
if not nemosis_cache.exists():
nemosis_cache.mkdir()
nemosis.cache_compiler(
analysis_start, analysis_end, "DISPATCHREGIONSUM", nemosis_cache, fformat="parquet"
)
Obtaining forecast demand data from NEMSEER
#
We will download REGIONSUM
to access the TOTALDEMAND
field in PREDISPATCH
forecasts.
We’ll first download the data we need and cache it so that it’s ready for computation.
download_raw_data(
"PREDISPATCH",
"REGIONSUM",
"nemseer_cache/",
forecasted_start=analysis_start,
forecasted_end=analysis_end,
)
Calculating regional forecast errors#
Below we calculate demand forecast error for PREDISPATCH
forecasts using forecast demand data and actual demand data.
Attention
The actual run time of PD is approximately 30 minutes before the nominal run time. We will adjust for this in this when calculating forecast ahead times. See the note in this section.
As data for the entire period is loaded into memory, adapt the length of the period you select to your machine specifications (e.g. a year’s worth of forecast data consumed ~10GB on the test machine).
Forecast error calculation functions#
The code below uses functionalities offered by NEMOSIS
, NEMSEER
and pandas
to calculate demand forecast error.
def calculate_predispatch_demand_forecast_error_vectorised(
analysis_start: str, analysis_end: str
) -> pd.DataFrame:
"""
Calculates PD demand forecast error (Actual - Forecast) for all forecasts
that are run for a given forecasted_time in a vectorised fashion.
Args:
forecasted_time: Datetime string in the form YYYY/mm/dd HH:MM:SS
Returns:
pandas DataFrame with forecast error in `TOTALDEMAND` columns, the ahead time
of the forecast run in `ahead_time`, and the forecasted time in
`forecasted_time`.
"""
def get_forecast_data(analysis_start: str, analysis_end: str) -> pd.DataFrame:
"""
Use NEMSEER to get PD forecast data. Also omits any intervention periods.
"""
# use NEMSEER functions to compile pre-cached data
forecasts_run_start, forecasts_run_end = generate_runtimes(
analysis_start, analysis_end, "PREDISPATCH"
)
forecast_df = compile_data(
forecasts_run_start,
forecasts_run_end,
analysis_start,
analysis_end,
"PREDISPATCH",
"REGIONSUM",
"nemseer_cache/",
)["REGIONSUM"]
# remove intervention periods
forecast_df = forecast_df.query("INTERVENTION == 0")
return forecast_df
def get_actual_data(analysis_start: str, analysis_end: str) -> pd.DataFrame:
"""
Use NEMOSIS to get actual data. Also omits any intervention periods
"""
# NEMOSIS start time must precede end of interval of interest by 5 minutes
nemosis_start = (
datetime.strptime(analysis_start, "%Y/%m/%d %H:%M:%S")
- timedelta(minutes=5)
).strftime("%Y/%m/%d %H:%M:%S")
# use NEMOSIS to compile pre-cached data and filter out interventions
actual_df = nemosis.dynamic_data_compiler(
nemosis_start,
analysis_end,
"DISPATCHREGIONSUM",
nemosis_cache,
filter_cols=["INTERVENTION"],
filter_values=([0],),
fformat="parquet",
)
return actual_df
def calculate_pd_forecast_demand_error(
actual_demand: pd.DataFrame, forecast_demand: pd.DataFrame
) -> pd.DataFrame:
"""
Calculate PD forecast demand error given actual and forecast demand
Ahead time calculation reflects the fact that PD actual run time is
30 minutes before the nominal run time.
"""
# merge the two types of demand
merged = pd.merge(
forecast_demand,
actual_demand,
on=["forecasted_time", "REGIONID"],
how="left",
)
if len(merged) > len(forecast_demand):
raise ValueError(
"Merge should return DataFrame with dimensions of forecast data"
)
# subtract 30 minutes from run time to get actual run time
merged["ahead_time"] = merged["forecasted_time"] - (
merged["run_time"] - timedelta(minutes=30)
)
# calculate forecast error
forecast_error = (
merged["TOTALDEMAND"] - merged["FORECAST_TOTALDEMAND"]
).rename("TOTALDEMAND")
# create the forecast error DataFrame
forecast_error = pd.concat(
[forecast_error, merged["ahead_time"], merged["REGIONID"]], axis=1
).set_index(merged["forecasted_time"])
return forecast_error
# get forecast data
forecast_df = get_forecast_data(analysis_start, analysis_end)
# rename columns in preparation for merge
forecast_df = forecast_df.rename(
columns={
"TOTALDEMAND": "FORECAST_TOTALDEMAND",
"DATETIME": "forecasted_time",
"PREDISPATCH_RUN_DATETIME": "run_time",
}
)
forecast_demand = forecast_df[
["run_time", "forecasted_time", "REGIONID", "FORECAST_TOTALDEMAND"]
]
# get actual data
actual_df = get_actual_data(analysis_start, analysis_end)
# rename columns in preparation for merge
actual_df = actual_df.rename(
columns={
"SETTLEMENTDATE": "forecasted_time",
"TOTALDEMAND": "TOTALDEMAND",
}
)
actual_demand = actual_df[["forecasted_time", "REGIONID", "TOTALDEMAND"]]
forecast_error = calculate_pd_forecast_demand_error(actual_demand, forecast_demand)
return forecast_error
forecast_error = calculate_predispatch_demand_forecast_error_vectorised(
analysis_start, analysis_end
)
Region-by-region error percentiles#
Below we plot regional error percentiles for all ahead times.
region_ahead_percentiles = {}
for region in (regions := ("QLD1", "NSW1", "VIC1", "SA1", "TAS1")):
percentile_data = []
region_error = forecast_error.query("REGIONID==@region")
for quantile in (0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95, 0.99):
quantile_result = region_error.groupby(
region_error["ahead_time"].dt.total_seconds() / (60**2)
)["TOTALDEMAND"].quantile(quantile)
percentile_result = pd.concat(
[
quantile_result,
pd.Series(
np.repeat(quantile * 100, len(quantile_result)),
index=quantile_result.index,
name="Percentile",
).astype(int),
],
axis=1,
)
percentile_data.append(percentile_result)
percentile_df = pd.concat(percentile_data, axis=0).reset_index()
region_ahead_percentiles[region] = percentile_df
figs = []
for region in regions:
fig = px.line(
region_ahead_percentiles[region],
x="ahead_time",
y="TOTALDEMAND",
color="Percentile",
title=f"PD {region} Demand Forecast Error, 2021<br><sup>Error = Actual - Forecast</sup>",
labels={
"TOTALDEMAND": "Demand Forecast Error (MW)",
"ahead_time": "Forecast Ahead Time (Hours, Actual Run Time)",
},
template=nemseer_template,
color_discrete_map={
1: "#E24A33",
5: "#348ABD",
10: "#988ED5",
25: "#777777",
50: "#FBC15E",
75: "#777777",
90: "#988ED5",
95: "#348ABD",
99: "#E24A33",
},
)
fig["layout"]["xaxis"]["autorange"] = "reversed"
figs.append(fig)
Why does the error drop off beyond ~24 hours?#
A limited number of periods during the day are actually forecasted beyond 24 hours out.
PREDISPATCH
is run until the end of the trading day for which bid price band submission has closed (1230 EST). So this means, for example:
The 1300 PD (nominal) run will forecast out til 4AM two days away (39 hours)
But the 1400 PD (nominal) run will still only forecast out til 4AM two days away (38 hours)
And the 0800 PD (nominal) run the next day will still only forecast out til 4AM the next day (20 hours)
So because of this, the number of error samples drops off beyond 16 hours ahead (see figure below).
In addition, the runs closer to ~35 hours will be forecasts for periods in the early hours of the morning. These periods tend to have more predictable demand.
sample_count = px.line(
forecast_error.groupby(forecast_error["ahead_time"].dt.total_seconds() / (60**2))[
"TOTALDEMAND"
]
.count()
.rename("Computed Errors"),
labels={"value": "Count of Samples"},
template=nemseer_template,
)
sample_count.update_layout(legend_title="", xaxis=dict(title="Ahead Time (hours)"));
NEM-wide Demand Forecast Error, less than 24 hours#
Because of the reasons above, we’ll focus on ahead times of up to 24 hours.
nem_error = (
forecast_error.reset_index()
.groupby(["forecasted_time", "ahead_time"])["TOTALDEMAND"]
.sum()
.reset_index()
)
nem_percentile_data = []
for quantile in (0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95, 0.99):
nem_quantile_result = nem_error.groupby(
nem_error["ahead_time"].dt.total_seconds() / (60**2)
)["TOTALDEMAND"].quantile(quantile)
nem_percentile_result = pd.concat(
[
nem_quantile_result,
pd.Series(
np.repeat(quantile * 100, len(nem_quantile_result)),
index=nem_quantile_result.index,
name="Percentile",
).astype(int),
],
axis=1,
)
nem_percentile_data.append(nem_percentile_result)
nem_percentile_df = pd.concat(nem_percentile_data, axis=0).reset_index()
nemwide = px.line(
nem_percentile_df.query("ahead_time < 24"),
x="ahead_time",
y="TOTALDEMAND",
color="Percentile",
title=f"Pre-dispatch NEM-wide Demand Forecast Error, 2021<br><sup>Error = Actual - Forecast</sup>",
labels={
"TOTALDEMAND": "Demand Forecast Error (MW)",
"ahead_time": "Forecast Ahead Time (Hours, Actual Run Time)",
},
template=nemseer_template,
color_discrete_map={
1: "#E24A33",
5: "#348ABD",
10: "#988ED5",
25: "#777777",
50: "#FBC15E",
75: "#777777",
90: "#988ED5",
95: "#348ABD",
99: "#E24A33",
},
)
nemwide["layout"]["xaxis"]["autorange"] = "reversed"
Distributions of Day-Ahead Demand Forecast Error by Region#
We can see that the TOTALDEMAND day-ahead demand forecast error distribution is long-tailed for every region.
day_ahead = forecast_error[
forecast_error["ahead_time"].dt.total_seconds() / (60**2) == 24
]
da_dists = px.histogram(
day_ahead,
x="TOTALDEMAND",
facet_row="REGIONID",
title="Pre-dispatch Demand Forecast Error, 2021<br><sup>Day-Ahead (24 hours ahead)</sup>",
template=nemseer_template,
)
da_dists.update_layout(xaxis=dict(title="Demand Forecast Error (MW)"));