{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Kapitel 19: Die moderne Portfoliotheorie nach Markowitz\n",
    "\n",
    "Begleitnotebook zu *Optimierte Entscheidungsfindung mit Python*. Die Codezellen sind identisch mit den im Buch abgedruckten Programmen.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Einmalig ausfuehren: installiert alle im Buch verwendeten Pakete.\n",
    "# Lokal in einer virtuellen Umgebung genauso gueltig wie in Google Colab.\n",
    "%pip install --quiet ortools highspy cvxpy scipy numpy pandas polars \\\n",
    "    scikit-learn matplotlib plotly pyomo linopy pymoo pydantic openpyxl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Warum Diversifikation funktioniert\n",
    "\n",
    "`Diversifikation_Demo.py`\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!/usr/bin/env python3\n",
    "\n",
    "# Diversifikation_Demo.py\n",
    "\"\"\"\n",
    "Kapitel Markowitz: Portfoliorisiko in Abhaengigkeit von der Korrelation - die\n",
    "Handrechnung zum Gratis-Mittagessen, geprueft ueber die analytische Formel und\n",
    "ueber simulierte, korrelierte Renditen.\n",
    "\"\"\"\n",
    "\n",
    "import numpy as np\n",
    "\n",
    "MU = 0.08\n",
    "SIGMA = 0.20\n",
    "\n",
    "\n",
    "def sigma_portfolio(rho: float, w1: float = 0.5) -> float:\n",
    "    \"\"\"Analytische Formel: sigma_p^2 = w1^2 s1^2 + w2^2 s2^2 + 2 w1 w2 rho s1 s2.\"\"\"\n",
    "    w2 = 1 - w1\n",
    "    varianz = w1**2 * SIGMA**2 + w2**2 * SIGMA**2 + 2 * w1 * w2 * rho * SIGMA**2\n",
    "    return float(np.sqrt(varianz))\n",
    "\n",
    "\n",
    "def simuliere(rho: float, n: int = 500_000, seed: int = 7) -> float:\n",
    "    \"\"\"Erzeugt n korrelierte Renditepaare und misst die Portfolio-Volatilitaet empirisch.\"\"\"\n",
    "    rng = np.random.default_rng(seed)\n",
    "    kovarianz = np.array([[SIGMA**2, rho * SIGMA**2],\n",
    "                           [rho * SIGMA**2, SIGMA**2]])\n",
    "    renditen = rng.multivariate_normal([MU, MU], kovarianz, size=n)\n",
    "    portfolio = 0.5 * renditen[:, 0] + 0.5 * renditen[:, 1]\n",
    "    return float(portfolio.std(ddof=1))\n",
    "\n",
    "\n",
    "if __name__ == \"__main__\":\n",
    "    print(\"=\" * 78)\n",
    "    print(\"  PORTFOLIORISIKO IN ABHAENGIGKEIT VON DER KORRELATION\")\n",
    "    print(\"=\" * 78)\n",
    "    print(f\"{'rho':>6} | {'sigma_p (Formel)':>18} | {'sigma_p (Simulation)':>20} | {'Risikoreduktion':>16}\")\n",
    "    print(\"-\" * 78)\n",
    "    for rho in [1.0, 0.5, 0.0, -0.5, -1.0]:\n",
    "        formel = sigma_portfolio(rho)\n",
    "        sim = simuliere(rho)\n",
    "        reduktion = (1 - formel / SIGMA) * 100\n",
    "        print(f\"{rho:6.1f} | {formel*100:16.2f} % | {sim*100:18.2f} % | {reduktion:14.0f} %\")\n",
    "\n",
    "    print(\"\\nDie simulierten Werte (500.000 gezogene Renditepaare je rho) bestaetigen\")\n",
    "    print(\"die Formel aus der Handrechnung bis auf statistisches Rauschen - und das,\")\n",
    "    print(\"obwohl Formel und Simulation nichts voneinander wissen.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Vollimplementierung mit CVXPY\n",
    "\n",
    "`Markowitz_CVXPY.py`\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!/usr/bin/env python3\n",
    "\n",
    "# Markowitz_CVXPY.py\n",
    "\"\"\"\n",
    "Kapitel Markowitz: Markowitz-Mean-Variance-Optimierung mit CVXPY.\n",
    "Enthaelt GMV, Maximum Sharpe (Korn-Transformation), Effizienzgrenze\n",
    "und institutionelle Restriktionen.\n",
    "\n",
    "Eigenschaften:\n",
    "  * Spaltenreihenfolge erzwungen; Sektoren ueber NAMEN statt Indizes\n",
    "  * Effizienzgrenze bis zur TATSAECHLICH erreichbaren Maximalrendite\n",
    "    (nicht etwa max(mu)*0.95 - das waere unter Restriktionen oft\n",
    "     unerreichbar und die Kurve wuerde stillschweigend vorzeitig abbrechen)\n",
    "  * Pruefung aller Restriktionen nach dem Loesen\n",
    "  * Vergleich mit Gleichgewichtung als Realitaetscheck\n",
    "\"\"\"\n",
    "\n",
    "import os\n",
    "\n",
    "import cvxpy as cp\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import yfinance as yf\n",
    "from sklearn.covariance import LedoitWolf\n",
    "import matplotlib\n",
    "matplotlib.use(\"Agg\")\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "OUTPUT_DIR = os.path.join(os.path.dirname(os.path.abspath(__file__)), \"output\")\n",
    "os.makedirs(OUTPUT_DIR, exist_ok=True)\n",
    "\n",
    "TICKER = [\"AAPL\", \"MSFT\", \"NVDA\", \"AMZN\", \"JNJ\", \"PFE\", \"JPM\", \"GS\", \"XOM\", \"CVX\"]\n",
    "\n",
    "# Sektoren ueber Namen definiert, nicht ueber Positionen\n",
    "SEKTOREN = {\n",
    "    \"Technologie\": [\"AAPL\", \"MSFT\", \"NVDA\", \"AMZN\"],\n",
    "    \"Gesundheit\":  [\"JNJ\", \"PFE\"],\n",
    "    \"Finanzen\":    [\"JPM\", \"GS\"],\n",
    "    \"Energie\":     [\"XOM\", \"CVX\"],\n",
    "}\n",
    "SEKTORGRENZEN = {\"Technologie\": 0.35, \"Gesundheit\": 0.40,\n",
    "                 \"Finanzen\": 0.40, \"Energie\": 0.40}\n",
    "\n",
    "MAX_GEWICHT = 0.20      # hoechstens 20 % je Einzeltitel\n",
    "RISIKOFREI = 0.03       # 3 % p.a.\n",
    "HANDELSTAGE = 252\n",
    "\n",
    "\n",
    "def lade_daten():\n",
    "    \"\"\"Laedt Kurse und garantiert die Spaltenreihenfolge.\"\"\"\n",
    "    ende = pd.Timestamp.today().normalize()\n",
    "    start = ende - pd.DateOffset(years=2)\n",
    "    roh = yf.download(TICKER, start=start, end=ende, auto_adjust=True, progress=False)\n",
    "    if roh.empty:\n",
    "        raise SystemExit(\"Download fehlgeschlagen (Netz, Ticker oder Rate-Limit pruefen).\")\n",
    "\n",
    "    if isinstance(roh.columns, pd.MultiIndex):\n",
    "        kurse = roh[\"Close\"][TICKER].dropna()      # [TICKER] erzwingt die Reihenfolge\n",
    "    else:\n",
    "        kurse = roh[[\"Close\"]].dropna()\n",
    "        kurse.columns = TICKER\n",
    "\n",
    "    assert list(kurse.columns) == TICKER, \"Spaltenreihenfolge weicht ab!\"\n",
    "    renditen = kurse.pct_change().dropna()\n",
    "    mu = renditen.mean().values * HANDELSTAGE\n",
    "    sigma = LedoitWolf().fit(renditen.values).covariance_ * HANDELSTAGE\n",
    "    return kurse, renditen, mu, sigma\n",
    "\n",
    "\n",
    "def sektor_indizes():\n",
    "    \"\"\"Uebersetzt Sektornamen einmalig in Positionsindizes - mit Pruefung.\"\"\"\n",
    "    ergebnis = {}\n",
    "    for sektor, titel in SEKTOREN.items():\n",
    "        fehlend = [t for t in titel if t not in TICKER]\n",
    "        if fehlend:\n",
    "            raise ValueError(f\"Sektor '{sektor}': Ticker {fehlend} nicht im Universum.\")\n",
    "        ergebnis[sektor] = [TICKER.index(t) for t in titel]\n",
    "    return ergebnis\n",
    "\n",
    "\n",
    "def basis_restriktionen(w, skala=None):\n",
    "    \"\"\"\n",
    "    Standardrestriktionen. Bei der Korn-Transformation muessen ALLE Grenzen\n",
    "    mit kappa mitskaliert werden - dafuer dient der Parameter 'skala'.\n",
    "    \"\"\"\n",
    "    eins = 1.0 if skala is None else skala\n",
    "    idx = sektor_indizes()\n",
    "    bedingungen = [w >= 0, w <= MAX_GEWICHT * eins]\n",
    "    for sektor, positionen in idx.items():\n",
    "        bedingungen.append(cp.sum(w[positionen]) <= SEKTORGRENZEN[sektor] * eins)\n",
    "    return bedingungen\n",
    "\n",
    "\n",
    "def loese_gmv(sigma):\n",
    "    \"\"\"Global Minimum Variance: minimiere Risiko, ignoriere Rendite.\"\"\"\n",
    "    w = cp.Variable(len(sigma))\n",
    "    problem = cp.Problem(cp.Minimize(0.5 * cp.quad_form(w, sigma)),\n",
    "                         [cp.sum(w) == 1] + basis_restriktionen(w))\n",
    "    problem.solve()\n",
    "    if problem.status not in (\"optimal\", \"optimal_inaccurate\"):\n",
    "        raise SystemExit(f\"GMV nicht loesbar: {problem.status}\")\n",
    "    return w.value\n",
    "\n",
    "\n",
    "def loese_max_sharpe(mu, sigma, r_f):\n",
    "    \"\"\"Maximum Sharpe Ratio ueber die Korn-Transformation.\"\"\"\n",
    "    n = len(mu)\n",
    "    ueberrendite = mu - r_f\n",
    "    if np.all(ueberrendite <= 0):\n",
    "        raise SystemExit(\"Kein Titel schlaegt den risikofreien Zins - \"\n",
    "                         \"Max-Sharpe-Portfolio existiert nicht.\")\n",
    "\n",
    "    y = cp.Variable(n)\n",
    "    kappa = cp.Variable(nonneg=True)\n",
    "    bedingungen = ([ueberrendite @ y == 1, cp.sum(y) == kappa]\n",
    "                   + basis_restriktionen(y, skala=kappa))\n",
    "    problem = cp.Problem(cp.Minimize(0.5 * cp.quad_form(y, sigma)), bedingungen)\n",
    "    problem.solve()\n",
    "    if problem.status not in (\"optimal\", \"optimal_inaccurate\"):\n",
    "        raise SystemExit(f\"Max-Sharpe nicht loesbar: {problem.status}\")\n",
    "    return y.value / kappa.value\n",
    "\n",
    "\n",
    "def max_erreichbare_rendite(mu, sigma):\n",
    "    \"\"\"\n",
    "    Achtung: max(mu)*0.95 als obere Grenze der Frontier anzunehmen, ist unter\n",
    "    Positions- und Sektorgrenzen oft unerreichbar; die Kurve wuerde dann\n",
    "    stillschweigend abbrechen. Hier wird die tatsaechliche Obergrenze berechnet.\n",
    "    \"\"\"\n",
    "    w = cp.Variable(len(mu))\n",
    "    problem = cp.Problem(cp.Maximize(mu @ w), [cp.sum(w) == 1] + basis_restriktionen(w))\n",
    "    problem.solve()\n",
    "    return float(problem.value)\n",
    "\n",
    "\n",
    "def berechne_frontier(mu, sigma, ret_min, ret_max, punkte=40):\n",
    "    \"\"\"Effizienzgrenze durch Variation der Mindestrendite.\"\"\"\n",
    "    n = len(mu)\n",
    "    w = cp.Variable(n)\n",
    "    ziel_rendite = cp.Parameter()\n",
    "    problem = cp.Problem(\n",
    "        cp.Minimize(0.5 * cp.quad_form(w, sigma)),\n",
    "        [cp.sum(w) == 1, mu @ w >= ziel_rendite] + basis_restriktionen(w))\n",
    "\n",
    "    volas, renditen, uebersprungen = [], [], 0\n",
    "    for ziel in np.linspace(ret_min, ret_max, punkte):\n",
    "        ziel_rendite.value = ziel\n",
    "        problem.solve()\n",
    "        if problem.status in (\"optimal\", \"optimal_inaccurate\"):\n",
    "            volas.append(float(np.sqrt(w.value @ sigma @ w.value)))\n",
    "            renditen.append(float(mu @ w.value))\n",
    "        else:\n",
    "            uebersprungen += 1\n",
    "    if uebersprungen:\n",
    "        print(f\"  Hinweis: {uebersprungen} Zielrenditen waren nicht erreichbar.\")\n",
    "    return np.array(volas), np.array(renditen)\n",
    "\n",
    "\n",
    "def pruefe_restriktionen(w, bezeichnung):\n",
    "    \"\"\"Nach dem Loesen: haelt die Loesung wirklich alle Regeln ein?\"\"\"\n",
    "    assert abs(w.sum() - 1) < 1e-6, f\"{bezeichnung}: Summe != 1\"\n",
    "    assert w.min() > -1e-6, f\"{bezeichnung}: negatives Gewicht\"\n",
    "    assert w.max() < MAX_GEWICHT + 1e-6, f\"{bezeichnung}: Positionsgrenze verletzt\"\n",
    "    for sektor, positionen in sektor_indizes().items():\n",
    "        anteil = w[positionen].sum()\n",
    "        assert anteil < SEKTORGRENZEN[sektor] + 1e-6, \\\n",
    "            f\"{bezeichnung}: Sektor {sektor} bei {anteil:.3f} ueber Grenze\"\n",
    "\n",
    "\n",
    "def kennzahlen(w, mu, sigma, r_f):\n",
    "    rendite = float(mu @ w)\n",
    "    vola = float(np.sqrt(w @ sigma @ w))\n",
    "    return rendite, vola, (rendite - r_f) / vola\n",
    "\n",
    "\n",
    "if __name__ == \"__main__\":\n",
    "    kurse, renditen, mu, sigma = lade_daten()\n",
    "    n = len(TICKER)\n",
    "\n",
    "    w_gmv = loese_gmv(sigma)\n",
    "    w_sharpe = loese_max_sharpe(mu, sigma, RISIKOFREI)\n",
    "    w_gleich = np.ones(n) / n\n",
    "\n",
    "    pruefe_restriktionen(w_gmv, \"GMV\")\n",
    "    pruefe_restriktionen(w_sharpe, \"Max Sharpe\")\n",
    "\n",
    "    print(\"=\" * 88)\n",
    "    print(\"        ERGEBNISSE DER MEAN-VARIANCE-OPTIMIERUNG\")\n",
    "    print(\"=\" * 88)\n",
    "    print(f\"Datenbasis: {len(renditen)} Handelstage, {n} Titel, \"\n",
    "          f\"risikofreier Zins {RISIKOFREI*100:.1f} %\")\n",
    "    print(f\"Restriktionen: max. {MAX_GEWICHT*100:.0f} % je Titel, \"\n",
    "          f\"Sektorgrenzen {SEKTORGRENZEN}\\n\")\n",
    "\n",
    "    print(f\"{'Portfolio':<26} {'Rendite':>10} {'Volatilitaet':>13} {'Sharpe':>9}\")\n",
    "    print(\"-\" * 88)\n",
    "    for name, w in [(\"Global Minimum Variance\", w_gmv),\n",
    "                    (\"Maximum Sharpe Ratio\", w_sharpe),\n",
    "                    (\"Gleichgewichtung (1/N)\", w_gleich)]:\n",
    "        r, v, sr = kennzahlen(w, mu, sigma, RISIKOFREI)\n",
    "        print(f\"{name:<26} {r*100:>9.2f} % {v*100:>12.2f} % {sr:>9.2f}\")\n",
    "\n",
    "    # --- Gewichte je Titel ------------------------------------------------\n",
    "    print(\"\\n--- Optimierte Portfoliogewichte ---\")\n",
    "    sektor_je_titel = {t: s for s, titel in SEKTOREN.items() for t in titel}\n",
    "    tabelle = pd.DataFrame({\n",
    "        \"Ticker\": TICKER,\n",
    "        \"Sektor\": [sektor_je_titel[t] for t in TICKER],\n",
    "        \"Rendite p.a.\": [f\"{r*100:+6.1f} %\" for r in mu],\n",
    "        \"Vola p.a.\": [f\"{np.sqrt(sigma[i, i])*100:5.1f} %\" for i in range(n)],\n",
    "        \"GMV\": [f\"{w*100:5.1f} %\" for w in w_gmv],\n",
    "        \"Max Sharpe\": [f\"{w*100:5.1f} %\" for w in w_sharpe],\n",
    "    })\n",
    "    print(tabelle.to_string(index=False))\n",
    "\n",
    "    print(\"\\n--- Sektoraufteilung (Kontrolle) ---\")\n",
    "    print(f\"{'Sektor':<14} {'Grenze':>8} {'GMV':>9} {'Max Sharpe':>12}\")\n",
    "    for sektor, positionen in sektor_indizes().items():\n",
    "        print(f\"{sektor:<14} {SEKTORGRENZEN[sektor]*100:>7.0f} % \"\n",
    "              f\"{w_gmv[positionen].sum()*100:>8.1f} % \"\n",
    "              f\"{w_sharpe[positionen].sum()*100:>11.1f} %\")\n",
    "\n",
    "    # --- Effizienzgrenze --------------------------------------------------\n",
    "    ret_gmv = float(mu @ w_gmv)\n",
    "    ret_max = max_erreichbare_rendite(mu, sigma)\n",
    "    print(f\"\\nEffizienzgrenze von {ret_gmv*100:.2f} % bis {ret_max*100:.2f} % \"\n",
    "          f\"(unter Restriktionen tatsaechlich erreichbar)\")\n",
    "    print(f\"  Zum Vergleich: bester Einzeltitel {mu.max()*100:.2f} % - \"\n",
    "          f\"durch die Grenzen nicht erreichbar.\")\n",
    "    volas, rets = berechne_frontier(mu, sigma, ret_gmv, ret_max)\n",
    "\n",
    "    # --- Diagramm ---------------------------------------------------------\n",
    "    plt.figure(figsize=(11, 6.5))\n",
    "    plt.plot(volas * 100, rets * 100, \"b-\", lw=2.5, label=\"Effizienzgrenze (restringiert)\")\n",
    "\n",
    "    for w, farbe, marker, groesse, name in [\n",
    "            (w_gmv, \"green\", \"o\", 150, \"GMV\"),\n",
    "            (w_sharpe, \"red\", \"*\", 260, \"Max Sharpe\"),\n",
    "            (w_gleich, \"purple\", \"D\", 110, \"Gleichgewichtung\")]:\n",
    "        r, v, sr = kennzahlen(w, mu, sigma, RISIKOFREI)\n",
    "        plt.scatter([v * 100], [r * 100], color=farbe, marker=marker, s=groesse,\n",
    "                    zorder=5, label=f\"{name} (SR={sr:.2f})\")\n",
    "\n",
    "    # Kapitalmarktlinie durch den risikofreien Zins und das Tangentialportfolio\n",
    "    r_s, v_s, _ = kennzahlen(w_sharpe, mu, sigma, RISIKOFREI)\n",
    "    x_linie = np.array([0, v_s * 1.25])\n",
    "    plt.plot(x_linie * 100, (RISIKOFREI + (r_s - RISIKOFREI) / v_s * x_linie) * 100,\n",
    "             color=\"orange\", ls=\":\", lw=2, label=\"Kapitalmarktlinie\")\n",
    "\n",
    "    for i, t in enumerate(TICKER):\n",
    "        plt.scatter(np.sqrt(sigma[i, i]) * 100, mu[i] * 100, color=\"gray\", alpha=0.5, s=40)\n",
    "        plt.annotate(t, (np.sqrt(sigma[i, i]) * 100 + 0.4, mu[i] * 100), fontsize=8)\n",
    "\n",
    "    plt.title(\"Markowitz-Effizienzgrenze mit Sektor- und Positionsgrenzen\", fontsize=12)\n",
    "    plt.xlabel(\"Annualisierte Volatilitaet [%]\")\n",
    "    plt.ylabel(\"Annualisierte erwartete Rendite [%]\")\n",
    "    plt.grid(True, linestyle=\":\", alpha=0.6)\n",
    "    plt.legend(loc=\"best\")\n",
    "    plt.tight_layout()\n",
    "    ziel = os.path.join(OUTPUT_DIR, \"markowitz_efficient_frontier.png\")\n",
    "    plt.savefig(ziel, dpi=150)\n",
    "    print(f\"\\nDiagramm gespeichert unter '{ziel}'\")\n",
    "    print(\"=\" * 88)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Finde den Denkfehler\n",
    "\n",
    "`Renditeschaetzung_Falle.py`\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!/usr/bin/env python3\n",
    "\n",
    "# Renditeschaetzung_Falle.py\n",
    "\"\"\"\n",
    "Kapitel Markowitz: Warum geschaetzte Renditen noch gefaehrlicher sind als geschaetzte\n",
    "Kovarianzen.\n",
    "\n",
    "Das Kapitel Finanzdaten hat gezeigt, was eine schlecht geschaetzte Kovarianzmatrix\n",
    "anrichtet. Bei den erwarteten RENDITEN ist es schlimmer - aus einem einfachen\n",
    "statistischen Grund: Um eine Rendite von 0,04 % pro Tag von null zu\n",
    "unterscheiden, braucht man bei 1,2 % Tagesschwankung ueber tausend\n",
    "Beobachtungen. Kovarianzen konvergieren dagegen deutlich schneller.\n",
    "\n",
    "Das Experiment ist bewusst so gebaut, dass die Antwort feststeht: ALLE zwoelf\n",
    "Anlagen haben exakt dieselbe wahre Rendite und dieselbe Schwankung. Ein\n",
    "korrektes Ergebnis waere also Gleichgewichtung. Was der Optimierer stattdessen\n",
    "tut, zeigt der Lauf.\n",
    "\n",
    "Wiederholt wird 200 mal, damit das Ergebnis nicht von einer Zufallsstichprobe\n",
    "abhaengt.\n",
    "\n",
    "Benoetigt: numpy, cvxpy\n",
    "\"\"\"\n",
    "\n",
    "from __future__ import annotations\n",
    "\n",
    "import numpy as np\n",
    "import cvxpy as cp\n",
    "\n",
    "ANZAHL_ANLAGEN = 12\n",
    "BEOBACHTUNGEN = 250              # ein Jahr Tagesdaten\n",
    "BEOBACHTUNGEN_SPAETER = 2000     # der \"spaetere Verlauf\"\n",
    "WIEDERHOLUNGEN = 200\n",
    "\n",
    "# Die Wahrheit, die der Optimierer nicht kennt: alle Anlagen sind gleich.\n",
    "WAHRE_RENDITE = 0.0004           # 0,04 % pro Tag\n",
    "WAHRE_SCHWANKUNG = 0.012         # 1,2 % pro Tag\n",
    "RISIKOAVERSION = 10.0\n",
    "\n",
    "RNG = np.random.default_rng(5)\n",
    "\n",
    "\n",
    "def baue_optimierer() -> tuple[cp.Problem, cp.Variable, cp.Parameter, cp.Parameter]:\n",
    "    \"\"\"Baut das Markowitz-Problem EINMAL mit Parametern.\n",
    "\n",
    "    Bei 200 Wiederholungen ist das der Unterschied zwischen Sekunden und\n",
    "    Minuten: cp.Parameter erlaubt es, nur die Daten zu tauschen, statt den\n",
    "    Ausdrucksbaum jedes Mal neu zu kompilieren (siehe Kapitel Oekosystem).\n",
    "    \"\"\"\n",
    "    gewichte = cp.Variable(ANZAHL_ANLAGEN, nonneg=True)\n",
    "    renditen = cp.Parameter(ANZAHL_ANLAGEN)\n",
    "    kovarianz = cp.Parameter((ANZAHL_ANLAGEN, ANZAHL_ANLAGEN), PSD=True)\n",
    "    problem = cp.Problem(\n",
    "        cp.Maximize(renditen @ gewichte\n",
    "                    - RISIKOAVERSION * cp.quad_form(gewichte, kovarianz)),\n",
    "        [cp.sum(gewichte) == 1])\n",
    "    return problem, gewichte, renditen, kovarianz\n",
    "\n",
    "\n",
    "def erzeuge(perioden: int) -> np.ndarray:\n",
    "    \"\"\"Renditen, bei denen alle Anlagen identisch verteilt und unabhaengig sind.\"\"\"\n",
    "    return RNG.normal(WAHRE_RENDITE, WAHRE_SCHWANKUNG,\n",
    "                      (perioden, ANZAHL_ANLAGEN))\n",
    "\n",
    "\n",
    "if __name__ == \"__main__\":\n",
    "    problem, gewichte, renditen, kovarianz = baue_optimierer()\n",
    "\n",
    "    im_zeitraum, spaeter, spaeter_gleich = [], [], []\n",
    "    groesstes_gewicht, schaetzspanne = [], []\n",
    "\n",
    "    for _ in range(WIEDERHOLUNGEN):\n",
    "        schaetzdaten = erzeuge(BEOBACHTUNGEN)\n",
    "        spaetere_daten = erzeuge(BEOBACHTUNGEN_SPAETER)\n",
    "\n",
    "        geschaetzte_rendite = schaetzdaten.mean(axis=0)\n",
    "        geschaetzte_kovarianz = np.cov(schaetzdaten, rowvar=False)\n",
    "\n",
    "        renditen.value = geschaetzte_rendite\n",
    "        # Symmetrisieren und minimal anheben: numerisches Rauschen kann die\n",
    "        # Matrix sonst knapp unter die PSD-Grenze druecken (Kapitel Fundament).\n",
    "        kovarianz.value = ((geschaetzte_kovarianz + geschaetzte_kovarianz.T) / 2\n",
    "                           + 1e-10 * np.eye(ANZAHL_ANLAGEN))\n",
    "        problem.solve()\n",
    "\n",
    "        w = np.array(gewichte.value).ravel()\n",
    "        im_zeitraum.append(float((schaetzdaten @ w).mean()))\n",
    "        spaeter.append(float((spaetere_daten @ w).mean()))\n",
    "        spaeter_gleich.append(float(spaetere_daten.mean()))\n",
    "        groesstes_gewicht.append(float(w.max()))\n",
    "        schaetzspanne.append(float(geschaetzte_rendite.max()\n",
    "                                   - geschaetzte_rendite.min()))\n",
    "\n",
    "    prozent = lambda werte: float(np.mean(werte)) * 100\n",
    "\n",
    "    print(\"=\" * 78)\n",
    "    print(\"  WENN DER OPTIMIERER GESCHAETZTE RENDITEN GLAUBT\")\n",
    "    print(\"=\" * 78)\n",
    "    print(f\"{ANZAHL_ANLAGEN} Anlagen, {WIEDERHOLUNGEN} Wiederholungen, \"\n",
    "          f\"je {BEOBACHTUNGEN} Beobachtungen zur Schaetzung.\")\n",
    "    print(f\"In Wahrheit haben ALLE dieselbe Tagesrendite von \"\n",
    "          f\"{WAHRE_RENDITE * 100:.4f} % und dieselbe Schwankung.\")\n",
    "    print(\"Die richtige Antwort waere also: gleichgewichten.\")\n",
    "    print()\n",
    "    print(f\"Mittlere Spanne der geschaetzten Renditen: \"\n",
    "          f\"{prozent(schaetzspanne):.4f} Prozentpunkte\")\n",
    "    print(f\"  -> Das ist das {prozent(schaetzspanne) / (WAHRE_RENDITE * 100):.1f}-Fache \"\n",
    "          f\"des wahren Werts. Aus reinem Rauschen.\")\n",
    "    print(f\"Mittleres groesstes Einzelgewicht: \"\n",
    "          f\"{prozent(groesstes_gewicht):.1f} % \"\n",
    "          f\"(bei Gleichgewichtung waeren es {100 / ANZAHL_ANLAGEN:.1f} %)\")\n",
    "    print()\n",
    "    print(f\"{'':<38} {'Tagesrendite':>14}\")\n",
    "    print(\"-\" * 78)\n",
    "    print(f\"{'Wahrheit (alle Anlagen)':<38} {WAHRE_RENDITE * 100:>13.4f} %\")\n",
    "    print(f\"{'Optimierer IM Schaetzzeitraum':<38} {prozent(im_zeitraum):>13.4f} %\"\n",
    "          f\"   <- Illusion\")\n",
    "    print(f\"{'Optimierer AUSSERHALB':<38} {prozent(spaeter):>13.4f} %\")\n",
    "    print(f\"{'Gleichgewichtung AUSSERHALB':<38} {prozent(spaeter_gleich):>13.4f} %\")\n",
    "    print(\"-\" * 78)\n",
    "\n",
    "    faktor = prozent(im_zeitraum) / (WAHRE_RENDITE * 100)\n",
    "    vorsprung = prozent(spaeter) - prozent(spaeter_gleich)\n",
    "    print(f\"\\nIm Schaetzzeitraum verspricht der Optimierer das \"\n",
    "          f\"{faktor:.1f}-Fache der wahren Rendite.\")\n",
    "    print(f\"Ausserhalb bleibt davon nichts: Er liegt {abs(vorsprung):.4f} \"\n",
    "          f\"Prozentpunkte\")\n",
    "    print(f\"{'schlechter' if vorsprung < 0 else 'besser'} als blosse Gleichgewichtung.\")\n",
    "    print()\n",
    "    print(\"Der Grund ist statistisch, nicht finanzwirtschaftlich: Um eine\")\n",
    "    print(f\"Rendite von {WAHRE_RENDITE*100:.2f} % taeglich von null zu unterscheiden,\")\n",
    "    print(f\"braucht man bei {WAHRE_SCHWANKUNG*100:.1f} % Schwankung ueber tausend\")\n",
    "    print(\"Beobachtungen. Mit 250 misst man ueberwiegend Rauschen - und der\")\n",
    "    print(\"Optimierer nimmt jedes Rauschen fuer bare Muenze.\")\n",
    "    print()\n",
    "    print(\"Drei praktische Konsequenzen:\")\n",
    "    print(\"  1. Verzichten Sie auf Renditeschaetzungen, wo es geht. Das\")\n",
    "    print(\"     Minimum-Varianz-Portfolio braucht gar keine.\")\n",
    "    print(\"  2. Wenn Sie welche brauchen: schrumpfen Sie sie stark zur Mitte\")\n",
    "    print(\"     (James-Stein, Black-Litterman) - viel staerker als Kovarianzen.\")\n",
    "    print(\"  3. Begrenzen Sie Einzelgewichte. Was der Optimierer nicht darf,\")\n",
    "    print(\"     kann er auch nicht auf Rauschen setzen.\")\n",
    "    print(\"=\" * 78)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
