.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/tutorials/plot_tuto_benchmark_TS.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_tutorials_plot_tuto_benchmark_TS.py: ========================= Benchmark for time series ========================= In this tutorial, we show how to use Qolmat to benchmark several imputation methods and a multivariate time series dataset. We use Beijing Multi-Site Air-Quality Data Set. It consists in hourly air pollutants data from 12 chinese nationally-controlled air-quality monitoring sites. .. GENERATED FROM PYTHON SOURCE LINES 13-14 First import some libraries .. GENERATED FROM PYTHON SOURCE LINES 14-33 .. code-block:: Python import numpy as np np.random.seed(1234) import matplotlib.ticker as plticker from matplotlib import pyplot as plt tab10 = plt.get_cmap("tab10") from sklearn.linear_model import LinearRegression from qolmat.benchmark import comparator, missing_patterns from qolmat.imputations import imputers from qolmat.utils import data, plot from sklearn import utils as sku seed = 1234 rng = sku.check_random_state(seed) .. GENERATED FROM PYTHON SOURCE LINES 34-46 1. Data --------------------------------------------------------------- We use the public Beijing Multi-Site Air-Quality Data Set. It consists in hourly air pollutants data from 12 chinese nationally-controlled air-quality monitoring sites. In this way, each column has missing values. We group the data by day and only consider 5 columns. For the purpose of this notebook, we corrupt the data, with the ``qolmat.utils.data.add_holes`` function on three variables: "TEMP", "PRES" and "WSPM" and the imputation methods will have access to two additional features: "DEWP" and "RAIN". .. GENERATED FROM PYTHON SOURCE LINES 46-53 .. code-block:: Python df_data = data.get_data("Beijing") df_data = df_data[["TEMP", "PRES", "DEWP", "RAIN", "WSPM"]] df_data = df_data.groupby(level=["station", "date"]).mean() cols_to_impute = ["TEMP", "PRES", "WSPM"] df = data.add_holes(df_data, ratio_masked=0.15, mean_size=50) df[["DEWP", "RAIN"]] = df_data[["DEWP", "RAIN"]] .. GENERATED FROM PYTHON SOURCE LINES 54-55 Let's take a look at one station, for instance "Aotizhongxin" .. GENERATED FROM PYTHON SOURCE LINES 55-66 .. code-block:: Python station = "Aotizhongxin" fig, ax = plt.subplots(len(cols_to_impute), 1, figsize=(13, 8)) for i, col in enumerate(cols_to_impute): ax[i].plot(df.loc[station, col]) ax[i].set_ylabel(col) fig.align_labels() ax[0].set_title(station, fontsize=14) plt.tight_layout() plt.show() .. image-sg:: /examples/tutorials/images/sphx_glr_plot_tuto_benchmark_TS_001.png :alt: Aotizhongxin :srcset: /examples/tutorials/images/sphx_glr_plot_tuto_benchmark_TS_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 67-81 2. Time series imputation methods --------------------------------------------------------------- All presented methods are group-wise: here each station is imputed independently. For example ImputerMean computes the mean of each variable in each station and uses the result for imputation; ImputerInterpolation interpolates temporal signals corresponding to each variable on each station. We consider five imputation methods: ``median`` for a baseline imputation; ``interpolation`` since this method is really simple and often works well for time series; ``residuals`` which is to be compared directly with interpolation since it is also a simple linear interpolation, but no longer on the series directly but on the residuals. This method works very well when the time series displays seasonal patterns; ``TSOU`` which assumes the time series follow a VAR(1) process and finally ``mice`` which is known to be very effective. We use mice with linear regressions. .. GENERATED FROM PYTHON SOURCE LINES 81-138 .. code-block:: Python ratio_masked = 0.1 imputer_median = imputers.ImputerSimple(groups=("station",), strategy="median") imputer_interpol = imputers.ImputerInterpolation( groups=("station",), method="linear" ) imputer_residuals = imputers.ImputerResiduals( groups=("station",), period=365, model_tsa="additive", extrapolate_trend="freq", method_interpolation="linear", ) imputer_tsou = imputers.ImputerEM( groups=("station",), model="VAR", method="sample", max_iter_em=30, n_iter_ou=15, dt=1e-3, p=1, random_state=rng ) imputer_mice = imputers.ImputerMICE( groups=("station",), estimator=LinearRegression(), sample_posterior=False, max_iter=100, ) generator_holes = missing_patterns.EmpiricalHoleGenerator( n_splits=4, groups=("station",), subset=cols_to_impute, ratio_masked=ratio_masked, random_state=rng ) dict_imputers = { "median": imputer_median, "interpolation": imputer_interpol, "residuals": imputer_residuals, "TSOU": imputer_tsou, "mice": imputer_mice, } n_imputers = len(dict_imputers) comparison = comparator.Comparator( dict_imputers, generator_holes=generator_holes, metrics=["mae", "wmape", "kl_columnwise", "wasserstein_columnwise"], max_evals=10, ) results = comparison.compare(df) results.style.highlight_min(color="lightsteelblue", axis=1) .. raw:: html
    median interpolation residuals TSOU mice
kl_columnwise PRES 20.641637 0.266796 0.180567 0.495054 0.491604
TEMP 25.288092 0.288897 0.085867 0.074901 0.046404
WSPM 15.925819 0.684028 0.227816 0.696523 0.592234
mae PRES 8.649781 5.041545 5.308747 4.094000 4.265280
TEMP 10.477251 2.990938 2.614373 2.028266 2.970154
WSPM 0.534205 0.593273 0.672117 0.428670 0.444895
wasserstein_columnwise PRES 6.517742 1.191430 1.213633 1.311710 1.181843
TEMP 7.781333 1.456686 0.580367 0.495105 0.840348
WSPM 0.435260 0.130819 0.104882 0.206771 0.172547
wmape PRES 0.008560 0.004989 0.005253 0.004051 0.004221
TEMP 0.697256 0.199104 0.174262 0.134933 0.197616
WSPM 0.317137 0.352568 0.399242 0.254728 0.264032


.. GENERATED FROM PYTHON SOURCE LINES 139-146 We have considered four metrics for comparison. ``mae`` and ``wmape`` are point-wise metrics, while ``kl_columnwise`` and ``wasserstein_columnwise`` are metrics that compare distributions. Since we treat time series with strong seasonal patterns, imputation on residuals works very well. From these results, users can choose the imputer that suits them best. .. GENERATED FROM PYTHON SOURCE LINES 148-153 3. Visualisation --------------------------------------------------------------- We will now proceed to the visualisation stage. TIn order to make the figures readable, we select a single station: Aotizhongxin .. GENERATED FROM PYTHON SOURCE LINES 153-182 .. code-block:: Python df_plot = df[cols_to_impute] dfs_imputed = { name: imp.fit_transform(df_plot) for name, imp in dict_imputers.items() } station = "Aotizhongxin" df_station = df_plot.loc[station] dfs_imputed_station = { name: df_plot.loc[station] for name, df_plot in dfs_imputed.items() } fig, axs = plt.subplots( 3, 1, sharex=True, figsize=(10, 3 * len(cols_to_impute)) ) for col, ax in zip(cols_to_impute, axs.flatten()): values_orig = df_station[col] ax.plot(values_orig, ".", color="black", label="original") for ind, (name, model) in enumerate(list(dict_imputers.items())): values_imp = dfs_imputed_station[name][col].copy() values_imp[values_orig.notna()] = np.nan ax.plot(values_imp, ".", color=tab10(ind), label=name, alpha=1) ax.set_ylabel(col, fontsize=12) ax.legend(fontsize=10, ncol=3) ax.tick_params(axis="both", which="major", labelsize=10) loc = plticker.MultipleLocator(base=365) ax.xaxis.set_major_locator(loc) fig.align_labels() plt.show() .. image-sg:: /examples/tutorials/images/sphx_glr_plot_tuto_benchmark_TS_002.png :alt: plot tuto benchmark TS :srcset: /examples/tutorials/images/sphx_glr_plot_tuto_benchmark_TS_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 183-188 We can also check the covariance. We simply plot one variable versus one another. One observes the methods provide similar visual results: it's difficult to compare them based on this criterion, except the median imputation that greatly differs. Black points and ellipses are original dataframes while colored ones are imputed dataframes. .. GENERATED FROM PYTHON SOURCE LINES 188-210 .. code-block:: Python n_columns = len(dfs_imputed_station) fig = plt.figure(figsize=(10, 10)) i_plot = 1 for i, col in enumerate(cols_to_impute[:-1]): for i_imputer, (name_imputer, df_imp) in enumerate( dfs_imputed_station.items() ): ax = fig.add_subplot(n_columns, n_imputers, i_plot) plot.compare_covariances( df_station, df_imp, col, cols_to_impute[i + 1], ax, color=tab10(i_imputer), label=name_imputer, ) ax.set_title(name_imputer, fontsize=10) i_plot += 1 plt.tight_layout() plt.show() .. image-sg:: /examples/tutorials/images/sphx_glr_plot_tuto_benchmark_TS_003.png :alt: median, interpolation, residuals, TSOU, mice, median, interpolation, residuals, TSOU, mice :srcset: /examples/tutorials/images/sphx_glr_plot_tuto_benchmark_TS_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 31.531 seconds) .. _sphx_glr_download_examples_tutorials_plot_tuto_benchmark_TS.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_tuto_benchmark_TS.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_tuto_benchmark_TS.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_tuto_benchmark_TS.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_