added files from exp_hdbscan_paper

examples/research/robustness-of-algorithms/delta_roam_vs_acc_exp.ipynb

@@ -0,0 +1,222 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "fc8fe351",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import pandas as pd\n",
+    "import geopandas as gpd\n",
+    "import shapely\n",
+    "import numpy as np\n",
+    "import json\n",
+    "from functools import partial\n",
+    "from tqdm import tqdm\n",
+    "\n",
+    "import matplotlib.pyplot as plt\n",
+    "import matplotlib.colors as mcolors\n",
+    "import matplotlib.patches as mpatches\n",
+    "from matplotlib.ticker import MaxNLocator\n",
+    "import seaborn as sns\n",
+    "from statsmodels.nonparametric.smoothers_lowess import lowess\n",
+    "\n",
+    "import nomad.io.base as loader\n",
+    "import nomad.stop_detection.utils as utils\n",
+    "import nomad.stop_detection.lachesis as LACHESIS\n",
+    "import nomad.stop_detection.dbscan as TADBSCAN\n",
+    "import nomad.stop_detection.grid_based as GRID_BASED # for oracle visits\n",
+    "import nomad.stop_detection.postprocessing as pp\n",
+    "\n",
+    "import nomad.visit_attribution.visit_attribution as visits\n",
+    "import nomad.filters as filters\n",
+    "import nomad.city_gen as cg\n",
+    "\n",
+    "from nomad.map_utils import blocks_to_mercator_gdf\n",
+    "from nomad.contact_estimation import compute_stop_detection_metrics"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "b7b4fc41",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "with open('config_2_stops.json', 'r', encoding='utf-8') as f:\n",
+    "    config3 = json.load(f)\n",
+    "\n",
+    "config=config3"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "78f6e650",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "poi_table = gpd.read_file(config[\"city_file\"], layer=\"buildings\").rename(columns={\"id\":\"location\"})\n",
+    "poi_table = blocks_to_mercator_gdf(poi_table, block_size=15, false_easting=-4265699.0, false_northing=4392976.0)\n",
+    "\n",
+    "sparse_path=config[\"output_files\"][\"sparse_path\"]\n",
+    "diaries_path=config[\"output_files\"][\"diaries_path\"]\n",
+    "\n",
+    "sparse_df = loader.from_file(sparse_path, format=\"parquet\")\n",
+    "diaries_df = loader.from_file(diaries_path, format=\"parquet\").rename(columns={\"identifier\":\"user_id\"})"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "ced917c8",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import os\n",
+    "os.makedirs('figures', exist_ok=True)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "8f199a72",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "delta_roam_values = np.linspace(1, 125, 100)\n",
+    "\n",
+    "config[\"algos\"] = {\n",
+    "    **{f\"lachesis_delta_{int(dr)}\": {\n",
+    "        \"func\": LACHESIS.lachesis,\n",
+    "        \"params\": {\n",
+    "            'dt_max': 60,\n",
+    "            'delta_roam': dr\n",
+    "        }\n",
+    "    } for dr in delta_roam_values}\n",
+    "}\n",
+    "\n",
+    "results_list = []\n",
+    "for user in tqdm(diaries_df.user_id.unique(), desc='Processing users'):\n",
+    "    sparse = sparse_df.query(\"user_id==@user\").copy()\n",
+    "    truth = diaries_df.query(\"user_id==@user\").copy()\n",
+    "    \n",
+    "    # Run all Lachesis variations\n",
+    "    for algo_name, algo_config in config[\"algos\"].items():\n",
+    "        # Run stop detection\n",
+    "        stops = algo_config[\"func\"](sparse, **algo_config[\"params\"], x=\"x\", y=\"y\", timestamp='timestamp')\n",
+    "        \n",
+    "        # Map stops to buildings\n",
+    "        stops[\"location\"] = visits.point_in_polygon(\n",
+    "            stops, \n",
+    "            poi_table=poi_table, \n",
+    "            data_crs='EPSG:3857',\n",
+    "            max_distance=10, \n",
+    "            location_id='location', \n",
+    "            x='x', \n",
+    "            y='y'\n",
+    "        )\n",
+    "        \n",
+    "        # Compute metrics\n",
+    "        metrics = compute_stop_detection_metrics(\n",
+    "            stops=stops,\n",
+    "            truth=truth,\n",
+    "            user_id=user,\n",
+    "            algorithm=algo_name,\n",
+    "            traj_cols={'location_id': 'location'},\n",
+    "            timestamp='timestamp'\n",
+    "        )\n",
+    "        \n",
+    "        # Add delta_roam value to metrics\n",
+    "        metrics['delta_roam'] = algo_config['params']['delta_roam']\n",
+    "        \n",
+    "        results_list.append(metrics)\n",
+    "\n",
+    "# Convert to DataFrame\n",
+    "results_df = pd.DataFrame(results_list)\n",
+    "\n",
+    "print(f\"Computed metrics for {len(results_df)} user-algorithm combinations\")\n",
+    "print(results_df.head())\n",
+    "\n",
+    "# Plotting function adapted for delta_roam\n",
+    "def plot_metric_delta_roam(metric, title, save_path='figures'):\n",
+    "    \"\"\"Plot a metric vs delta_roam for Lachesis.\"\"\"\n",
+    "    chart_df = results_df.groupby(['delta_roam'])[metric].agg(['mean', 'sem']).reset_index()\n",
+    "    \n",
+    "    fig, ax = plt.subplots(figsize=(10, 6))\n",
+    "    \n",
+    "    # Add scatter plot background (individual data points)\n",
+    "    sns.scatterplot(data=results_df, x='delta_roam', y=metric,\n",
+    "                   alpha=0.2, s=30, ax=ax, color='steelblue')\n",
+    "    \n",
+    "    # Add line plot (mean values)\n",
+    "    ax.plot(chart_df['delta_roam'], chart_df['mean'], \n",
+    "            color='darkblue', linewidth=3, label='Mean')\n",
+    "    \n",
+    "    # Add confidence band (mean ± SEM)\n",
+    "    ax.fill_between(chart_df['delta_roam'], \n",
+    "                    chart_df['mean'] - chart_df['sem'],\n",
+    "                    chart_df['mean'] + chart_df['sem'],\n",
+    "                    alpha=0.3, color='steelblue')\n",
+    "    \n",
+    "    # Styling\n",
+    "    ax.set_xlabel('Delta Roam (meters)', fontsize=11, labelpad=10)\n",
+    "    ax.set_ylabel(title, fontsize=11, labelpad=10)\n",
+    "    ax.set_ylim(0, 1)\n",
+    "    ax.set_title(f'{title} vs Delta Roam (Sequential)', fontsize=13, pad=15, fontweight='bold')\n",
+    "    ax.grid(True, linestyle=':', alpha=0.6, linewidth=0.7)\n",
+    "    ax.tick_params(axis='both', which='major', labelsize=10, length=6, width=1.2)\n",
+    "    ax.minorticks_on()\n",
+    "    \n",
+    "    plt.tight_layout()\n",
+    "    plt.savefig(f\"{save_path}/lachesis_delta_roam_{metric}.svg\", bbox_inches='tight')\n",
+    "    plt.savefig(f\"{save_path}/lachesis_delta_roam_{metric}.png\", dpi=600, bbox_inches='tight')\n",
+    "    plt.show(block=False)\n",
+    "    plt.close()\n",
+    "\n",
+    "# Plot all metrics\n",
+    "os.makedirs('figures', exist_ok=True)\n",
+    "\n",
+    "# metrics = {\n",
+    "#     'precision': 'Precision',\n",
+    "#     'recall': 'Accuracy',\n",
+    "#     'f1': 'F1 Score',\n",
+    "#     'missed_fraction': 'Proportion of Stops Missed',\n",
+    "#     'merged_fraction': 'Proportion of Stops Merged',\n",
+    "#     'split_fraction': 'Proportion of Stops Split'\n",
+    "# }\n",
+    "\n",
+    "metrics = {\n",
+    "    'recall': 'Accuracy',\n",
+    "}\n",
+    "\n",
+    "for metric, title in metrics.items():\n",
+    "    plot_metric_delta_roam(metric, title)\n",
+    "\n",
+    "print(f\"Number of unique delta_roam values: {results_df['delta_roam'].nunique()}\")"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "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.10.0"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}