{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "---\n",
        "title: \"Italian Tourism Occupancy, 1956–2025: From Mass Growth to a New Record\"\n",
        "description: \"A reproducible analysis of ISTAT arrivals and nights spent in Italian accommodation establishments, from the long historical series to the post-pandemic recovery.\"\n",
        "author: \"Federico Viscioletti\"\n",
        "image: \"images/italian-tourism-occupancy-analysis.png\"\n",
        "date: \"2026-08-10\"\n",
        "lang: en\n",
        "categories: [data science, python, tourism, data visualization]\n",
        "pillar: data-science-insights\n",
        "cluster: data-cleaning\n",
        "pillar-stage: case-study\n",
        "jupyter: python3\n",
        "format:\n",
        "  html:\n",
        "    toc: true\n",
        "    code-fold: false\n",
        "    code-tools: true\n",
        "    code-overflow: wrap\n",
        "    embed-resources: false\n",
        "    fig-responsive: true\n",
        "execute:\n",
        "  warning: false\n",
        "  message: false\n",
        "  echo: true\n",
        "  cache: false\n",
        "\n",
        "fig-width: 7\n",
        "fig-height: 4.2\n",
        "fig-align: center\n",
        "out-width: 100%\n",
        "---\n",
        "\n",
        "<style>\n",
        "img.cover {\n",
        "  width: 100%;\n",
        "  aspect-ratio: 16 / 9;\n",
        "  object-fit: cover;\n",
        "  object-position: center 45%;\n",
        "}\n",
        ".image-credit {\n",
        "  margin-top: -0.5rem;\n",
        "  color: #7a7974;\n",
        "  font-size: 0.8rem;\n",
        "  text-align: right;\n",
        "}\n",
        "</style>\n",
        "\n",
        "<script src=\"https://cdn.plot.ly/plotly-3.6.0.min.js\"></script>\n",
        "\n",
        "<img src=\"images/italian-tourism-occupancy-analysis.png\" title=\"Aerial view of blue and white beach umbrellas\" class=\"cover\" />\n",
        "\n",
        "<p class=\"image-credit\">Original photo by <a href=\"https://unsplash.com/it/@rgaleriacom?utm_source=unsplash&utm_medium=referral&utm_content=creditCopyText\">Ricardo Gomez Angel</a> on <a href=\"https://unsplash.com/it/foto/fotografia-aerea-di-ombrelloni-da-giardino-blu-e-bianchi-hlmuzhcpCkY?utm_source=unsplash&utm_medium=referral&utm_content=creditCopyText\">Unsplash</a>.</p>\n",
        "\n",
        "# Introduction\n",
        "\n",
        "It is summer, it is hot, and I am thinking about going on holiday. More specifically, I am thinking about the seaside: crowded beaches, full hotels, and the usual question of where everyone is going.\n",
        "\n",
        "That made me want to look at tourism data. Not just to find out how many people visit Italy, but to understand how tourism has changed over time, which places receive the most visitors, where those visitors come from, and how strongly the whole system still depends on a few summer months.\n",
        "\n",
        "The data comes from three ISTAT workbooks. One contains the national series from 1956 to 2025; another describes recent tourism by province and guests’ country of residence; the third breaks the data down to Italian municipalities and individual months. They are not exactly designed to tell a story, so the first part of the work is making the different levels and definitions fit together.\n",
        "\n",
        "The main question I wanted to answer was:\n",
        "\n",
        "> **How has Italian tourism changed over time, and what does the latest recovery look like once we separate arrivals, nights, destinations, guests’ residence, and seasonality?**\n",
        "\n",
        "## Quick answers about Italian tourism\n",
        "\n",
        "**Did Italian tourism grow over the long run?** Yes. Nights spent in accommodation establishments increased from 52.6 million in 1956 to 436.7 million in 2019. The 2025 provisional figure is higher still, at 535.5 million nights.\n",
        "\n",
        "**What did the pandemic do?** Nights spent fell to 208.4 million in 2020 — about 52% below 2019 — before recovering to 466.2 million in 2024 and 535.5 million in 2025.\n",
        "\n",
        "**Where is tourism concentrated?** In 2025, Rome, Venice, Bolzano, Milan, and Trento were the leading provinces by nights spent. The municipal table tells a similar story, with Rome, Milan, Venice, and Florence at the top.\n",
        "\n",
        "**Who are the international visitors?** In the provincial-by-residence table, Germany is the largest foreign market in 2025, followed by the United States, France, and the United Kingdom.\n",
        "\n",
        "**Is Italian tourism seasonal?** Very much so. August accounts for about 17.6% of the nights recorded in the 2025 monthly municipal file, with July and September also prominent.\n",
        "\n",
        "The point of this article is not to produce one “tourism score”. It is to show how a few careful reshapes and checks turn official workbooks into a coherent, auditable analysis.\n",
        "\n",
        "# Loading the ISTAT workbooks\n",
        "\n",
        "To answer these questions, I used three ISTAT workbooks. They describe Italian tourism at different levels of detail, and each one is useful for a different part of the analysis:\n",
        "\n",
        "- `1_serie_storica.xlsx`: Italy-wide arrivals and nights spent from 1956 to 2025\n",
        "- `3_dati_provinciali_per_provenienza.xlsx`: tourism by province and guests’ country of residence from 2019 to 2025\n",
        "- `2_dati_comunali.xlsx`: annual data by municipality from 2014 to 2025, plus monthly municipal data from 2022 to 2025\n",
        "\n",
        "The files do not use exactly the same structure or unit of measurement. The historical workbook reports values in thousands, while the provincial and municipal workbooks report individual arrivals and nights. I kept this distinction explicit throughout the analysis rather than converting everything immediately and risking confusion between the series.\n",
        "\n",
        "Before looking at the results, I needed to clean the files, name the columns clearly, and reshape some of the wider tables into a format that was easier to analyse."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the imports and visual style\"\n",
        "from pathlib import Path\n",
        "\n",
        "import matplotlib.pyplot as plt\n",
        "from matplotlib.ticker import FuncFormatter, MaxNLocator\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "import plotly.graph_objects as go\n",
        "\n",
        "HISTORICAL_FILE = Path(\"1_serie_storica.xlsx\")\n",
        "PROVINCIAL_FILE = Path(\"3_dati_provinciali_per_provenienza.xlsx\")\n",
        "MUNICIPAL_FILE = Path(\"2_dati_comunali.xlsx\")\n",
        "\n",
        "BG = \"#f7f6f2\"\n",
        "PANEL = \"#f9f8f5\"\n",
        "TEXT = \"#28251d\"\n",
        "MUTED = \"#7a7974\"\n",
        "GRID = \"#dcd9d5\"\n",
        "TEAL = \"#01696f\"\n",
        "CORAL = \"#e76f51\"\n",
        "\n",
        "COUNTRY_EN = {\n",
        "    \"Germania\": \"Germany\", \"Stati Uniti d'America\": \"United States\",\n",
        "    \"Francia\": \"France\", \"Regno Unito\": \"United Kingdom\",\n",
        "    \"Svizzera e Liechtenstein\": \"Switzerland and Liechtenstein\",\n",
        "    \"Paesi Bassi\": \"Netherlands\", \"Polonia\": \"Poland\", \"Austria\": \"Austria\",\n",
        "    \"Spagna\": \"Spain\", \"Belgio\": \"Belgium\", \"Canada\": \"Canada\",\n",
        "    \"Australia\": \"Australia\", \"Brasile\": \"Brazil\", \"Israele\": \"Israel\",\n",
        "    \"India\": \"India\", \"Cina\": \"China\", \"Giappone\": \"Japan\",\n",
        "    \"Nuova Zelanda\": \"New Zealand\", \"Turchia\": \"Turkey\", \"Russia\": \"Russia\",\n",
        "    \"Danimarca\": \"Denmark\", \"Svezia\": \"Sweden\", \"Norvegia\": \"Norway\",\n",
        "    \"Finlandia\": \"Finland\", \"Irlanda\": \"Ireland\", \"Portogallo\": \"Portugal\",\n",
        "    \"Grecia\": \"Greece\", \"Repubblica Ceca\": \"Czechia\", \"Romania\": \"Romania\",\n",
        "}\n",
        "\n",
        "PROVINCE_EN = {\n",
        "    \"Bolzano-Bozen\": \"Bolzano\", \"Venezia\": \"Venice\", \"Milano\": \"Milan\",\n",
        "    \"Roma\": \"Rome\", \"Firenze\": \"Florence\", \"Napoli\": \"Naples\",\n",
        "    \"Torino\": \"Turin\", \"Sassari\": \"Sassari\", \"Livorno\": \"Livorno\",\n",
        "    \"Brescia\": \"Brescia\", \"Verona\": \"Verona\", \"Trento\": \"Trento\",\n",
        "    \"Rimini\": \"Rimini\", \"Como\": \"Como\", \"Siena\": \"Siena\", \"Palermo\": \"Palermo\",\n",
        "    \"Udine\": \"Udine\",\n",
        "}\n",
        "\n",
        "plt.rcParams.update({\n",
        "    \"figure.facecolor\": BG,\n",
        "    \"axes.facecolor\": PANEL,\n",
        "    \"axes.edgecolor\": GRID,\n",
        "    \"axes.labelcolor\": TEXT,\n",
        "    \"axes.titlecolor\": TEXT,\n",
        "    \"xtick.color\": MUTED,\n",
        "    \"ytick.color\": MUTED,\n",
        "    \"text.color\": TEXT,\n",
        "    \"font.size\": 11,\n",
        "    \"axes.titlesize\": 18,\n",
        "    \"axes.titleweight\": \"bold\",\n",
        "    \"legend.frameon\": False,\n",
        "})\n",
        "pd.options.display.float_format = \"{:,.2f}\".format"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## From wide Excel layouts to tidy tables\n",
        "\n",
        "The first challenge - as always in EDA - was not analysing the data, but reading it properly. The three workbooks were built for publication rather than for analysis: they have different header structures, different levels of detail, and several columns that need to be interpreted before they can be used.\n",
        "\n",
        "For example, the historical series has two header rows before the actual data begins. The provincial and municipal files instead contain wide blocks of columns split by accommodation type, guest residence, or month. I wrote a small parsing function for each workbook so that these choices were explicit rather than hidden in one long cleaning step."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the workbook parsing functions\"\n",
        "def read_historical(path):\n",
        "    raw = pd.read_excel(path, header=None).iloc[6:, :19].copy()\n",
        "    raw.columns = [\n",
        "        \"year\", \"arrivals_residents\", \"arrivals_non_residents\", \"arrivals_total\",\n",
        "        \"arrivals_hotels_residents\", \"arrivals_hotels_non_residents\", \"arrivals_hotels_total\",\n",
        "        \"arrivals_other_residents\", \"arrivals_other_non_residents\", \"arrivals_other_total\",\n",
        "        \"nights_residents\", \"nights_non_residents\", \"nights_total\",\n",
        "        \"nights_hotels_residents\", \"nights_hotels_non_residents\", \"nights_hotels_total\",\n",
        "        \"nights_other_residents\", \"nights_other_non_residents\", \"nights_other_total\",\n",
        "    ]\n",
        "    raw[\"year\"] = pd.to_numeric(raw[\"year\"].astype(str).str.extract(r\"(\\d{4})\")[0], errors=\"coerce\")\n",
        "    for col in raw.columns[1:]:\n",
        "        raw[col] = pd.to_numeric(raw[col], errors=\"coerce\")\n",
        "    return raw.dropna(subset=[\"year\", \"arrivals_total\", \"nights_total\"]).astype({\"year\": int})\n",
        "\n",
        "\n",
        "def read_provincial(path):\n",
        "    raw = pd.read_excel(path, header=None).iloc[2:, :21].copy()\n",
        "    raw.columns = [\n",
        "        \"year\", \"region_code\", \"region\", \"province_code\", \"province\",\n",
        "        \"country_code\", \"country\", *[f\"accommodation_{i}\" for i in range(14)],\n",
        "    ]\n",
        "    data = raw[[\"year\", \"region\", \"province_code\", \"province\", \"country_code\", \"country\",\n",
        "                \"accommodation_12\", \"accommodation_13\"]].copy()\n",
        "    data.columns = [\"year\", \"region\", \"province_code\", \"province\", \"country_code\", \"country\",\n",
        "                    \"arrivals\", \"nights\"]\n",
        "    data[\"year\"] = pd.to_numeric(data[\"year\"], errors=\"coerce\")\n",
        "    data[\"province_code\"] = (data[\"province_code\"].astype(\"string\")\n",
        "                              .str.replace(r\"\\.0$\", \"\", regex=True)\n",
        "                              .str.zfill(3))\n",
        "    data[\"country_code\"] = pd.to_numeric(data[\"country_code\"], errors=\"coerce\")\n",
        "    data[\"arrivals\"] = pd.to_numeric(data[\"arrivals\"], errors=\"coerce\")\n",
        "    data[\"nights\"] = pd.to_numeric(data[\"nights\"], errors=\"coerce\")\n",
        "    data = data.dropna(subset=[\"year\", \"country_code\", \"province\", \"country\", \"arrivals\", \"nights\"])\n",
        "    data[\"year\"] = data[\"year\"].astype(int)\n",
        "    data[[\"arrivals\", \"nights\"]] = data[[\"arrivals\", \"nights\"]].round().astype(\"Int64\")\n",
        "    return data\n",
        "\n",
        "\n",
        "def read_municipal_annual(path):\n",
        "    raw = pd.read_excel(path, sheet_name=0, header=None).iloc[6:, :].copy()\n",
        "    data = pd.DataFrame({\n",
        "        \"year\": pd.to_numeric(raw.iloc[:, 0], errors=\"coerce\"),\n",
        "        \"region\": raw.iloc[:, 2].astype(\"string\").str.strip(),\n",
        "        \"province\": raw.iloc[:, 4].astype(\"string\").str.strip(),\n",
        "        \"municipality\": raw.iloc[:, 5].astype(\"string\").str.strip(),\n",
        "        \"arrivals_residents\": pd.to_numeric(raw.iloc[:, 8], errors=\"coerce\"),\n",
        "        \"arrivals_non_residents\": pd.to_numeric(raw.iloc[:, 9], errors=\"coerce\"),\n",
        "        \"arrivals\": pd.to_numeric(raw.iloc[:, 10], errors=\"coerce\"),\n",
        "        \"nights_residents\": pd.to_numeric(raw.iloc[:, 17], errors=\"coerce\"),\n",
        "        \"nights_non_residents\": pd.to_numeric(raw.iloc[:, 18], errors=\"coerce\"),\n",
        "        \"nights\": pd.to_numeric(raw.iloc[:, 19], errors=\"coerce\"),\n",
        "    })\n",
        "    data = data.dropna(subset=[\"year\", \"municipality\", \"arrivals\", \"nights\"]).astype({\"year\": int})\n",
        "    count_columns = [column for column in data if column.startswith((\"arrivals\", \"nights\"))]\n",
        "    data[count_columns] = data[count_columns].round().astype(\"Int64\")\n",
        "    return data\n",
        "\n",
        "\n",
        "def read_municipal_monthly(path):\n",
        "    raw = pd.read_excel(path, sheet_name=1, header=None).iloc[4:, :].copy()\n",
        "    data = pd.DataFrame({\n",
        "        \"year\": pd.to_numeric(raw.iloc[:, 0], errors=\"coerce\"),\n",
        "        \"province\": raw.iloc[:, 4].astype(\"string\").str.strip(),\n",
        "        \"municipality\": raw.iloc[:, 5].astype(\"string\").str.strip(),\n",
        "        **{f\"arrivals_{m:02d}\": pd.to_numeric(raw.iloc[:, 8 + m - 1], errors=\"coerce\") for m in range(1, 13)},\n",
        "        **{f\"nights_{m:02d}\": pd.to_numeric(raw.iloc[:, 20 + m - 1], errors=\"coerce\") for m in range(1, 13)},\n",
        "    })\n",
        "    data = data.dropna(subset=[\"year\", \"municipality\"]).astype({\"year\": int})\n",
        "    count_columns = [column for column in data if column.startswith((\"arrivals\", \"nights\"))]\n",
        "    data[count_columns] = data[count_columns].fillna(0).round().astype(\"Int64\")\n",
        "    return data\n",
        "\n",
        "\n",
        "historical = read_historical(HISTORICAL_FILE)\n",
        "provincial = read_provincial(PROVINCIAL_FILE)\n",
        "municipal_annual = read_municipal_annual(MUNICIPAL_FILE)\n",
        "municipal_monthly = read_municipal_monthly(MUNICIPAL_FILE)\n",
        "\n",
        "historical.agg(start_year=(\"year\", \"min\"), end_year=(\"year\", \"max\"), rows=(\"year\", \"size\"))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The main historical table contains 66 reported annual observations from 1956 to 2025. The 2025 row is marked provisional in the source and also introduces a definition change: other privately rented accommodation is included in the tourism-flow data from that year. I therefore use 2024 when comparing the long series with the pre-pandemic period, and label 2025 as provisional wherever it appears.\n",
        "\n",
        "# The long-term picture: tourism became much bigger\n",
        "\n",
        "I started with nights spent rather than arrivals. Arrivals tell us how many trips took place, but nights tell us something closer to the actual pressure and scale of tourism: a person staying for one night and a family spending two weeks somewhere should not count in the same way.\n",
        "\n",
        "The chart below shows the total number of nights spent in Italian accommodation establishments between 1956 and 2025."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the long-run chart code\"\n",
        "history_plot = historical.copy()\n",
        "history_plot[\"nights_millions\"] = history_plot[\"nights_total\"] / 1_000\n",
        "\n",
        "fig, ax = plt.subplots(figsize=(10.5, 5.8))\n",
        "ax.plot(history_plot[\"year\"], history_plot[\"nights_millions\"], color=TEAL, linewidth=2.8, marker=\"o\", markersize=4.5, zorder=2)\n",
        "ax.scatter([2019, 2020, 2025], history_plot.set_index(\"year\").loc[[2019, 2020, 2025], \"nights_millions\"],\n",
        "           color=[CORAL, CORAL, TEAL], s=42, zorder=3)\n",
        "ax.axvspan(2020, 2021, color=CORAL, alpha=0.08, label=\"Pandemic shock\")\n",
        "ax.set_title(\"Nights spent in Italian accommodation establishments\", loc=\"left\", pad=38)\n",
        "ax.text(0, 1.02, \"Italy-wide historical series; 2025 is provisional and includes a definition change\", transform=ax.transAxes, fontsize=10.5, color=MUTED)\n",
        "ax.set_xlabel(\"Year\")\n",
        "ax.set_ylabel(\"Nights spent (millions)\")\n",
        "ax.grid(axis=\"y\", color=GRID, linewidth=0.8, alpha=0.85)\n",
        "ax.grid(axis=\"x\", visible=False)\n",
        "ax.xaxis.set_major_locator(MaxNLocator(integer=True, nbins=9))\n",
        "ax.legend(loc=\"upper left\")\n",
        "for spine in [\"top\", \"right\"]:\n",
        "    ax.spines[spine].set_visible(False)\n",
        "ax.spines[\"left\"].set_color(GRID)\n",
        "ax.spines[\"bottom\"].set_color(GRID)\n",
        "plt.tight_layout()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The overall direction is clear: Italian tourism became much bigger over the last seventy years. Nights spent rose from 52.6 million in 1956 to 436.7 million in 2019. The provisional figure for 2025 reaches 535.5 million.\n",
        "\n",
        "But this is not a simple growth story. The series includes changes in how accommodation is classified, a visible discontinuity in the late 1980s, and the sudden collapse caused by the pandemic in 2020.\n",
        "\n",
        "## What happened around 1987?\n",
        "\n",
        "The drop around 1987 immediately caught my attention. Total nights spent fall from 342.3 million in 1986 to 249.7 million in 1987, which looks like a 27.1% collapse in a single year.\n",
        "\n",
        "But arrivals fall by only 5.3% over the same period. The average stay also declines, from 5.94 to 4.57 nights, but not enough to fully explain such a large fall in total nights. That made me suspect that the break was at least partly about definitions and classification rather than tourists suddenly disappearing."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the late-1980s chart code\"\n",
        "late_eighties = historical.loc[historical[\"year\"].between(1980, 2005)].copy()\n",
        "for column in [\"nights_total\", \"nights_hotels_total\", \"nights_other_total\"]:\n",
        "    late_eighties[column + \"_millions\"] = late_eighties[column] / 1_000\n",
        "\n",
        "fig, ax = plt.subplots(figsize=(10.5, 5.8))\n",
        "ax.plot(late_eighties[\"year\"], late_eighties[\"nights_total_millions\"], color=TEXT, linewidth=2.8, marker=\"o\", markersize=4.5, label=\"Total nights\", zorder=3)\n",
        "ax.plot(late_eighties[\"year\"], late_eighties[\"nights_hotels_total_millions\"], color=TEAL, linewidth=2.4, label=\"Hotels and similar establishments\")\n",
        "ax.plot(late_eighties[\"year\"], late_eighties[\"nights_other_total_millions\"], color=CORAL, linewidth=2.4, label=\"Other accommodation establishments\")\n",
        "ax.axvline(1987, color=CORAL, linestyle=(0, (4, 4)), linewidth=1.4, alpha=0.85)\n",
        "ax.annotate(\"1987: break concentrated in\\nother accommodation\", xy=(1987, 249.7), xytext=(1988.2, 315),\n",
        "            arrowprops=dict(arrowstyle=\"->\", color=CORAL, linewidth=1.1), color=TEXT, fontsize=10.5,\n",
        "            bbox=dict(boxstyle=\"round,pad=0.35\", facecolor=BG, edgecolor=\"none\", alpha=0.92))\n",
        "ax.set_title(\"The late-1980s trough is mostly an accommodation-series break\", loc=\"left\", pad=38)\n",
        "ax.text(0, 1.02, \"Nights spent by accommodation group, Italy, 1980–2005\", transform=ax.transAxes, fontsize=10.5, color=MUTED)\n",
        "ax.set_xlabel(\"Year\")\n",
        "ax.set_ylabel(\"Nights spent (millions)\")\n",
        "ax.set_xticks(late_eighties[\"year\"][::2])\n",
        "ax.grid(axis=\"y\", color=GRID, linewidth=0.8, alpha=0.85)\n",
        "ax.grid(axis=\"x\", visible=False)\n",
        "for spine in [\"top\", \"right\"]:\n",
        "    ax.spines[spine].set_visible(False)\n",
        "ax.spines[\"left\"].set_color(GRID)\n",
        "ax.spines[\"bottom\"].set_color(GRID)\n",
        "ax.legend(loc=\"upper left\", ncol=3)\n",
        "plt.tight_layout()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The accommodation breakdown makes the issue clearer. Hotel nights increase slightly, from 176.7 million to 183.1 million, while nights in “other accommodation” fall from 165.6 million to 66.5 million.\n",
        "\n",
        "ISTAT notes a classification change around this point: tourist residences move from extra-hotel to hotel accommodation from 1986. However, that note does not completely explain the fall in the total series, because the total should still include both accommodation groups.\n",
        "\n",
        "My reading is that the 1987 break combines a genuine reduction in average length of stay with a break in comparability in the accommodation data. I would therefore avoid describing it as a 27% collapse in tourism demand. It is better treated as a point where the historical series becomes less directly comparable.\n",
        "\n",
        "After 1987, total nights remain around 250–300 million through much of the 1990s and return to the 1986 level only around 2000–2001. Arrivals recover sooner, which suggests that the slower recovery in nights was mainly about shorter stays and changes in the accommodation mix, rather than simply fewer people travelling to Italy.\n",
        "\n",
        "## The pandemic was a shock, not a new normal\n",
        "\n",
        "The long historical series shows the scale of the 2020 collapse, but the provincial workbook gives a more detailed view of what happened during the recovery. Unlike the historical series, these data are reported as individual counts rather than thousands.\n",
        "\n",
        "To calculate national totals from this file, I selected the provincial total rows (`country_code == 0`) and summed them across provinces. I kept individual foreign countries separate, so that aggregate rows would not be counted twice."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "prov_total = (\n",
        "    provincial[provincial[\"country_code\"].eq(0)]\n",
        "    .groupby(\"year\", as_index=False)[[\"arrivals\", \"nights\"]].sum()\n",
        ")\n",
        "prov_total[\"average_stay\"] = prov_total[\"nights\"] / prov_total[\"arrivals\"]\n",
        "prov_total[\"nights_millions\"] = prov_total[\"nights\"] / 1_000_000\n",
        "prov_total[\"change_vs_2019\"] = prov_total[\"nights\"] / prov_total.loc[prov_total.year.eq(2019), \"nights\"].iloc[0] - 1\n",
        "prov_total"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the recovery chart code\"\n",
        "fig, ax = plt.subplots(figsize=(10.5, 5.8))\n",
        "ax.plot(prov_total[\"year\"], prov_total[\"nights_millions\"], color=TEAL, linewidth=2.8, marker=\"o\", markersize=5.5)\n",
        "ax.axhline(prov_total.loc[prov_total.year.eq(2019), \"nights_millions\"].iloc[0], color=MUTED, linestyle=(0, (4, 4)), linewidth=1.4, label=\"2019 level\")\n",
        "ax.set_title(\"The post-pandemic recovery moved beyond the old peak\", loc=\"left\", pad=38)\n",
        "ax.text(0, 1.02, \"Nights spent by year, summed across Italian provinces\", transform=ax.transAxes, fontsize=10.5, color=MUTED)\n",
        "ax.set_xlabel(\"Year\")\n",
        "ax.set_ylabel(\"Nights spent (millions)\")\n",
        "ax.set_xticks(prov_total[\"year\"])\n",
        "ax.grid(axis=\"y\", color=GRID, linewidth=0.8, alpha=0.85)\n",
        "ax.grid(axis=\"x\", visible=False)\n",
        "ax.legend(loc=\"upper left\")\n",
        "for spine in [\"top\", \"right\"]:\n",
        "    ax.spines[spine].set_visible(False)\n",
        "ax.spines[\"left\"].set_color(GRID)\n",
        "ax.spines[\"bottom\"].set_color(GRID)\n",
        "plt.tight_layout()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The recovery also changed the mix of visitors. Foreign residents accounted for about 50.5% of nights in 2019, fell to 31.4% in 2020, and reached 54.5% in 2024 in the historical series. That makes the rebound look less like a simple return of domestic demand and more like the reopening of Italy's international tourism circuit.\n",
        "\n",
        "# Where the nights are spent\n",
        "\n",
        "The provincial totals show a concentrated geography: a small group of destinations captures a large share of the national market. The ranking is not just a ranking of large cities. It also includes resort systems and provinces whose accommodation capacity serves a whole landscape — the Dolomites, the Venetian lagoon, Lake Garda, or the Roman metropolitan area."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "latest_provinces = (\n",
        "    provincial.query(\"year == 2025 and country_code == 0\")\n",
        "    .nlargest(10, \"nights\")[[\"province\", \"nights\", \"arrivals\"]]\n",
        "    .assign(nights_millions=lambda d: d[\"nights\"] / 1_000_000)\n",
        ")\n",
        "latest_provinces[\"province_en\"] = (\n",
        "    latest_provinces[\"province\"].astype(\"string\").str.strip().str.title()\n",
        "    .map(PROVINCE_EN).fillna(latest_provinces[\"province\"])\n",
        ")\n",
        "latest_provinces = latest_provinces.drop(columns=\"province\")\n",
        "latest_provinces = latest_provinces.rename(columns={\"province_en\": \"province\"})\n",
        "latest_provinces = latest_provinces[[\"province\"] + [col for col in latest_provinces.columns if col != \"province\"]]\n",
        "latest_provinces"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "fig, ax = plt.subplots(figsize=(10.5, 5.8))\n",
        "plot_data = latest_provinces.sort_values(\"nights_millions\")\n",
        "ax.barh(plot_data[\"province\"], plot_data[\"nights_millions\"], color=TEAL, alpha=0.8, edgecolor=\"none\", height=0.65)\n",
        "ax.set_title(\"The ten leading provinces by nights spent in 2025\", loc=\"left\", pad=38)\n",
        "ax.text(0, 1.02, \"Provisional provincial totals; English labels used for readability\", transform=ax.transAxes, fontsize=10.5, color=MUTED)\n",
        "ax.set_xlabel(\"Nights spent (millions)\", color=MUTED, fontsize=11, labelpad=10)\n",
        "ax.grid(axis=\"x\", color=GRID, linewidth=0.8, alpha=0.85)\n",
        "ax.grid(axis=\"y\", visible=False)\n",
        "for spine in [\"top\", \"right\", \"left\"]:\n",
        "    ax.spines[spine].set_visible(False)\n",
        "ax.spines[\"bottom\"].set_color(GRID)\n",
        "plt.tight_layout()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Rome led the 2025 provincial table with 56.0 million nights, followed by Venice and Bolzano. The municipal table provides a useful cross-check: it puts Rome, Milan, Venice, and Florence at the top when the analysis is performed at destination-municipality level. The small difference between the provincial and municipal totals is a reminder that official tables can have different coverage and disclosure rules; they should be compared as complementary views, not mechanically stacked together.\n",
        "\n",
        "# Who is visiting Italy?\n",
        "\n",
        "The residence table is especially useful because it turns “international tourism” into a set of visible flows. The chart below excludes `TOTALE`, Italian regions, and other aggregate codes, then ranks individual foreign residences by nights spent in 2025."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "foreign_2025 = provincial.query(\"year == 2025 and country_code.between(1, 887)\")\n",
        "top_markets = (\n",
        "    foreign_2025.groupby(\"country\", as_index=False)[[\"arrivals\", \"nights\"]].sum()\n",
        "    .nlargest(12, \"nights\")\n",
        "    .assign(nights_millions=lambda d: d[\"nights\"] / 1_000_000)\n",
        ")\n",
        "top_markets[\"country_en\"] = top_markets[\"country\"].map(COUNTRY_EN).fillna(top_markets[\"country\"])\n",
        "top_markets = top_markets.drop(columns=\"country\")\n",
        "top_markets = top_markets.rename(columns={\"country_en\": \"country\"})\n",
        "top_markets = top_markets[[\"country\"] + [col for col in top_markets.columns if col != \"country\"]]\n",
        "top_markets"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "fig, ax = plt.subplots(figsize=(10.5, 5.8))\n",
        "plot_data = top_markets.sort_values(\"nights_millions\")\n",
        "ax.barh(plot_data[\"country\"], plot_data[\"nights_millions\"], color=CORAL, alpha=0.8, edgecolor=\"none\", height=0.65)\n",
        "ax.set_title(\"Germany is Italy's largest foreign tourism market\", loc=\"left\", pad=38)\n",
        "ax.text(0, 1.02, \"Top foreign residences by nights spent, 2025 provisional provincial data\", transform=ax.transAxes, fontsize=10.5, color=MUTED)\n",
        "ax.set_xlabel(\"Nights spent (millions)\", color=MUTED, fontsize=11, labelpad=10)\n",
        "ax.grid(axis=\"x\", color=GRID, linewidth=0.8, alpha=0.85)\n",
        "ax.grid(axis=\"y\", visible=False)\n",
        "for spine in [\"top\", \"right\", \"left\"]:\n",
        "    ax.spines[spine].set_visible(False)\n",
        "ax.spines[\"bottom\"].set_color(GRID)\n",
        "plt.tight_layout()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Germany contributed 71.0 million nights in the 2025 provincial file, more than twice the United States total. The ranking is a useful descriptive view of demand, not a measure of the value or profitability of each market: average stay, accommodation type, spending, and geography differ across countries.\n",
        "\n",
        "# August still owns the calendar\n",
        "\n",
        "Annual totals hide the operational rhythm of tourism. The monthly municipal file covers 2022–2025; filtering it to 2025 and summing the twelve monthly columns shows how concentrated demand is in the summer."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "monthly_2025 = municipal_monthly.query(\"year == 2025\").copy()\n",
        "monthly_nights = pd.Series({\n",
        "    month: monthly_2025[f\"nights_{month:02d}\"].sum()\n",
        "    for month in range(1, 13)\n",
        "})\n",
        "monthly_nights.index = [\"Jan\", \"Feb\", \"Mar\", \"Apr\", \"May\", \"Jun\", \"Jul\", \"Aug\", \"Sep\", \"Oct\", \"Nov\", \"Dec\"]\n",
        "monthly_nights = monthly_nights / 1_000_000\n",
        "monthly_nights.to_frame(\"nights_millions\")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the monthly seasonality chart code\"\n",
        "fig, ax = plt.subplots(figsize=(10.5, 5.8))\n",
        "colors = [CORAL if month == \"Aug\" else TEAL for month in monthly_nights.index]\n",
        "ax.bar(monthly_nights.index, monthly_nights, color=colors, alpha=0.8, edgecolor=\"none\")\n",
        "ax.set_title(\"August is the single busiest month\", loc=\"left\", pad=38)\n",
        "ax.text(0, 1.02, \"Nights spent across municipalities, 2025 monthly data\", transform=ax.transAxes, fontsize=10.5, color=MUTED)\n",
        "ax.set_xlabel(\"Month\", color=MUTED, fontsize=11, labelpad=10)\n",
        "ax.set_ylabel(\"Nights spent (millions)\")\n",
        "ax.grid(axis=\"y\", color=GRID, linewidth=0.8, alpha=0.85)\n",
        "ax.grid(axis=\"x\", visible=False)\n",
        "for spine in [\"top\", \"right\"]:\n",
        "    ax.spines[spine].set_visible(False)\n",
        "ax.spines[\"left\"].set_color(GRID)\n",
        "ax.spines[\"bottom\"].set_color(GRID)\n",
        "plt.tight_layout()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "August accounts for roughly 17.6% of the nights in the 2025 monthly municipal file. July and September add another large block around it, which matters for transport, staffing, water use, housing pressure, and the environmental footprint of destinations. A national annual total is therefore only the beginning of the planning question.\n",
        "\n",
        "# Reading the flows: from residence to destination\n",
        "\n",
        "The previous charts rank places and markets independently. The two views below keep the relationship between them visible. First, a Sankey diagram shows how the leading foreign residences connect to the leading Italian provinces in 2025. The colours follow the article palette: coral identifies origin markets and teal identifies Italian destinations."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the Sankey chart code\"\n",
        "SANKey_YEAR = 2025\n",
        "TOP_COUNTRIES = 12\n",
        "TOP_PROVINCES = 15\n",
        "MIN_PRESENCES = 1_000\n",
        "\n",
        "flow_data = provincial.query(\n",
        "    \"year == @SANKey_YEAR and country_code.between(1, 887) and nights >= @MIN_PRESENCES\"\n",
        ").copy()\n",
        "flow_data[\"country\"] = flow_data[\"country\"].astype(str).str.strip()\n",
        "flow_data[\"province\"] = flow_data[\"province\"].astype(str).str.strip().str.title()\n",
        "flow_data = (\n",
        "    flow_data.groupby([\"country_code\", \"country\", \"province\"], as_index=False)\n",
        "    .agg(arrivals=(\"arrivals\", \"sum\"), nights=(\"nights\", \"sum\"))\n",
        ")\n",
        "\n",
        "top_countries_flow = (\n",
        "    flow_data.groupby([\"country_code\", \"country\"], as_index=False)[\"nights\"]\n",
        "    .sum()\n",
        "    .nlargest(TOP_COUNTRIES, \"nights\")\n",
        ")\n",
        "top_provinces_flow = (\n",
        "    flow_data.groupby(\"province\", as_index=False)[\"nights\"]\n",
        "    .sum()\n",
        "    .nlargest(TOP_PROVINCES, \"nights\")\n",
        ")\n",
        "flow_data = flow_data.merge(top_countries_flow[[\"country_code\", \"country\"]], on=[\"country_code\", \"country\"])\n",
        "flow_data = flow_data.merge(top_provinces_flow[[\"province\"]], on=\"province\")\n",
        "\n",
        "province_coordinates = {\n",
        "    \"Bolzano-Bozen\": (46.4983, 11.3548), \"Verona\": (45.4385, 10.9938),\n",
        "    \"Livorno\": (43.5443, 10.3262), \"Brescia\": (45.5356, 10.2147),\n",
        "    \"Udine\": (46.0711, 13.2346),\n",
        "    \"Torino\": (45.0703, 7.6869),\n",
        "    \"Venezia\": (45.4333, 12.3500), \"Trento\": (46.0667, 11.1333),\n",
        "    \"Rimini\": (44.0667, 12.5667), \"Sassari\": (40.7167, 8.5667),\n",
        "    \"Como\": (45.8000, 9.0833), \"Siena\": (43.3167, 11.3000),\n",
        "    \"Palermo\": (38.1166, 13.3636), \"Milano\": (45.4643, 9.1895),\n",
        "    \"Napoli\": (40.8522, 14.2681), \"Roma\": (41.9000, 12.4833),\n",
        "    \"Firenze\": (43.7667, 11.2500),\n",
        "}\n",
        "\n",
        "province_totals_flow = flow_data.groupby(\"province\", as_index=False)[\"nights\"].sum()\n",
        "missing_provinces = sorted(set(province_totals_flow[\"province\"]) - set(province_coordinates))\n",
        "if missing_provinces:\n",
        "    raise ValueError(f\"Missing coordinates for Sankey provinces: {missing_provinces}\")\n",
        "\n",
        "country_totals_flow = (\n",
        "    flow_data.groupby([\"country_code\", \"country\"], as_index=False)[\"nights\"]\n",
        "    .sum()\n",
        "    .sort_values(\"nights\", ascending=False)\n",
        "    .reset_index(drop=True)\n",
        ")\n",
        "province_totals_flow = province_totals_flow.assign(\n",
        "    latitude=province_totals_flow[\"province\"].map(lambda p: province_coordinates[p][0])\n",
        ").sort_values([\"latitude\", \"province\"], ascending=[False, True]).reset_index(drop=True)\n",
        "\n",
        "ISO2_BY_NUMERIC = {\n",
        "    1: \"FR\", 3: \"NL\", 4: \"DE\", 6: \"GB\", 7: \"IE\", 8: \"DK\", 9: \"GR\", 10: \"PT\", 11: \"ES\",\n",
        "    17: \"BE\", 18: \"LU\", 24: \"IS\", 28: \"NO\", 30: \"SE\", 32: \"FI\", 36: \"CH\", 38: \"AT\",\n",
        "    52: \"TR\", 55: \"LT\", 60: \"PL\", 61: \"CZ\", 63: \"SK\", 64: \"HU\", 66: \"RO\", 68: \"BG\",\n",
        "    72: \"UA\", 75: \"RU\", 400: \"US\", 404: \"CA\", 508: \"BR\", 624: \"IL\", 664: \"IN\",\n",
        "    720: \"CN\", 732: \"JP\", 800: \"AU\", 804: \"NZ\",\n",
        "}\n",
        "\n",
        "def flag(iso2):\n",
        "    if not iso2 or len(iso2) != 2:\n",
        "        return \"🌍\"\n",
        "    return \"\".join(chr(127397 + ord(c)) for c in iso2.upper())\n",
        "\n",
        "\n",
        "country_nodes = [\n",
        "    f\"{flag(ISO2_BY_NUMERIC.get(int(row.country_code)))} {COUNTRY_EN.get(row.country, row.country)}\"\n",
        "    for row in country_totals_flow.itertuples(index=False)\n",
        "]\n",
        "province_nodes = [f\"🇮🇹 {PROVINCE_EN.get(p, p)}\" for p in province_totals_flow[\"province\"]]\n",
        "labels = country_nodes + province_nodes\n",
        "country_y = np.linspace(0.03, 0.97, len(country_nodes))\n",
        "province_y = np.linspace(0.03, 0.97, len(province_nodes))\n",
        "node_x = [0.01] * len(country_nodes) + [0.99] * len(province_nodes)\n",
        "node_y = country_y.tolist() + province_y.tolist()\n",
        "\n",
        "country_index = {\n",
        "    (int(row.country_code), row.country): i\n",
        "    for i, row in enumerate(country_totals_flow.itertuples(index=False))\n",
        "}\n",
        "province_index = {\n",
        "    p: len(country_nodes) + i\n",
        "    for i, p in enumerate(province_totals_flow[\"province\"])\n",
        "}\n",
        "sources = [country_index[(int(r.country_code), r.country)] for r in flow_data.itertuples(index=False)]\n",
        "targets = [province_index[r.province] for r in flow_data.itertuples(index=False)]\n",
        "link_codes = flow_data[\"country_code\"].astype(int).tolist()\n",
        "\n",
        "default_links = [\"rgba(1, 105, 111, 0.28)\"] * len(flow_data)\n",
        "default_nodes = ([\"rgba(231, 111, 81, 0.88)\"] * len(country_nodes)\n",
        "                 + [\"rgba(1, 105, 111, 0.78)\"] * len(province_nodes))\n",
        "hover = [\n",
        "    f\"<b>{flag(ISO2_BY_NUMERIC.get(int(r.country_code)))} {COUNTRY_EN.get(r.country, r.country)} → 🇮🇹 {PROVINCE_EN.get(r.province, r.province)}</b><br>\"\n",
        "    f\"Nights spent: {r.nights:,.0f}<br>Arrivals: {r.arrivals:,.0f}<br>\"\n",
        "    f\"Average stay: {r.nights / r.arrivals:.2f} nights\"\n",
        "    if r.arrivals else f\"<b>{COUNTRY_EN.get(r.country, r.country)} → {PROVINCE_EN.get(r.province, r.province)}</b><br>Nights spent: {r.nights:,.0f}\"\n",
        "    for r in flow_data.itertuples(index=False)\n",
        "]\n",
        "\n",
        "sankey_fig = go.Figure(go.Sankey(\n",
        "    arrangement=\"snap\",\n",
        "    node=dict(\n",
        "        pad=24, thickness=18, x=node_x, y=node_y,\n",
        "        line=dict(color=\"#ffffff\", width=0.8), label=labels, color=default_nodes,\n",
        "        hovertemplate=\"%{label}<br>%{value:,.0f} nights spent<extra></extra>\",\n",
        "    ),\n",
        "    link=dict(\n",
        "        source=sources, target=targets, value=flow_data[\"nights\"], color=default_links,\n",
        "        customdata=hover, hovertemplate=\"%{customdata}<extra></extra>\",\n",
        "    ),\n",
        "))\n",
        "\n",
        "sankey_buttons = [dict(\n",
        "    label=\"All countries\", method=\"restyle\",\n",
        "    args=[{\"link.color\": [default_links], \"node.color\": [default_nodes], \"node.label\": [labels]}],\n",
        ")]\n",
        "for row in country_totals_flow.itertuples(index=False):\n",
        "    selected = int(row.country_code)\n",
        "    selected_provinces = set(flow_data.loc[flow_data[\"country_code\"].eq(selected), \"province\"])\n",
        "    selected_links = [\n",
        "        \"rgba(231, 111, 81, 0.90)\" if code == selected else \"rgba(108, 117, 125, 0.035)\"\n",
        "        for code in link_codes\n",
        "    ]\n",
        "    selected_nodes = [\n",
        "        \"rgba(231, 111, 81, 0.95)\" if int(code) == selected else \"rgba(108, 117, 125, 0)\"\n",
        "        for code in country_totals_flow[\"country_code\"]\n",
        "    ] + [\n",
        "        \"rgba(1, 105, 111, 0.95)\" if province in selected_provinces else \"rgba(108, 117, 125, 0)\"\n",
        "        for province in province_totals_flow[\"province\"]\n",
        "    ]\n",
        "    selected_labels = [\n",
        "        label if int(code) == selected else \"\"\n",
        "        for label, code in zip(country_nodes, country_totals_flow[\"country_code\"])\n",
        "    ] + [\n",
        "        label if province in selected_provinces else \"\"\n",
        "        for label, province in zip(province_nodes, province_totals_flow[\"province\"])\n",
        "    ]\n",
        "    sankey_buttons.append(dict(\n",
        "        label=f\"{flag(ISO2_BY_NUMERIC.get(selected))} {COUNTRY_EN.get(row.country, row.country)}\", method=\"restyle\",\n",
        "        args=[{\"link.color\": [selected_links], \"node.color\": [selected_nodes], \"node.label\": [selected_labels]}],\n",
        "    ))\n",
        "\n",
        "sankey_fig.update_layout(\n",
        "    title=dict(text=f\"From foreign residences to Italian provinces: tourism flows {SANKey_YEAR}\", x=0.02, xanchor=\"left\"),\n",
        "    font=dict(family=\"Arial, sans-serif\", size=13, color=TEXT),\n",
        "    paper_bgcolor=BG, plot_bgcolor=BG, margin=dict(l=20, r=20, t=145, b=20), height=900,\n",
        "    updatemenus=[dict(\n",
        "        buttons=sankey_buttons, direction=\"down\", x=0.02, y=1.16, xanchor=\"left\", yanchor=\"top\",\n",
        "        bgcolor=PANEL, bordercolor=GRID, borderwidth=1, font=dict(size=13),\n",
        "    )],\n",
        "    annotations=[dict(\n",
        "        text=\"Select a country to highlight its destination flows.\", x=0.02, y=1.08,\n",
        "        xref=\"paper\", yref=\"paper\", showarrow=False, xanchor=\"left\", font=dict(size=12, color=MUTED),\n",
        "    )],\n",
        ")\n",
        "sankey_fig"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The Sankey is deliberately filtered to the leading markets, destinations, and links above a minimum size. That keeps the visual legible while preserving the same country-to-province logic as the original standalone analysis. Select a country from the dropdown to isolate its network.\n",
        "\n",
        "## Following the flows across years\n",
        "\n",
        "The second view animates the same idea over time. To keep the article self-contained and reproducible, it uses a small local coordinate table for the leading foreign markets and province capitals rather than downloading geographic files at render time. The result retains the original map's useful interaction: press play or move the year slider to see destinations and international origins change from 2019 to 2025."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {},
      "source": [
        "#| code-fold: true\n",
        "#| code-summary: \"Show the animated map code\"\n",
        "from IPython.display import HTML, display\n",
        "\n",
        "TOP_COUNTRIES_MAP = 10\n",
        "TOP_DESTINATIONS_PER_COUNTRY = 3\n",
        "\n",
        "country_coordinates = {\n",
        "    \"Germania\": (51.16, 10.45), \"Stati Uniti d'America\": (39.83, -98.58),\n",
        "    \"Francia\": (46.23, 2.21), \"Regno Unito\": (55.38, -3.44),\n",
        "    \"Svizzera e Liechtenstein\": (46.82, 8.23), \"Paesi Bassi\": (52.13, 5.29),\n",
        "    \"Polonia\": (51.92, 19.15), \"Austria\": (47.52, 14.55),\n",
        "    \"Spagna\": (40.46, -3.75), \"Belgio\": (50.50, 4.47),\n",
        "    \"Canada\": (56.13, -106.35), \"Australia\": (-25.27, 133.78),\n",
        "    \"Brasile\": (-14.24, -51.93), \"Israele\": (31.05, 34.85),\n",
        "}\n",
        "\n",
        "def marker_size(values, low, high):\n",
        "    values = np.asarray(values, dtype=float)\n",
        "    if values.max() == values.min():\n",
        "        return np.repeat((low + high) / 2, len(values))\n",
        "    roots = np.sqrt(values)\n",
        "    return low + (roots - roots.min()) / (roots.max() - roots.min()) * (high - low)\n",
        "\n",
        "\n",
        "province_map_coordinates = {\n",
        "    p: coordinates for p, coordinates in province_coordinates.items()\n",
        "}\n",
        "\n",
        "def map_frame_data(year):\n",
        "    yearly = provincial.query(\"year == @year and country_code.between(1, 887)\").copy()\n",
        "    top_countries = (\n",
        "        yearly.groupby(\"country\", as_index=False)[\"nights\"]\n",
        "        .sum().nlargest(TOP_COUNTRIES_MAP, \"nights\")\n",
        "    )\n",
        "    flows = yearly.merge(top_countries[[\"country\"]], on=\"country\")\n",
        "    flows[\"rank\"] = flows.groupby(\"country\")[\"nights\"].rank(method=\"first\", ascending=False)\n",
        "    flows = flows[flows[\"rank\"] <= TOP_DESTINATIONS_PER_COUNTRY].copy()\n",
        "    flows[\"province\"] = flows[\"province\"].astype(str).str.strip().str.title()\n",
        "    flows[\"country_lat\"] = flows[\"country\"].map(lambda value: country_coordinates.get(value, (np.nan, np.nan))[0])\n",
        "    flows[\"country_lon\"] = flows[\"country\"].map(lambda value: country_coordinates.get(value, (np.nan, np.nan))[1])\n",
        "    flows[\"province_lat\"] = flows[\"province\"].map(lambda value: province_map_coordinates.get(value, (np.nan, np.nan))[0])\n",
        "    flows[\"province_lon\"] = flows[\"province\"].map(lambda value: province_map_coordinates.get(value, (np.nan, np.nan))[1])\n",
        "    flows = flows.dropna(subset=[\"country_lon\", \"country_lat\", \"province_lon\", \"province_lat\"])\n",
        "    destinations = flows.groupby([\"province\", \"province_lon\", \"province_lat\"], as_index=False).agg(\n",
        "        arrivals=(\"arrivals\", \"sum\"), nights=(\"nights\", \"sum\")\n",
        "    )\n",
        "    destinations[\"average_stay\"] = destinations[\"nights\"] / destinations[\"arrivals\"].replace(0, np.nan)\n",
        "    origins = flows.groupby([\"country\", \"country_lon\", \"country_lat\"], as_index=False).agg(nights=(\"nights\", \"sum\"))\n",
        "    return flows, destinations, origins\n",
        "\n",
        "\n",
        "map_years = sorted(provincial[\"year\"].unique())\n",
        "map_country_order = sorted({\n",
        "    country\n",
        "    for year in map_years\n",
        "    for country in (\n",
        "        provincial.query(\"year == @year and country_code.between(1, 887)\")\n",
        "        .groupby(\"country\")[\"nights\"]\n",
        "        .sum()\n",
        "        .nlargest(TOP_COUNTRIES_MAP)\n",
        "        .index\n",
        "    )\n",
        "})\n",
        "\n",
        "\n",
        "def map_traces(flows, destinations, origins):\n",
        "    traces = []\n",
        "    for country_index, country in enumerate(map_country_order):\n",
        "        country_flows = flows[flows[\"country\"].eq(country)].copy()\n",
        "        if country_flows.empty:\n",
        "            traces.extend([\n",
        "                go.Scattergeo(lon=[], lat=[], mode=\"lines\", showlegend=False),\n",
        "                go.Scattergeo(lon=[], lat=[], mode=\"markers\", showlegend=False),\n",
        "                go.Scattergeo(lon=[], lat=[], mode=\"markers\", showlegend=False),\n",
        "            ])\n",
        "            continue\n",
        "\n",
        "        line_lon, line_lat = [], []\n",
        "        for row in country_flows.itertuples():\n",
        "            line_lon += [row.country_lon, row.province_lon, None]\n",
        "            line_lat += [row.country_lat, row.province_lat, None]\n",
        "        traces.append(go.Scattergeo(\n",
        "            lon=line_lon, lat=line_lat, mode=\"lines\", hoverinfo=\"skip\", showlegend=False,\n",
        "            line=dict(color=\"rgba(1, 105, 111, 0.27)\", width=1.2),\n",
        "        ))\n",
        "\n",
        "        country_destinations = country_flows.groupby(\n",
        "            [\"province\", \"province_lon\", \"province_lat\"], as_index=False\n",
        "        ).agg(arrivals=(\"arrivals\", \"sum\"), nights=(\"nights\", \"sum\"))\n",
        "        country_destinations[\"average_stay\"] = (\n",
        "            country_destinations[\"nights\"] / country_destinations[\"arrivals\"].replace(0, np.nan)\n",
        "        )\n",
        "        country_destinations[\"province_display\"] = country_destinations[\"province\"].map(PROVINCE_EN).fillna(country_destinations[\"province\"])\n",
        "        traces.append(go.Scattergeo(\n",
        "            lon=country_destinations.province_lon, lat=country_destinations.province_lat, mode=\"markers\",\n",
        "            marker=dict(size=marker_size(country_destinations.nights, 8, 34), color=TEAL, opacity=0.82,\n",
        "                        line=dict(color=\"white\", width=0.8)),\n",
        "            customdata=np.c_[country_destinations.province_display, country_destinations.arrivals,\n",
        "                             country_destinations.nights, country_destinations.average_stay],\n",
        "            hovertemplate=(\"<b>%{customdata[0]}</b><br>Arrivals: %{customdata[1]:,.0f}\"\n",
        "                           \"<br>Nights spent: %{customdata[2]:,.0f}<br>Average stay: %{customdata[3]:.2f} nights<extra></extra>\"),\n",
        "            name=\"Destination provinces\", showlegend=country_index == 0,\n",
        "        ))\n",
        "\n",
        "        country_origin = country_flows.iloc[0]\n",
        "        country_nights = country_flows[\"nights\"].sum()\n",
        "        country_display = COUNTRY_EN.get(country, country)\n",
        "        traces.append(go.Scattergeo(\n",
        "            lon=[country_origin.country_lon], lat=[country_origin.country_lat], mode=\"markers+text\",\n",
        "            text=[country_display], textposition=\"top center\", textfont=dict(size=10, color=TEXT),\n",
        "            marker=dict(size=marker_size([country_nights], 7, 20)[0], color=CORAL, opacity=0.9,\n",
        "                        line=dict(color=\"white\", width=0.8)),\n",
        "            customdata=[[country_display, country_nights]],\n",
        "            hovertemplate=\"<b>%{customdata[0]}</b><br>Nights in the top three destinations: %{customdata[1]:,.0f}<extra></extra>\",\n",
        "            name=\"Foreign residences\", showlegend=country_index == 0,\n",
        "        ))\n",
        "    return traces\n",
        "\n",
        "\n",
        "map_frames = []\n",
        "for year in map_years:\n",
        "    map_flows, map_destinations, map_origins = map_frame_data(year)\n",
        "    map_frames.append(go.Frame(name=str(year), data=map_traces(map_flows, map_destinations, map_origins)))\n",
        "\n",
        "first_flows, first_destinations, first_origins = map_frame_data(map_years[0])\n",
        "flow_map_fig = go.Figure(\n",
        "    data=map_traces(first_flows, first_destinations, first_origins),\n",
        "    frames=map_frames,\n",
        ")\n",
        "flow_map_fig.update_layout(\n",
        "    title=dict(text=\"International tourism flows to Italian provinces\", x=0.02, xanchor=\"left\"),\n",
        "    margin=dict(l=15, r=15, t=105, b=90), height=1000, paper_bgcolor=BG,\n",
        "    font=dict(color=TEXT), legend=dict(orientation=\"h\", y=0.16, x=0.02),\n",
        "    geo=dict(\n",
        "        scope=\"world\", domain=dict(x=[0, 1], y=[0.12, 0.91]),\n",
        "        projection=dict(type=\"natural earth\", scale=1.6),\n",
        "        center=dict(lon=-35, lat=47), showland=True, landcolor=PANEL,\n",
        "        showcountries=True, countrycolor=GRID, showocean=True, oceancolor=\"#e8f1f2\",\n",
        "        lonaxis=dict(range=[-110, 35]), lataxis=dict(range=[0, 80]), bgcolor=BG,\n",
        "    ),\n",
        "    updatemenus=[dict(\n",
        "        type=\"buttons\", direction=\"left\", x=0.02, y=-0.1, xanchor=\"left\", yanchor=\"middle\",\n",
        "        showactive=False, pad=dict(t=4, r=8, b=4, l=8),\n",
        "        buttons=[\n",
        "        dict(label=\"▶ Play\", method=\"animate\", args=[None, {\"frame\": {\"duration\": 900, \"redraw\": True}, \"fromcurrent\": True}]),\n",
        "        dict(label=\"❚❚ Pause\", method=\"animate\", args=[[None], {\"frame\": {\"duration\": 0, \"redraw\": False}, \"mode\": \"immediate\"}]),\n",
        "    ]), dict(\n",
        "        type=\"buttons\", direction=\"left\", x=0.98, y=1.02, xanchor=\"right\", yanchor=\"top\",\n",
        "        showactive=False, pad=dict(t=4, r=8, b=4, l=8),\n",
        "        buttons=[dict(\n",
        "            label=\"↺ All countries\", method=\"restyle\",\n",
        "            args=[{\"visible\": [True] * (len(map_country_order) * 3)}],\n",
        "        )],\n",
        "    )],\n",
        "    sliders=[dict(\n",
        "        active=0, x=0.22, y=0.045, len=0.62, pad=dict(t=8, b=0),\n",
        "        currentvalue={\"prefix\": \"Year: \", \"font\": {\"size\": 18}},\n",
        "        steps=[dict(label=str(year), method=\"animate\", args=[[str(year)], {\"mode\": \"immediate\", \"frame\": {\"duration\": 0, \"redraw\": True}}]) for year in map_years],\n",
        "    )],\n",
        "    annotations=[dict(\n",
        "        text=f\"Click a foreign market to show its top {TOP_DESTINATIONS_PER_COUNTRY} destinations; use All countries to reset.\",\n",
        "        x=0.02, y=0.96, xref=\"paper\", yref=\"paper\", showarrow=False,\n",
        "        font=dict(size=12, color=MUTED), align=\"left\",\n",
        "    )],\n",
        ")\n",
        "flow_map_fig.update_layout(height=1000)\n",
        "\n",
        "map_click_script = \"\"\"\n",
        "const plot = document.getElementById('{plot_id}');\n",
        "const groupCount = %d;\n",
        "const totalTraces = groupCount * 3;\n",
        "plot.on('plotly_click', (event) => {\n",
        "  const point = event.points[0];\n",
        "  if (point.curveNumber %% 3 !== 2) return;\n",
        "  const selectedGroup = Math.floor(point.curveNumber / 3);\n",
        "  const visible = Array(totalTraces).fill(false);\n",
        "  visible[selectedGroup * 3] = true;\n",
        "  visible[selectedGroup * 3 + 1] = true;\n",
        "  visible[selectedGroup * 3 + 2] = true;\n",
        "  Plotly.restyle(plot, {visible});\n",
        "});\n",
        "\"\"\" % len(map_country_order)\n",
        "\n",
        "display(HTML(flow_map_fig.to_html(\n",
        "    full_html=False, include_plotlyjs=False, div_id=\"tourism-flow-map\", post_script=map_click_script\n",
        ")))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The animated map highlights a different aspect of the data than the Sankey. The Sankey is best for comparing the structure of the 2025 network; the map is better for seeing how the leading international markets and destinations move through the recovery period.\n",
        "\n",
        "# What to take away\n",
        "\n",
        "The three workbooks describe the same tourism system at different resolutions:\n",
        "\n",
        "- **Time:** Italian accommodation grew enormously over the long run, suffered a historic interruption in 2020, and reached a new provisional high in 2025.\n",
        "- **Place:** demand is geographically concentrated in a relatively small group of urban, cultural, mountain, and coastal destinations.\n",
        "- **Residence:** international tourism is central to the rebound, with Germany the largest foreign market in 2025.\n",
        "- **Calendar:** August remains the peak month, so annual growth does not distribute evenly through the year.\n",
        "\n",
        "There are also two methodological lessons. First, arrivals and nights answer different questions: arrivals measure the volume of trips, while nights capture duration. Second, 2025 should be handled carefully. The source marks it provisional and notes a change in coverage, so a clean analysis should flag it rather than present it as perfectly comparable with every prior year.\n",
        "\n",
        "The next natural extension would be to combine these occupancy tables with accommodation supply — beds, establishments, and occupancy rates — or to build a destination-level typology that separates city tourism from beach, mountain, and lake destinations.\n",
        "\n",
        "## Reproduce the analysis\n",
        "\n",
        "The source workbooks used here are the official ISTAT tourism tables included with this article. For definitions and future releases, consult the [ISTAT tourism statistics portal](https://www.istat.it/en/statistics-by-topic/industry-and-services/tourism/).\n",
        "\n",
        "The original files came from the `DCSC_Occupancy_in_collective_accommodation` analysis folder. The article code deliberately reads the workbooks with `header=None`, names the fields explicitly, and filters aggregate rows before calculating rankings. That makes the cleaning decisions visible and easy to adapt.\n",
        "\n",
        "*This article is descriptive. It does not estimate causal effects, forecast tourism, or treat provisional 2025 values as final.*"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 4
}
