<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 tsam - Animation of pareto-optimal aggregation
Date
Step1:  Input data
Read in time series from testdata.csv with pandas
Step2:  Setup the hyperparameter instance
Step3:  Load the resulting combination
Step4:  Create the animated aggregations
Drop all results with timesteps below 1% of the original data set since they are not meaningful.
Step5:  Append the original time series
Step6:  And reverse the order
Step7:  And create a dictionary with all aggregations we want to show in the animation
Step8:  And then append a last aggregation with the novel duration/distribution represenation
Step9:  Create the animation
Let animation warp - slow in the beginning and slow in the end
Step10:  Create the plot and the animation loop
Step12:  And save as animation parelllized with ffmpeg since the default matplotlib implemenation takes too long. Faster implemntation than matplotib from here
Step13:  You can also show it inline but it takes quite long.
<ASSISTANT_TASK:>
Python Code:
%load_ext autoreload
%autoreload 2
import copy
import os
import pandas as pd
import matplotlib.pyplot as plt
import tsam.timeseriesaggregation as tsam
import tsam.hyperparametertuning as tune
import tqdm
import numpy as np
from matplotlib.animation import FuncAnimation
%matplotlib inline
import matplotlib
import subprocess
import multiprocessing
import math
matplotlib.rcParams['animation.embed_limit'] = 1000 
Explanation: tsam - Animation of pareto-optimal aggregation
Date: 02.05.2022
Author: Leander Kotzur
Import pandas and the relevant time series aggregation class
End of explanation
raw = pd.read_csv('testdata.csv', index_col = 0)
raw=raw.rename(columns={'T': 'Temperature [°C]', 'Load':'Demand [kW]', 'Wind':'Wind [m/s]', 'GHI': 'Solar [W/m²]'})
Explanation: Input data
Read in time series from testdata.csv with pandas
End of explanation
tunedAggregations = tune.HyperTunedAggregations(
    tsam.TimeSeriesAggregation(
        raw,
        hoursPerPeriod=24,
        clusterMethod="hierarchical",
        representationMethod="medoidRepresentation",
        rescaleClusterPeriods=False,
        segmentation=True,
    )
)
Explanation: Setup the hyperparameter instance
End of explanation
results = pd.read_csv(os.path.join("results","paretoOptimalAggregation.csv"),index_col=0)
results["time_steps"] = results["segments"] * results["periods"]
Explanation: Load the resulting combination
End of explanation
results = results[results["time_steps"]>80]
Explanation: Create the animated aggregations
Drop all results with timesteps below 1% of the original data set since they are not meaningful.
End of explanation
results=results.append({"segments":24, "periods":365, "time_steps":len(raw)}, ignore_index=True)
Explanation: Append the original time series
End of explanation
results=results.iloc[::-1]
Explanation: And reverse the order
End of explanation
animation_list = []
for i, index in enumerate(tqdm.tqdm(results.index)):
    segments = results.loc[index,:].to_dict()["segments"]
    periods = results.loc[index,:].to_dict()["periods"]
    # aggregate to the selected set
    tunedAggregations._testAggregation(noTypicalPeriods=periods, noSegments=segments)
    # and repredict the data
    prediction = tunedAggregations.aggregationHistory[-1].predictOriginalData()
    # relative reduction of time steps
    reduction = 1 -  ( float(tunedAggregations._segmentHistory[-1] * tunedAggregations._periodHistory[-1]) /len(raw))
    # add a change layer which shows the difference of the latest aggregation to the previous
    if i > 0:
        # difference
        diff_val = previouspredictedPeriods - prediction
        # all fields that changed
        diff_bool = abs(diff_val)>1e-10
        # make sure that when any change is there it gets set to Nan
        prediction_diff = copy.deepcopy(prediction)
        prediction_diff[diff_bool.max(axis=1)] = np.nan
        # what changes ? segments or periods
        if segments == animation_list[-1]["Segments"]:
            misc = "Clustering periods"
        else:
            misc = "Clustering segments"
        animation_list.append({
            "Prediction": prediction_diff,
            "Segments": segments,
            "Periods": periods,
            "Reduction": reduction,
            "Misc": "Medoid representation",#misc,
        })
    animation_list.append({
        "Prediction": prediction,
        "Segments": segments,
        "Periods": periods,
        "Reduction": reduction,
        "Misc": "Medoid representation",
    })
    # and set previous prediction preiods
    previouspredictedPeriods = prediction
Explanation: And create a dictionary with all aggregations we want to show in the animation
End of explanation
aggregation=tsam.TimeSeriesAggregation(
        raw,
        hoursPerPeriod=24,
        noSegments=segments,
        noTypicalPeriods=periods,
        clusterMethod="hierarchical",
        rescaleClusterPeriods=False,
        segmentation=True,
        representationMethod="durationRepresentation",
        distributionPeriodWise=False
    )
    
animation_list.append({
        "Prediction": aggregation.predictOriginalData(),
        "Segments": segments,
        "Periods": periods,
        "Reduction": reduction,
        "Misc": "Distribution representation",
    })
Explanation: And then append a last aggregation with the novel duration/distribution represenation
End of explanation
iterator = []
for i in range(len(animation_list )):
    if i < 1:
        iterator+=[i]*100
    elif i < 3:
        iterator+=[i]*50
    elif i < 6:
        iterator+=[i]*30
    elif i < 20:
        iterator+=[i]*10
    elif i >= len(animation_list )-1:
        iterator+=[i]*150
    elif i > len(animation_list )-3:
        iterator+=[i]*50
    elif i > len(animation_list )-6:
        iterator+=[i]*30
    else:
        iterator+=[i]
Explanation: Create the animation
Let animation warp - slow in the beginning and slow in the end
End of explanation
import matplotlib.ticker as tick
fig, axes = plt.subplots(figsize = [7, 5], dpi = 300, nrows = raw.shape[1], ncols = 1)
cmap = plt.cm.get_cmap("Spectral_r").copy()
cmap.set_bad((.7, .7, .7, 1))
for ii, column in enumerate(raw.columns):
    data = raw[column]
    stacked, timeindex = tsam.unstackToPeriods(copy.deepcopy(data), tunedAggregations.base_aggregation.hoursPerPeriod)
    cax = axes[ii].imshow(stacked.values.T, interpolation = 'nearest', vmin = raw[column].min(), vmax = raw[column].max(), origin='lower', cmap = cmap)
    axes[ii].set_aspect('auto')  
    axes[ii].set_ylabel('Hour')
    plt.xlabel('Day in the year')
    cbar=plt.colorbar(cax, ax=axes[ii], pad=0.01, aspect=7)
    cbar.set_label(column, fontsize="small")
fig.suptitle('Time series aggregation by xx.xx %',
    y=0.97,
    x=0.3,
    horizontalalignment="left",)
text=fig.text(
        0.27,
        0.91,
        "with xxx periods and xx segments - Medoid representation",
        horizontalalignment="left",
        fontsize="small",)
fig.subplots_adjust(right = 1.1, hspace = 0.05)
def animate(iter):
    i=iterator[iter]
    predictedPeriods = animation_list[i]["Prediction"]
    fig.suptitle(
        'Time series aggregation by ' + str(round(animation_list[i]["Reduction"]*100,2)) + " % ",
        y=0.97,
        x=0.3,
        horizontalalignment="left",)
    text.set_text(
        "with " + 
        str(animation_list[i]["Periods"]) + " periods and " + str(animation_list[i]["Segments"]) + " segments - " + animation_list[i]["Misc"],
        )
    for ii, column in enumerate(raw.columns):
        data = predictedPeriods[column]        
        stacked, timeindex = tsam.unstackToPeriods(copy.deepcopy(data), tunedAggregations.base_aggregation.hoursPerPeriod)
        cax = axes[ii].imshow(stacked.values.T, interpolation = 'nearest', vmin = raw[column].min(), vmax = raw[column].max(), origin='lower', cmap = cmap)
        axes[ii].set_aspect('auto')  
    fig.subplots_adjust(right = 1.0, hspace = 0.05)
ani = FuncAnimation(fig, animate, repeat_delay=600, interval=20, frames=len(iterator) )
Explanation: Create the plot and the animation loop
End of explanation
def chunks(lst, n):
    Yield successive n-sized chunks from lst.
    for i in range(0, len(lst), n):
        yield lst[i:i + n]
threads = multiprocessing.cpu_count()
frames=[i for i in range(len(iterator))]
# divide the frame equally
i_length=math.ceil(len(frames)/(threads))
frame_sets=list(chunks(frames,i_length))
filenames=[]
for i in range(len(frame_sets)):
    filenames.append("_temp_video_chunk_" + str(i) + ".mp4")
def ani_to_mp4(frame_range, filename):
    canvas_width, canvas_height = fig.canvas.get_width_height()
    # Open an ffmpeg process
    outf = os.path.join("results", filename)
    cmdstring = ('ffmpeg', 
                    '-y', '-r', '100', # fps
                    '-s', '%dx%d' % (canvas_width, canvas_height), # size of image string
                    '-pix_fmt', 'argb', # formats
                    '-f', 'rawvideo',  '-i', '-', # tell ffmpeg to expect raw video from the pipe
                    '-vcodec', 'mpeg4', outf) # output encoding
    p = subprocess.Popen(cmdstring, stdin=subprocess.PIPE)
    # Draw frames and write to the pipe
    for frame in frame_range:
        # draw the frame
        animate(frame)
        fig.canvas.draw()
        # extract the image as an ARGB string
        string = fig.canvas.tostring_argb()
        # write to pipe
        p.stdin.write(string)
    # Finish up
    p.communicate()
with multiprocessing.Pool() as pool:
    pool.starmap(ani_to_mp4, zip(frame_sets, filenames))
filename_list=os.path.join("results", "filenames.txt")
with open(filename_list, "w") as textfile:
    for filename in filenames:
        textfile.write("file '" + filename + "'\n")
cmdstring = ('ffmpeg', '-y',
                '-f', 'concat', 
                '-safe', '0', 
                '-i', filename_list,
                '-c', 'copy', 
                os.path.join("results", "animation.mp4") # output encoding
p = subprocess.Popen(cmdstring, stdin=subprocess.PIPE)
Explanation: And save as animation parelllized with ffmpeg since the default matplotlib implemenation takes too long. Faster implemntation than matplotib from here: https://stackoverflow.com/a/31315362/3253411 
Parallelize animation to video
End of explanation
from IPython.display import HTML
HTML(ani.to_jshtml())
Explanation: You can also show it inline but it takes quite long.
End of explanation
<END_TASK>