Maps (Interactive)¶

Next to static map plotting, PyPSA allows for exploring networks on a map, interactively. With n.explore(), you can explore the location of all components, including buses, lines, links, transformers, their component attributes and map results or other properties to the bus sizes, branch widths, colors, etc. Calling the method returns a standard pydeck.Deck object than can be layered on top of other pydeck.Deck objects (see https://deckgl.readthedocs.io/en/latest/layer.html). They can also be exported in self-contained HTML files for sharing. In this notebook, we demonstrate the features of n.explore() using the SciGRID example.

Input data¶

In [2]:

import geopandas as gpd

import pypsa

n = pypsa.examples.scigrid_de()

import geopandas as gpd import pypsa n = pypsa.examples.scigrid_de()

INFO:pypsa.network.io:Imported network 'SciGrid-DE' has buses, carriers, generators, lines, loads, storage_units, transformers

Preparation¶

For illustrative purposes, we cluster the network based on federal states. For a more detailed guide on clustering, please go to Network Clustering. For the scope of this guide, you can ignore the following blocks.

In [3]:

n.calculate_dependent_values()
n.lines = n.lines.reindex(columns=n.components["Line"]["defaults"].index[1:])
n.lines["type"] = "Al/St 240/40 2-bundle 220.0"
n.buses = n.buses.reindex(columns=n.components["Bus"]["defaults"].index[1:])
n.buses["frequency"] = 50

url = "https://media.githubusercontent.com/media/wmgeolab/geoBoundaries/9469f09592ced973a3448cf66b6100b741b64c0d/releaseData/gbOpen/DEU/ADM1/geoBoundaries-DEU-ADM1-all.zip"
states = gpd.read_file(url, layer="geoBoundaries-DEU-ADM1_simplified")
states["shapeName"] = states["shapeName"].apply(
    lambda x: x.encode("latin1").decode("utf-8")
)  # fix encoding issue

bus_coords = gpd.GeoDataFrame(
    geometry=gpd.points_from_xy(n.buses.x, n.buses.y, crs=4326), index=n.buses.index
)
busmap = bus_coords.to_crs(3035).sjoin_nearest(states.to_crs(3035), how="left").shapeISO
nc = n.cluster.spatial.cluster_by_busmap(busmap)

n.calculate_dependent_values() n.lines = n.lines.reindex(columns=n.components["Line"]["defaults"].index[1:]) n.lines["type"] = "Al/St 240/40 2-bundle 220.0" n.buses = n.buses.reindex(columns=n.components["Bus"]["defaults"].index[1:]) n.buses["frequency"] = 50 url = "https://media.githubusercontent.com/media/wmgeolab/geoBoundaries/9469f09592ced973a3448cf66b6100b741b64c0d/releaseData/gbOpen/DEU/ADM1/geoBoundaries-DEU-ADM1-all.zip" states = gpd.read_file(url, layer="geoBoundaries-DEU-ADM1_simplified") states["shapeName"] = states["shapeName"].apply( lambda x: x.encode("latin1").decode("utf-8") ) # fix encoding issue bus_coords = gpd.GeoDataFrame( geometry=gpd.points_from_xy(n.buses.x, n.buses.y, crs=4326), index=n.buses.index ) busmap = bus_coords.to_crs(3035).sjoin_nearest(states.to_crs(3035), how="left").shapeISO nc = n.cluster.spatial.cluster_by_busmap(busmap)

Let’s take an initial look at the network. By default, n.explore() displays information in the tooltip for each component type, including component names and their default sizes or widths. You can hover over the elements in the map below to inspect the data interactively. Optionally, you can disable tooltip by setting tooltip=False - this can help reduce processing time and decrease the file size when exporting to HTML.

In [4]:

nc.explore()

nc.explore()

Out[4]:

With help(n.explore), we can check what parameters the method accepts:

Docstring:
Create an interactive map of the PyPSA network using Pydeck.

Parameters
----------
branch_width_factor : float, default 1.0
    Branch widths are scaled by this factor.
bus_size : float/dict/pandas.Series
    Sizes of bus points in km² (corresponds to circle area), defaults to 25 km².
bus_size_factor : float, default 1.0
    Bus sizes are scaled by this factor.
bus_split_circle : bool, default False
    Draw half circles if bus_size is a pandas.Series with a Multiindex.
    If set to true, the upper half circle per bus then includes all positive values
    of the series, the lower half circle all negative values. Defaults to False.
bus_color : str/dict/pandas.Series/None
    Colors for the buses, defaults to "cadetblue". If bus_size is a
    pandas.Series with a Multiindex, bus_color defaults to the
    n.c.carriers.static['color'] column.
bus_cmap : mcolors.Colormap/str, default 'Reds'
    If bus_color are floats, this color map will assign the colors.
bus_cmap_norm : mcolors.Normalize/None
    Normalization for bus_cmap, defaults to None.
bus_alpha : float/dict/pandas.Series
    Add alpha channel to buses, defaults to 0.9.
line_flow : float/dict/pandas.Series, default 0
    Series of line flows indexed by line names, defaults to 0. If 0, no arrows will be created.
    If a float is provided, it will be used as a constant flow for all lines.
line_color : str/dict/pandas.Series
    Colors for the lines, defaults to 'rosybrown'.
line_alpha : float/dict/pandas.Series
    Add alpha channel to lines, defaults to 0.9.
line_width : float/dict/pandas.Series, default 2
    Widths of line component in km.
link_flow : float/dict/pandas.Series, default 0
    Series of link flows indexed by link names, defaults to 0. If 0, no arrows will be created.
    If a float is provided, it will be used as a constant flow for all links.
link_color : str/dict/pandas.Series
    Colors for the links, defaults to 'darkseagreen'.
link_alpha : float/dict/pandas.Series
    Add alpha channel to links, defaults to 0.9.
link_width : float/dict/pandas.Series, default 2
    Widths of link component in km.
tooltip : bool, default True
    Whether to add a tooltip to the bus layer.

Other Parameters
----------------
branch_components : list, set, optional, default ['Line', 'Link', 'Transformer']
    Branch components to be plotted.
branch_width_max : float, default 10
    Maximum width of branch component in km when `auto_scale` is True.
bus_size_max : float, default 10000
    Maximum area size of bus component in km² when `auto_scale` is True.
line_cmap : mcolors.Colormap/str, default 'viridis'
    If line_color are floats, this color map will assign the colors.
line_cmap_norm : mcolors.Normalize
    The norm applied to the line_cmap.
link_cmap : mcolors.Colormap/str, default 'viridis'
    If link_color are floats, this color map will assign the colors.
link_cmap_norm : mcolors.Normalize|matplotlib.colors.*Norm
    The norm applied to the link_cmap.
transformer_flow : float/dict/pandas.Series, default 0
    Series of transformer flows indexed by transformer names, defaults to 0. If 0, no arrows will be created.
    If a float is provided, it will be used as a constant flow for all transformers.
transformer_color : str/dict/pandas.Series
    Colors for the transformers, defaults to 'orange'.
transformer_cmap : mcolors.Colormap/str, default 'viridis'
    If transformer_color are floats, this color map will assign the colors.
transformer_cmap_norm : matplotlib.colors.Normalize|matplotlib.colors.*Norm
    The norm applied to the transformer_cmap.
transformer_alpha : float/dict/pandas.Series
    Add alpha channel to transformers, defaults to 0.9.
transformer_width : float/dict/pandas.Series, default 2
    Widths of transformer in km.
arrow_size_factor : float, default 1.5
    Multiplier on branch flows to scale the arrow size.
arrow_color : str/dict/pandas.Series | None, default None
    Colors for the arrows. If not specified, defaults to the same colors as the respective branch component.
arrow_cmap : str/matplotlib.colors.Colormap, default 'viridis'
    Colormap to use if arrow_color is a numeric pandas.Series.
arrow_cmap_norm : matplotlib.colors.Normalize, optional
    Normalization to use if arrow_color is a numeric pandas.Series.
arrow_alpha : float/dict/pandas.Series, default 0.9
    Add alpha channel to arrows, defaults to 0.9.
bus_columns : list, default None
    List of bus columns to include.
    Specify additional columns to include in the tooltip.
line_columns : list, default None
    List of line columns to include. If None, only the bus0 and bus1 columns are used.
    Specify additional columns to include in the tooltip.
link_columns : list, default None
    List of link columns to include. If None, only the bus0 and bus1 columns are used.
    Specify additional columns to include in the tooltip.
transformer_columns : list, default None
    List of transformer columns to include. If None, only the bus0 and bus1 columns are used.
    Specify additional columns to include in the tooltip.
geomap : bool, default False
    Whether to add a geomap layer to the plot.
geomap_alpha : float, default 0.9
    Alpha transparency for the geomap features.
geomap_color : dict | None, default None
    Dictionary specifying colors for different geomap features. If None, default colors will be used: `{'land': 'whitesmoke', 'ocean': 'lightblue'}
geomap_resolution : {'110m', '50m', '10m'}, default '50m'
    Resolution of the geomap features. One of '110m', '50m', or '10m'.
geometry : bool, default False
    Whether to use the geometry column of the branch components.

Returns
-------
PydeckPlotter
    The PydeckPlotter instance with the created layers.

In the docstrings, you can see that the method allows for passing the same arguments as its static counterpart n.plot() and even a few more. Please be aware that due to how pydeck and matplotlib.pyplot handle numeric values, passing a value for e.g. line_width will achieve a different scalings in each method. In n.explore(), all widths and flows passed are translated into kilometers on the interactive map. Values passed for bus_size translate into km², accordingly. By default, all branches are rendered at a width of 2 km, buses at a size of 25 km².

Retrieving Results Data¶

To map result to parameters of the interactive map, we first solve the network and then use n.statistics() to calculate relevant metrics.

In [5]:

# We reduce logging output for clarity
import logging

logging.getLogger("pypsa").setLevel(logging.ERROR)
logging.getLogger("linopy").setLevel(logging.ERROR)

nc.optimize()

# We reduce logging output for clarity import logging logging.getLogger("pypsa").setLevel(logging.ERROR) logging.getLogger("linopy").setLevel(logging.ERROR) nc.optimize()

/tmp/ipykernel_9166/319203369.py:7: FutureWarning: The default value of `include_objective_constant` will change from True to False in version 2.0. Set `include_objective_constant` explicitly to suppress this warning. Using False improves LP numerical conditioning by not including the objective constant as a variable.
  nc.optimize()

Writing constraints.:   0%|          | 0/16 [00:00<?, ?it/s]

Writing constraints.: 100%|██████████| 16/16 [00:00<00:00, 171.46it/s]

Writing continuous variables.:   0%|          | 0/5 [00:00<?, ?it/s]

Writing continuous variables.: 100%|██████████| 5/5 [00:00<00:00, 299.32it/s]

Out[5]:

('ok', 'optimal')

From above we learned that bus_size accepts parameters of type float, dict, and pd.Series. When passing a multi-index pd.Series, its values will be mapped to pie chart slices.

In [6]:

eb = (
    nc.statistics.energy_balance(
        groupby=["bus", "carrier"],
        components=["Generator", "Load", "StorageUnit"],
    )
    .groupby(["bus", "carrier"])
    .sum()
)

eb = ( nc.statistics.energy_balance( groupby=["bus", "carrier"], components=["Generator", "Load", "StorageUnit"], ) .groupby(["bus", "carrier"]) .sum() )

We also extract branch results, e.g., line and link flows in this example.

In [7]:

line_flow = nc.lines_t.p0.sum(axis=0)
link_flow = nc.links_t.p0.sum(axis=0)

line_flow = nc.lines_t.p0.sum(axis=0) link_flow = nc.links_t.p0.sum(axis=0)

Note that for the pie slices to be plotted and colored correctly, passing a multi-index pd.Series requires all carrier colors to exist. Colors can be specified by their hex code representation or from the list of matplotlib names. In n.statistics.energy_balance() load is also included, so we also need to include a color for the load carrier.

In [8]:

colors = {
    "Multiple": "pink",
    "AC": "black",
    "Brown Coal": "saddlebrown",
    "Gas": "darkorange",
    "Geothermal": "firebrick",
    "Hard Coal": "darkslategray",
    "Nuclear": "mediumorchid",
    "Oil": "peru",
    "Other": "dimgray",
    "Pumped Hydro": "cornflowerblue",
    "Run of River": "royalblue",
    "Solar": "gold",
    "Storage Hydro": "navy",
    "Waste": "olive",
    "Wind Offshore": "teal",
    "Wind Onshore": "turquoise",
}

nc.carriers.color = nc.carriers.index.map(colors)

colors = { "Multiple": "pink", "AC": "black", "Brown Coal": "saddlebrown", "Gas": "darkorange", "Geothermal": "firebrick", "Hard Coal": "darkslategray", "Nuclear": "mediumorchid", "Oil": "peru", "Other": "dimgray", "Pumped Hydro": "cornflowerblue", "Run of River": "royalblue", "Solar": "gold", "Storage Hydro": "navy", "Waste": "olive", "Wind Offshore": "teal", "Wind Onshore": "turquoise", } nc.carriers.color = nc.carriers.index.map(colors)

As the carriers for loads are missing, we need to add them, manually.

In [9]:

nc.carriers.loc["", "color"] = "darkred"
nc.carriers.loc["-", "color"] = "darkred"

nc.carriers.loc["", "color"] = "darkred" nc.carriers.loc["-", "color"] = "darkred"

Applying n.explore()¶

Finally, let's pass the results to n.explore(). By default the map_style='road' is used, we pass dark for illustrative purposes. Setting bus_split_circle=True maps negative values to the bottom half and positive values to the positive half. If set to False, bottom half circles are not used and negative values will automatically be omitted. As the values scale proportionally to the bus area, they are directly translated into km². This may not achieve the outcome we want, this is why we set auto_scale=True. This scales the maximum value to bus_size_max. The same applies to branch_width and branch_flow (branch_width_max). By defaults, arrows are scaled by 1.5, so that arrowheads are visible. If you want those to be less or more prominent, scale accordingly. We can pass additional columns for each component type that we want to include in the tooltip. Note that we disable the tooltip in the documentation due to file size limitations of our documentation.

Optionally, you can pass a pdk.ViewState object or a dict. By default, PyPSA will set zoom level of 4 and calculate the initial view based on all coordinates in n.buses. For details on how to use view_state, we refer to the pydeck documentation.

In [10]:

view_state = {}
view_state["zoom"] = 6
view_state["pitch"] = 35  # Up/down angle relative to the map's plane

map = nc.explore(
    view_state=view_state,
    map_style="dark",
    bus_size=eb,  # MWh -> km²
    bus_split_circle=True,
    bus_size_max=7000,  # km²
    line_color="yellow",
    line_width=line_flow,  # MWh -> km
    link_width=link_flow,  # MWh -> km
    line_flow=line_flow,  # MWh -> km
    link_flow=link_flow,  # MWh -> km
    branch_width_max=16,  # km
    auto_scale=True,
    bus_columns=["v_nom"],
    line_columns=["s_nom"],
    link_columns=["p_nom"],
    arrow_size_factor=2,
    tooltip=True,  # disabled here for technical limits of mkdocs-jupyter plugin
)

view_state = {} view_state["zoom"] = 6 view_state["pitch"] = 35 # Up/down angle relative to the map's plane map = nc.explore( view_state=view_state, map_style="dark", bus_size=eb, # MWh -> km² bus_split_circle=True, bus_size_max=7000, # km² line_color="yellow", line_width=line_flow, # MWh -> km link_width=link_flow, # MWh -> km line_flow=line_flow, # MWh -> km link_flow=link_flow, # MWh -> km branch_width_max=16, # km auto_scale=True, bus_columns=["v_nom"], line_columns=["s_nom"], link_columns=["p_nom"], arrow_size_factor=2, tooltip=True, # disabled here for technical limits of mkdocs-jupyter plugin )

Map Export¶

To export the interactive map, we use pydeck's built-in features (pdk.to_html()). Passing offline=True embeds deck.gls JavaScript library. Note that you require an internet connection if you want to be able to view the map tiles in the background. If you want country shapes also be included in the self-contained HTML, set geomap=True. This will however increase the file size noticeably.

In [11]:

map.to_html("exploring.html")

map.to_html("exploring.html")

Out[11]:

Static Equivalent¶

This would be its static equivalent. As mentioned above, scaling is handled differently in pydeck and matplotlib.

In [12]:

nc.plot(
    bus_size=eb / 3e6,
    bus_split_circle=True,
    line_width=line_flow / 1e4,
    link_width=link_flow / 1e4,
    line_flow=line_flow / 5e4,
)

nc.plot( bus_size=eb / 3e6, bus_split_circle=True, line_width=line_flow / 1e4, link_width=link_flow / 1e4, line_flow=line_flow / 5e4, )

Out[12]:

{'nodes': {'Bus': <matplotlib.collections.PatchCollection at 0x7e8550f06510>},
 'branches': {'Line': <matplotlib.collections.LineCollection at 0x7e8550f06660>},
 'flows': {'Line': <matplotlib.collections.PatchCollection at 0x7e8550fa3ed0>}}

Stacking Pydeck Layers¶

We can use pydecks built-in layer functionalities to extend the interactive map with features that are completely unrelated to pypsa. For example, we can use the GeoDataFrame containing Polygon and MultiPolygon geometries to color them by average prices (€/MWh).

In [13]:

avg_prices = nc.statistics.prices()
avg_prices.head()

avg_prices = nc.statistics.prices() avg_prices.head()

Out[13]:

name
DE-BB     9.10477
DE-BE     9.21549
DE-BW    10.49886
DE-BY    10.41145
DE-HB     7.80802
Name: objective, dtype: float64

In [14]:

# Map average prices by shapeISO
states["avg_price"] = states["shapeISO"].map(avg_prices).round(2)

# Map average prices by shapeISO states["avg_price"] = states["shapeISO"].map(avg_prices).round(2)

Now we can map colors to the shapes using a colormap and store them as in RGBA formatted lists in the states GeoDataFrame.

In [15]:

import matplotlib.colors as mcolors
import matplotlib.pyplot as plt

values = states["avg_price"]
cmap = plt.get_cmap("Reds")
norm = mcolors.Normalize(vmin=values.min(), vmax=values.max())

import matplotlib.colors as mcolors import matplotlib.pyplot as plt values = states["avg_price"] cmap = plt.get_cmap("Reds") norm = mcolors.Normalize(vmin=values.min(), vmax=values.max())

We need a small helper to convert mcolors to pydeck ready RGBA lists

In [16]:

def price_to_color(price, alpha=0.7):
    color = cmap(norm(price))  # RGBA in 0-1
    rgb = [round(c * 255) for c in color[:3]]  # only RGB
    a = round(alpha * 255)
    return rgb + [a]

def price_to_color(price, alpha=0.7): color = cmap(norm(price)) # RGBA in 0-1 rgb = [round(c * 255) for c in color[:3]] # only RGB a = round(alpha * 255) return rgb + [a]

In [17]:

states["color"] = states["avg_price"].apply(price_to_color)
states.head()

states["color"] = states["avg_price"].apply(price_to_color) states.head()

Out[17]:

	shapeName	shapeISO	shapeID	shapeGroup	shapeType	geometry	avg_price	color
0	Baden-Württemberg	DE-BW	10402087B60055985875400	DEU	ADM1	MULTIPOLYGON (((9.12593 47.66864, 9.12068 47.6...	10.50	[148, 11, 19, 178]
1	Bayern	DE-BY	10402087B60477050509260	DEU	ADM1	POLYGON ((9.60208 47.58434, 9.60589 47.5857, 9...	10.41	[163, 15, 21, 178]
2	Berlin	DE-BE	10402087B20892132820961	DEU	ADM1	POLYGON ((13.48006 52.67465, 13.47601 52.67039...	9.22	[250, 104, 73, 178]
3	Brandenburg	DE-BB	10402087B40185768535592	DEU	ADM1	MULTIPOLYGON (((13.05103 51.64768, 13.15453 51...	9.10	[251, 115, 83, 178]
4	Bremen	DE-HB	10402087B44391416804171	DEU	ADM1	MULTIPOLYGON (((8.6164 53.19703, 8.61522 53.19...	7.81	[254, 234, 224, 178]

We now create a new layer based on our states GeoDataFrame and insert it to map.layers. If we append it, the states layer would lie above the pie charts, which would obtrude our previous map.

In [18]:

import pydeck as pdk

# Add a custom tooltip column (HTML or plain text)
states["tooltip_html"] = (
    "<b>State:</b> "
    + states["shapeName"]
    + "<br>"
    + "<b>GID:</b> "
    + states["shapeISO"]
    + "<br><b>Avg. Price:</b> "
    + states["avg_price"].astype(str)
    + " €/MWh"
)

# Create layer
states_layer = pdk.Layer(
    "GeoJsonLayer",
    states,
    stroked=True,
    filled=True,
    get_fill_color="color",
    get_line_color=[255, 255, 255, 255],
    line_width_min_pixels=1,
    pickable=True,
    auto_highlight=True,
)

import pydeck as pdk # Add a custom tooltip column (HTML or plain text) states["tooltip_html"] = ( "<b>State:</b> " + states["shapeName"] + "<br>" + "<b>GID:</b> " + states["shapeISO"] + "<br><b>Avg. Price:</b> " + states["avg_price"].astype(str) + " €/MWh" ) # Create layer states_layer = pdk.Layer( "GeoJsonLayer", states, stroked=True, filled=True, get_fill_color="color", get_line_color=[255, 255, 255, 255], line_width_min_pixels=1, pickable=True, auto_highlight=True, )

In [19]:

map.layers.insert(0, states_layer)

map.layers.insert(0, states_layer)

In [20]:

map.show()

map.show()

Out[20]: