scripts/COMPARE/3a_extrinsic_analysis_metrics.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<img src=\"../../images/BikeDNA_logo.svg\" width=\"250\"  alt=\"BikeDNA logo\" style=\"display:block; margin-left: auto; margin-right: auto;\">\n",
    "<a href=\"https://github.com/anerv/BikeDNA\">Github</a>"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3a. Extrinsic Analysis: Comparison of OSM & Reference Data\n",
    "\n",
    "This notebook compares the provided reference bicycle infrastructure data set with OSM data in the same area in a so-called extrinsic quality assessment. To run this part of the analysis, a reference data set thus must be available for comparison. \n",
    "\n",
    "This analysis is based on comparing the reference data set to OSM and highlighting how and where they differ, both in terms of *how much* bicycle infrastructure is mapped in the two data sets, and of *how* the infrastructure is mapped, pinpointing differences in network structure. \n",
    "\n",
    "**All differences are computed for the reference data in relation to OSM, taking the OSM data as the base line.** For example, the difference in network density is computed by calculating reference density minus OSM density. Hence, positive difference values (over 0) indicate how much higher the reference value is; negative difference values (below 0) indicate how much lower the reference value is. Accordingly, if differences are given in percent, the OSM value is taken to be the total value (100%).\n",
    "\n",
    "While the analysis is based on a comparison, it makes no a priori assumptions about which data set is better. The same goes for the identified differences: BikeDNA does not allow an automatic conclusion as to which data set is of better quality, but instead requires the user to interpret the meaning of the differences found, e.g., whether differing features are results of errors of omission or commission, and which data set is more correct. However, many low values can be an indication that the reference data is of lower completeness than the OSM data.\n",
    "\n",
    "The goal is that the identified differences can be used to both assess the quality of the reference and OSM data sets, and to support the decision of which data set should be used for further analysis."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Prerequisites &amp; Input/Output</b>\n",
    "<p>\n",
    "The extrinsic comparison is based almost exclusively on results from notebooks 2a and 2b and expects results and plots to be stored using the file paths used in the previous notebooks. Notebooks 1a, 2a, 1b, and 2b must be run successfully before this notebook can be used. Detailed explanations of the applied concepts and metrics can be found in the intrinsic analysis notebooks (2a and 2b).\n",
    "</p>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Familiarity required</b>\n",
    "<p>\n",
    "For a correct interpretation of some of the metrics for spatial data quality, some familiarity with the area is necessary.\n",
    "</p>\n",
    "</div>\n",
    "\n",
    "**Sections**\n",
    "* [Data completeness](#Data-completeness)\n",
    "    * [Network length](#Network-length)\n",
    "    * [Network density](#Network-density)\n",
    "* [Network topology](#Network-topology)\n",
    "    * [Simplification outcomes](#Simplification-outcomes)\n",
    "    * [Alpha, beta, and gamma indices](#Alpha,-beta,-and-gamma-indices)\n",
    "    * [Dangling nodes](#Dangling-nodes)\n",
    "    * [Under/overshoots](#Under/overshoots)\n",
    "* [Network components](#Network-components)\n",
    "    * [Disconnected components](#Disconnected-components)\n",
    "    * [Component size distribution](#Component-size-distribution)\n",
    "    * [Largest connected component](#Largest-connected-component)\n",
    "    * [Missing links](#Missing-links)\n",
    "    * [Components per grid cell](#Components-per-grid-cell)\n",
    "    * [Component connectivity](#Component-connectivity)\n",
    "* [Summary](#Summary)\n",
    "\n",
    "<br />"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Load libraries, settings and data\n",
    "\n",
    "import json\n",
    "import pickle\n",
    "import warnings\n",
    "from collections import Counter\n",
    "\n",
    "import contextily as cx\n",
    "import folium\n",
    "import geopandas as gpd\n",
    "import matplotlib as mpl\n",
    "import matplotlib.colors as colors\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import osmnx as ox\n",
    "import pandas as pd\n",
    "import yaml\n",
    "from matplotlib import cm\n",
    "from collections import Counter\n",
    "\n",
    "from bikedna import evaluation_functions as eval_func\n",
    "from bikedna import graph_functions as graph_func\n",
    "from bikedna import plotting_functions as plot_func\n",
    "\n",
    "# Read in dictionaries with settings\n",
    "%run ../settings/yaml_variables.py\n",
    "%run ../settings/plotting.py\n",
    "%run ../settings/tiledict.py\n",
    "%run ../settings/df_styler.py\n",
    "\n",
    "# Load data\n",
    "%run ../settings/load_osmdata.py\n",
    "%run ../settings/load_refdata.py\n",
    "%run ../settings/load_results.py\n",
    "\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**OSM versus reference network**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_static_maps_fp + \"area_network_osm\",\n",
    "        ref_results_static_maps_fp + \"area_network_reference\",\n",
    "    ]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data completeness\n",
    "\n",
    "This section compares the OSM and reference data sets in terms of data completeness. The goal is to identify whether one data set has more bicycle infrastructure mapped than the other, and if so, whether those differences are concentrated in some areas.\n",
    "\n",
    "The section starts with a comparison of the total length of the infrastructure in both data sets. Then, infrastructure, node and dangling node densities (i.e., the length of infrastructure/nodes per km2) is compared first at a global (study area) and at local (grid cell) level. Finally, density differences for protected and unprotected bicycle infrastructure are compared separately.\n",
    "\n",
    "Computing gridded local density differences as a measure of data quality has also been applied by e.g. [Haklay (2010)](https://journals-sagepub-com.ep.ituproxy.kb.dk/doi/10.1068/b35097).\n",
    "\n",
    "**Method**\n",
    "\n",
    "To account for differences in how bicycle infrastructure has been mapped, the computation of network length and density is based on the infrastructure length, not the geometric length of the network edges. For example, a 100 meter long **bi**directional path (geometric length: 100m) contributes with 200 meters of bicycle infrastructure (infrastructure length: 200m).\n",
    "\n",
    "**Interpretation**\n",
    "\n",
    "Density differences can point to incomplete data. For instance, if a grid cell has a significantly higher edge density in the OSM than in the reference data set, this can indicate unmapped, missing features in the reference data set, or that a street mistakenly has been tagged as bicycle infrastructure in OSM.\n",
    "\n",
    "### Network length"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.compare_print_network_length(\n",
    "    osm_edges_simplified.infrastructure_length.sum(),\n",
    "    ref_edges_simplified.infrastructure_length.sum(),\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot length comparison\n",
    "\n",
    "set_renderer(renderer_plot)\n",
    "\n",
    "bar_labels = [\"OSM\", reference_name]\n",
    "x_positions = [1, 2]\n",
    "bar_colors = [pdict[\"osm_base\"], pdict[\"ref_base\"]]\n",
    "\n",
    "# Infrastructure length density\n",
    "data = [osm_edges_simplified.infrastructure_length.sum()/1000,\n",
    "        ref_edges_simplified.infrastructure_length.sum()/1000\n",
    "    ]\n",
    "y_label = \"Network edge lengths [km]\"\n",
    "filepath = compare_results_plots_fp + \"network_length_compare\"\n",
    "title = area_name\n",
    "\n",
    "plot = plot_func.make_bar_plot(\n",
    "    data=data,\n",
    "    bar_labels=bar_labels,\n",
    "    y_label=y_label,\n",
    "    x_positions=x_positions,\n",
    "    title=title,\n",
    "    bar_colors=bar_colors,\n",
    "    filepath=filepath,\n",
    "    figsize=pdict[\"fsbar_small\"]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Network Density\n",
    "\n",
    "#### Global network densities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.print_network_densities(osm_intrinsic_results, \"OSM\")\n",
    "plot_func.print_network_densities(ref_intrinsic_results, \"reference\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Global network densities (per km2)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot global difference\n",
    "\n",
    "set_renderer(renderer_plot)\n",
    "\n",
    "# Infrastructure length density\n",
    "subplotdata = [\n",
    "    (\n",
    "        osm_intrinsic_results[\"network_density\"][\"edge_density_m_sqkm\"],\n",
    "        ref_intrinsic_results[\"network_density\"][\"edge_density_m_sqkm\"],\n",
    "    ),\n",
    "    (\n",
    "        osm_intrinsic_results[\"network_density\"][\"node_density_count_sqkm\"],\n",
    "        ref_intrinsic_results[\"network_density\"][\"node_density_count_sqkm\"],\n",
    "    ),\n",
    "    (\n",
    "        osm_intrinsic_results[\"network_density\"][\"dangling_node_density_count_sqkm\"],\n",
    "        ref_intrinsic_results[\"network_density\"][\"dangling_node_density_count_sqkm\"],\n",
    "    ),\n",
    "]\n",
    "subplotx_positions = [[1,2] for j in range(3)]\n",
    "subplotbar_labels = [\n",
    "    \"Infrastructure length (m/km2)\",\n",
    "    \"Total node density (1/km2)\",\n",
    "    \"Dangling node density (1/km2)\",\n",
    "]\n",
    "filepath = compare_results_plots_fp + \"network_densities_compare\"\n",
    "subplottitle = [\n",
    "    \"Avg infrastructure length per km2\",\n",
    "    \"Avg nodes per km2\",\n",
    "    \"Avg dangling nodes per km2\",\n",
    "]\n",
    "\n",
    "plot = plot_func.make_bar_subplots(\n",
    "    subplot_data=subplotdata,\n",
    "    nrows=1,\n",
    "    ncols=3,\n",
    "    bar_labels=[[\"OSM\", reference_name] for j in range(3)],\n",
    "    y_label=subplotbar_labels,\n",
    "    x_positions=subplotx_positions,\n",
    "    title=subplottitle,\n",
    "    bar_colors=[pdict[\"osm_base\"], pdict[\"ref_base\"]],\n",
    "    filepath=filepath,\n",
    "    wspace=0.4\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Local network densities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Plot density comparisons\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "\n",
    "# Edge density\n",
    "plot_cols = [\"osm_edge_density\", \"ref_edge_density\"]\n",
    "plot_func.plot_multiple_grid_results(\n",
    "    grid=grid,\n",
    "    plot_cols=plot_cols,\n",
    "    plot_titles=[\"OSM edge density\", f\"{reference_name} edge density\"],\n",
    "    filepath=compare_results_static_maps_fp + \"density_edge_compare\",\n",
    "    cmap=pdict[\"pos\"],\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    cx_tile=cx_tile_2,\n",
    "    no_data_cols=[\"count_osm_edges\", \"count_ref_edges\"],\n",
    "    use_norm=True,\n",
    "    norm_min=0,\n",
    "    norm_max = max(grid[plot_cols].max()),\n",
    "    figsize = pdict[\"fsmap\"],\n",
    "    legend = False,\n",
    ")\n",
    "\n",
    "# Node density\n",
    "plot_cols = [\"osm_node_density\", \"ref_node_density\"]\n",
    "plot_func.plot_multiple_grid_results(\n",
    "    grid=grid,\n",
    "    plot_cols=plot_cols,\n",
    "    plot_titles=[\"OSM node density\", f\"{reference_name} node density\"],\n",
    "    filepath=compare_results_static_maps_fp + \"density_node_compare\",\n",
    "    cmap=pdict[\"pos\"],\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    cx_tile=cx_tile_2,\n",
    "    no_data_cols=[\"count_osm_nodes\", \"count_ref_nodes\"],\n",
    "    use_norm=True,\n",
    "    norm_min=0,\n",
    "    norm_max = max(grid[plot_cols].max()),\n",
    "    figsize = pdict[\"fsmap\"],\n",
    "    legend = False,\n",
    ")\n",
    "\n",
    "# Dangling node density\n",
    "plot_cols = [\"osm_dangling_node_density\", \"ref_dangling_node_density\"]\n",
    "plot_func.plot_multiple_grid_results(\n",
    "    grid=grid,\n",
    "    plot_cols=plot_cols,\n",
    "    plot_titles=[\"OSM dangling node density\", f\"{reference_name} dangling node density\"],\n",
    "    filepath=compare_results_static_maps_fp + \"density_danglingnode_compare\",\n",
    "    cmap=pdict[\"pos\"],\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    cx_tile=cx_tile_2,\n",
    "    no_data_cols=[\"count_osm_nodes\", \"count_ref_nodes\"],\n",
    "    use_norm=True,\n",
    "    norm_min=0,\n",
    "    norm_max = max(grid[plot_cols].max()),\n",
    "    figsize = pdict[\"fsmap\"],\n",
    "    legend = False,\n",
    ")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Local differences in network densities**\n",
    "\n",
    "The densities in the OSM data are taken as base line for comparison, with absolute differences computed as `reference value` - `OSM value`. Hence, positive values indicate that the reference density of the infrastructure type is higher than the OSM density; negative values indicate that the reference density is lower than the OSM density."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "grid[\"edge_density_diff\"] = grid.ref_edge_density.fillna(\n",
    "    0\n",
    ") - grid.osm_edge_density.fillna(0)\n",
    "\n",
    "grid[\"node_density_diff\"] = grid.ref_node_density.fillna(\n",
    "    0\n",
    ") - grid.osm_node_density.fillna(0)\n",
    "\n",
    "grid[\"dangling_node_density_diff\"] = grid.ref_dangling_node_density.fillna(\n",
    "    0\n",
    ") - grid.osm_dangling_node_density.fillna(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Network density grid plots\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "\n",
    "plot_cols = [\"edge_density_diff\", \"node_density_diff\", \"dangling_node_density_diff\"]\n",
    "plot_titles = [\n",
    "    area_name + f\": {reference_name} edge density differences to OSM (m/km2)\",\n",
    "    area_name + f\": {reference_name} node density differences to OSM (m/km2)\",\n",
    "    area_name + f\": {reference_name} dangling node density differences to OSM (m/km2)\"\n",
    "]\n",
    "filepaths = [\n",
    "    compare_results_static_maps_fp + \"edge_density_compare\",\n",
    "    compare_results_static_maps_fp + \"node_density_compare\",\n",
    "    compare_results_static_maps_fp + \"dangling_node_density_compare\",\n",
    "]\n",
    "\n",
    "cmaps = [pdict[\"diff\"]] * 3\n",
    "\n",
    "# Cols for no-data plots\n",
    "no_data_cols = [\n",
    "    (\"osm_edge_density\", \"ref_edge_density\"),\n",
    "    (\"osm_node_density\", \"ref_node_density\"),\n",
    "    (\"osm_dangling_node_density\", \"ref_dangling_node_density\"),\n",
    "]\n",
    "\n",
    "cblim_edge = max(\n",
    "    abs(min(grid[\"edge_density_diff\"].fillna(value=0))),\n",
    "    max(grid[\"edge_density_diff\"].fillna(value=0)),\n",
    ")\n",
    "\n",
    "cblim_node = max(\n",
    "    abs(min(grid[\"node_density_diff\"].fillna(value=0))),\n",
    "    max(grid[\"node_density_diff\"].fillna(value=0)),\n",
    ")\n",
    "\n",
    "cblim_d_node = max(\n",
    "    abs(min(grid[\"dangling_node_density_diff\"].fillna(value=0))),\n",
    "    max(grid[\"dangling_node_density_diff\"].fillna(value=0)),\n",
    ")\n",
    "\n",
    "norm_min = [-cblim_edge, -cblim_node, -cblim_d_node]\n",
    "norm_max = [cblim_edge, cblim_node, cblim_d_node]\n",
    "\n",
    "plot_func.plot_grid_results(\n",
    "    grid=grid,\n",
    "    plot_cols=plot_cols,\n",
    "    plot_titles=plot_titles,\n",
    "    filepaths=filepaths,\n",
    "    cmaps=cmaps,\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    cx_tile=cx_tile_2,\n",
    "    no_data_cols=no_data_cols,\n",
    "    use_norm=True,\n",
    "    norm_min=norm_min,\n",
    "    norm_max=norm_max,\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Densities of protected and unprotected bicycle infrastructure\n",
    "\n",
    "**Global network densities for protected/unprotected infrastructure**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot global differences\n",
    "\n",
    "set_renderer(renderer_plot)\n",
    "\n",
    "data_labels = [\"protected_density\", \"unprotected_density\", \"mixed_density\"]\n",
    "legend_labels = [\"OSM\", reference_name]\n",
    "data_osm = [\n",
    "    osm_intrinsic_results[\"network_density\"][label + \"_m_sqkm\"] for label in data_labels\n",
    "]\n",
    "data_ref = [\n",
    "    ref_intrinsic_results[\"network_density\"][label + \"_m_sqkm\"] for label in data_labels\n",
    "]\n",
    "bar_colors = [pdict[\"osm_base\"], pdict[\"ref_base\"]]\n",
    "title = f\"{area_name}\"\n",
    "\n",
    "x_labels = [\"Protected\", \"Unprotected\", \"Mixed\"]\n",
    "x_axis = [np.arange(len(x_labels)) * 2 - 0.25, np.arange(len(x_labels)) * 2 + 0.25]\n",
    "x_ticks = np.arange(len(x_labels)) * 2\n",
    "\n",
    "y_label = \"Infrastructure densities (m/km2)\"\n",
    "\n",
    "filepath = compare_results_plots_fp + \"infrastructure_type_density_diff_compare\"\n",
    "\n",
    "fig = plot_func.make_bar_plot_side(\n",
    "    x_axis=x_axis,\n",
    "    data_osm=data_osm,\n",
    "    data_ref=data_ref,\n",
    "    bar_colors=bar_colors,\n",
    "    legend_labels=legend_labels,\n",
    "    title=title,\n",
    "    x_ticks=x_ticks,\n",
    "    x_labels=x_labels,\n",
    "    x_label=None,\n",
    "    y_label=y_label,\n",
    "    filepath=filepath,\n",
    "    figsize=pdict[\"fsbar_small\"]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Local network densities for protected/unprotected infrastructure**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if \"ref_protected_density\" in grid.columns:\n",
    "\n",
    "    plot_func.plot_saved_maps(\n",
    "        [\n",
    "            osm_results_static_maps_fp + \"density_protected_OSM\",\n",
    "            ref_results_static_maps_fp + \"density_protected_reference\",\n",
    "        ]\n",
    "    )\n",
    "\n",
    "else:\n",
    "    plot_func.plot_saved_maps(\n",
    "        [osm_results_static_maps_fp + \"density_protected_OSM\"] * 2, alpha=[1, 0]\n",
    "    )\n",
    "    print(f\"No infrastructure is mapped as protected in the {reference_name} data.\")\n",
    "\n",
    "\n",
    "if \"ref_unprotected_density\" in grid.columns:\n",
    "\n",
    "    plot_func.plot_saved_maps(\n",
    "        [\n",
    "            osm_results_static_maps_fp + \"density_unprotected_OSM\",\n",
    "            ref_results_static_maps_fp + \"density_unprotected_reference\",\n",
    "        ]\n",
    "    )\n",
    "\n",
    "else:\n",
    "    plot_func.plot_saved_maps(\n",
    "        [osm_results_static_maps_fp + \"density_unprotected_OSM\"] * 2, alpha=[1, 0]\n",
    "    )\n",
    "    print(f\"No infrastructure is mapped as unprotected in the {reference_name} data.\")\n",
    "\n",
    "if \"ref_mixed_density\" in grid.columns:\n",
    "\n",
    "    plot_func.plot_saved_maps(\n",
    "        [\n",
    "            osm_results_static_maps_fp + \"density_mixed_OSM\",\n",
    "            ref_results_static_maps_fp + \"density_mixed_reference\",\n",
    "        ]\n",
    "    )\n",
    "\n",
    "else:\n",
    "    plot_func.plot_saved_maps(\n",
    "        [osm_results_static_maps_fp + \"density_mixed_OSM\"] * 2, alpha=[1, 0]\n",
    "    )\n",
    "    print(f\"No infrastructure is mapped as mixed protected/unprotected in the {reference_name} data.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Differences in infrastructure type density**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Computing difference in infrastructure type density\n",
    "\n",
    "# In case no infrastructure with mixed protected/unprotected exist\n",
    "if \"osm_mixed_density\" not in grid.columns:\n",
    "    grid[\"osm_mixed_density\"] = 0\n",
    "\n",
    "if \"ref_mixed_density\" not in grid.columns:\n",
    "    grid[\"ref_mixed_density\"] = 0\n",
    "\n",
    "if \"ref_protected_density\" not in grid.columns:\n",
    "    grid[\"ref_protected_density\"] = 0\n",
    "\n",
    "if \"ref_unprotected_density\" not in grid.columns:\n",
    "    grid[\"ref_unprotected_density\"] = 0\n",
    "\n",
    "grid[\"protected_density_diff\"] = grid.ref_protected_density.fillna(\n",
    "    0\n",
    ") - grid.osm_protected_density.fillna(0)\n",
    "grid[\"unprotected_density_diff\"] = grid.ref_unprotected_density.fillna(\n",
    "    0\n",
    ") - grid.osm_unprotected_density.fillna(0)\n",
    "grid[\"mixed_density_diff\"] = grid.ref_mixed_density.fillna(\n",
    "    0\n",
    ") - grid.osm_mixed_density.fillna(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Infrastructure type density grid plots\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "\n",
    "plot_cols = [\"protected_density_diff\", \"unprotected_density_diff\", \"mixed_density_diff\"]\n",
    "plot_cols = [c for c in plot_cols if c in grid.columns]\n",
    "plot_titles = [\n",
    "    area_name + f\": {reference_name} protected infrastructure density differences to OSM (m/km2)\",\n",
    "    area_name + f\": {reference_name} unprotected infrastructure density differences to OSM (1/km2)\",\n",
    "    area_name + f\": {reference_name} mixed infrastructure density differences to OSM (1/km2)\",\n",
    "]\n",
    "\n",
    "filepaths = [\n",
    "    compare_results_static_maps_fp + \"protected_density_compare\",\n",
    "    compare_results_static_maps_fp + \"unprotected_density_compare\",\n",
    "    compare_results_static_maps_fp + \"mixed_density_compare\",\n",
    "]\n",
    "\n",
    "cmaps = [pdict[\"diff\"]] * len(plot_cols)\n",
    "\n",
    "# Cols for no-data plots\n",
    "no_data_cols = [\n",
    "    (\"osm_edge_density\", \"ref_edge_density\"),\n",
    "    (\"osm_edge_density\", \"ref_edge_density\"),\n",
    "    (\"osm_edge_density\", \"ref_edge_density\"),\n",
    "]\n",
    "\n",
    "# Create symmetrical color range around zero\n",
    "cblim = max(\n",
    "    abs(min(grid[plot_cols].fillna(value=0).min())),\n",
    "    max(grid[plot_cols].fillna(value=0).max()),\n",
    ")\n",
    "\n",
    "norm_min = [-cblim] * len(plot_cols)\n",
    "norm_max = [cblim] * len(plot_cols)\n",
    "\n",
    "plot_func.plot_grid_results(\n",
    "    grid=grid,\n",
    "    plot_cols=plot_cols,\n",
    "    plot_titles=plot_titles,\n",
    "    filepaths=filepaths,\n",
    "    cmaps=cmaps,\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    cx_tile=cx_tile_2,\n",
    "    no_data_cols=no_data_cols,\n",
    "    use_norm=True,\n",
    "    norm_min=norm_min,\n",
    "    norm_max=norm_max,\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Network topology\n",
    "\n",
    "After having compared data completeness, i.e. *how much* infrastructure is mapped, here we focus on differences in network *topology*, which gives information about *how* the infrastructure is mapped in both data sets. Here we also analyze the extent to which network edges are connected to one or more other edges, or if they end in a dangling node. The extent to which edges are properly connected to adjacent edges are important for, for example, analyzes of accessibility and routing. \n",
    "\n",
    "When working with data on bicycle networks, a data set without gaps between actually connected network elements is preferred - while of course reflecting the real conditions.\n",
    "Identifying the dangling nodes in a network is a quick and easy way to identify edges that end in a 'dead end'. Under- and overshoots offer a more precise picture of respectively network gaps and overextended edges, that give a misleading count of dangling nodes.\n",
    "\n",
    "**Method**\n",
    "\n",
    "To identify potential gaps or missing links in the data, first the dangling nodes in both data sets are plotted. Then, the local percentage of dangling nodes out of all nodes in each data set is plotted separately. Finally, we show the local difference in the percent of dangling nodes.\n",
    "\n",
    "Under and overshoots in both OSM and reference data are finally plotted together in an interactive plot for further inspection.\n",
    "\n",
    "**Interpretation**\n",
    "\n",
    "If an edge ends in a dangling node in one data set but not the other, this indicates a problem with the data quality. There either is a missing connection in the data, or two edges have been connected erroneously. Similarly, different local rates in the share of dangling nodes indicates differences in how the bicycle networks have been mapped - although differences in data completeness of course should be considered in the interpretation.\n",
    "\n",
    "Undershoots are clear indications of misleading gaps in network data - although they might also represent actual gaps in bicycle infrastructure. Comparing undershoots in one data set with another data set can help identify whether it is a question of data quality or the quality of the actual infrastructure. Systematic differences in the presence of undershoots or gaps across intersections might be an indication in differing digitizing strategies, since some approaches will map a bike lane crossing a street as a connected stretch, while others will introduce a gap in the width of the crossing street. While both approaches are valid, data sets created with the former method are more suited for routing-based analysis.\n",
    "\n",
    "Overshoots will often be less consequential for analysis, but a high number of overshoots will introduce false dangling nodes and distort measures for network structure based on e.g., node degree or the ratio between nodes and edges.\n",
    "\n",
    "### Simplification outcomes\n",
    "\n",
    "<br />"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\n",
    "    f\"Simplifying the OSM network decreased the number of edges by {osm_intrinsic_results['simplification_outcome']['edge_percent_diff']:.1f}%.\"\n",
    ")\n",
    "print(\n",
    "    f\"Simplifying the OSM network decreased the number of nodes by {osm_intrinsic_results['simplification_outcome']['node_percent_diff']:.1f}%.\"\n",
    ")\n",
    "print(\"\\n\")\n",
    "print(\n",
    "    f\"Simplifying the {reference_name} network decreased the number of edges by {ref_intrinsic_results['simplification_outcome']['edge_percent_diff']:.1f}%.\"\n",
    ")\n",
    "print(\n",
    "    f\"Simplifying the {reference_name} network decreased the number of nodes by {ref_intrinsic_results['simplification_outcome']['node_percent_diff']:.1f}%.\"\n",
    ")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Node degree distribution**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "degree_sequence_before_osm = sorted((d for n, d in osm_graph.degree()), reverse=True)\n",
    "degree_sequence_after_osm = sorted(\n",
    "    (d for n, d in osm_graph_simplified.degree()), reverse=True\n",
    ")\n",
    "\n",
    "degree_sequence_before_ref = sorted((d for n, d in ref_graph.degree()), reverse=True)\n",
    "degree_sequence_after_ref = sorted(\n",
    "    (d for n, d in ref_graph_simplified.degree()), reverse=True\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Display already saved degree distribution plots\n",
    "\n",
    "print(\"Note that the two figures below have different y-axis scales.\")\n",
    "\n",
    "plot_func.plot_saved_maps(\n",
    "    [osm_results_plots_fp + \"degree_dist_osm\",\n",
    "     ref_results_plots_fp + \"degree_dist_reference\"]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Alpha, beta, and gamma indices\n",
    "\n",
    "In this subsection, we compute and contrast the three aggregated network metrics alpha, beta, and gamma. These metrics are often used to describe network structure, but as measures of data quality, they are only meaningful when compared to the values of a corresponding data set. For this reason, alpha, beta, and gamma are only part of the extrinsic analysis and not included in the intrinsic notebooks.\n",
    "\n",
    "While no conclusion can be drawn about data quality based on any of the three metrics by itself, a comparison of the metrics for the two data sets can indicate differences in network topology, and hence differences in how the infrastructure has been mapped.\n",
    "\n",
    "#### Method\n",
    "All three indices are computed with `eval_func.compute_alpha_beta_gamma`.\n",
    "\n",
    "The **alpha** value is the ratio of actual to possible *cycles* in the network. A network cycle is defined as a closed loop - i.e. a path that ends on the same node that it started from. The value of alpha ranges from 0 to 1. An alpha value of 0 means that the network has no cycles at all, i.e. it is a tree. An alpha value of 1 means that the network is fully connected, which is very rarely the case.\n",
    "\n",
    "The **beta** value is the ratio of existing edges to existing nodes in the network. The value of beta ranges from 0 to N-1, where N is the number of existing nodes. A beta value of 0 means that the network has no edges; a beta value of N-1 means that the network is fully connected (see also gamma value of 1). The higher the beta value, the more different paths (on average) can be chosen between any pair of nodes.\n",
    "\n",
    "The **gamma** value is the ratio of existing to *possible* edges in the network. Any edge that connects two of the existing network nodes is defined as \"possible\". Hence, the value of gamma ranges from 0 to 1. A gamma value of 0 means that the network has no edges; a gamma value of 1 means that every node of the network is connected to every other node.\n",
    "\n",
    "For all three indices, see [Ducruet and Rodrigue, 2020](https://transportgeography.org/contents/methods/graph-theory-measures-indices/). All three indices can be interpreted in respect to network connectivity: The higher the alpha value, the more cycles are present in the network; the higher the beta value, the higher the number of paths and thus the higher the complexity of the network; and the higher the gamma value, the fewer edges lie between any pair of nodes.\n",
    "\n",
    "#### Interpretation\n",
    "These metrics do not say much about the data quality itself, nor are they useful for a topological comparison of networks of similar size. However, some conclusions can be drawn through a comparison. For example, if the indices are very similar for the two networks, despite the networks e.g. having very different geometric lengths, this suggests that the data sets have been mapped in roughly the same way, but that one simply includes more features than the other. However, if the networks have roughly the same total geometric length, but the values from alpha, beta and gamma differ, this can be an indication that the structure and topology of the two data sets are fundamentally different.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "osm_alpha, osm_beta, osm_gamma = eval_func.compute_alpha_beta_gamma(\n",
    "    edges=osm_edges_simplified,\n",
    "    nodes=osm_nodes_simplified,\n",
    "    G=osm_graph_simplified,\n",
    "    planar=True,\n",
    ")  # We assume network to be planar or approximately planar\n",
    "\n",
    "print(f\"Alpha for the simplified OSM network: {osm_alpha:.2f}\")\n",
    "print(f\"Beta for the simplified OSM network: {osm_beta:.2f}\")\n",
    "print(f\"Gamma for the simplified OSM network: {osm_gamma:.2f}\")\n",
    "\n",
    "print(\"\\n\")\n",
    "\n",
    "ref_alpha, ref_beta, ref_gamma = eval_func.compute_alpha_beta_gamma(\n",
    "    ref_edges_simplified, ref_nodes_simplified, ref_graph_simplified, planar=True\n",
    ")  # We assume network to be planar or approximately planar\n",
    "\n",
    "print(f\"Alpha for the simplified {reference_name} network: {ref_alpha:.2f}\")\n",
    "print(f\"Beta for the simplified {reference_name} network: {ref_beta:.2f}\")\n",
    "print(f\"Gamma for the simplified {reference_name} network: {ref_gamma:.2f}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot alpha, beta, gamma\n",
    "\n",
    "set_renderer(renderer_plot)\n",
    "\n",
    "bar_labels = [\"OSM\", reference_name]\n",
    "\n",
    "bar_colors = [pdict[\"osm_base\"], pdict[\"ref_base\"]]\n",
    "\n",
    "subplot_data = [\n",
    "    (len(osm_nodes_simplified), len(ref_nodes_simplified)),\n",
    "    (len(osm_edges_simplified), len(ref_edges_simplified)),\n",
    "    (\n",
    "        osm_edges_simplified.geometry.length.sum()/1000,\n",
    "        ref_edges_simplified.geometry.length.sum()/1000,\n",
    "    ),\n",
    "    (osm_alpha, ref_alpha),\n",
    "    (osm_beta, ref_beta),\n",
    "    (osm_gamma, ref_gamma),\n",
    "]\n",
    "y_label = [\n",
    "    \"\",\n",
    "    \"\",\n",
    "    \"\",\n",
    "    \"\",\n",
    "    \"\",\n",
    "    \"\",\n",
    "]\n",
    "subplottitle = [\"Total nodes\", \"Total edges\", \"Network length (km)\", \"Alpha\", \"Beta\", \"Gamma\"]\n",
    "filepath = compare_results_plots_fp + \"alpha_beta_gamma\"\n",
    "\n",
    "plot = plot_func.make_bar_subplots(\n",
    "    subplot_data=subplot_data,\n",
    "    nrows=2,\n",
    "    ncols=3,\n",
    "    bar_labels=[[\"OSM\", reference_name] for j in range(6)],\n",
    "    y_label=y_label,\n",
    "    x_positions=[[1,2] for j in range(6)],\n",
    "    title=subplottitle,\n",
    "    bar_colors=bar_colors,\n",
    "    filepath=filepath,\n",
    "    wspace=0.4\n",
    ");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Dangling nodes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "osm_dangling_nodes = gpd.read_file(osm_results_data_fp + \"dangling_nodes.gpkg\")\n",
    "ref_dangling_nodes = gpd.read_file(ref_results_data_fp + \"dangling_nodes.gpkg\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Dangling nodes in OSM & reference networks**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "interactive"
    ]
   },
   "outputs": [],
   "source": [
    "# Interactive plot of dangling nodes\n",
    "\n",
    "osm_edges_simplified_folium = plot_func.make_edgefeaturegroup(\n",
    "    gdf=osm_edges_simplified,\n",
    "    mycolor=pdict[\"osm_base\"],\n",
    "    myweight=pdict[\"line_base\"],\n",
    "    nametag=\"OSM edges\",\n",
    "    show_edges=True,\n",
    ")\n",
    "\n",
    "osm_nodes_simplified_folium = plot_func.make_nodefeaturegroup(\n",
    "    gdf=osm_nodes_simplified,\n",
    "    mysize=pdict[\"mark_base\"],\n",
    "    mycolor=pdict[\"osm_base\"],\n",
    "    nametag=\"OSM all nodes\",\n",
    "    show_nodes=True,\n",
    ")\n",
    "\n",
    "osm_dangling_nodes_folium = plot_func.make_nodefeaturegroup(\n",
    "    gdf=osm_dangling_nodes,\n",
    "    mysize=pdict[\"mark_emp\"],\n",
    "    mycolor=pdict[\"osm_contrast\"],\n",
    "    nametag=\"OSM dangling nodes\",\n",
    "    show_nodes=True,\n",
    ")\n",
    "\n",
    "ref_edges_simplified_folium = plot_func.make_edgefeaturegroup(\n",
    "    gdf=ref_edges_simplified,\n",
    "    mycolor=pdict[\"ref_base\"],\n",
    "    myweight=pdict[\"line_base\"],\n",
    "    nametag=\"Reference edges\",\n",
    "    show_edges=True,\n",
    ")\n",
    "\n",
    "ref_nodes_simplified_folium = plot_func.make_nodefeaturegroup(\n",
    "    gdf=ref_nodes_simplified,\n",
    "    mysize=pdict[\"mark_base\"],\n",
    "    mycolor=pdict[\"ref_base\"],\n",
    "    nametag=\"Reference all nodes\",\n",
    "    show_nodes=True,\n",
    ")\n",
    "\n",
    "ref_dangling_nodes_folium = plot_func.make_nodefeaturegroup(\n",
    "    gdf=ref_dangling_nodes,\n",
    "    mysize=pdict[\"mark_emp\"],\n",
    "    mycolor=pdict[\"ref_contrast2\"],\n",
    "    nametag=\"Reference dangling nodes\",\n",
    "    show_nodes=True,\n",
    ")\n",
    "\n",
    "m = plot_func.make_foliumplot(\n",
    "    feature_groups=[\n",
    "        osm_edges_simplified_folium,\n",
    "        osm_nodes_simplified_folium,\n",
    "        osm_dangling_nodes_folium,\n",
    "        ref_edges_simplified_folium,\n",
    "        ref_nodes_simplified_folium,\n",
    "        ref_dangling_nodes_folium,\n",
    "    ],\n",
    "    layers_dict=folium_layers,\n",
    "    center_gdf=osm_nodes_simplified,\n",
    "    center_crs=osm_nodes_simplified.crs,\n",
    ")\n",
    "\n",
    "bounds = plot_func.compute_folium_bounds(osm_nodes_simplified)\n",
    "m.fit_bounds(bounds)\n",
    "m.save(compare_results_inter_maps_fp + \"danglingmap_compare.html\")\n",
    "\n",
    "display(m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Interactive map saved at \" + compare_results_inter_maps_fp.lstrip(\"../\") + \"danglingmap_compare.html\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Local values for dangling nodes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute pct difference relative to OSM\n",
    "\n",
    "grid[\"dangling_nodes_diff_pct\"] = np.round(\n",
    "    100\n",
    "    * (grid.count_ref_dangling_nodes - grid.count_osm_dangling_nodes)\n",
    "    / grid.count_osm_dangling_nodes,\n",
    "    2,\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Dangling nodes as percentage of all nodes**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_static_maps_fp + \"pct_dangling_nodes_osm\",\n",
    "        ref_results_static_maps_fp + \"pct_dangling_nodes_reference\",\n",
    "    ]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Local differences in dangling nodes percentages**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "\n",
    "# norm color bar\n",
    "cbnorm_dang_diff = colors.Normalize(vmin=-100, vmax=100)  # from -max to +max\n",
    "\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "from mpl_toolkits.axes_grid1 import make_axes_locatable\n",
    "divider = make_axes_locatable(ax)\n",
    "cax = divider.append_axes(\"right\", size=\"3.5%\", pad=\"1%\")\n",
    "\n",
    "grid.plot(\n",
    "    cax=cax,\n",
    "    ax=ax,\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    column=\"dangling_nodes_diff_pct\",\n",
    "    cmap=pdict[\"diff\"],\n",
    "    legend=True,\n",
    "    norm=cbnorm_dang_diff,\n",
    ")\n",
    "\n",
    "# Add no data patches\n",
    "grid[grid[\"dangling_nodes_diff_pct\"].isnull()].plot(\n",
    "    cax=cax,\n",
    "    ax=ax,\n",
    "    facecolor=pdict[\"nodata_face\"],\n",
    "    edgecolor=pdict[\"nodata_edge\"],\n",
    "    linewidth= pdict[\"line_nodata\"],\n",
    "    hatch=pdict[\"nodata_hatch\"],\n",
    "    alpha=pdict[\"alpha_nodata\"],\n",
    ")\n",
    "\n",
    "ax.legend(handles=[nodata_patch], loc=\"upper right\")\n",
    "\n",
    "ax.set_title(\n",
    "    area_name + f\": {reference_name} percent difference to OSM in dangling nodes\"\n",
    ")\n",
    "ax.set_axis_off()\n",
    "cx.add_basemap(ax=ax, crs=study_crs, source=cx_tile_2)\n",
    "\n",
    "plot_func.save_fig(fig, compare_results_static_maps_fp + \"dangling_nodes_pct_diff_compare\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Under/overshoots"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# USER INPUT: LENGTH TOLERANCE FOR OVER- AND UNDERSHOOTS\n",
    "length_tolerance_over = 3\n",
    "length_tolerance_under = 3\n",
    "\n",
    "for s in [length_tolerance_over, length_tolerance_under]:\n",
    "    assert isinstance(s, int) or isinstance(s, float), print(\n",
    "        \"Settings must be integer or float values!\"\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "osm_overshoot_ids = pd.read_csv(\n",
    "    osm_results_data_fp + f\"overshoot_edges_{length_tolerance_over}.csv\"\n",
    ")[\"edge_id\"].to_list()\n",
    "osm_undershoot_ids = pd.read_csv(\n",
    "    osm_results_data_fp + f\"undershoot_nodes_{length_tolerance_under}.csv\"\n",
    ")[\"node_id\"].to_list()\n",
    "\n",
    "ref_overshoot_ids = pd.read_csv(\n",
    "    ref_results_data_fp + f\"overshoot_edges_{length_tolerance_over}.csv\"\n",
    ")[\"edge_id\"].to_list()\n",
    "ref_undershoot_ids = pd.read_csv(\n",
    "    ref_results_data_fp + f\"undershoot_nodes_{length_tolerance_under}.csv\"\n",
    ")[\"node_id\"].to_list()\n",
    "\n",
    "osm_overshoots = osm_edges_simplified.loc[\n",
    "    osm_edges_simplified.edge_id.isin(osm_overshoot_ids)\n",
    "]\n",
    "ref_overshoots = ref_edges_simplified.loc[\n",
    "    ref_edges_simplified.edge_id.isin(ref_overshoot_ids)\n",
    "]\n",
    "ref_undershoots = ref_nodes_simplified.loc[\n",
    "    ref_nodes_simplified.nodeID.isin(ref_undershoot_ids)\n",
    "]\n",
    "osm_undershoots = osm_nodes_simplified.loc[\n",
    "    osm_nodes_simplified.osmid.isin(osm_undershoot_ids)\n",
    "]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Over and undershoots in OSM and reference networks**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "interactive"
    ]
   },
   "outputs": [],
   "source": [
    "# Interactive plot of over/undershoots\n",
    "\n",
    "feature_groups = []\n",
    "\n",
    "# OSM feature groups\n",
    "osm_edges_simplified_folium = plot_func.make_edgefeaturegroup(\n",
    "    gdf=osm_edges_simplified,\n",
    "    mycolor=pdict[\"osm_base\"],\n",
    "    myweight=pdict[\"line_base\"],\n",
    "    nametag=\"OSM network\",\n",
    "    show_edges=True,\n",
    ")\n",
    "\n",
    "feature_groups.append(osm_edges_simplified_folium)\n",
    "\n",
    "\n",
    "if len(osm_overshoots) > 0:\n",
    "    osm_overshoots_folium = plot_func.make_edgefeaturegroup(\n",
    "        gdf=osm_overshoots,\n",
    "        mycolor=pdict[\"osm_contrast\"],\n",
    "        myweight=pdict[\"line_emp\"],\n",
    "        nametag=\"OSM overshoots\",\n",
    "        show_edges=True,\n",
    "    )\n",
    "\n",
    "    feature_groups.append(osm_overshoots_folium)\n",
    "\n",
    "if len(osm_undershoots) > 0:\n",
    "    osm_undershoot_nodes_folium = plot_func.make_nodefeaturegroup(\n",
    "        gdf=osm_undershoots,\n",
    "        mysize=pdict[\"mark_emp\"],\n",
    "        mycolor=pdict[\"osm_contrast2\"],\n",
    "        nametag=\"OSM undershoot nodes\",\n",
    "        show_nodes=True,\n",
    "    )\n",
    "\n",
    "    feature_groups.append(osm_undershoot_nodes_folium)\n",
    "\n",
    "\n",
    "# Reference feature groups\n",
    "ref_edges_simplified_folium = plot_func.make_edgefeaturegroup(\n",
    "    gdf=ref_edges_simplified,\n",
    "    mycolor=pdict[\"ref_base\"],\n",
    "    myweight=pdict[\"line_base\"],\n",
    "    nametag=f\"{reference_name} network\",\n",
    "    show_edges=True,\n",
    ")\n",
    "\n",
    "feature_groups.append(ref_edges_simplified_folium)\n",
    "\n",
    "if len(ref_overshoots) > 0:\n",
    "    ref_overshoots_folium = plot_func.make_edgefeaturegroup(\n",
    "        gdf=ref_overshoots,\n",
    "        mycolor=pdict[\"ref_contrast\"],\n",
    "        myweight=pdict[\"line_emp\"],\n",
    "        nametag=f\"{reference_name} overshoots\",\n",
    "        show_edges=True,\n",
    "    )\n",
    "\n",
    "    feature_groups.append(ref_overshoots_folium)\n",
    "\n",
    "if len(ref_undershoots) > 0:\n",
    "\n",
    "    ref_undershoot_nodes_folium = plot_func.make_nodefeaturegroup(\n",
    "        gdf=ref_undershoots,\n",
    "        mysize=pdict[\"mark_emp\"],\n",
    "        mycolor=pdict[\"ref_contrast2\"],\n",
    "        nametag=f\"{reference_name} undershoot nodes\",\n",
    "        show_nodes=True,\n",
    "    )\n",
    "\n",
    "    feature_groups.append(ref_undershoot_nodes_folium)\n",
    "\n",
    "m = plot_func.make_foliumplot(\n",
    "    feature_groups=feature_groups,\n",
    "    layers_dict=folium_layers,\n",
    "    center_gdf=osm_nodes_simplified,\n",
    "    center_crs=osm_nodes_simplified.crs,\n",
    ")\n",
    "\n",
    "bounds = plot_func.compute_folium_bounds(osm_nodes_simplified)\n",
    "m.fit_bounds(bounds)\n",
    "\n",
    "m.save(\n",
    "    compare_results_inter_maps_fp\n",
    "    + f\"overundershoots_{length_tolerance_over}_{length_tolerance_under}_compare.html\"\n",
    ")\n",
    "\n",
    "display(m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Interactive map saved at \" + compare_results_inter_maps_fp.lstrip(\"../\")\n",
    "    + f\"overundershoots_{length_tolerance_over}_{length_tolerance_under}_compare.html\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Network components\n",
    "\n",
    "This section takes a close look at the network component characteristics for the two data sets.\n",
    "\n",
    "Disconnected components do not share any elements (nodes/edges). In other words, there is no network path that could lead from one disconnected component to the other. As mentioned above, most real-world networks of bicycle infrastructure do consist of many disconnected components ([Natera Orozco et al., 2020](https://onlinelibrary.wiley.com/doi/pdfdirect/10.1111/gean.12324)). However, when two disconnected components are very close to each other, it might be a sign of a missing edge or another digitizing error.\n",
    "\n",
    "**Method**\n",
    "\n",
    "To compare the number and pattern of disconnected components in OSM and reference data, all component results from the intrinsic analyses are juxtaposed and two new plots showing respectively components gaps for OSM and reference data and the difference in component connectivity are produced.\n",
    "\n",
    "**Interpretation**\n",
    "\n",
    "The fragmented nature of many bicycle networks make it hard to assess whether disconnected components are a question of a lack of data quality of a lack of properly connected bicycle infrastructure. Comparing disconnected components in two data sets enables a more accurate assessment of whether a disconnected component is a data or a planning issue.\n",
    "\n",
    "<br />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Disconnected components"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\n",
    "    f\"The OSM network in the study area consists of {osm_intrinsic_results['component_analysis']['component_count']} disconnected components.\"\n",
    ")\n",
    "print(\n",
    "    f\"The {reference_name} network in the study area consists of {ref_intrinsic_results['component_analysis']['component_count']} disconnected components.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_static_maps_fp + \"all_components_osm\",\n",
    "        ref_results_static_maps_fp + \"all_components_reference\",\n",
    "    ]\n",
    ")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Component length distribution\n",
    "\n",
    "The distribution of all network component lengths can be visualized in a so-called *Zipf plot*, which orders the lengths of each component by rank, showing the largest component's length on the left, then the second largest component's length, etc., until the smallest component's length on the right. When a Zipf plot follows a straight line in [log-log scale](https://en.wikipedia.org/wiki/Logarithmic_scale), it means that there is a much higher chance to find small disconnected components than expected from traditional distributions [(Clauset et al., 2009)]( https://epubs.siam.org/doi/abs/10.1137/070710111). This can mean that there has been no consolidation of the network, only piece-wise or random additions [(Szell et al., 2022)](https://www.nature.com/articles/s41598-022-10783-y), or that the data itself suffers from many gaps and topology errors resulting in small disconnected components.\n",
    "\n",
    "However, it can also happen that the largest connected component (the leftmost marker in the plot at rank $10^0$) is a clear outlier, while the rest of the plot follows a different shape. This can mean that at the infrastructure level, most of the infrastructure has been connected to one large component, and that the data reflects this - i.e. the data is not suffering from gaps and missing links to a large extent. Bicycle networks might also be somewhere inbetween, with several large components as outliers.\n",
    "\n",
    "In case of a comparison over the same region, as shown below, if one data set shows a clear outlier in its largest connected component while the other does not, and if it is also at least as large, it can in general be interpreted as being more complete."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_plots_fp + \"component_length_distribution_osm\",\n",
    "        ref_results_plots_fp + \"component_length_distribution_reference\",\n",
    "    ],\n",
    "    figsize=pdict[\"fsmap\"]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Largest connected component"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Read largest cc\n",
    "osm_largest_cc = gpd.read_file(osm_results_data_fp + \"largest_connected_component.gpkg\")\n",
    "ref_largest_cc = gpd.read_file(ref_results_data_fp + \"largest_connected_component.gpkg\")\n",
    "\n",
    "print(\n",
    "    f\"The largest connected component in the OSM network contains {osm_intrinsic_results['component_analysis']['largest_cc_pct_size']:.2f}% of the network length.\"\n",
    ")\n",
    "print(\n",
    "    f\"The largest connected component in the {reference_name} network contains {ref_intrinsic_results['component_analysis']['largest_cc_pct_size']:.2f}% of the network length.\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_static_maps_fp + \"largest_conn_comp_osm\",\n",
    "        ref_results_static_maps_fp + \"largest_conn_comp_reference\",\n",
    "    ]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Overlay of largest connected component in OSM and reference networks**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "\n",
    "osm_largest_cc.plot(ax=ax, linewidth=3.5, color=pdict[\"osm_base\"], label=\"OSM\")\n",
    "ref_largest_cc.plot(ax=ax, linewidth=1.25, color=pdict[\"ref_base\"], label=reference_name)\n",
    "\n",
    "ax.set_title(f\" {area_name}: largest connected components\")\n",
    "ax.set_axis_off()\n",
    "ax.legend()\n",
    "cx.add_basemap(ax=ax, crs=study_crs, source=cx_tile_2)\n",
    "\n",
    "plot_func.save_fig(fig, compare_results_static_maps_fp + \"largest_cc_overlay_compare\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "# Plot again for potential report titlepage\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "\n",
    "osm_largest_cc.plot(ax=ax, linewidth=3, color=pdict[\"osm_base\"], label=\"OSM\")\n",
    "ref_largest_cc.plot(ax=ax, linewidth=1, color=pdict[\"ref_base\"], label=reference_name)\n",
    "ax.set_axis_off()\n",
    "\n",
    "plot_func.save_fig(fig, compare_results_static_maps_fp + \"titleimage\",plot_res=\"high\")\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Missing links\n",
    "In the plot of potential missing links between components, all edges that are within the specified distance of an edge on another component are plotted. The gaps between disconnected edges are highlighted with a marker. The map thus highlights edges which, despite being in close proximity of each other, are disconnected and where it thus would not be possible to bike on cycling infrastructure between the edges."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# DEFINE MAX BUFFER DISTANCE BETWEEN COMPONENTS CONSIDERED A GAP/MISSING LINK\n",
    "component_min_distance = 10\n",
    "\n",
    "assert isinstance(component_min_distance, int) or isinstance(\n",
    "    component_min_distance, float\n",
    "), print(\"Setting must be integer or float value!\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Read results with component gaps\n",
    "\n",
    "osm_cg_edge_ids = pd.read_csv(\n",
    "    osm_results_data_fp + f\"component_gaps_edges_{component_min_distance}.csv\"\n",
    ")[\"edge_id\"].to_list()\n",
    "osm_component_gaps_edges = osm_edges_simplified.loc[\n",
    "    osm_edges_simplified.edge_id.isin(osm_cg_edge_ids)\n",
    "]\n",
    "\n",
    "ref_cg_edge_ids = pd.read_csv(\n",
    "    ref_results_data_fp + f\"component_gaps_edges_{component_min_distance}.csv\"\n",
    ")[\"edge_id\"].to_list()\n",
    "ref_component_gaps_edges = ref_edges_simplified.loc[\n",
    "    ref_edges_simplified.edge_id.isin(ref_cg_edge_ids)\n",
    "]\n",
    "\n",
    "osm_component_gaps = gpd.read_file(\n",
    "    osm_results_data_fp + f\"component_gaps_centroids_{component_min_distance}.gpkg\"\n",
    ")\n",
    "ref_component_gaps = gpd.read_file(\n",
    "    ref_results_data_fp + f\"component_gaps_centroids_{component_min_distance}.gpkg\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "interactive"
    ]
   },
   "outputs": [],
   "source": [
    "# Interactive plot of adjacent components\n",
    "\n",
    "feature_groups = []\n",
    "\n",
    "if len(osm_component_gaps_edges) > 0:\n",
    "\n",
    "    # Feature groups for OSM\n",
    "    osm_edges_simplified_folium = plot_func.make_edgefeaturegroup(\n",
    "        gdf=osm_edges_simplified,\n",
    "        mycolor=pdict[\"osm_base\"],\n",
    "        myweight=pdict[\"line_base\"],\n",
    "        nametag=\"OSM network\",\n",
    "        show_edges=True,\n",
    "    )\n",
    "\n",
    "    osm_component_gaps_edges_folium = plot_func.make_edgefeaturegroup(\n",
    "        gdf=osm_component_gaps_edges,\n",
    "        mycolor=pdict[\"osm_emp\"],\n",
    "        myweight=pdict[\"line_emp\"],\n",
    "        nametag=\"OSM: Adjacent disconnected edges\",\n",
    "        show_edges=True,\n",
    "    )\n",
    "\n",
    "    osm_component_gaps_folium = plot_func.make_markerfeaturegroup(\n",
    "        gdf=osm_component_gaps, nametag=\"OSM: Component gaps\", show_markers=True\n",
    "    )\n",
    "\n",
    "    feature_groups.extend(\n",
    "        [\n",
    "            osm_edges_simplified_folium,\n",
    "            osm_component_gaps_edges_folium,\n",
    "            osm_component_gaps_folium,\n",
    "        ]\n",
    "    )\n",
    "\n",
    "# Feature groups for reference\n",
    "if len(ref_component_gaps_edges) > 0:\n",
    "\n",
    "    ref_edges_simplified_folium = plot_func.make_edgefeaturegroup(\n",
    "        gdf=ref_edges_simplified,\n",
    "        mycolor=pdict[\"ref_base\"],\n",
    "        myweight=pdict[\"line_base\"],\n",
    "        nametag=f\"{reference_name} network\",\n",
    "        show_edges=True,\n",
    "    )\n",
    "\n",
    "    ref_component_gaps_edges_folium = plot_func.make_edgefeaturegroup(\n",
    "        gdf=ref_component_gaps_edges,\n",
    "        mycolor=pdict[\"ref_emp\"],\n",
    "        myweight=pdict[\"line_emp\"],\n",
    "        nametag=f\"{reference_name}: Adjacent disconnected edges\",\n",
    "        show_edges=True,\n",
    "    )\n",
    "\n",
    "    ref_component_gaps_folium = plot_func.make_markerfeaturegroup(\n",
    "        gdf=ref_component_gaps, nametag=f\"{reference_name}: Component gaps\", show_markers=True\n",
    "    )\n",
    "\n",
    "    feature_groups.extend(\n",
    "        [\n",
    "            ref_edges_simplified_folium,\n",
    "            ref_component_gaps_edges_folium,\n",
    "            ref_component_gaps_folium,\n",
    "        ]\n",
    "    )\n",
    "\n",
    "m = plot_func.make_foliumplot(\n",
    "    feature_groups=feature_groups,\n",
    "    layers_dict=folium_layers,\n",
    "    center_gdf=osm_nodes_simplified,\n",
    "    center_crs=osm_nodes_simplified.crs,\n",
    ")\n",
    "\n",
    "bounds = plot_func.compute_folium_bounds(osm_nodes_simplified)\n",
    "m.fit_bounds(bounds)\n",
    "m.save(compare_results_inter_maps_fp + \"component_gaps_compare.html\")\n",
    "\n",
    "display(m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Interactive map saved at \" + compare_results_inter_maps_fp.lstrip(\"../\") + \"component_gaps_compare.html\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Components per grid cell\n",
    "\n",
    "The plots below show the number of components intersecting a grid cell. A high number of components in a grid cell is generally an indication of poor network connectivity - either due to fragmented infrastructure or because of problems with the data quality."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_static_maps_fp + \"number_of_components_in_grid_cells_osm\",\n",
    "        ref_results_static_maps_fp + \"number_of_components_in_grid_cells_reference\",\n",
    "    ]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Component connectivity\n",
    "\n",
    "Here we visualize differences between how many cells can be reached from each cell. The metric is a crude measure for network connectivity but has the benefit of being computationally cheap and thus able to quickly highlight stark differences in network connectivity.\n",
    "\n",
    "In the plot showing the difference in percent cells reached, positive values indicate a higher connectivity using the reference data set, while negative values means that more cells can be reached from a particular cell in the OSM data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_func.plot_saved_maps(\n",
    "    [\n",
    "        osm_results_static_maps_fp + \"percent_cells_reachable_grid_osm\",\n",
    "        ref_results_static_maps_fp + \"percent_cells_reachable_grid_reference\",\n",
    "    ]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute difference in cell reach percentage (where data for both OSM and REF is available)\n",
    "\n",
    "grid[\"cell_reach_pct_diff\"] = (\n",
    "    grid[\"cells_reached_ref_pct\"] - grid[\"cells_reached_osm_pct\"]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot\n",
    "\n",
    "set_renderer(renderer_map)\n",
    "\n",
    "# norm color bar\n",
    "cbnorm_diff = colors.Normalize(vmin=-100, vmax=100)\n",
    "\n",
    "fig, ax = plt.subplots(1, figsize=pdict[\"fsmap\"])\n",
    "from mpl_toolkits.axes_grid1 import make_axes_locatable\n",
    "divider = make_axes_locatable(ax)\n",
    "cax = divider.append_axes(\"right\", size=\"3.5%\", pad=\"1%\")\n",
    "\n",
    "grid.plot(\n",
    "    cax=cax,\n",
    "    ax=ax,\n",
    "    alpha=pdict[\"alpha_grid\"],\n",
    "    column=\"cell_reach_pct_diff\",\n",
    "    cmap=pdict[\"diff\"],\n",
    "    legend=True,\n",
    "    norm=cbnorm_diff,\n",
    ")\n",
    "\n",
    "# Add no data patches\n",
    "grid[grid[\"cell_reach_pct_diff\"].isnull()].plot(\n",
    "    cax=cax,\n",
    "    ax=ax,\n",
    "    facecolor=pdict[\"nodata_face\"],\n",
    "    edgecolor=pdict[\"nodata_edge\"],\n",
    "    linewidth= pdict[\"line_nodata\"],\n",
    "    hatch=pdict[\"nodata_hatch\"],\n",
    "    alpha=pdict[\"alpha_nodata\"],\n",
    ")\n",
    "\n",
    "# osm_edges_simplified.plot(ax=ax, color=pdict[\"osm_base\"], alpha=1,linewidth=2)\n",
    "# ref_edges_simplified.plot(ax=ax, color=pdict[\"ref_base\"], alpha=1,linewidth=1)\n",
    "\n",
    "ax.legend(handles=[nodata_patch], loc=\"upper right\")\n",
    "\n",
    "ax.set_title(\n",
    "    f\"{area_name}: {reference_name} difference to OSM in percent of cells reached\"\n",
    ")\n",
    "ax.set_axis_off()\n",
    "cx.add_basemap(ax=ax, crs=study_crs, source=cx_tile_2)\n",
    "\n",
    "plot_func.save_fig(fig, compare_results_static_maps_fp + \"percent_cell_reached_diff_compare\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load results from intrinsic\n",
    "osm_intrinsic_df = pd.read_csv(\n",
    "    osm_results_data_fp + \"intrinsic_summary_results.csv\",\n",
    "    index_col=0,\n",
    "    names=[\"OSM\"],\n",
    "    header=0,\n",
    ")\n",
    "\n",
    "ref_intrinsic_df = pd.read_csv(\n",
    "    ref_results_data_fp + \"intrinsic_summary_results.csv\",\n",
    "    index_col=0,\n",
    "    names=[reference_name],\n",
    "    header=0,\n",
    ")\n",
    "\n",
    "# Drop rows from OSM results not available for reference\n",
    "osm_intrinsic_df.drop(\n",
    "    [\"Incompatible tag combinations\", \"Missing intersection nodes\"],\n",
    "    axis=0,\n",
    "    inplace=True,\n",
    ")\n",
    "\n",
    "# Save new results\n",
    "osm_intrinsic_df.at[\"Alpha\", \"OSM\"] = osm_alpha\n",
    "osm_intrinsic_df.at[\"Beta\", \"OSM\"] = osm_beta\n",
    "osm_intrinsic_df.at[\"Gamma\", \"OSM\"] = osm_gamma\n",
    "\n",
    "ref_intrinsic_df.at[\"Alpha\", reference_name] = ref_alpha\n",
    "ref_intrinsic_df.at[\"Beta\", reference_name] = ref_beta\n",
    "ref_intrinsic_df.at[\"Gamma\", reference_name] = ref_gamma\n",
    "\n",
    "# Combine\n",
    "extrinsic_df = osm_intrinsic_df.join(ref_intrinsic_df)\n",
    "assert len(extrinsic_df) == len(osm_intrinsic_df) == len(ref_intrinsic_df)\n",
    "extrinsic_df.style.pipe(format_extrinsic_style)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "## Save results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "extrinsic_df.to_csv(\n",
    "    compare_results_data_fp + \"extrinsic_summary_results.csv\", index=True\n",
    ")\n",
    "\n",
    "with open(\n",
    "    compare_results_data_fp + f\"grid_results_extrinsic.pickle\", \"wb\"\n",
    ") as f:\n",
    "    pickle.dump(grid, f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "source": [
    "***"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "noex"
    ]
   },
   "outputs": [],
   "source": [
    "from time import strftime\n",
    "print(\"Time of analysis: \" + strftime(\"%a, %d %b %Y %H:%M:%S\"))"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "kernelspec": {
   "display_name": "bikedna",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.3"
  },
  "vscode": {
   "interpreter": {
    "hash": "25d89d047ff73a4f08e1f58df0505313a978ba3f16552d4fb4e98a48d36b765b"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}