Segmentation#

How to use time step segmentation to reduce the number of timesteps per period.

Author: Maximilian Hoffmann

Import pandas and the relevant time series aggregation class

[1]:
%load_ext autoreload
%autoreload 2

import pandas as pd
import plotly.express as px
import plotly.io as pio

import tsam
from tsam import ClusterConfig, SegmentConfig

pio.renderers.default = "notebook"

Input data#

Read in time series from testdata.csv with pandas

[2]:
raw = pd.read_csv("testdata.csv", index_col=0)

Create a plot function for the temperature for a visual comparison of the time series

[3]:
# Use tsam.unstack_to_periods() with plotly for heatmap visualization
# px.imshow(unstacked["column"].values.T) creates interactive heatmaps

Hierarchical aggregation with medoid representation and 10 typical days with 24 hourly segments#

Initialize an aggregation class object with hierarchical as method for eight typical days

[4]:
result = tsam.aggregate(
    raw,
    n_clusters=10,
    period_duration=24,
    cluster=ClusterConfig(method="hierarchical"),
)

Create the typical periods

[5]:
cluster_representatives = result.cluster_representatives
cluster_representatives
[5]:
GHI Load T Wind
timestep
0 0 0.000000 383.795921 4.421453 5.043609
1 0.000000 371.689230 4.026704 4.933965
2 0.000000 374.310544 3.829330 4.714678
3 0.000000 369.720991 3.730643 4.605034
4 0.000000 376.053584 3.335894 5.591827
... ... ... ... ... ...
9 19 2.199516 478.910919 10.145308 2.192873
20 0.000000 457.575577 10.342682 2.192873
21 0.000000 444.225788 10.540056 2.192873
22 0.000000 430.132842 10.441369 2.192873
23 0.000000 401.874887 10.441369 2.192873

240 rows × 4 columns

Predict original data

[6]:
reconstructed = result.reconstructed
reconstructed
[6]:
GHI Load T Wind
2009-12-31 23:30:00 0.0 383.795921 4.421453 5.043609
2010-01-01 00:30:00 0.0 371.689230 4.026704 4.933965
2010-01-01 01:30:00 0.0 374.310544 3.829330 4.714678
2010-01-01 02:30:00 0.0 369.720991 3.730643 4.605034
2010-01-01 03:30:00 0.0 376.053584 3.335894 5.591827
... ... ... ... ...
2010-12-31 18:30:00 0.0 524.855994 3.335894 3.289310
2010-12-31 19:30:00 0.0 496.026033 2.842458 4.385747
2010-12-31 20:30:00 0.0 480.154017 2.645084 4.385747
2010-12-31 21:30:00 0.0 456.008193 2.546397 4.495390
2010-12-31 22:30:00 0.0 417.990119 2.151648 2.412161

8760 rows × 4 columns

Get accuracy indicators

[7]:
result.accuracy
[7]:
AccuracyMetrics(
  rmse=0.1018 (mean),
  mae=0.0715 (mean),
  rmse_duration=0.0305 (mean)
)

Hierarchical aggregation with medoid representation and 20 typical days with 12 irregular segments#

[8]:
result_segmented = tsam.aggregate(
    raw,
    n_clusters=20,
    period_duration=24,
    cluster=ClusterConfig(method="hierarchical"),
    segments=SegmentConfig(n_segments=12),
)

Create the typical periods

[9]:
cluster_representatives_segmented = result_segmented.cluster_representatives
cluster_representatives_segmented
[9]:
GHI Load T Wind
Segment Step Segment Duration
0 0 4 0.000000 403.310708 1.025025 2.226777
1 2 0.000000 424.393351 1.122830 3.340165
2 1 0.000000 495.118712 1.171733 3.340165
3 2 11.182790 541.665766 1.562953 3.340165
4 4 74.551932 547.697631 2.272040 2.226777
... ... ... ... ... ... ...
19 7 3 76.681988 471.099097 -2.447057 3.340165
8 4 34.080883 469.022225 -1.664616 2.226777
9 1 0.000000 499.029025 -2.055836 2.226777
10 4 0.000000 498.349573 -2.178093 1.113388
11 2 0.000000 474.908494 -2.495959 2.226777

240 rows × 4 columns

Predict original data

[10]:
reconstructed_segmented = result_segmented.reconstructed
reconstructed_segmented
[10]:
GHI Load T Wind
2009-12-31 23:30:00 0.000000 409.327552 1.660758 5.566942
2010-01-01 00:30:00 0.000000 396.579747 1.660758 4.453553
2010-01-01 01:30:00 0.000000 386.327886 1.758564 6.680330
2010-01-01 02:30:00 0.000000 386.327886 1.758564 6.680330
2010-01-01 03:30:00 3.834099 378.749459 1.621636 5.344264
... ... ... ... ...
2010-12-31 18:30:00 0.000000 497.614656 2.051979 3.340165
2010-12-31 19:30:00 0.000000 497.614656 2.051979 3.340165
2010-12-31 20:30:00 0.000000 477.379017 2.003076 2.783471
2010-12-31 21:30:00 0.000000 477.379017 2.003076 2.783471
2010-12-31 22:30:00 0.000000 450.506015 1.856369 2.226777

8760 rows × 4 columns

Get accuracy indicators

[11]:
result_segmented.accuracy
[11]:
AccuracyMetrics(
  rmse=0.0934 (mean),
  mae=0.0658 (mean),
  rmse_duration=0.0190 (mean)
)

Comparison of the aggregations#

It was shown for the temperature, but both times all four time series have been aggregated. Therefore, we compare here also the duration curves of the electrical load for the original time series, the aggregation with k-mean, and the hierarchical aggregation including peak periods.

[12]:
# Compare duration curves using plotly express
comparison_data = {
    "Original": raw,
    "10 with 24 hours": reconstructed,
    "20 with 12 Seg": reconstructed_segmented,
}

frames = []
for name, df in comparison_data.items():
    sorted_vals = df["Load"].sort_values(ascending=False).reset_index(drop=True)
    frames.append(
        pd.DataFrame(
            {"Hour": range(len(sorted_vals)), "Load": sorted_vals, "Method": name}
        )
    )
long_df = pd.concat(frames, ignore_index=True)

px.line(
    long_df,
    x="Hour",
    y="Load",
    color="Method",
    title="Duration Curve Comparison - Load",
)
[13]:
param = "GHI"
[14]:
# Original data heatmap
unstacked_orig = tsam.unstack_to_periods(raw, period_duration=24)
px.imshow(
    unstacked_orig[param].values.T,
    labels={"x": "Day", "y": "Hour", "color": param},
    title=f"Original {param}",
    aspect="auto",
)
[15]:
# 10 periods with 24 hours heatmap
unstacked_10 = tsam.unstack_to_periods(reconstructed, period_duration=24)
px.imshow(
    unstacked_10[param].values.T,
    labels={"x": "Day", "y": "Hour", "color": param},
    title=f"10 with 24 hours - {param}",
    aspect="auto",
)
[16]:
# 20 periods with 12 segments heatmap
unstacked_20 = tsam.unstack_to_periods(reconstructed_segmented, period_duration=24)
px.imshow(
    unstacked_20[param].values.T,
    labels={"x": "Day", "y": "Hour", "color": param},
    title=f"20 with 12 Seg - {param}",
    aspect="auto",
)
[17]:
# Time slice comparison using plotly express
frames = []
for name, df in comparison_data.items():
    sliced = df.loc["20100210":"20100218", ["Load"]].copy()
    sliced["Method"] = name
    frames.append(sliced)
long_df = pd.concat(frames).reset_index(names="Time")

px.line(
    long_df,
    x="Time",
    y="Load",
    color="Method",
    title="Time Slice Comparison - Load (Feb 10-18)",
)

Validation#

Check that the means of the original time series and the predicted ones are the same.

[18]:
raw.mean()
[18]:
GHI     110.990183
Load    450.260335
T         7.790616
Wind      3.057306
dtype: float64
[19]:
reconstructed.mean()
[19]:
GHI     110.990183
Load    450.260335
T         7.790616
Wind      3.057306
dtype: float64
[20]:
reconstructed_segmented.mean()
[20]:
GHI     110.990179
Load    450.260335
T         7.790616
Wind      3.057306
dtype: float64

Check that a segmented period has the same column-wise means as a non-segmented period for if the periods are the same.

[21]:
# Mean of first period with non-segmented aggregation
result.cluster_representatives.loc[0, :].mean()
[21]:
GHI      58.149701
Load    465.032475
T         3.710083
Wind      4.979650
dtype: float64
[22]:
# Segmented aggregation with same number of periods for comparison
result_segmented_test = tsam.aggregate(
    raw,
    n_clusters=10,
    period_duration=24,
    cluster=ClusterConfig(method="hierarchical"),
    segments=SegmentConfig(n_segments=12),
)
[23]:
# Get segment durations
segment_durations = result_segmented_test.segment_durations
print("Segment durations:", segment_durations)
Segment durations: ((6, 1, 2, 2, 1, 2, 2, 2, 2, 2, 1, 1), (3, 2, 2, 2, 1, 1, 3, 1, 2, 1, 4, 2), (1, 4, 1, 1, 3, 2, 2, 3, 3, 2, 1, 1), (1, 4, 1, 1, 1, 1, 2, 3, 2, 2, 4, 2), (6, 1, 1, 1, 1, 3, 1, 2, 2, 2, 2, 2), (3, 2, 3, 2, 2, 2, 1, 1, 1, 3, 3, 1), (6, 1, 1, 1, 1, 2, 3, 1, 1, 1, 4, 2), (4, 2, 1, 2, 4, 2, 1, 1, 3, 2, 1, 1), (3, 1, 1, 2, 2, 1, 1, 3, 1, 5, 2, 2), (2, 4, 1, 2, 1, 2, 3, 2, 1, 2, 3, 1))
[24]:
# Weighted mean of first period (should match non-segmented period mean)
period_0 = result_segmented_test.cluster_representatives.loc[0, :].reset_index(
    0, drop=True
)
# Convert segment_durations tuple-of-tuples to DataFrame and get values for period 0
segment_durations_df = pd.DataFrame(segment_durations)
segment_durations_values = segment_durations_df.loc[0, :].values.flatten().tolist()
weighted_mean = period_0.mul(segment_durations_values, axis=0).sum() / sum(
    segment_durations_values
)
weighted_mean
[24]:
GHI      58.149701
Load    465.032475
T         3.710083
Wind      4.979650
dtype: float64

Print out the (segmented) typical periods.

[25]:
# Display segmented typical periods
result_segmented.cluster_representatives
[25]:
GHI Load T Wind
Segment Step Segment Duration
0 0 4 0.000000 403.310708 1.025025 2.226777
1 2 0.000000 424.393351 1.122830 3.340165
2 1 0.000000 495.118712 1.171733 3.340165
3 2 11.182790 541.665766 1.562953 3.340165
4 4 74.551932 547.697631 2.272040 2.226777
... ... ... ... ... ... ...
19 7 3 76.681988 471.099097 -2.447057 3.340165
8 4 34.080883 469.022225 -1.664616 2.226777
9 1 0.000000 499.029025 -2.055836 2.226777
10 4 0.000000 498.349573 -2.178093 1.113388
11 2 0.000000 474.908494 -2.495959 2.226777

240 rows × 4 columns

[26]:
# Display non-segmented typical periods
result.cluster_representatives
[26]:
GHI Load T Wind
timestep
0 0 0.000000 383.795921 4.421453 5.043609
1 0.000000 371.689230 4.026704 4.933965
2 0.000000 374.310544 3.829330 4.714678
3 0.000000 369.720991 3.730643 4.605034
4 0.000000 376.053584 3.335894 5.591827
... ... ... ... ... ...
9 19 2.199516 478.910919 10.145308 2.192873
20 0.000000 457.575577 10.342682 2.192873
21 0.000000 444.225788 10.540056 2.192873
22 0.000000 430.132842 10.441369 2.192873
23 0.000000 401.874887 10.441369 2.192873

240 rows × 4 columns