diff --git a/README.md b/README.md index a6600cba2..c8323fad9 100644 --- a/README.md +++ b/README.md @@ -28,6 +28,10 @@ Or the development version: `> pip install https://github.com/DHI/modelskill/archive/main.zip` +Interactive plots (`backend="plotly"`) require an extra dependency: + +`> pip install "modelskill[plotly]"` + ## Example notebooks @@ -101,10 +105,22 @@ ms.plotting.spatial_overview([HKNA, EPL, c2], mr, figsize=(7,7)) ### Timeseries plot -Timeseries plots can either be static and report-friendly ([matplotlib](https://matplotlib.org/)) or interactive with zoom functionality ([plotly](https://plotly.com/python/)). +Plots can either be static and report-friendly ([matplotlib](https://matplotlib.org/), the default) or interactive with zoom functionality ([plotly](https://plotly.com/python/)). ```python -cc["HKNA"].plot.timeseries(width=1000, backend="plotly") +cc["HKNA"].plot.timeseries(figsize=(10, 4), backend="plotly") ``` ![timeseries](https://raw.githubusercontent.com/DHI/modelskill/main/images/plotly_timeseries.png) + +Every plot takes a `backend` argument -- `scatter`, `hist`, `kde`, `qq`, `box`, +`residual_hist` and `taylor` on both `Comparer` and `ComparerCollection`, +`Comparer.plot.timeseries`, `ComparerCollection.plot.spatial_overview` and +`.temporal_coverage`, the `timeseries` and `hist` plots on observations and model results, +and the standalone functions in `ms.plotting`. + +The same arguments (`title`, `figsize` in inches, `xlim`, `ylim`, ...) work with both +backends. The matplotlib backend returns a `matplotlib.axes.Axes` (or `Figure` for +`taylor`), the plotly backend a `plotly.graph_objects.Figure`. Extra `**kwargs` go to the +underlying matplotlib call, or to +[`Figure.update_layout`](https://plotly.com/python/reference/layout/) respectively. diff --git a/docs/user-guide/plotting.qmd b/docs/user-guide/plotting.qmd index 5fdce145d..ff3114d5a 100644 --- a/docs/user-guide/plotting.qmd +++ b/docs/user-guide/plotting.qmd @@ -1,5 +1,22 @@ # Plotting +## Backends + +Every plot can be rendered by either of two backends, selected with the `backend` argument: + +* `"matplotlib"` (the default) - static, report-friendly figures, returns a `matplotlib.axes.Axes` (or a `Figure` for `taylor`) +* `"plotly"` - [interactive](https://plotly.com/python/) figures with zoom and hover, returns a `plotly.graph_objects.Figure` + +The plotly backend requires an optional dependency: + +```bash +pip install "modelskill[plotly]" +``` + +Both backends take the same arguments — `title`, `figsize` (in inches for both), `xlim`, `ylim` and so on. Extra `**kwargs` go to the underlying matplotlib call, or to plotly's [`Figure.update_layout`](https://plotly.com/python/reference/layout/) respectively, so a matplotlib-only argument such as `cmap` is rejected by the plotly backend. `ax` is matplotlib-only; the plotly backend always returns a new figure. + +Both backends return the figure rather than showing it, so in a script you need `.show()` (or `fig.write_html(...)`). + ## Plotting observations and model results [](`~modelskill.PointObservation`)s and [](`~modelskill.PointModelResult`)s can be plotted using their `plot` accessor: @@ -17,6 +34,10 @@ mr = ms.PointModelResult('../data/SW/ts_storm_4.dfs0', item=0) # TODO coords o.plot.timeseries(); ``` +```{python} +o.plot.timeseries(backend="plotly") +``` + ```{python} mr.plot.timeseries(); ``` @@ -71,6 +92,12 @@ cmp.plot.timeseries(); cmp.plot.scatter(); ``` +The same plot with the interactive backend: + +```{python} +cmp.plot.scatter(backend="plotly") +``` + ## Taylor diagrams diff --git a/pyproject.toml b/pyproject.toml index 890fc7574..29767610a 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -44,11 +44,12 @@ classifiers = [ [project.optional-dependencies] networks = ["mikeio1d", "networkx"] +plotly = ["plotly >= 4.5"] [dependency-groups] dev = ["pytest", "plotly >= 4.5", "ruff==0.6.2", "netCDF4", "dask"] -docs = ["quartodoc==0.11.1", "nbformat", "nbconvert", "ipykernel", "griffe<2"] +docs = ["quartodoc==0.11.1", "nbformat", "nbconvert", "ipykernel", "griffe<2", "plotly >= 4.5"] test = [ "pytest", @@ -58,6 +59,7 @@ test = [ "mypy==1.19.1", "types-PyYAML", "geopandas", + "plotly >= 4.5", ] notebooks = ["nbformat", "nbconvert", "jupyter", "plotly", "shapely", "seaborn"] diff --git a/src/modelskill/comparison/_collection_plotter.py b/src/modelskill/comparison/_collection_plotter.py index b1df0cd2f..3f431ed51 100644 --- a/src/modelskill/comparison/_collection_plotter.py +++ b/src/modelskill/comparison/_collection_plotter.py @@ -4,7 +4,6 @@ TYPE_CHECKING, Any, List, - Literal, Mapping, Sequence, Tuple, @@ -17,13 +16,24 @@ if TYPE_CHECKING: from ._collection import ComparerCollection -from matplotlib.figure import Figure import numpy as np import pandas as pd from .. import metrics as mtr -from ..plotting import TaylorPoint, scatter, taylor_diagram -from ..plotting._misc import _get_fig_ax, _xtick_directional, _ytick_directional +from ..plotting import TaylorPoint, scatter, taylor_diagram, _plotly +from ..plotting._backend import ( + Backend, + FigureResult, + PlotResult, + reject_matplotlib_axes, + validate_backend, +) +from ..plotting._misc import ( + RESIDUAL_COLOR, + _get_fig_ax, + _xtick_directional, + _ytick_directional, +) from ..settings import options from ..utils import _get_idx from ._comparer_plotter import quantiles_xy @@ -49,7 +59,7 @@ def __init__(self, cc: ComparerCollection) -> None: self.cc = cc self.is_directional = False - def __call__(self, *args: Any, **kwds: Any) -> Axes | list[Axes]: + def __call__(self, *args: Any, **kwds: Any) -> PlotResult | list[PlotResult]: return self.scatter(*args, **kwds) def scatter( @@ -62,7 +72,7 @@ def scatter( show_hist: bool | None = None, show_density: bool | None = None, norm: colors.Normalize | None = None, - backend: Literal["matplotlib", "plotly"] = "matplotlib", + backend: Backend = "matplotlib", figsize: Tuple[float, float] = (8, 8), xlim: Tuple[float, float] | None = None, ylim: Tuple[float, float] | None = None, @@ -73,7 +83,7 @@ def scatter( skill_table: Union[str, List[str], Mapping[str, str], bool] | None = None, ax: Axes | None = None, **kwargs, - ) -> Axes | list[Axes]: + ) -> PlotResult | list[PlotResult]: """Scatter plot tailored for comparing model output with observations. Optionally, with density histogram. @@ -152,6 +162,9 @@ def scatter( >>> cc.sel(observations=['c2','HKNA']).plot.scatter() """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cc = self.cc mod_names = cc.mod_names @@ -191,7 +204,7 @@ def _scatter_one_model( show_points: bool | int | float | None, show_hist: bool | None, show_density: bool | None, - backend: Literal["matplotlib", "plotly"], + backend: Backend, figsize: Tuple[float, float], xlim: Tuple[float, float] | None, ylim: Tuple[float, float] | None, @@ -202,7 +215,7 @@ def _scatter_one_model( skill_table: Union[str, List[str], Mapping[str, str], bool] | None, ax, **kwargs, - ): + ) -> PlotResult: assert ( mod_name in self.cc.mod_names ), f"Model {mod_name} not found in collection {self.cc.mod_names}" @@ -267,33 +280,42 @@ def _scatter_one_model( skill_scores=skill_scores, skill_score_unit=skill_score_unit, ax=ax, + directional=self.is_directional, **kwargs, ) - if backend == "matplotlib" and self.is_directional: - _xtick_directional(ax, xlim) - _ytick_directional(ax, ylim) - return ax - def kde(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: + def kde( + self, + *, + ax=None, + figsize=None, + title=None, + backend: Backend = "matplotlib", + **kwargs, + ) -> PlotResult: """Plot kernel density estimate of observation and model data. Parameters ---------- ax : Axes, optional - matplotlib axes, by default None + matplotlib axes (matplotlib backend only), by default None figsize : tuple, optional - width and height of the figure, by default None + width and height of the figure in inches, by default None title : str, optional plot title, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - passed to pandas.DataFrame.plot.kde() + passed to pandas.DataFrame.plot.kde() (matplotlib backend) or + fig.update_layout() (plotly backend); `bw_method` is passed to + the kernel density estimate by both backends Returns ------- - Axes - matplotlib axes + Axes or plotly.graph_objects.Figure Examples -------- @@ -302,9 +324,32 @@ def kde(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: >>> cc.plot.kde(bw_method='silverman') """ - _, ax = _get_fig_ax(ax, figsize) + validate_backend(backend) + reject_matplotlib_axes(ax, backend) df = self.cc._to_long_dataframe() + title = ( + _default_univarate_title("Density plot", self.cc) + if title is None + else title + ) + + if backend == "plotly": + series = {"Observation": df.obs_val.values} + series.update( + {m: df[df.model == m].mod_val.values for m in self.cc.mod_names} + ) + return _plotly.kde( + series=series, + title=title, + xlabel=self.cc._unit_text, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) + + _, ax = _get_fig_ax(ax, figsize) + ax = df.obs_val.plot.kde( ax=ax, linestyle="dashed", label="Observation", **kwargs ) @@ -315,11 +360,6 @@ def kde(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: ax.set_xlabel(f"{self.cc._unit_text}") - title = ( - _default_univarate_title("Density plot", self.cc) - if title is None - else title - ) ax.set_title(title) ax.legend() @@ -348,12 +388,11 @@ def hist( alpha: float = 0.5, ax=None, figsize: Tuple[float, float] | None = None, + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult | list[PlotResult]: """Plot histogram of specific model and all observations. - Wraps pandas.DataFrame hist() method. - Parameters ---------- bins : int, optional @@ -365,15 +404,20 @@ def hist( alpha : float, optional alpha transparency fraction, by default 0.5 ax : matplotlib axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None figsize : tuple, optional - width and height of the figure, by default None + width and height of the figure in inches, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to df.hist() + other keyword arguments to df.hist() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib axes + Axes or plotly.graph_objects.Figure + one per model, or a list if the collection has multiple models Examples -------- @@ -385,23 +429,25 @@ def hist( pandas.Series.hist matplotlib.axes.Axes.hist """ - - mod_names = self.cc.mod_names - - axes = [] - for mod_name in mod_names: - ax_mod = self._hist_one_model( - mod_name=mod_name, - bins=bins, - title=title, - density=density, - alpha=alpha, - ax=ax, - figsize=figsize, - **kwargs, + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + + figs = [] + for mod_name in self.cc.mod_names: + figs.append( + self._hist_one_model( + mod_name=mod_name, + bins=bins, + title=title, + density=density, + alpha=alpha, + ax=ax, + figsize=figsize, + backend=backend, + **kwargs, + ) ) - axes.append(ax_mod) - return axes[0] if len(axes) == 1 else axes + return figs[0] if len(figs) == 1 else figs def _hist_one_model( self, @@ -413,12 +459,11 @@ def _hist_one_model( alpha: float, ax, figsize: Tuple[float, float] | None, + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult: from ._comparison import MOD_COLORS - _, ax = _get_fig_ax(ax, figsize) - assert ( mod_name in self.cc.mod_names ), f"Model {mod_name} not found in collection" @@ -428,21 +473,37 @@ def _hist_one_model( _default_univarate_title("Histogram", self.cc) if title is None else title ) - cmp = self.cc - df = cmp._to_long_dataframe() + df = self.cc._to_long_dataframe() + obs_color = self.cc[0].data["Observation"].attrs["color"] + xlabel = f"{self.cc[df.observation.iloc[0]]._unit_text}" + + if backend == "plotly": + return _plotly.histogram( + series={ + mod_name: df[df.model == mod_name].mod_val.values, + "observations": df.obs_val.values, + }, + colors=[MOD_COLORS[mod_idx], obs_color], + bins=bins, + density=density, + alpha=alpha, + title=title, + xlabel=xlabel, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) + + _, ax = _get_fig_ax(ax, figsize) + kwargs["alpha"] = alpha kwargs["density"] = density df.mod_val.hist(bins=bins, color=MOD_COLORS[mod_idx], ax=ax, **kwargs) - df.obs_val.hist( - bins=bins, - color=self.cc[0].data["Observation"].attrs["color"], - ax=ax, - **kwargs, - ) + df.obs_val.hist(bins=bins, color=obs_color, ax=ax, **kwargs) ax.legend([mod_name, "observations"]) ax.set_title(title) - ax.set_xlabel(f"{self.cc[df.observation.iloc[0]]._unit_text}") + ax.set_xlabel(xlabel) if density: ax.set_ylabel("density") @@ -463,7 +524,8 @@ def taylor( marker: str = "o", marker_size: float = 6.0, title: str = "Taylor diagram", - ) -> Figure | None: + backend: Backend = "matplotlib", + ) -> FigureResult | None: """Taylor diagram for model skill comparison. Taylor diagram showing model std and correlation to observation @@ -484,10 +546,13 @@ def taylor( size of the marker, by default 6 title : str, optional title of the plot, by default "Taylor diagram" + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" Returns ------- - matplotlib.figure.Figure + matplotlib.figure.Figure or plotly.graph_objects.Figure Examples ------ @@ -539,26 +604,38 @@ def taylor( figsize=figsize, normalize_std=normalize_std, title=title, + backend=backend, ) - def box(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: + def box( + self, + *, + ax=None, + figsize=None, + title=None, + backend: Backend = "matplotlib", + **kwargs, + ) -> PlotResult: """Plot box plot of observations and model data. Parameters ---------- ax : Axes, optional - matplotlib axes, by default None + matplotlib axes (matplotlib backend only), by default None figsize : tuple, optional - width and height of the figure, by default None + width and height of the figure in inches, by default None title : str, optional plot title, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - passed to pandas.DataFrame.plot.box() + passed to pandas.DataFrame.plot.box() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - Axes - matplotlib axes + Axes or plotly.graph_objects.Figure Examples -------- @@ -566,7 +643,8 @@ def box(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: >>> cc.plot.box(showmeans=True) >>> cc.plot.box(ax=ax, title="Box plot") """ - _, ax = _get_fig_ax(ax, figsize) + validate_backend(backend) + reject_matplotlib_axes(ax, backend) df = self.cc._to_long_dataframe() @@ -579,6 +657,22 @@ def box(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: df_model = df[df.model == model] data[model] = df_model.mod_val.values + title = ( + _default_univarate_title("Box plot", self.cc) if title is None else title + ) + + if backend == "plotly": + return _plotly.box( + series=data, + title=title, + ylabel=f"{self.cc._unit_text}", + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) + + _, ax = _get_fig_ax(ax, figsize) + data = {k: pd.Series(v) for k, v in data.items()} df = pd.DataFrame(data) @@ -588,10 +682,6 @@ def box(self, *, ax=None, figsize=None, title=None, **kwargs) -> Axes: ax = df.plot.box(ax=ax, **kwargs) ax.set_ylabel(f"{self.cc._unit_text}") - - title = ( - _default_univarate_title("Box plot", self.cc) if title is None else title - ) ax.set_title(title) if self.is_directional: @@ -606,8 +696,9 @@ def qq( title=None, ax=None, figsize=None, + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult: """Make quantile-quantile (q-q) plot of model data and observations. Primarily used to compare multiple models. @@ -621,26 +712,54 @@ def qq( title : str, optional plot title, default: "Q-Q plot for [observation name]" ax : matplotlib.axes.Axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None figsize : tuple, optional - figure size, by default None + figure size in inches, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to plt.plot() + other keyword arguments to plt.plot() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib axes + Axes or plotly.graph_objects.Figure Examples -------- >>> cc.plot.qq() """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cc = self.cc + df = cc._to_long_dataframe() + title = ( + _default_univarate_title("Q-Q plot for ", self.cc) + if title is None + else title + ) - _, ax = _get_fig_ax(ax, figsize) + if backend == "plotly": + quantile_pairs = {} + for model in cc.mod_names: + df_model = df[df.model == model] + quantile_pairs[model] = quantiles_xy( + df_model.obs_val.values, df_model.mod_val.values, quantiles + ) + return _plotly.qq( + quantiles=quantile_pairs, + title=title, + xlabel="Observation, " + cc._unit_text, + ylabel="Model, " + cc._unit_text, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) - df = cc._to_long_dataframe() + _, ax = _get_fig_ax(ax, figsize) xmin, xmax, ymin, ymax = np.inf, -np.inf, np.inf, -np.inf @@ -676,11 +795,6 @@ def qq( ax.legend() ax.set_xlabel("Observation, " + cc._unit_text) ax.set_ylabel("Model, " + cc._unit_text) - title = ( - _default_univarate_title("Q-Q plot for ", self.cc) - if title is None - else title - ) ax.set_title(title) if self.is_directional: @@ -690,8 +804,15 @@ def qq( return ax def residual_hist( - self, bins=100, title=None, color=None, figsize=None, ax=None, **kwargs - ) -> Axes | list[Axes]: + self, + bins=100, + title=None, + color=None, + figsize=None, + ax=None, + backend: Backend = "matplotlib", + **kwargs, + ) -> PlotResult | list[PlotResult]: """plot histogram of residual values Parameters @@ -703,16 +824,24 @@ def residual_hist( color : str, optional residual color, by default "#8B8D8E" figsize : tuple, optional - figure size, by default None + figure size in inches, by default None ax : Axes | list[Axes], optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to plt.hist() + other keyword arguments to plt.hist() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - Axes | list[Axes] + Axes or plotly.graph_objects.Figure + one per model, or a list if the collection has multiple models """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cc = self.cc if cc.n_models == 1: @@ -723,6 +852,7 @@ def residual_hist( figsize=figsize, ax=ax, mod_name=cc.mod_names[0], + backend=backend, **kwargs, ) @@ -739,6 +869,7 @@ def residual_hist( color=color, figsize=figsize, ax=axs[i], + backend=backend, **kwargs, ) axs[i] = ax_mod @@ -753,24 +884,38 @@ def _residual_hist_one_model( figsize=None, ax=None, mod_name=None, + backend: Backend = "matplotlib", **kwargs, - ) -> Axes: + ) -> PlotResult: """Residual histogram for one model only""" - _, ax = _get_fig_ax(ax, figsize) - df = self.cc.sel(model=mod_name)._to_long_dataframe() residuals = df.mod_val.values - df.obs_val.values - default_color = "#8B8D8E" - color = default_color if color is None else color title = ( _default_univarate_title(f"Residuals, Model {mod_name}", self.cc) if title is None else title ) + xlabel = f"Residuals of {self.cc._unit_text}" + + if backend == "plotly": + return _plotly.residual_hist( + residuals=residuals, + bins=bins, + color=color, + title=title, + xlabel=xlabel, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) + + _, ax = _get_fig_ax(ax, figsize) + + color = RESIDUAL_COLOR if color is None else color ax.hist(residuals, bins=bins, color=color, **kwargs) ax.set_title(title) - ax.set_xlabel(f"Residuals of {self.cc._unit_text}") + ax.set_xlabel(xlabel) if self.is_directional: ticks = np.linspace(-180, 180, 9) @@ -784,7 +929,8 @@ def spatial_overview( ax=None, figsize: Tuple | None = None, title: str | None = None, - ) -> Axes: + backend: Backend = "matplotlib", + ) -> PlotResult: """Plot observation points on a map showing the model domain Parameters @@ -795,18 +941,23 @@ def spatial_overview( figure size, by default None title: str, optional plot title, default empty + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" Returns ------- - matplotlib.axes.Axes - The matplotlib axes object + matplotlib.axes.Axes or plotly.graph_objects.Figure + The axes (matplotlib backend) or figure (plotly backend) """ from ..plotting import spatial_overview obs = [cmp._to_observation() for cmp in self.cc] # TODO how to add model domain(s) - return spatial_overview(obs, ax=ax, figsize=figsize, title=title) + return spatial_overview( + obs, ax=ax, figsize=figsize, title=title, backend=backend + ) def temporal_coverage( self, @@ -815,7 +966,8 @@ def temporal_coverage( ax: Any | None = None, figsize: Any | None = None, title: Any | None = None, - ) -> Axes: + backend: Backend = "matplotlib", + ) -> PlotResult: """Plot graph showing temporal coverage for all observations and models Parameters @@ -831,6 +983,14 @@ def temporal_coverage( size of figure, by default (7, 0.45*n_lines) title: str, optional plot title, default empty + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" + + Returns + ------- + matplotlib.axes.Axes or plotly.graph_objects.Figure + The axes (matplotlib backend) or figure (plotly backend) """ from ..plotting import temporal_coverage @@ -845,4 +1005,5 @@ def temporal_coverage( ax=ax, figsize=figsize, title=title, + backend=backend, ) diff --git a/src/modelskill/comparison/_comparer_plotter.py b/src/modelskill/comparison/_comparer_plotter.py index d226252e2..f11bca680 100644 --- a/src/modelskill/comparison/_comparer_plotter.py +++ b/src/modelskill/comparison/_comparer_plotter.py @@ -1,6 +1,5 @@ from __future__ import annotations from typing import ( - Literal, Union, List, Tuple, @@ -16,11 +15,21 @@ from ._comparison import Comparer import numpy as np # type: ignore +import pandas as pd from .. import metrics as mtr from ..utils import _get_idx import matplotlib.colors as colors +from ..plotting import _plotly +from ..plotting._backend import ( + Backend, + FigureResult, + PlotResult, + reject_matplotlib_axes, + validate_backend, +) from ..plotting._misc import ( + RESIDUAL_COLOR, _get_fig_ax, _xtick_directional, _ytick_directional, @@ -47,9 +56,7 @@ def __init__(self, comparer: Comparer) -> None: self.comparer = comparer self.is_directional = comparer.quantity.is_directional - def __call__( - self, *args, **kwargs - ) -> matplotlib.axes.Axes | list[matplotlib.axes.Axes]: + def __call__(self, *args, **kwargs) -> PlotResult | list[PlotResult]: """Plot scatter plot of modelled vs observed data""" return self.scatter(*args, **kwargs) @@ -60,9 +67,9 @@ def timeseries( ylim: Tuple[float, float] | None = None, ax=None, figsize: Tuple[float, float] | None = None, - backend: str = "matplotlib", + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult: """Timeseries plot showing compared data: observation vs modelled Parameters @@ -72,82 +79,74 @@ def timeseries( ylim : (float, float), optional plot range for the model (ymin, ymax), by default None ax : matplotlib.axes.Axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None figsize : (float, float), optional - figure size, by default None + figure size in inches, by default None backend : str, optional - use "plotly" (interactive) or "matplotlib" backend, + "matplotlib" (static) or "plotly" (interactive), by default "matplotlib" **kwargs - other keyword arguments to fig.update_layout (plotly backend) + other keyword arguments to pandas.Series.plot() (matplotlib + backend) or fig.update_layout() (plotly backend) Returns ------- matplotlib.axes.Axes or plotly.graph_objects.Figure + + Examples + -------- + >>> cmp.plot.timeseries() + >>> cmp.plot.timeseries(backend="plotly") """ from ._comparison import MOD_COLORS - cmp = self.comparer + validate_backend(backend) + reject_matplotlib_axes(ax, backend) - if title is None: - title = cmp.name - - if backend == "matplotlib": - fig, ax = _get_fig_ax(ax, figsize) - for j in range(cmp.n_models): - key = cmp.mod_names[j] - mod = cmp.raw_mod_data[key]._values_as_series - mod.plot(ax=ax, color=MOD_COLORS[j]) - - ax.scatter( - cmp.time, - cmp.data[cmp._obs_name].values, - marker=".", - color=cmp.data[cmp._obs_name].attrs["color"], - ) - ax.set_ylabel(cmp._unit_text) - ax.legend([*cmp.mod_names, cmp._obs_name]) - ax.set_ylim(ylim) - if self.is_directional: - _ytick_directional(ax, ylim) - ax.set_title(title) - return ax - - elif backend == "plotly": # pragma: no cover - import plotly.graph_objects as go # type: ignore - - mod_scatter_list = [] - for j in range(cmp.n_models): - key = cmp.mod_names[j] - mod = cmp.raw_mod_data[key]._values_as_series - mod_scatter_list.append( - go.Scatter( - x=mod.index, - y=mod.values, - name=key, - line=dict(color=MOD_COLORS[j]), - ) - ) - - fig = go.Figure( - [ - *mod_scatter_list, - go.Scatter( - x=cmp.time, - y=cmp.data[cmp._obs_name].values, - name=cmp._obs_name, - mode="markers", - marker=dict(color=cmp.data[cmp._obs_name].attrs["color"]), - ), - ] + cmp = self.comparer + title = cmp.name if title is None else title + + if backend == "plotly": + return _plotly.timeseries( + obs=self._obs_series, + obs_color=cmp.data[cmp._obs_name].attrs["color"], + mods={k: cmp.raw_mod_data[k]._values_as_series for k in cmp.mod_names}, + mod_colors=MOD_COLORS, + title=title, + ylabel=cmp._unit_text, + ylim=ylim, + figsize=figsize, + directional=self.is_directional, + **kwargs, ) - fig.update_layout(title=title, yaxis_title=cmp._unit_text, **kwargs) - fig.update_yaxes(range=ylim) + _, ax = _get_fig_ax(ax, figsize) + for j in range(cmp.n_models): + key = cmp.mod_names[j] + mod = cmp.raw_mod_data[key]._values_as_series + mod.plot(ax=ax, color=MOD_COLORS[j], **kwargs) + + ax.scatter( + cmp.time, + cmp.data[cmp._obs_name].values, + marker=".", + color=cmp.data[cmp._obs_name].attrs["color"], + ) + ax.set_ylabel(cmp._unit_text) + ax.legend([*cmp.mod_names, cmp._obs_name]) + ax.set_ylim(ylim) + if self.is_directional: + _ytick_directional(ax, ylim) + ax.set_title(title) + return ax - return fig - else: - raise ValueError(f"Plotting backend: {backend} not supported") + @property + def _obs_series(self) -> pd.Series: + """Observation values as a named, time-indexed series""" + cmp = self.comparer + return pd.Series( + cmp.data[cmp._obs_name].values, index=cmp.time, name=cmp._obs_name + ) def hist( self, @@ -158,12 +157,11 @@ def hist( figsize: Tuple[float, float] | None = None, density: bool = True, alpha: float = 0.5, + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult | list[PlotResult]: """Plot histogram of model data and observations. - Wraps pandas.DataFrame hist() method. - Parameters ---------- bins : int, optional @@ -171,44 +169,52 @@ def hist( title : str, optional plot title, default: [model name] vs [observation name] ax : matplotlib.axes.Axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None figsize : tuple, optional - figure size, by default None + figure size in inches, by default None density: bool, optional If True, draw and return a probability density alpha : float, optional alpha transparency fraction, by default 0.5 + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to df.plot.hist() + other keyword arguments to df.plot.hist() (matplotlib backend) + or fig.update_layout() (plotly backend) Returns ------- - matplotlib axes + matplotlib.axes.Axes or plotly.graph_objects.Figure + one per model, or a list if the comparer has multiple models See also -------- pandas.Series.plot.hist matplotlib.axes.Axes.hist """ - cmp = self.comparer + validate_backend(backend) + reject_matplotlib_axes(ax, backend) - mod_names = cmp.mod_names + cmp = self.comparer - axes = [] - for mod_name in mod_names: - ax_mod = self._hist_one_model( - mod_name=mod_name, - bins=bins, - title=title, - ax=ax, - figsize=figsize, - density=density, - alpha=alpha, - **kwargs, + figs = [] + for mod_name in cmp.mod_names: + figs.append( + self._hist_one_model( + mod_name=mod_name, + bins=bins, + title=title, + ax=ax, + figsize=figsize, + density=density, + alpha=alpha, + backend=backend, + **kwargs, + ) ) - axes.append(ax_mod) - return axes[0] if len(axes) == 1 else axes + return figs[0] if len(figs) == 1 else figs def _hist_one_model( self, @@ -220,8 +226,9 @@ def _hist_one_model( figsize: Tuple[float, float] | None, density: bool | None, alpha: float | None, + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult: from ._comparison import MOD_COLORS # TODO move to here cmp = self.comparer @@ -229,6 +236,24 @@ def _hist_one_model( mod_idx = _get_idx(mod_name, cmp.mod_names) title = f"{mod_name} vs {cmp.name}" if title is None else title + obs_color = cmp.data[cmp._obs_name].attrs["color"] + + if backend == "plotly": + return _plotly.histogram( + series={ + mod_name: cmp.data[mod_name].values, + cmp._obs_name: cmp.data[cmp._obs_name].values, + }, + colors=[MOD_COLORS[mod_idx], obs_color], + bins=bins if bins is not None else 100, + density=bool(density), + alpha=alpha if alpha is not None else 0.5, + title=title, + xlabel=cmp._unit_text, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) _, ax = _get_fig_ax(ax, figsize) @@ -242,9 +267,7 @@ def _hist_one_model( .hist(bins=bins, color=MOD_COLORS[mod_idx], **kwargs) ) - cmp.data[cmp._obs_name].to_series().hist( - bins=bins, color=cmp.data[cmp._obs_name].attrs["color"], **kwargs - ) + cmp.data[cmp._obs_name].to_series().hist(bins=bins, color=obs_color, **kwargs) ax.legend([mod_name, cmp._obs_name]) ax.set_title(title) ax.set_xlabel(f"{cmp._unit_text}") @@ -258,25 +281,35 @@ def _hist_one_model( return ax - def kde(self, ax=None, title=None, figsize=None, **kwargs) -> matplotlib.axes.Axes: + def kde( + self, + ax=None, + title=None, + figsize=None, + backend: Backend = "matplotlib", + **kwargs, + ) -> PlotResult: """Plot kde (kernel density estimates of distributions) of model data and observations. - Wraps pandas.DataFrame kde() method. - Parameters ---------- ax : matplotlib.axes.Axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None title : str, optional plot title, default: "KDE plot for [observation name]" figsize : tuple, optional - figure size, by default None + figure size in inches, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to df.plot.kde() + other keyword arguments to df.plot.kde() (matplotlib backend) + or fig.update_layout() (plotly backend); `bw_method` is passed + to the kernel density estimate by both backends Returns ------- - matplotlib.axes.Axes + matplotlib.axes.Axes or plotly.graph_objects.Figure Examples -------- @@ -289,7 +322,23 @@ def kde(self, ax=None, title=None, figsize=None, **kwargs) -> matplotlib.axes.Ax -------- pandas.Series.plot.kde """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cmp = self.comparer + title = f"KDE plot for {cmp.name}" if title is None else title + + if backend == "plotly": + series = {"Observation": cmp.data.Observation.values} + series.update({m: cmp.data[m].values for m in cmp.mod_names}) + return _plotly.kde( + series=series, + title=title, + xlabel=cmp._unit_text, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) _, ax = _get_fig_ax(ax, figsize) @@ -308,7 +357,6 @@ def kde(self, ax=None, title=None, figsize=None, **kwargs) -> matplotlib.axes.Ax ax.yaxis.set_visible(False) ax.tick_params(axis="y", which="both", length=0) ax.set_ylabel("") - title = f"KDE plot for {cmp.name}" if title is None else title ax.set_title(title) # remove box around plot @@ -328,8 +376,9 @@ def qq( title=None, ax=None, figsize=None, + backend: Backend = "matplotlib", **kwargs, - ): + ) -> PlotResult: """Make quantile-quantile (q-q) plot of model data and observations. Primarily used to compare multiple models. @@ -343,26 +392,48 @@ def qq( title : str, optional plot title, default: "Q-Q plot for [observation name]" ax : matplotlib.axes.Axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None figsize : tuple, optional - figure size, by default None + figure size in inches, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to plt.plot() + other keyword arguments to plt.plot() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib axes + matplotlib.axes.Axes or plotly.graph_objects.Figure Examples -------- >>> cmp.plot.qq() """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cmp = self.comparer + title = f"Q-Q plot for {cmp.name}" if title is None else title + x = cmp.data.Observation.values + + if backend == "plotly": + return _plotly.qq( + quantiles={ + m: quantiles_xy(x, cmp.data[m].values, quantiles) + for m in cmp.mod_names + }, + title=title, + xlabel="Observation, " + cmp._unit_text, + ylabel="Model, " + cmp._unit_text, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) _, ax = _get_fig_ax(ax, figsize) - x = cmp.data.Observation.values xmin, xmax = x.min(), x.max() ymin, ymax = np.inf, -np.inf @@ -401,7 +472,7 @@ def qq( ax.legend() ax.set_xlabel("Observation, " + cmp._unit_text) ax.set_ylabel("Model, " + cmp._unit_text) - ax.set_title(title or f"Q-Q plot for {cmp.name}") + ax.set_title(title) if self.is_directional: _xtick_directional(ax) @@ -409,25 +480,35 @@ def qq( return ax - def box(self, *, ax=None, title=None, figsize=None, **kwargs): + def box( + self, + *, + ax=None, + title=None, + figsize=None, + backend: Backend = "matplotlib", + **kwargs, + ) -> PlotResult: """Make a box plot of model data and observations. - Wraps pandas.DataFrame boxplot() method. - Parameters ---------- ax : matplotlib.axes.Axes, optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None title : str, optional plot title, default: [observation name] figsize : tuple, optional - figure size, by default None + figure size in inches, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to df.boxplot() + other keyword arguments to df.boxplot() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib axes + matplotlib.axes.Axes or plotly.graph_objects.Figure Examples -------- @@ -440,15 +521,29 @@ def box(self, *, ax=None, title=None, figsize=None, **kwargs): pandas.DataFrame.boxplot matplotlib.pyplot.boxplot """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cmp = self.comparer + title = cmp.name if title is None else title + cols = ["Observation"] + cmp.mod_names + + if backend == "plotly": + return _plotly.box( + series={c: cmp.data[c].values for c in cols}, + title=title, + ylabel=cmp._unit_text, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) _, ax = _get_fig_ax(ax, figsize) - cols = ["Observation"] + cmp.mod_names df = cmp.data[cols].to_dataframe()[cols] df.boxplot(ax=ax, **kwargs) ax.set_ylabel(cmp._unit_text) - ax.set_title(title or cmp.name) + ax.set_title(title) if self.is_directional: _ytick_directional(ax) @@ -466,7 +561,7 @@ def scatter( show_hist: bool | None = None, show_density: bool | None = None, norm: colors.Normalize | None = None, - backend: Literal["matplotlib", "plotly"] = "matplotlib", + backend: Backend = "matplotlib", figsize: Tuple[float, float] = (8, 8), xlim: Tuple[float, float] | None = None, ylim: Tuple[float, float] | None = None, @@ -477,7 +572,7 @@ def scatter( skill_table: Union[str, List[str], Mapping[str, str], bool] | None = None, ax: matplotlib.axes.Axes | None = None, **kwargs, - ) -> matplotlib.axes.Axes | list[matplotlib.axes.Axes]: + ) -> PlotResult | list[PlotResult]: """Scatter plot tailored for model-observation comparison. Optionally, with density histogram. @@ -596,7 +691,7 @@ def _scatter_one_model( show_hist: bool | None, show_density: bool | None, norm: colors.Normalize | None, - backend: Literal["matplotlib", "plotly"], + backend: Backend, figsize: Tuple[float, float], xlim: Tuple[float, float] | None, ylim: Tuple[float, float] | None, @@ -606,7 +701,7 @@ def _scatter_one_model( ylabel: str | None, skill_table: Union[str, List[str], Mapping[str, str], bool] | None, **kwargs, - ): + ) -> PlotResult: """Scatter plot for one model only""" cmp = self.comparer @@ -665,13 +760,10 @@ def _scatter_one_model( ylabel=ylabel, skill_scores=skill_scores, skill_score_unit=skill_score_unit, + directional=self.is_directional, **kwargs, ) - if backend == "matplotlib" and self.is_directional: - _xtick_directional(ax, xlim) - _ytick_directional(ax, ylim) - return ax def taylor( @@ -682,7 +774,8 @@ def taylor( marker: str = "o", marker_size: float = 6.0, title: str = "Taylor diagram", - ): + backend: Backend = "matplotlib", + ) -> FigureResult | None: """Taylor diagram for model skill comparison. Taylor diagram showing model std and correlation to observation @@ -700,10 +793,13 @@ def taylor( size of the marker, by default 6 title : str, optional title of the plot, by default "Taylor diagram" + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" Returns ------- - matplotlib.figure.Figure + matplotlib.figure.Figure or plotly.graph_objects.Figure Examples -------- @@ -764,11 +860,19 @@ def taylor( obs_text=f"Obs: {cmp.name}", normalize_std=normalize_std, title=title, + backend=backend, ) def residual_hist( - self, bins=100, title=None, color=None, figsize=None, ax=None, **kwargs - ) -> matplotlib.axes.Axes | list[matplotlib.axes.Axes]: + self, + bins=100, + title=None, + color=None, + figsize=None, + ax=None, + backend: Backend = "matplotlib", + **kwargs, + ) -> PlotResult | list[PlotResult]: """plot histogram of residual values Parameters @@ -780,16 +884,24 @@ def residual_hist( color : str, optional residual color, by default "#8B8D8E" figsize : tuple, optional - figure size, by default None + figure size in inches, by default None ax : matplotlib.axes.Axes | list[matplotlib.axes.Axes], optional - axes to plot on, by default None + axes to plot on (matplotlib backend only), by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to plt.hist() + other keyword arguments to plt.hist() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib.axes.Axes | list[matplotlib.axes.Axes] + matplotlib.axes.Axes or plotly.graph_objects.Figure + one per model, or a list if the comparer has multiple models """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + cmp = self.comparer if cmp.n_models == 1: @@ -800,6 +912,7 @@ def residual_hist( figsize=figsize, ax=ax, mod_name=cmp.mod_names[0], + backend=backend, **kwargs, ) @@ -816,6 +929,7 @@ def residual_hist( color=color, figsize=figsize, ax=axs[i], + backend=backend, **kwargs, ) axs[i] = ax_mod @@ -830,21 +944,36 @@ def _residual_hist_one_model( figsize=None, ax=None, mod_name=None, + backend: Backend = "matplotlib", **kwargs, - ) -> matplotlib.axes.Axes: + ) -> PlotResult: """Residual histogram for one model only""" - _, ax = _get_fig_ax(ax, figsize) - - default_color = "#8B8D8E" - color = default_color if color is None else color + cmp = self.comparer title = ( - f"Residuals, Observation: {self.comparer.name}, Model: {mod_name}" + f"Residuals, Observation: {cmp.name}, Model: {mod_name}" if title is None else title ) - ax.hist(self.comparer._residual, bins=bins, color=color, **kwargs) + xlabel = f"Residuals of {cmp._unit_text}" + + if backend == "plotly": + return _plotly.residual_hist( + residuals=cmp._residual, + bins=bins, + color=color, + title=title, + xlabel=xlabel, + figsize=figsize, + directional=self.is_directional, + **kwargs, + ) + + _, ax = _get_fig_ax(ax, figsize) + + color = RESIDUAL_COLOR if color is None else color + ax.hist(cmp._residual, bins=bins, color=color, **kwargs) ax.set_title(title) - ax.set_xlabel(f"Residuals of {self.comparer._unit_text}") + ax.set_xlabel(xlabel) if self.is_directional: ticks = np.linspace(-180, 180, 9) diff --git a/src/modelskill/plotting/_backend.py b/src/modelskill/plotting/_backend.py new file mode 100644 index 000000000..e2cb4b19d --- /dev/null +++ b/src/modelskill/plotting/_backend.py @@ -0,0 +1,110 @@ +"""Plotting backend selection. + +modelskill plots can be rendered by either of two backends: + +* ``"matplotlib"`` - static, report-friendly figures, returns + :class:`matplotlib.axes.Axes` +* ``"plotly"`` - interactive figures with zoom/hover, returns + :class:`plotly.graph_objects.Figure` + +Both backends accept the same plot arguments (``title``, ``figsize``, +``xlim``, ``ylim``, ...). Backend-specific extras are passed via +``**kwargs``: to the underlying matplotlib/pandas call for the +matplotlib backend, and to :meth:`plotly.graph_objects.Figure.update_layout` +for the plotly backend. + +This module holds the backend vocabulary only; the plotly renderers and the +plotly layout interop live in `_plotly.py`. plotly is an optional dependency, +install it with ``pip install "modelskill[plotly]"``. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING, Any, Literal, Tuple + +import numpy as np +from typing_extensions import TypeAlias + +if TYPE_CHECKING: + from matplotlib.axes import Axes + from matplotlib.figure import Figure + import plotly.graph_objects as go + +Backend = Literal["matplotlib", "plotly"] + +# What a plot returns depends on the backend: axes for matplotlib, a figure for +# plotly. A few plots (taylor) return a matplotlib figure rather than axes. +PlotResult: TypeAlias = "Axes | go.Figure" +FigureResult: TypeAlias = "Figure | go.Figure" + +BACKENDS: Tuple[Backend, ...] = ("matplotlib", "plotly") + + +def validate_backend(backend: str) -> Backend: + """Check that a backend name is supported. + + Parameters + ---------- + backend : str + name of the plotting backend + + Returns + ------- + str + the validated backend name + + Raises + ------ + ValueError + if the backend is not one of the supported backends + """ + if backend not in BACKENDS: + raise ValueError( + f"Invalid backend '{backend}'. Valid options are: {list(BACKENDS)}" + ) + return backend # type: ignore[return-value] + + +def reject_matplotlib_axes(ax: Any, backend: str) -> None: + """Raise if matplotlib axes are passed to a non-matplotlib backend. + + Parameters + ---------- + ax : matplotlib.axes.Axes or None + the axes argument given by the user + backend : str + the selected backend + + Raises + ------ + ValueError + if ``ax`` is not None and the backend is not matplotlib + """ + if ax is not None and backend != "matplotlib": + raise ValueError( + f"Cannot pass matplotlib axes to the '{backend}' backend. " + f"The '{backend}' backend returns a new figure." + ) + + +def directional_ticks( + lim: Tuple[float, float] | None = None, n_sectors: int = 8 +) -> np.ndarray: + """Tick values for a directional (0-360 degrees) axis. + + Parameters + ---------- + lim : (float, float), optional + axis limits to clip the ticks to, by default None + n_sectors : int, optional + number of sectors, by default 8 + + Returns + ------- + np.ndarray + tick values + """ + ticks = np.linspace(0, 360, n_sectors + 1) + if lim is not None: + ticks = ticks[(ticks >= lim[0]) & (ticks <= lim[1])] + return ticks diff --git a/src/modelskill/plotting/_misc.py b/src/modelskill/plotting/_misc.py index e41a8e214..1b18bdb7a 100644 --- a/src/modelskill/plotting/_misc.py +++ b/src/modelskill/plotting/_misc.py @@ -10,6 +10,51 @@ from ..metrics import metric_has_units, defined_metrics, get_display_name from ..obs import unit_display_name +# grey used for residual histograms, shared by both plotting backends +RESIDUAL_COLOR = "#8B8D8E" + + +def series_range(series: Sequence) -> Tuple[float, float]: + """Combined min/max across a sequence of arrays, ignoring NaN. + + Parameters + ---------- + series : Sequence + arrays to take the range over + + Returns + ------- + (float, float) + overall minimum and maximum + """ + values = np.concatenate([np.asarray(s, dtype=float).ravel() for s in series]) + return float(np.nanmin(values)), float(np.nanmax(values)) + + +def reglabel(slope: float, intercept: float, fit_to_quantiles: bool) -> str: + """Legend label for a regression line. + + Parameters + ---------- + slope : float + slope of the fitted line + intercept : float + intercept of the fitted line + fit_to_quantiles : bool + whether the line was fitted to the quantiles rather than to all data + + Returns + ------- + str + label text + """ + sign = "" if intercept < 0 else "+" + if fit_to_quantiles: + fit = "QQ fit" + else: + fit = "Fit" + return f"{fit}: y={slope:.2f}x{sign}{intercept:.2f}" + def _get_ax(ax=None, figsize=None): if ax is None: diff --git a/src/modelskill/plotting/_plotly.py b/src/modelskill/plotting/_plotly.py new file mode 100644 index 000000000..82c3ae412 --- /dev/null +++ b/src/modelskill/plotting/_plotly.py @@ -0,0 +1,980 @@ +"""plotly renderers for the modelskill plots. + +Every function here takes plain data (arrays, series, labels) and returns a +:class:`plotly.graph_objects.Figure`. They are the plotly counterparts of the +matplotlib code in the plotter classes and in `_scatter.py`, and are selected +by the ``backend="plotly"`` argument on the plot methods. +""" + +from __future__ import annotations + +from types import ModuleType +from typing import TYPE_CHECKING, Any, Dict, List, Mapping, Sequence, Tuple + +import numpy as np +import pandas as pd + +if TYPE_CHECKING: + import plotly.graph_objects as go + +from ..metrics import _linear_regression +from ..settings import options +from ._backend import directional_ticks +from ._misc import RESIDUAL_COLOR, format_skill_table, reglabel, series_range + +# plotly sizes are in pixels, matplotlib figsize is in inches +PIXELS_PER_INCH = 100 + +_PLOTLY_INSTALL_HINT = ( + "The 'plotly' backend requires the optional plotly package. " + 'Install it with `pip install "modelskill[plotly]"`.' +) + + +def import_plotly_go() -> ModuleType: + """Import plotly.graph_objects with an actionable error if it is missing. + + Returns + ------- + module + the ``plotly.graph_objects`` module + + Raises + ------ + ImportError + if plotly is not installed + """ + try: + import plotly.graph_objects as go # type: ignore + except ImportError as e: + raise ImportError(_PLOTLY_INSTALL_HINT) from e + return go + + +def figsize_to_layout(figsize: Tuple[float, float] | None) -> Dict[str, float]: + """Translate a matplotlib figsize (inches) to plotly width/height (pixels). + + Parameters + ---------- + figsize : (float, float), optional + width and height in inches, by default None + + Returns + ------- + dict + ``{"width": ..., "height": ...}``, empty if figsize is None + """ + if figsize is None: + return {} + width, height = figsize + return {"width": width * PIXELS_PER_INCH, "height": height * PIXELS_PER_INCH} + + +def apply_layout( + fig: Any, + *, + figsize: Tuple[float, float] | None = None, + **kwargs: Any, +) -> Any: + """Apply modelskill and user layout arguments to a plotly figure. + + ``figsize`` is translated to plotly's ``width``/``height``; an explicit + ``width``/``height`` in ``kwargs`` wins. Remaining ``kwargs`` are passed + to ``fig.update_layout``. + + Parameters + ---------- + fig : plotly.graph_objects.Figure + figure to update + figsize : (float, float), optional + width and height in inches, by default None + **kwargs + keyword arguments for ``plotly.graph_objects.Figure.update_layout`` + + Returns + ------- + plotly.graph_objects.Figure + the updated figure + + Raises + ------ + ValueError + if a keyword argument is not a valid plotly layout property + """ + layout = {**figsize_to_layout(figsize), **kwargs} + layout = {k: v for k, v in layout.items() if v is not None} + try: + fig.update_layout(**layout) + except ValueError as e: + raise ValueError(_layout_error_message(layout, e)) from e + return fig + + +def _layout_error_message(layout: Dict[str, Any], error: ValueError) -> str: + invalid = _invalid_layout_keys(layout) + named = ", ".join(f"'{k}'" for k in invalid) if invalid else "argument" + return ( + f"Invalid plotly layout argument: {named}. The plotly backend passes " + "keyword arguments to plotly.graph_objects.Figure.update_layout, so " + "matplotlib-only arguments are not accepted. Valid layout properties " + "are documented at https://plotly.com/python/reference/layout/.\n" + f"Original plotly error: {error}" + ) + + +def _invalid_layout_keys(layout: Dict[str, Any]) -> List[str]: + go = import_plotly_go() + invalid = [] + for key in layout: + try: + go.Layout(**{key: layout[key]}) + except ValueError: + invalid.append(key) + return invalid + + +def directional_axis( + fig: Any, axis: str, lim: Tuple[float, float] | None = None +) -> None: + """Make a plotly axis directional (0-360 degrees with sector ticks). + + Parameters + ---------- + fig : plotly.graph_objects.Figure + figure to update + axis : str + "x" or "y" + lim : (float, float), optional + axis range, by default None which means (0, 360) + """ + ticks = directional_ticks(lim) + update = fig.update_xaxes if axis == "x" else fig.update_yaxes + if len(ticks) > 2: + update(tickmode="array", tickvals=ticks) + update(range=lim if lim is not None else (0, 360)) + + +def timeseries( + *, + obs: pd.Series, + obs_color: str, + mods: Mapping[str, pd.Series], + mod_colors: Sequence[str], + title: str | None = None, + ylabel: str | None = None, + ylim: Tuple[float, float] | None = None, + figsize: Tuple[float, float] | None = None, + directional: bool = False, + **kwargs: Any, +) -> go.Figure: + """Timeseries of observation and model data.""" + go = import_plotly_go() + + traces = [ + go.Scatter( + x=mod.index, + y=mod.values, + name=name, + line=dict(color=mod_colors[j]), + ) + for j, (name, mod) in enumerate(mods.items()) + ] + traces.append( + go.Scatter( + x=obs.index, + y=obs.values, + name=str(obs.name), + mode="markers", + marker=dict(color=obs_color), + ) + ) + + fig = go.Figure(traces) + apply_layout(fig, figsize=figsize, title=title, yaxis_title=ylabel, **kwargs) + if directional: + directional_axis(fig, "y", ylim) + else: + fig.update_yaxes(range=ylim) + return fig + + +def line( + *, + series: pd.Series, + color: str | None = None, + title: str | None = None, + ylabel: str | None = None, + figsize: Tuple[float, float] | None = None, + **kwargs: Any, +) -> go.Figure: + """Line plot of a single time series.""" + go = import_plotly_go() + + fig = go.Figure( + go.Scatter( + x=series.index, + y=series.values, + name=str(series.name), + line=dict(color=color), + ) + ) + apply_layout(fig, figsize=figsize, title=title, yaxis_title=ylabel, **kwargs) + return fig + + +def histogram( + *, + series: Mapping[str, np.ndarray], + colors: Sequence[str], + bins: int | Sequence = 100, + density: bool = True, + alpha: float = 0.5, + title: str | None = None, + xlabel: str | None = None, + figsize: Tuple[float, float] | None = None, + directional: bool = False, + **kwargs: Any, +) -> go.Figure: + """Overlaid histograms of the given named data series.""" + go = import_plotly_go() + + nbins, bin_edges = _hist_bins(bins) + + traces = [] + for i, (name, values) in enumerate(series.items()): + traces.append( + go.Histogram( + x=values, + name=name, + nbinsx=nbins, + xbins=bin_edges, + histnorm="probability density" if density else None, + opacity=alpha, + marker=dict(color=colors[i]), + ) + ) + + fig = go.Figure(traces) + apply_layout( + fig, + figsize=figsize, + title=title, + xaxis_title=xlabel, + yaxis_title="density" if density else "count", + barmode="overlay", + **kwargs, + ) + if directional: + directional_axis(fig, "x") + return fig + + +def _hist_bins(bins: int | Sequence) -> Tuple[int | None, Any]: + """Translate a matplotlib `bins` argument to plotly nbinsx/xbins. + + An int becomes plotly's `nbinsx`, which is an upper bound rather than an + exact bin count. A sequence of edges becomes `xbins`, which can only + express uniformly spaced bins. + """ + if isinstance(bins, (int, np.integer)): + return int(bins), None + edges = np.asarray(bins, dtype=float) + if edges.size < 2: + raise ValueError("`bins` must be an int or a sequence of at least two edges") + widths = np.diff(edges) + if not np.allclose(widths, widths[0]): + raise ValueError( + "the plotly backend supports only uniformly spaced bin edges, " + f"got widths {widths}" + ) + return None, dict(start=edges[0], end=edges[-1], size=widths[0]) + + +def kde( + *, + series: Mapping[str, np.ndarray], + title: str | None = None, + xlabel: str | None = None, + figsize: Tuple[float, float] | None = None, + directional: bool = False, + bw_method: Any = None, + n_points: int = 200, + **kwargs: Any, +) -> go.Figure: + """Kernel density estimates of the given named data series. + + The first series is drawn dashed, matching the matplotlib backend where + the observation is dashed and the models are solid. + """ + go = import_plotly_go() + from scipy.stats import gaussian_kde + + xmin, xmax = series_range(list(series.values())) + span = xmax - xmin + grid = np.linspace(xmin - 0.1 * span, xmax + 0.1 * span, n_points) + + traces = [] + for i, (name, values) in enumerate(series.items()): + density = gaussian_kde(np.asarray(values, dtype=float), bw_method=bw_method) + traces.append( + go.Scatter( + x=grid, + y=density(grid), + name=name, + mode="lines", + line=dict(dash="dash" if i == 0 else "solid"), + ) + ) + + fig = go.Figure(traces) + apply_layout(fig, figsize=figsize, title=title, xaxis_title=xlabel, **kwargs) + # the density scale carries no information the user needs, as in matplotlib + fig.update_yaxes(visible=False) + if directional: + directional_axis(fig, "x") + return fig + + +def qq( + *, + quantiles: Mapping[str, Tuple[np.ndarray, np.ndarray]], + title: str | None = None, + xlabel: str | None = None, + ylabel: str | None = None, + figsize: Tuple[float, float] | None = None, + directional: bool = False, + **kwargs: Any, +) -> go.Figure: + """Quantile-quantile plot with a 1:1 line, one trace per model.""" + go = import_plotly_go() + + all_values = [v for pair in quantiles.values() for v in pair] + xymin, xymax = series_range(all_values) + + traces = [ + go.Scatter( + x=[xymin, xymax], + y=[xymin, xymax], + name=options.plot.scatter.oneone_line.label, + mode="lines", + line=dict(color=options.plot.scatter.oneone_line.color), + ) + ] + for name, (xq, yq) in quantiles.items(): + traces.append( + go.Scatter(x=xq, y=yq, name=name, mode="lines+markers", marker=dict(size=4)) + ) + + fig = go.Figure(traces) + apply_layout( + fig, + figsize=figsize, + title=title, + xaxis_title=xlabel, + yaxis_title=ylabel, + yaxis=dict(scaleanchor="x", scaleratio=1), + **kwargs, + ) + if directional: + directional_axis(fig, "x") + directional_axis(fig, "y") + else: + fig.update_xaxes(range=(xymin, xymax)) + fig.update_yaxes(range=(xymin, xymax)) + return fig + + +def box( + *, + series: Mapping[str, np.ndarray], + title: str | None = None, + ylabel: str | None = None, + figsize: Tuple[float, float] | None = None, + directional: bool = False, + **kwargs: Any, +) -> go.Figure: + """Box plot with one box per named data series.""" + go = import_plotly_go() + + traces = [ + go.Box(y=np.asarray(values, dtype=float), name=name) + for name, values in series.items() + ] + + fig = go.Figure(traces) + apply_layout( + fig, + figsize=figsize, + title=title, + yaxis_title=ylabel, + showlegend=False, + **kwargs, + ) + if directional: + directional_axis(fig, "y") + return fig + + +def residual_hist( + *, + residuals: np.ndarray, + bins: int | Sequence = 100, + color: str | None = None, + title: str | None = None, + xlabel: str | None = None, + figsize: Tuple[float, float] | None = None, + directional: bool = False, + **kwargs: Any, +) -> go.Figure: + """Histogram of model residuals.""" + go = import_plotly_go() + + nbins, bin_edges = _hist_bins(bins) + + fig = go.Figure( + go.Histogram( + x=residuals, + nbinsx=nbins, + xbins=bin_edges, + marker=dict(color=color or RESIDUAL_COLOR), + ) + ) + apply_layout( + fig, + figsize=figsize, + title=title, + xaxis_title=xlabel, + yaxis_title="count", + showlegend=False, + **kwargs, + ) + if directional: + fig.update_xaxes( + tickmode="array", tickvals=np.linspace(-180, 180, 9), range=(-180, 180) + ) + return fig + + +def scatter( + *, + x, + y, + x_sample, + y_sample, + z, + xq, + yq, + x_trend, + show_density, + show_points, + norm, # matplotlib-only, plotly scales its own colorbar + show_hist, + nbins_hist, + reg_method, + xlabel, + ylabel, + figsize, + xlim, + ylim, + title, + skill_scores, + skill_score_unit, + fit_to_quantiles, + directional=False, + **kwargs, +) -> go.Figure: + """Scatter plot of observation vs model, with 1:1 line and regression.""" + go = import_plotly_go() + + data = [ + go.Scatter(x=xlim, y=xlim, name="1:1", mode="lines", line=dict(color="blue")), + ] + + if reg_method: + if fit_to_quantiles: + slope, intercept = _linear_regression( + obs=xq, model=yq, reg_method=reg_method + ) + else: + slope, intercept = _linear_regression(obs=x, model=y, reg_method=reg_method) + + regression_line = go.Scatter( + x=x_trend, + y=intercept + slope * x_trend, + name=reglabel( + slope=slope, intercept=intercept, fit_to_quantiles=fit_to_quantiles + ), + mode="lines", + line=dict(color="red"), + ) + data.append(regression_line) + + if show_hist: + data.append( + go.Histogram2d( + x=x, + y=y, + nbinsx=nbins_hist, + nbinsy=nbins_hist, + colorscale=[ + [0.0, "rgba(0,0,0,0)"], + [0.1, "purple"], + [0.5, "green"], + [1.0, "yellow"], + ], + colorbar=dict(title="# of points"), + ) + ) + + if show_points is None or show_points: + if show_density: + c = z + cbar = dict(thickness=20, title="# of points") + else: + c = "black" + cbar = None + data.append( + go.Scatter( + x=x_sample, + y=y_sample, + mode="markers", + name="Data", + marker=dict(color=c, opacity=0.5, size=3.0, colorbar=cbar), + ) + ) + if len(xq) > 0: + data.append( + go.Scatter( + x=xq, + y=yq, + name=options.plot.scatter.quantiles.label, + mode="markers", + marker_symbol="x", + marker_color=options.plot.scatter.quantiles.color, + marker_line_color="midnightblue", + marker_line_width=0.6, + ) + ) + + fig = go.Figure(data=data) + apply_layout( + fig, + figsize=figsize, + legend=dict(x=0.01, y=0.99), + yaxis=dict(scaleanchor="x", scaleratio=1), + title=dict(text=title, xanchor="center", yanchor="top", x=0.5, y=0.9), + yaxis_title=ylabel, + xaxis_title=xlabel, + **kwargs, + ) + if directional: + directional_axis(fig, "x", xlim) + directional_axis(fig, "y", ylim) + else: + fig.update_xaxes(range=xlim, nticks=10) + fig.update_yaxes(range=ylim, nticks=10) + + if skill_scores is not None: + _add_skill_table(fig, skill_scores=skill_scores, unit=skill_score_unit) + + return fig + + +def _add_skill_table(fig: Any, *, skill_scores: Mapping[str, float], unit: str) -> None: + table = format_skill_table(skill_scores=skill_scores, unit=unit) + lines = [ + f"{row['name']:<6} {row['sep']} {row['value']:<6}" + for _, row in table.iterrows() + ] + fig.add_annotation( + x=0.99, + y=0.01, + xref="paper", + yref="paper", + text="
".join(lines), + showarrow=False, + align="left", + bordercolor="black", + borderwidth=1, + borderpad=4, + bgcolor="white", + font=dict(family="Consolas, 'Liberation Mono', monospace"), + ) + + +def taylor( + *, + points: Sequence[Any], + obs_std: float, + obs_text: str = "Observations", + normalize_std: bool = False, + title: str = "Taylor diagram", + figsize: Tuple[float, float] | None = None, + n_rms_contours: int = 5, + **kwargs: Any, +) -> go.Figure: + """Taylor diagram in a single-quadrant polar plot, r=std and theta=arccos(cc). + + Parameters + ---------- + points : Sequence[TaylorPoint] + the model points to show + obs_std : float + standard deviation of the observations (the reference radius) + obs_text : str, optional + label of the reference point, by default "Observations" + normalize_std : bool, optional + model std is normalized with observation std, by default False + title : str, optional + plot title, by default "Taylor diagram" + figsize : (float, float), optional + figure size in inches, by default None + n_rms_contours : int, optional + number of dotted centered-RMS-difference contours, by default 5 + **kwargs + keyword arguments for fig.update_layout + """ + go = import_plotly_go() + + stds = [p.std / p.obs_std if normalize_std else p.std for p in points] + rmax = max([obs_std, *stds]) * 1.4 + + traces = [ + _rms_contour( + go, obs_std=obs_std, radius=r, name=f"RMSD={r:.2g}", showlegend=i == 0 + ) + for i, r in enumerate(_rms_contour_radii(obs_std, rmax, n_rms_contours)) + ] + + traces.append( + go.Scatterpolar( + r=[obs_std], + theta=[0.0], + name=obs_text, + mode="markers", + marker=dict(symbol="star", size=12, color="black"), + ) + ) + for p, std in zip(points, stds): + traces.append( + go.Scatterpolar( + r=[std], + theta=[np.degrees(np.arccos(np.clip(p.cc, -1.0, 1.0)))], + name=p.name, + mode="markers", + marker=dict(size=2 * p.marker_size), + hovertemplate=f"{p.name}
std=%{{r:.3g}}
cc={p.cc:.3f}", + ) + ) + + cc_ticks = np.array([0.0, 0.2, 0.4, 0.6, 0.7, 0.8, 0.9, 0.95, 0.99, 1.0]) + fig = go.Figure(traces) + fig.update_layout( + polar=dict( + sector=[0, 90], + radialaxis=dict( + range=[0, rmax], + title=dict( + text="Std. dev." + (" (normalized)" if normalize_std else "") + ), + angle=45, + ), + angularaxis=dict( + direction="counterclockwise", + tickmode="array", + tickvals=np.degrees(np.arccos(cc_ticks)), + ticktext=[f"{c:g}" for c in cc_ticks], + ), + ) + ) + apply_layout(fig, figsize=figsize, title=title, **kwargs) + return fig + + +def _rms_contour_radii(obs_std: float, rmax: float, n: int) -> list[float]: + """Radii of the centered-RMS-difference contours to draw""" + step = rmax / (n + 1) + return [step * (i + 1) for i in range(n)] + + +def _rms_contour( + go: Any, *, obs_std: float, radius: float, name: str, showlegend: bool +) -> Any: + """A circle of constant centered RMS difference, in Taylor diagram polar coordinates""" + t = np.linspace(0, 2 * np.pi, 180) + x = obs_std + radius * np.cos(t) + y = radius * np.sin(t) + keep = y >= 0 + x, y = x[keep], y[keep] + return go.Scatterpolar( + r=np.hypot(x, y), + theta=np.degrees(np.arctan2(y, x)), + mode="lines", + line=dict(color="lightgray", dash="dot", width=1), + name=name, + legendgroup="rmsd", + showlegend=showlegend, + hoverinfo="skip", + ) + + +def temporal_coverage( + *, + lines: Sequence[Tuple[str, Any, bool]], + xlim: Tuple[Any, Any] | None = None, + title: str | None = None, + figsize: Tuple[float, float] | None = None, + **kwargs: Any, +) -> go.Figure: + """Temporal coverage of observations and models, one row per data source. + + Parameters + ---------- + lines : Sequence of (name, times, is_model) + the rows to draw; models are drawn as a line from first to last time, + observations as markers at every time + xlim : (datetime, datetime), optional + limit the time axis, by default None + title : str, optional + plot title, by default None + figsize : (float, float), optional + figure size in inches, by default None + **kwargs + keyword arguments for fig.update_layout + """ + go = import_plotly_go() + + traces = [] + for name, times, is_model in lines: + if is_model: + x, mode = [times[0], times[-1]], "lines" + else: + x, mode = list(times), "markers" + traces.append( + go.Scatter( + x=x, + y=[name] * len(x), + name=name, + mode=mode, + marker=dict(symbol="line-ns", size=8, line=dict(width=1)), + ) + ) + + fig = go.Figure(traces) + apply_layout( + fig, + figsize=figsize, + title=title, + showlegend=False, + yaxis=dict(type="category"), + **kwargs, + ) + if xlim is not None: + fig.update_xaxes(range=list(xlim)) + return fig + + +def spatial_overview( + *, + outlines: Sequence[np.ndarray], + points: Sequence[Tuple[str, float, float]], + tracks: Sequence[Tuple[str, np.ndarray, np.ndarray]], + title: str | None = None, + figsize: Tuple[float, float] | None = None, + **kwargs: Any, +) -> go.Figure: + """Map of observation positions on the model domain outline. + + Parameters + ---------- + outlines : Sequence of (n, 2) arrays + model domain boundary polygons + points : Sequence of (name, x, y) + point observations, labelled on the map + tracks : Sequence of (name, x, y) + track observations + title : str, optional + plot title, by default "Spatial coverage" + figsize : (float, float), optional + figure size in inches, by default None + **kwargs + keyword arguments for fig.update_layout + """ + go = import_plotly_go() + + traces = [] + for i, xy in enumerate(outlines): + traces.append( + go.Scatter( + x=xy[:, 0], + y=xy[:, 1], + mode="lines", + line=dict(color="black", width=1), + name="Domain", + showlegend=i == 0, + hoverinfo="skip", + ) + ) + for name, x, y in tracks: + traces.append( + go.Scatter(x=x, y=y, mode="markers", name=name, marker=dict(size=3)) + ) + for name, px, py in points: + traces.append( + go.Scatter( + x=[px], + y=[py], + mode="markers+text", + name=name, + text=[name], + textposition="middle right", + marker=dict(symbol="x", size=8), + ) + ) + + fig = go.Figure(traces) + apply_layout( + fig, + figsize=figsize, + title=title if title else "Spatial coverage", + yaxis=dict(scaleanchor="x", scaleratio=1), + **kwargs, + ) + return fig + + +def wind_rose( + *, + dir_centers: np.ndarray, + dir_step: float, + densities: Sequence[np.ndarray], + mag_bins: np.ndarray, + mag_max: float, + labels: Sequence[str], + colorscales: Sequence[str], + dir_labels: Sequence[str], + dir_label_positions: np.ndarray, + calm: float, + calm_text: str = "Calm", + rmax: float, + r_ticks: np.ndarray, + title: str | None = None, + figsize: Tuple[float, float] | None = None, + secondary_dir_step_factor: float = 2.0, + legend: bool = True, + **kwargs: Any, +) -> go.Figure: + """Dual wind rose as stacked polar bars, with a calm hole in the centre. + + Parameters + ---------- + dir_centers : np.ndarray + centre of each directional sector, in degrees + dir_step : float + width of a directional sector, in degrees + densities : Sequence of (n_mag, n_dir) arrays + one array per dataset, fraction of data in each magnitude/direction bin + mag_bins : np.ndarray + magnitude bin edges, used for the legend labels; the last edge is an + open-ended catch-all + mag_max : float + magnitude the colorscale is normalized to, as in the matplotlib backend + labels : Sequence[str] + dataset names + colorscales : Sequence[str] + one plotly/matplotlib colorscale name per dataset + dir_labels : Sequence[str] + compass labels for the angular axis + dir_label_positions : np.ndarray + angles of the compass labels, in degrees + calm : float + radius of the calm hole + calm_text : str, optional + label of the calm hole, by default "Calm" + rmax : float + maximum radius beyond the calm hole + r_ticks : np.ndarray + radial tick positions (fractions), excluding the calm offset + title : str, optional + plot title, by default None + figsize : (float, float), optional + figure size in inches, by default None + secondary_dir_step_factor : float, optional + the secondary dataset is drawn with sectors this much narrower, + by default 2.0 + legend : bool, optional + show the magnitude legend, by default True + **kwargs + keyword arguments for fig.update_layout + """ + go = import_plotly_go() + + traces = [] + n_mag = len(densities[0]) + # a bin is colored by its upper magnitude edge, as in the matplotlib backend + color_positions = np.asarray(mag_bins[1 : n_mag + 1], dtype=float) / mag_max + for i, density in enumerate(densities): + width = dir_step if i == 0 else dir_step / secondary_dir_step_factor + colors = _sample_colorscale(colorscales[i], color_positions) + for j in range(n_mag): + # the last bin edge is an open-ended catch-all + is_last = j == n_mag - 1 + name = ( + f">= {mag_bins[j]:.3g}" + if is_last + else f"{mag_bins[j]:.3g} - {mag_bins[j + 1]:.3g}" + ) + traces.append( + go.Barpolar( + r=density[j, :], + theta=dir_centers, + width=[width] * len(dir_centers), + name=name, + legendgroup=labels[i], + legendgrouptitle_text=labels[i] if j == 0 else None, + marker=dict(color=colors[j], line=dict(width=0)), + hovertemplate=( + f"{labels[i]}
{name}
" + "%{theta}°
%{r:.1%}" + ), + ) + ) + + fig = go.Figure(traces) + fig.update_layout( + barmode="stack", + polar=dict( + hole=calm / (calm + rmax) if (calm + rmax) > 0 else 0, + radialaxis=dict( + range=[0, rmax], + tickmode="array", + tickvals=r_ticks, + ticktext=[f"{t * 100:.0f}%" for t in r_ticks], + angle=5, + ), + angularaxis=dict( + direction="clockwise", + rotation=90, + tickmode="array", + tickvals=dir_label_positions, + ticktext=list(dir_labels), + ), + ), + ) + if calm > 0: + fig.add_annotation( + x=0.5, y=0.5, xref="paper", yref="paper", text=calm_text, showarrow=False + ) + apply_layout(fig, figsize=figsize, title=title, showlegend=legend, **kwargs) + return fig + + +def _sample_colorscale(cmap: str, values: np.ndarray) -> list[str]: + """Colors at the given positions of a matplotlib colormap, as plotly rgb strings""" + import matplotlib as mpl + + colormap = mpl.colormaps[cmap] if isinstance(cmap, str) else cmap + return [ + "rgb({:.0f},{:.0f},{:.0f})".format(*(np.array(colormap(v)[:3]) * 255)) + for v in values + ] diff --git a/src/modelskill/plotting/_scatter.py b/src/modelskill/plotting/_scatter.py index 013aa1f84..692eba67d 100644 --- a/src/modelskill/plotting/_scatter.py +++ b/src/modelskill/plotting/_scatter.py @@ -1,5 +1,5 @@ from __future__ import annotations -from typing import Literal, Sequence, Tuple, Callable, TYPE_CHECKING, Mapping +from typing import Sequence, Tuple, Callable, TYPE_CHECKING, Mapping if TYPE_CHECKING: import matplotlib.axes @@ -17,7 +17,22 @@ from modelskill.settings import options from ..metrics import _linear_regression -from ._misc import quantiles_xy, sample_points, format_skill_table, _get_fig_ax +from ._backend import ( + Backend, + PlotResult, + reject_matplotlib_axes, + validate_backend, +) +from ._plotly import scatter as _scatter_plotly +from ._misc import ( + quantiles_xy, + reglabel, + sample_points, + format_skill_table, + _get_fig_ax, + _xtick_directional, + _ytick_directional, +) def scatter( @@ -31,7 +46,7 @@ def scatter( show_hist: bool | None = None, show_density: bool | None = None, norm: colors.Normalize | None = None, - backend: Literal["matplotlib", "plotly"] = "matplotlib", + backend: Backend = "matplotlib", figsize: Tuple[float, float] = (8, 8), xlim: Tuple[float, float] | None = None, ylim: Tuple[float, float] | None = None, @@ -43,8 +58,9 @@ def scatter( skill_scores: Mapping[str, float] | None = None, skill_score_unit: str | None = "", ax: Axes | None = None, + directional: bool = False, **kwargs, -) -> Axes: +) -> PlotResult: """Scatter plot tailored for model skill comparison. Scatter plot showing compared data: observation vs modelled @@ -86,9 +102,9 @@ def scatter( colormap normalization If None, defaults to matplotlib.colors.PowerNorm(vmin=1,gamma=0.5) backend : str, optional - use "plotly" (interactive) or "matplotlib" backend, by default "matplotlib" + "matplotlib" (static) or "plotly" (interactive), by default "matplotlib" figsize : tuple, optional - width and height of the figure, by default (8, 8) + width and height of the figure in inches, by default (8, 8) xlim : tuple, optional plot range for the observation (xmin, xmax), by default None ylim : tuple, optional @@ -119,12 +135,18 @@ def scatter( unit for skill_scores, by default None ax : matplotlib.axes.Axes, optional axes to plot on, by default None + directional : bool, optional + draw the axes as a 0-360 degree compass, for directional quantities, + by default False **kwargs + other keyword arguments to plt.scatter() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib.axes.Axes - The axes on which the scatter plot was drawn. + matplotlib.axes.Axes or plotly.graph_objects.Figure + The axes the scatter plot was drawn on (matplotlib backend) or the + figure (plotly backend). Examples -------- @@ -169,6 +191,11 @@ def scatter( nbins_hist, binsize = _get_bins(bins, xymin=xymin, xymax=xymax) + if directional: + # a directional axis spans the full compass unless the user says otherwise + xlim = (0.0, 360.0) if xlim is None else xlim + ylim = (0.0, 360.0) if ylim is None else ylim + if xlim is None: xlim = (xymin - binsize, xymax + binsize) @@ -208,8 +235,8 @@ def scatter( "plotly": _scatter_plotly, } - if backend not in PLOTTING_BACKENDS: - raise ValueError(f"backend must be one of {list(PLOTTING_BACKENDS.keys())}") + validate_backend(backend) + reject_matplotlib_axes(ax, backend) if skill_table: from modelskill import from_matched @@ -224,6 +251,10 @@ def scatter( skill = cmp.skill(metrics=metrics) skill_scores = skill.to_dict("records")[0] + backend_kwargs = dict(kwargs) + if backend == "matplotlib": + backend_kwargs["ax"] = ax + return PLOTTING_BACKENDS[backend]( x=x, y=y, @@ -248,8 +279,8 @@ def scatter( skill_scores=skill_scores, skill_score_unit=skill_score_unit, fit_to_quantiles=fit_to_quantiles, - ax=ax, - **kwargs, + directional=directional, + **backend_kwargs, ) @@ -278,6 +309,7 @@ def _scatter_matplotlib( skill_scores, skill_score_unit, fit_to_quantiles, + directional, ax, cmap=None, **kwargs, @@ -343,7 +375,7 @@ def _scatter_matplotlib( x_trend, intercept + slope * x_trend, **settings.get_option("plot.scatter.reg_line.kwargs"), - label=_reglabel( + label=reglabel( slope=slope, intercept=intercept, fit_to_quantiles=fit_to_quantiles ), zorder=2, @@ -412,165 +444,11 @@ def _scatter_matplotlib( ax.set_title(title) - return ax + if directional: + _xtick_directional(ax, xlim) + _ytick_directional(ax, ylim) - -def _scatter_plotly( - *, - x, - y, - x_sample, - y_sample, - z, - xq, - yq, - x_trend, - show_density, - show_points, - norm, # TODO not used by plotly, remove or keep for consistency? - show_hist, - nbins_hist, - reg_method, - xlabel, - ylabel, - figsize, # TODO not used by plotly, remove or keep for consistency? - xlim, - ylim, - title, - skill_scores, - skill_score_unit, - fit_to_quantiles, - **kwargs, -): - import plotly.graph_objects as go - - if "ax" in kwargs: - ax = kwargs.pop("ax") - if ax is not None: - raise ValueError("Cannot pass matplotlib axes to plotly backend.") - - data = [ - go.Scatter(x=xlim, y=xlim, name="1:1", mode="lines", line=dict(color="blue")), - ] - - if reg_method: - if fit_to_quantiles: - slope, intercept = _linear_regression( - obs=xq, model=yq, reg_method=reg_method - ) - else: - slope, intercept = _linear_regression(obs=x, model=y, reg_method=reg_method) - - regression_line = go.Scatter( - x=x_trend, - y=intercept + slope * x_trend, - name=_reglabel( - slope=slope, intercept=intercept, fit_to_quantiles=fit_to_quantiles - ), - mode="lines", - line=dict(color="red"), - ) - data.append(regression_line) - - if show_hist: - data.append( - go.Histogram2d( - x=x, - y=y, - nbinsx=nbins_hist, - nbinsy=nbins_hist, - colorscale=[ - [0.0, "rgba(0,0,0,0)"], - [0.1, "purple"], - [0.5, "green"], - [1.0, "yellow"], - ], - colorbar=dict(title="# of points"), - ) - ) - - if show_points is None or show_points: - if show_density: - c = z - cbar = dict(thickness=20, title="# of points") - else: - c = "black" - cbar = None - data.append( - go.Scatter( - x=x_sample, - y=y_sample, - mode="markers", - name="Data", - marker=dict(color=c, opacity=0.5, size=3.0, colorbar=cbar), - ) - ) - if len(xq) > 0: - data.append( - go.Scatter( - x=xq, - y=yq, - name=options.plot.scatter.quantiles.label, - mode="markers", - marker_symbol="x", - marker_color=options.plot.scatter.quantiles.color, - marker_line_color="midnightblue", - marker_line_width=0.6, - ) - ) - - defaults = {"width": 600, "height": 600} - defaults = {**defaults, **kwargs} - - layout = layout = go.Layout( - legend=dict(x=0.01, y=0.99), - yaxis=dict(scaleanchor="x", scaleratio=1), - title=dict(text=title, xanchor="center", yanchor="top", x=0.5, y=0.9), - yaxis_title=ylabel, - xaxis_title=xlabel, - **defaults, - ) - - fig = go.Figure(data=data, layout=layout) - fig.update_xaxes(range=xlim, nticks=10) - fig.update_yaxes(range=ylim, nticks=10) - - if skill_scores is not None: - table = format_skill_table( - skill_scores=skill_scores, - unit=skill_score_unit, - ) - lines = [ - f"{row['name']:<6} {row['sep']} {row['value']:<6}" - for _, row in table.iterrows() - ] - - # add text box - fig.add_annotation( - x=0.99, - y=0.01, - xref="paper", - yref="paper", - text="
".join(lines), - showarrow=False, - align="left", - bordercolor="black", - borderwidth=1, - borderpad=4, - bgcolor="white", - font=dict(family="Consolas, 'Liberation Mono', monospace"), - ) - - fig.show() # Should this be here - - -def _reglabel(slope: float, intercept: float, fit_to_quantiles: bool) -> str: - sign = "" if intercept < 0 else "+" - if fit_to_quantiles: - fit = "QQ fit" - else: - fit = "Fit" - return f"{fit}: y={slope:.2f}x{sign}{intercept:.2f}" + return ax def _get_bins(bins: int | float, xymin, xymax) -> Tuple[int, float]: diff --git a/src/modelskill/plotting/_spatial_overview.py b/src/modelskill/plotting/_spatial_overview.py index 739369147..dfdbf1b1a 100644 --- a/src/modelskill/plotting/_spatial_overview.py +++ b/src/modelskill/plotting/_spatial_overview.py @@ -2,7 +2,6 @@ from typing import Iterable, Tuple, TYPE_CHECKING if TYPE_CHECKING: - import matplotlib.axes from ..model import DfsuModelResult from mikeio import GeometryFM2D @@ -10,6 +9,12 @@ from ..model.track import TrackModelResult from ..model.vertical import VerticalModelResult from ..obs import Observation, PointObservation, TrackObservation, VerticalObservation +from ._backend import ( + Backend, + PlotResult, + reject_matplotlib_axes, + validate_backend, +) from ._misc import _get_ax @@ -25,7 +30,8 @@ def spatial_overview( ax=None, figsize: Tuple | None = None, title: str | None = None, -) -> matplotlib.axes.Axes: + backend: Backend = "matplotlib", +) -> PlotResult: """Plot observation points on a map showing the model domain Parameters @@ -40,6 +46,8 @@ def spatial_overview( figure size, by default None title: str, optional plot title, default empty + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), by default "matplotlib" See Also -------- @@ -47,8 +55,8 @@ def spatial_overview( Returns ------- - matplotlib.axes.Axes - The matplotlib axes object + matplotlib.axes.Axes or plotly.graph_objects.Figure + The axes (matplotlib backend) or figure (plotly backend) Examples -------- @@ -63,43 +71,44 @@ def spatial_overview( ms.plotting.spatial_overview([o1, o2], mr) ``` """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + obs = [] if obs is None else list(obs) if isinstance(obs, Iterable) else [obs] # type: ignore mods = [] if mod is None else list(mod) if isinstance(mod, Iterable) else [mod] # type: ignore - ax = _get_ax(ax=ax, figsize=figsize) + geometries = [_model_geometry(m) for m in mods] + points, tracks = _classify_observations(obs) - # TODO: support Gridded ModelResults - for m in mods: - if isinstance(m, (PointModelResult, TrackModelResult, VerticalModelResult)): - raise ValueError( - f"Model type {type(m)} not supported. Only DfsuModelResult and mikeio.GeometryFM supported!" - ) - if hasattr(m, "data") and hasattr(m.data, "geometry"): - # mod_name = m.name # TODO: better support for multiple models - g = m.data.geometry - else: - g = m + if backend == "plotly": + from . import _plotly - # mikeio's 3D geometries (GeometryFM3D) cannot be plotted directly - if hasattr(g, "to_2d_geometry"): - g = g.to_2d_geometry() + return _plotly.spatial_overview( + outlines=[ + polygon.xy + for g in geometries + for polygon in g.boundary_polygons.exteriors + ], + points=points, + tracks=tracks, + title=title, + figsize=figsize, + ) + + ax = _get_ax(ax=ax, figsize=figsize) - # TODO this is not supported for all model types + for g in geometries: g.plot.outline(ax=ax) # type: ignore - for o in obs: - if isinstance(o, (PointObservation, VerticalObservation)): - ax.scatter(x=o.x, y=o.y, marker="x") - elif isinstance(o, TrackObservation): - if o.n_points < 10000: - ax.scatter(x=o.x, y=o.y, marker=".") - else: - print(f"{o.name}: Too many points to plot") - # TODO: group by lonlat bin or sample randomly + for _, x, y in points: + ax.scatter(x=x, y=y, marker="x") + + for name, x, y in tracks: + if len(x) < 10000: + ax.scatter(x=x, y=y, marker=".") else: - raise ValueError( - f"Could not show observation {o}. Only PointObservation and TrackObservation supported." - ) + print(f"{name}: Too many points to plot") + # TODO: group by lonlat bin or sample randomly xlim = ax.get_xlim() offset_x = 0.02 * (xlim[1] - xlim[0]) @@ -114,3 +123,43 @@ def spatial_overview( ax.set_title(title) return ax + + +def _classify_observations(obs): + """Split observations into labelled points and tracks, for either backend + + Raises + ------ + ValueError + if an observation is neither a point nor a track observation + """ + points, tracks = [], [] + for o in obs: + if isinstance(o, (PointObservation, VerticalObservation)): + points.append((o.name, o.x, o.y)) + elif isinstance(o, TrackObservation): + tracks.append((o.name, o.x, o.y)) + else: + raise ValueError( + f"Could not show observation {o}. Only PointObservation and TrackObservation supported." + ) + return points, tracks + + +def _model_geometry(m): + """The 2D flexible mesh geometry of a model result or geometry""" + # TODO: support Gridded ModelResults + if isinstance(m, (PointModelResult, TrackModelResult, VerticalModelResult)): + raise ValueError( + f"Model type {type(m)} not supported. Only DfsuModelResult and mikeio.GeometryFM supported!" + ) + if hasattr(m, "data") and hasattr(m.data, "geometry"): + # TODO: better support for multiple models + g = m.data.geometry + else: + g = m + + # mikeio's 3D geometries (GeometryFM3D) cannot be plotted directly + if hasattr(g, "to_2d_geometry"): + g = g.to_2d_geometry() + return g diff --git a/src/modelskill/plotting/_taylor_diagram.py b/src/modelskill/plotting/_taylor_diagram.py index 3ac913cc6..4a54d795a 100644 --- a/src/modelskill/plotting/_taylor_diagram.py +++ b/src/modelskill/plotting/_taylor_diagram.py @@ -1,14 +1,17 @@ from __future__ import annotations from dataclasses import dataclass import warnings -from typing import TYPE_CHECKING, Collection - -if TYPE_CHECKING: - import matplotlib.figure +from typing import Collection from matplotlib.axes import Axes import matplotlib.pyplot as plt +from ._backend import ( + Backend, + FigureResult, + reject_matplotlib_axes, + validate_backend, +) from ._taylor_diagram_external import TaylorDiagram @@ -30,7 +33,8 @@ def taylor_diagram( normalize_std: bool = False, ax: Axes | None = None, title: str = "Taylor diagram", -) -> matplotlib.figure.Figure: + backend: Backend = "matplotlib", +) -> FigureResult: """ Plot a Taylor diagram using the given observations and points. @@ -48,12 +52,31 @@ def taylor_diagram( Whether to normalize the standard deviation of the points by the standard deviation of the observations. Default is False. title : str, optional Title of the plot. Default is "Taylor diagram". + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), by default "matplotlib" Returns -------- - matplotlib.figure.Figure - The matplotlib figure object + matplotlib.figure.Figure or plotly.graph_objects.Figure + The figure object """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + + if isinstance(points, TaylorPoint): + points = [points] + + if backend == "plotly": + from . import _plotly + + return _plotly.taylor( + points=list(points), + obs_std=obs_std, + obs_text=obs_text, + normalize_std=normalize_std, + title=title, + figsize=figsize, + ) if figsize[0] != figsize[1]: warnings.warn( @@ -71,8 +94,6 @@ def taylor_diagram( contours = td.add_contours(levels=8, colors="0.5", linestyles="dotted") plt.clabel(contours, inline=1, fontsize=10, fmt="%.2f") - if isinstance(points, TaylorPoint): - points = [points] for p in points: assert isinstance(p, TaylorPoint) m = "o" if p.marker is None else p.marker diff --git a/src/modelskill/plotting/_temporal_coverage.py b/src/modelskill/plotting/_temporal_coverage.py index 5052ac8cc..458411cca 100644 --- a/src/modelskill/plotting/_temporal_coverage.py +++ b/src/modelskill/plotting/_temporal_coverage.py @@ -1,12 +1,15 @@ from __future__ import annotations -from typing import Sequence, TYPE_CHECKING - -if TYPE_CHECKING: - import matplotlib.axes +from typing import Sequence import matplotlib.pyplot as plt import numpy as np +from ._backend import ( + Backend, + PlotResult, + reject_matplotlib_axes, + validate_backend, +) from ._misc import _get_fig_ax @@ -19,7 +22,8 @@ def temporal_coverage( ax=None, figsize=None, title=None, -) -> matplotlib.axes.Axes: + backend: Backend = "matplotlib", +) -> PlotResult: """Plot graph showing temporal coverage for all observations and models Parameters @@ -32,13 +36,15 @@ def temporal_coverage( Show temporal coverage only for period covered by the model, by default True marker : str, optional - plot marker for observations, by default "_" + plot marker for observations (matplotlib backend only), by default "_" ax: matplotlib.axes, optional Adding to existing axis, instead of creating new fig figsize : Tuple(float, float), optional size of figure, by default (7, 0.45*n_lines) title: str, optional plot title, default empty + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), by default "matplotlib" See Also -------- @@ -46,8 +52,8 @@ def temporal_coverage( Returns ------- - matplotlib.axes.Axes - The matplotlib axes object + matplotlib.axes.Axes or plotly.graph_objects.Figure + The axes (matplotlib backend) or figure (plotly backend) Examples -------- @@ -71,13 +77,30 @@ def temporal_coverage( ms.plotting.temporal_coverage(mod=[mr1, mr2], figsize=(5,3)) ``` """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + obs = [] if obs is None else list(obs) if isinstance(obs, Sequence) else [obs] mod = [] if mod is None else list(mod) if isinstance(mod, Sequence) else [mod] n_lines = len(obs) + len(mod) if figsize is None: - ysize = max(2.0, 0.45 * n_lines) - figsize = (7, ysize) + figsize = (7, max(2.0, 0.45 * n_lines)) + + if backend == "plotly": + from . import _plotly + + # models first, so that the rows match the matplotlib backend + lines = [(mr.name, mr.time, True) for mr in mod] + lines += [(o.name, o.time, False) for o in obs] + xlim = ( + (mod[0].time[0], mod[0].time[-1]) + if (len(mod) > 0 and limit_to_model_period) + else None + ) + return _plotly.temporal_coverage( + lines=lines, xlim=xlim, title=title, figsize=figsize + ) fig, ax = _get_fig_ax(ax=ax, figsize=figsize) y = np.repeat(0.0, 2) diff --git a/src/modelskill/plotting/_wind_rose.py b/src/modelskill/plotting/_wind_rose.py index 76180bd88..d78b31620 100644 --- a/src/modelskill/plotting/_wind_rose.py +++ b/src/modelskill/plotting/_wind_rose.py @@ -1,9 +1,6 @@ from __future__ import annotations from dataclasses import dataclass -from typing import List, Tuple, Union, TYPE_CHECKING - -if TYPE_CHECKING: - import matplotlib.axes +from typing import List, Tuple, Union import matplotlib as mpl import matplotlib.pyplot as plt @@ -13,6 +10,13 @@ from matplotlib.legend import Legend from matplotlib.patches import Polygon, Rectangle +from ._backend import ( + Backend, + PlotResult, + reject_matplotlib_axes, + validate_backend, +) + @dataclass class DirectionalHistogram: @@ -148,7 +152,8 @@ def wind_rose( figsize: tuple[float, float] = (8, 8), ax=None, title=None, -) -> matplotlib.axes.Axes: + backend: Backend = "matplotlib", +) -> PlotResult: """Plots a (dual) wind (wave or current) roses with calms. The size of the calm is determined by the primary (measurement) data. @@ -191,11 +196,13 @@ def wind_rose( Matplotlib axis to plot on defined as polar, it can be done using "subplot_kw = dict(projection = 'polar')". Default = None, new axis created. title: str Default= None title of the plot + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), by default "matplotlib" Returns ------- - matplotlib.axes.Axes - Matplotlib axis with the plot + matplotlib.axes.Axes or plotly.graph_objects.Figure + The axes (matplotlib backend) or figure (plotly backend) Examples -------- @@ -207,6 +214,9 @@ def wind_rose( ms.plotting.wind_rose(df) ``` """ + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + if hasattr(data, "to_numpy"): data = data.to_numpy() @@ -258,6 +268,29 @@ def wind_rose( if calm_size is not None: calm = calm_size + if backend == "plotly": + from . import _plotly + + return _plotly.wind_rose( + dir_centers=dh.dir_centers, + dir_step=dir_step, + densities=[dh.density, dh2.density] if dual else [dh.density], + mag_bins=ui, + mag_max=vmax, + labels=labels if dual else labels[:1], + colorscales=[cmap1, cmap2], + dir_labels=directional_labels(n_dir_labels), + dir_label_positions=np.linspace(0, 360, n_dir_labels + 1)[:-1], + calm=calm, + calm_text=calm_text, + rmax=rmax, + r_ticks=ri, + title=title, + figsize=figsize, + secondary_dir_step_factor=secondary_dir_step_factor, + legend=legend, + ) + cmap = _get_cmap(cmap1) if ax is None: diff --git a/src/modelskill/timeseries/_plotter.py b/src/modelskill/timeseries/_plotter.py index e66374a48..84d57fa83 100644 --- a/src/modelskill/timeseries/_plotter.py +++ b/src/modelskill/timeseries/_plotter.py @@ -1,35 +1,45 @@ -from typing import Protocol +from __future__ import annotations +from typing import TYPE_CHECKING, Any, Tuple -class TimeSeriesPlotter(Protocol): - def __init__(self, ts) -> None: - pass +# modelskill.plotting depends on obs/model, which in turn depend on this +# module, so plotting is imported inside the methods to avoid a circular import +if TYPE_CHECKING: + from ..plotting._backend import Backend, PlotResult - def __call__(self): - pass - def timeseries(self): - pass +class TimeSeriesPlotter: + """Plotter for TimeSeries (observations and model results) - def hist(self): - pass + Both plots are available with the "matplotlib" (static) and the + "plotly" (interactive) backend. + Examples + -------- + >>> obs.plot.timeseries() + >>> obs.plot.hist(backend="plotly") + """ -class MatplotlibTimeSeriesPlotter(TimeSeriesPlotter): def __init__(self, ts) -> None: self._ts = ts - def __call__(self, **kwargs): + def __call__(self, **kwargs) -> PlotResult: # default to timeseries plot - self.timeseries(**kwargs) + return self.timeseries(**kwargs) def timeseries( - self, title=None, color=None, marker=".", linestyle="None", **kwargs - ): + self, + title: str | None = None, + color: str | None = None, + marker: str = ".", + linestyle: str = "None", + ax=None, + figsize: Tuple[float, float] | None = None, + backend: Backend = "matplotlib", + **kwargs: Any, + ) -> PlotResult: """Plot timeseries - Wraps pandas.DataFrame plot() method. - Parameters ---------- title : str, optional @@ -37,28 +47,69 @@ def timeseries( color : str, optional plot color, by default '#d62728' marker : str, optional - plot marker, by default '.' + plot marker (matplotlib backend only), by default '.' linestyle : str, optional - line style, by default None + line style (matplotlib backend only), by default None + ax : matplotlib.axes.Axes, optional + axes to plot on (matplotlib backend only), by default None + figsize : (float, float), optional + figure size in inches, by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to df.plot() + other keyword arguments to df.plot() (matplotlib backend) or + fig.update_layout() (plotly backend) + + Returns + ------- + matplotlib.axes.Axes or plotly.graph_objects.Figure """ - kwargs["color"] = self._ts._color if color is None else color - ax = self._ts._values_as_series.plot( - marker=marker, linestyle=linestyle, **kwargs + from ..plotting._backend import reject_matplotlib_axes, validate_backend + + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + + ts = self._ts + title = ts.name if title is None else title + color = ts._color if color is None else color + + if backend == "plotly": + from ..plotting import _plotly + + return _plotly.line( + series=ts._values_as_series, + color=color, + title=title, + ylabel=str(ts.quantity), + figsize=figsize, + **kwargs, + ) + + if figsize is not None: + kwargs["figsize"] = figsize + if ax is not None: + kwargs["ax"] = ax + + ax = ts._values_as_series.plot( + marker=marker, linestyle=linestyle, color=color, **kwargs ) - - title = self._ts.name if title is None else title ax.set_title(title) - - ax.set_ylabel(str(self._ts.quantity)) + ax.set_ylabel(str(ts.quantity)) return ax - def hist(self, bins=100, title=None, color=None, **kwargs): + def hist( + self, + bins: int = 100, + title: str | None = None, + color: str | None = None, + figsize: Tuple[float, float] | None = None, + ax=None, + backend: Backend = "matplotlib", + **kwargs: Any, + ) -> PlotResult: """Plot histogram of timeseries values - Wraps pandas.DataFrame hist() method. - Parameters ---------- bins : int, optional @@ -67,61 +118,56 @@ def hist(self, bins=100, title=None, color=None, **kwargs): plot title, default: observation name color : str, optional plot color, by default "#d62728" + figsize : (float, float), optional + figure size in inches, by default None + ax : matplotlib.axes.Axes, optional + axes to plot on (matplotlib backend only), by default None + backend : str, optional + "matplotlib" (static) or "plotly" (interactive), + by default "matplotlib" **kwargs - other keyword arguments to df.hist() + other keyword arguments to df.hist() (matplotlib backend) or + fig.update_layout() (plotly backend) Returns ------- - matplotlib axes + matplotlib.axes.Axes or plotly.graph_objects.Figure """ - title = self._ts.name if title is None else title - - kwargs["color"] = self._ts._color if color is None else color - - ax = self._ts._values_as_series.hist(bins=bins, **kwargs) + from ..plotting._backend import reject_matplotlib_axes, validate_backend + + validate_backend(backend) + reject_matplotlib_axes(ax, backend) + + ts = self._ts + title = ts.name if title is None else title + color = ts._color if color is None else color + + if backend == "plotly": + from ..plotting import _plotly + + return _plotly.histogram( + series={ts.name: ts._values_as_series.values}, + colors=[color], + bins=bins, + density=False, + alpha=1.0, + title=title, + xlabel=str(ts.quantity), + figsize=figsize, + **kwargs, + ) + + if figsize is not None: + kwargs["figsize"] = figsize + if ax is not None: + kwargs["ax"] = ax + + ax = ts._values_as_series.hist(bins=bins, color=color, **kwargs) ax.set_title(title) - ax.set_xlabel(str(self._ts.quantity)) + ax.set_xlabel(str(ts.quantity)) return ax -class PlotlyTimeSeriesPlotter(TimeSeriesPlotter): - def __init__(self, ts) -> None: - self._ts = ts - - def __call__(self): - # default to timeseries plot - self.timeseries() - - def timeseries(self): - """Plot timeseries - - Wraps plotly.express.line() function. - """ - import plotly.express as px # type: ignore - - fig = px.line( - self._ts._values_as_series, color_discrete_sequence=[self._ts._color] - ) - fig.show() - - def hist(self, bins=100, **kwargs): - """Plot histogram of timeseries values - - Wraps plotly.express.histogram() function. - - Parameters - ---------- - bins : int, optional - specification of bins, by default 100 - **kwargs - other keyword arguments to df.hist() - """ - import plotly.express as px # type: ignore - - fig = px.histogram( - self._ts._values_as_series, - nbins=bins, - color_discrete_sequence=[self._ts._color], - **kwargs, - ) - fig.show() +# kept as an alias: the plotter used to be selected by class, it is now +# selected by the `backend` argument on each plot method +MatplotlibTimeSeriesPlotter = TimeSeriesPlotter diff --git a/src/modelskill/timeseries/_timeseries.py b/src/modelskill/timeseries/_timeseries.py index bba5d78f7..c284e1b32 100644 --- a/src/modelskill/timeseries/_timeseries.py +++ b/src/modelskill/timeseries/_timeseries.py @@ -10,7 +10,7 @@ from ..types import GeometryType from ..quantity import Quantity -from ._plotter import TimeSeriesPlotter, MatplotlibTimeSeriesPlotter +from ._plotter import TimeSeriesPlotter from .. import __version__ T = TypeVar("T", bound="TimeSeries") @@ -170,7 +170,7 @@ class TimeSeries: """Time series data""" data: xr.Dataset - plotter: ClassVar = MatplotlibTimeSeriesPlotter # TODO is this the best option to choose a plotter? Can we use the settings module? + plotter: ClassVar = TimeSeriesPlotter def __init__(self, data: xr.Dataset) -> None: self.data = data if self._is_input_validated(data) else _validate_dataset(data) diff --git a/tests/plot/test_plotly_backend.py b/tests/plot/test_plotly_backend.py new file mode 100644 index 000000000..f44621e25 --- /dev/null +++ b/tests/plot/test_plotly_backend.py @@ -0,0 +1,369 @@ +"""Tests that plotly is a peer backend to matplotlib. + +Same plots, same arguments, and a figure returned rather than shown. +""" + +import sys + +import matplotlib +import matplotlib.figure +import numpy as np +import plotly.graph_objects as go +import pytest +from matplotlib.axes import Axes + +import modelskill as ms + +# every plot method that takes a `backend` argument +PLOT_KINDS = ["timeseries", "scatter", "hist", "kde", "qq", "box", "residual_hist"] +COLLECTION_PLOT_KINDS = ["scatter", "hist", "kde", "qq", "box", "residual_hist"] + + +@pytest.fixture(autouse=True) +def use_non_interactive_matplotlib(): + matplotlib.use("Agg") + + +@pytest.fixture +def o1(): + fn = "tests/testdata/SW/HKNA_Hm0.dfs0" + return ms.PointObservation(fn, item=0, x=4.2420, y=52.6887, name="HKNA") + + +@pytest.fixture +def o2(): + fn = "tests/testdata/SW/eur_Hm0.dfs0" + return ms.PointObservation(fn, item=0, x=3.2760, y=51.9990, name="EPL") + + +@pytest.fixture +def mr1(): + fn = "tests/testdata/SW/HKZN_local_2017_DutchCoast.dfsu" + return ms.model_result(fn, item=0, name="SW_1") + + +@pytest.fixture +def mr2(): + fn = "tests/testdata/SW/HKZN_local_2017_DutchCoast_v2.dfsu" + return ms.model_result(fn, item=0, name="SW_2") + + +@pytest.fixture +def cmp(o1, mr1): + return ms.match(obs=o1, mod=mr1) + + +@pytest.fixture +def cmp_two_models(o1, mr1, mr2): + return ms.match(obs=o1, mod=[mr1, mr2]) + + +@pytest.fixture +def cc(o1, o2, mr1): + return ms.match([o1, o2], mr1) + + +@pytest.fixture +def directional_cmp(): + """Comparer of a directional quantity, where axes are 0-360 degrees""" + import pandas as pd + + time = pd.date_range("2017-01-01", periods=100, freq="h") + rng = np.random.default_rng(42) + df = pd.DataFrame( + {"obs": rng.uniform(0, 360, 100), "model": rng.uniform(0, 360, 100)}, + index=time, + ) + return ms.from_matched( + df, + obs_item="obs", + mod_items=["model"], + quantity=ms.Quantity("Wave direction", "degree", is_directional=True), + ) + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_comparer_plot_returns_plotly_figure(cmp, kind): + fig = getattr(cmp.plot, kind)(backend="plotly") + + assert isinstance(fig, go.Figure) + assert len(fig.data) > 0 + + +@pytest.mark.parametrize("kind", COLLECTION_PLOT_KINDS) +def test_collection_plot_returns_plotly_figure(cc, kind): + fig = getattr(cc.plot, kind)(backend="plotly") + + assert isinstance(fig, go.Figure) + assert len(fig.data) > 0 + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_matplotlib_is_still_the_default(cmp, kind): + assert isinstance(getattr(cmp.plot, kind)(), Axes) + + +@pytest.mark.parametrize("kind", ["timeseries", "hist"]) +def test_observation_plot_returns_plotly_figure(o1, kind): + fig = getattr(o1.plot, kind)(backend="plotly") + + assert isinstance(fig, go.Figure) + assert isinstance(getattr(o1.plot, kind)(), Axes) + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_figsize_sets_the_plotly_figure_size(cmp, kind): + """figsize is given in inches for both backends""" + fig = getattr(cmp.plot, kind)(backend="plotly", figsize=(4, 3)) + + assert fig.layout.width == 400 + assert fig.layout.height == 300 + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_title_is_set_for_both_backends(cmp, kind): + fig = getattr(cmp.plot, kind)(backend="plotly", title="my title") + ax = getattr(cmp.plot, kind)(title="my title") + + assert fig.layout.title.text == "my title" + assert ax.get_title() == "my title" + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_unknown_backend_is_rejected(cmp, kind): + with pytest.raises(ValueError, match="Invalid backend 'plotLY'"): + getattr(cmp.plot, kind)(backend="plotLY") + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_matplotlib_axes_are_rejected_by_the_plotly_backend(cmp, kind): + _, ax = matplotlib.pyplot.subplots() + + with pytest.raises(ValueError, match="Cannot pass matplotlib axes"): + getattr(cmp.plot, kind)(backend="plotly", ax=ax) + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_matplotlib_only_argument_gives_an_actionable_error(cmp, kind): + with pytest.raises(ValueError, match="Invalid plotly layout argument: 'cmap'"): + getattr(cmp.plot, kind)(backend="plotly", cmap="OrRd") + + +@pytest.mark.parametrize("kind", PLOT_KINDS) +def test_plotly_layout_arguments_are_forwarded(cmp, kind): + fig = getattr(cmp.plot, kind)(backend="plotly", width=1234) + + assert fig.layout.width == 1234 + + +def test_plot_per_model_for_multiple_models(cmp_two_models): + figs = cmp_two_models.plot.hist(backend="plotly") + + assert len(figs) == 2 + assert all(isinstance(f, go.Figure) for f in figs) + + +@pytest.mark.parametrize("kind", ["timeseries", "qq", "box"]) +def test_directional_quantity_gets_a_compass_axis(directional_cmp, kind): + fig = getattr(directional_cmp.plot, kind)(backend="plotly") + + axis = fig.layout.xaxis if kind == "qq" else fig.layout.yaxis + assert axis.range == (0, 360) + assert list(axis.tickvals) == list(np.linspace(0, 360, 9)) + + +def test_directional_scatter_gets_a_compass_axis_on_both_backends(directional_cmp): + compass = list(np.linspace(0, 360, 9)) + fig = directional_cmp.plot.scatter(backend="plotly") + ax = directional_cmp.plot.scatter() + + assert list(fig.layout.xaxis.tickvals) == compass + assert list(fig.layout.yaxis.tickvals) == compass + assert fig.layout.xaxis.range == (0, 360) + assert fig.layout.yaxis.range == (0, 360) + assert list(ax.get_xticks()) == compass + assert list(ax.get_yticks()) == compass + assert ax.get_xlim() == (0.0, 360.0) + assert ax.get_ylim() == (0.0, 360.0) + + +def test_scatter_skill_table_is_shown_in_plotly(cmp): + fig = cmp.plot.scatter(backend="plotly", skill_table=True) + + assert len(fig.layout.annotations) == 1 + assert "BIAS" in fig.layout.annotations[0].text + + +def test_hist_density_switches_the_plotly_normalisation(cmp): + density = cmp.plot.hist(backend="plotly", density=True) + counts = cmp.plot.hist(backend="plotly", density=False) + + assert density.data[0].histnorm == "probability density" + assert counts.data[0].histnorm is None + assert density.layout.yaxis.title.text == "density" + assert counts.layout.yaxis.title.text == "count" + + +def test_bin_edges_are_translated_to_plotly_bins(cmp): + fig = cmp.plot.hist(bins=[0.0, 0.5, 1.0, 1.5], backend="plotly") + + assert fig.data[0].xbins.start == 0.0 + assert fig.data[0].xbins.end == 1.5 + assert fig.data[0].xbins.size == 0.5 + + +def test_non_uniform_bin_edges_are_rejected_by_the_plotly_backend(cmp): + with pytest.raises(ValueError, match="uniformly spaced bin edges"): + cmp.plot.hist(bins=[0.0, 0.5, 1.0, 3.0], backend="plotly") + + +def test_plotly_scatter_traces_cover_1to1_regression_points_and_quantiles(cmp): + fig = cmp.plot.scatter(backend="plotly") + + names = [t.name for t in fig.data] + assert "1:1" in names + assert "Data" in names + assert "Q-Q" in names + assert any(n.startswith("Fit:") for n in names) + + +@pytest.mark.parametrize("normalize_std", [False, True]) +def test_taylor_returns_a_polar_plotly_figure(cc, normalize_std): + fig = cc.plot.taylor(backend="plotly", normalize_std=normalize_std) + + assert isinstance(fig, go.Figure) + assert fig.layout.polar.sector == (0, 90) + assert all(t.type == "scatterpolar" for t in fig.data) + + +def std_mod(obs, model): + """Standard deviation of the model, as a custom metric""" + return model.std() + + +def test_taylor_places_the_models_at_arccos_of_the_correlation(cmp): + sk = cmp.skill(metrics=[ms.metrics.cc, std_mod]).to_dataframe() + fig = cmp.plot.taylor(backend="plotly") + + # a single-model Comparer labels its taylor point "model", not the model name + model = fig.data[-1] + assert model.theta[0] == pytest.approx( + np.degrees(np.arccos(sk["cc"].iloc[0])), abs=1e-6 + ) + assert model.r[0] == pytest.approx(sk["std_mod"].iloc[0], rel=1e-6) + + +def test_taylor_matplotlib_still_returns_a_matplotlib_figure(cc): + assert isinstance(cc.plot.taylor(), matplotlib.figure.Figure) + + +def test_temporal_coverage_has_one_row_per_data_source(o1, o2, mr1): + fig = ms.plotting.temporal_coverage([o1, o2], mr1, backend="plotly") + + assert isinstance(fig, go.Figure) + assert [t.name for t in fig.data] == ["SW_1", "HKNA", "EPL"] + + +def test_temporal_coverage_limits_the_time_axis_to_the_model_period(o1, mr1): + limited = ms.plotting.temporal_coverage(o1, mr1, backend="plotly") + unlimited = ms.plotting.temporal_coverage( + o1, mr1, limit_to_model_period=False, backend="plotly" + ) + + assert limited.layout.xaxis.range is not None + assert unlimited.layout.xaxis.range is None + + +def test_spatial_overview_shows_the_domain_and_the_observations(o1, o2, mr1): + fig = ms.plotting.spatial_overview([o1, o2], mr1, backend="plotly") + + assert isinstance(fig, go.Figure) + assert [t.name for t in fig.data] == ["Domain", "HKNA", "EPL"] + # equal aspect ratio, as for a map + assert fig.layout.yaxis.scaleanchor == "x" + + +def test_spatial_overview_track_observations_are_drawn_as_points(o1, mr1): + track = ms.TrackObservation( + "tests/testdata/SW/Alti_c2_Dutch.dfs0", item=3, name="c2" + ) + fig = ms.plotting.spatial_overview([o1, track], mr1, backend="plotly") + + c2 = [t for t in fig.data if t.name == "c2"][0] + assert c2.mode == "markers" + assert len(c2.x) == track.n_points + + +@pytest.mark.parametrize("backend", ["matplotlib", "plotly"]) +def test_spatial_overview_rejects_an_unsupported_observation(o1, mr1, backend): + class NotAnObservation: + name = "nope" + + with pytest.raises(ValueError, match="Could not show observation"): + ms.plotting.spatial_overview([o1, NotAnObservation()], mr1, backend=backend) + + +@pytest.fixture +def wave_dir_dataframe(): + import mikeio + + ds = mikeio.read("tests/testdata/wave_dir.dfs0") + return ds[[0, 2, 1, 3]].to_dataframe() + + +def test_wind_rose_is_a_stacked_polar_bar_chart(wave_dir_dataframe): + fig = ms.plotting.wind_rose(wave_dir_dataframe, backend="plotly") + + assert isinstance(fig, go.Figure) + assert fig.layout.barmode == "stack" + assert all(t.type == "barpolar" for t in fig.data) + # north up, clockwise, as in the matplotlib version + assert fig.layout.polar.angularaxis.direction == "clockwise" + assert fig.layout.polar.angularaxis.rotation == 90 + + +def test_wind_rose_dual_has_a_legend_group_per_dataset(wave_dir_dataframe): + dual = ms.plotting.wind_rose(wave_dir_dataframe, backend="plotly") + single = ms.plotting.wind_rose(wave_dir_dataframe.iloc[:, :2], backend="plotly") + + assert set(t.legendgroup for t in dual.data) == {"Measurement", "Model"} + assert set(t.legendgroup for t in single.data) == {"Measurement"} + assert len(dual.data) == 2 * len(single.data) + + +def test_wind_rose_has_compass_labels_like_matplotlib(wave_dir_dataframe): + fig = ms.plotting.wind_rose(wave_dir_dataframe, backend="plotly") + + labels = list(fig.layout.polar.angularaxis.ticktext) + assert labels[:5] == ["N", "NNE", "NE", "ENE", "E"] + assert list(fig.layout.polar.angularaxis.tickvals)[:3] == [0.0, 22.5, 45.0] + + +def test_wind_rose_calm_becomes_the_polar_hole(wave_dir_dataframe): + fig = ms.plotting.wind_rose(wave_dir_dataframe, backend="plotly") + + assert 0 < fig.layout.polar.hole < 1 + + +def test_wind_rose_matplotlib_is_unchanged(wave_dir_dataframe): + ax = ms.plotting.wind_rose(wave_dir_dataframe) + + assert ax.name == "polar" + + +# --- plotly layout interop --- + + +def test_explicit_width_wins_over_figsize(cmp): + fig = cmp.plot.scatter(backend="plotly", figsize=(3, 4), width=1000) + + assert fig.layout.width == 1000 + assert fig.layout.height == 400 + + +def test_missing_plotly_gives_an_actionable_error(cmp, monkeypatch): + # setting a sys.modules entry to None makes the import raise ImportError + monkeypatch.setitem(sys.modules, "plotly.graph_objects", None) + + with pytest.raises(ImportError, match=r'pip install "modelskill\[plotly\]"'): + cmp.plot.hist(backend="plotly") diff --git a/tests/test_multimodelcompare.py b/tests/test_multimodelcompare.py index 729d7af75..0873c11f7 100644 --- a/tests/test_multimodelcompare.py +++ b/tests/test_multimodelcompare.py @@ -304,11 +304,20 @@ def test_mm_scatter(cc): cc.sel(model="SW_2").plot.scatter(show_points=0.75, show_density=True) cc.sel(model="SW_2", observation="HKNA").plot.scatter(skill_table=True) cc.sel(model="SW_2").plot.scatter(fit_to_quantiles=True) - # cc.sel(model="SW_2").plot.scatter(binsize=0.5, backend="plotly") assert True plt.close("all") +def test_mm_scatter_plotly_backend(cc): + import plotly.graph_objects as go + + fig = cc.sel(model="SW_2").plot.scatter(bins=0.5, backend="plotly") + assert isinstance(fig, go.Figure) + + with pytest.raises(ValueError, match="plotLY"): + cc.sel(model="SW_2").plot.scatter(backend="plotLY") + + def cm_1(obs, model): """Custom metric #1""" return np.mean(obs / model) @@ -378,8 +387,7 @@ def test_mm_plot_timeseries(cc): cc["EPL"].plot.timeseries() cc["EPL"].plot.timeseries(title="t", figsize=(3, 3)) - # cc["EPL"].plot_timeseries(backend="plotly") - with pytest.raises(ValueError): + with pytest.raises(ValueError, match="mpl"): cc["EPL"].plot.timeseries(backend="mpl") ax = cc["EPL"].plot.timeseries() @@ -388,6 +396,16 @@ def test_mm_plot_timeseries(cc): plt.close("all") +def test_mm_plot_timeseries_plotly_backend(cc): + import plotly.graph_objects as go + + fig = cc["EPL"].plot.timeseries(backend="plotly") + assert isinstance(fig, go.Figure) + + with pytest.raises(ValueError, match="plotLY"): + cc["EPL"].plot.timeseries(backend="plotLY") + + def test_match_including_dummy(mr1, mr2, o1, o2, o3): mr3 = ms.DummyModelResult(strategy="constant", data=0.0) cc = ms.match([o1, o2, o3], [mr3, mr1, mr2])