{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Kapitel 20: Tail-Risiko, CVaR und Transaktionskosten{idx:Transaktionskosten}\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": [ "## Die zwei Schwächen des Markowitz-Modells\n", "\n", "`VaR_CVaR_Demo.py`\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!/usr/bin/env python3\n", "\n", "# VaR_CVaR_Demo.py\n", "\"\"\"\n", "Kapitel CVaR: Fat Tails, VaR und CVaR anschaulich.\n", "\n", "Teil 1: Wie oft treten \"unmoegliche\" Tage wirklich auf?\n", "Teil 2: Warum ist der VaR nicht subadditiv - ein Gegenbeispiel zum Nachrechnen.\n", "\"\"\"\n", "\n", "import numpy as np\n", "from scipy import stats\n", "\n", "\n", "def var_quantil(verluste: np.ndarray, alpha: float = 0.95) -> float:\n", " \"\"\"VaR = Quantil der Verlustverteilung (Verluste positiv, Gewinne negativ).\"\"\"\n", " return float(np.quantile(verluste, alpha))\n", "\n", "\n", "def cvar_rockafellar(verluste: np.ndarray, alpha: float = 0.95) -> float:\n", " \"\"\"\n", " CVaR ueber die Rockafellar-Uryasev-Formel:\n", " CVaR = min_gamma { gamma + 1/(1-alpha) * E[max(Verlust - gamma, 0)] }\n", "\n", " WICHTIG: Der naheliegende Weg \"Mittelwert aller Werte >= VaR\" ist FALSCH,\n", " sobald die Verteilung Atome hat (z. B. genau zwei moegliche Verluste).\n", " Dann liegt der VaR selbst auf einem Atom, und der Vergleich '>=' erfasst\n", " zu viel Wahrscheinlichkeitsmasse. Die Formel unten behandelt das korrekt -\n", " und ist zugleich genau der Ausdruck, den wir im Abschnitt 'Value at Risk\n", " und Conditional Value at Risk' optimieren.\n", " \"\"\"\n", " kandidaten = np.unique(verluste) # Optimum liegt immer auf einem Datenpunkt\n", " return float(min(g + np.mean(np.maximum(verluste - g, 0.0)) / (1.0 - alpha)\n", " for g in kandidaten))\n", "\n", "\n", "if __name__ == \"__main__\":\n", " rng = np.random.default_rng(2026)\n", "\n", " # --- Teil 1: Fat Tails (analytisch, nicht simuliert) ------------------\n", " print(\"=\" * 88)\n", " print(\" TEIL 1: WIE OFT TRITT DAS 'UNMOEGLICHE' EIN?\")\n", " print(\"=\" * 88)\n", " print(\"Vergleich: Normalverteilung gegen t-Verteilung mit 3 Freiheitsgraden\")\n", " print(\"(beide auf Standardabweichung 1 normiert).\\n\")\n", "\n", " t_verteilung = stats.t(df=3)\n", " skalierung = t_verteilung.std() # auf Varianz 1 bringen\n", "\n", " print(f\"{'Ereignis':<22} {'Normal':>14} {'t (df=3)':>14} {'Faktor':>11} \"\n", " f\"{'Normal: 1 Tag in':>18}\")\n", " print(\"-\" * 88)\n", " for k in [3, 4, 5, 6]:\n", " p_normal = 2 * stats.norm.sf(k) # beidseitig\n", " p_t = 2 * t_verteilung.sf(k * skalierung)\n", " print(f\"Abweichung > {k} Sigma {p_normal*100:>13.6f} % {p_t*100:>13.6f} % \"\n", " f\"{p_t/p_normal:>10.1f}x {1/p_normal/252:>15,.0f} Jahre\")\n", "\n", " print(\"\\nDeutung: Ein 5-Sigma-Tag ist unter Normalverteilung ein Ereignis von\")\n", " print(\"etwa einmal in 6.900 Jahren. Reale Aktienmaerkte liefern solche Tage\")\n", " print(\"mehrfach pro Jahrzehnt. Wer allein mit Varianz steuert, plant fuer\")\n", " print(\"eine Welt, in der Crashs praktisch nicht vorkommen.\")\n", "\n", " # --- Teil 2: VaR ist nicht subadditiv ---------------------------------\n", " print(\"\\n\" + \"=\" * 88)\n", " print(\" TEIL 2: WARUM DER VaR KEIN KOHAERENTES RISIKOMASS IST\")\n", " print(\"=\" * 88)\n", " print(\"Zwei unabhaengige Anleihen, je 100 EUR Nominal.\")\n", " print(\"Jede faellt mit 4 % Wahrscheinlichkeit aus (Verlust 100),\")\n", " print(\"sonst zahlt sie 2 EUR Kupon (Verlust -2).\\n\")\n", "\n", " ziehungen = 2_000_000\n", " verlust_a = np.where(rng.random(ziehungen) < 0.04, 100.0, -2.0)\n", " verlust_b = np.where(rng.random(ziehungen) < 0.04, 100.0, -2.0)\n", " verlust_ab = verlust_a + verlust_b\n", "\n", " print(f\"{'':<28} {'VaR 95%':>12} {'CVaR 95%':>12}\")\n", " print(\"-\" * 88)\n", " werte = {}\n", " for name, v in [(\"Anleihe A allein\", verlust_a),\n", " (\"Anleihe B allein\", verlust_b),\n", " (\"Portfolio A+B\", verlust_ab)]:\n", " werte[name] = (var_quantil(v), cvar_rockafellar(v))\n", " print(f\"{name:<28} {werte[name][0]:>12.2f} {werte[name][1]:>12.2f}\")\n", "\n", " var_summe = werte[\"Anleihe A allein\"][0] + werte[\"Anleihe B allein\"][0]\n", " cvar_summe = werte[\"Anleihe A allein\"][1] + werte[\"Anleihe B allein\"][1]\n", " print(f\"{'Summe der Einzelwerte':<28} {var_summe:>12.2f} {cvar_summe:>12.2f}\")\n", "\n", " var_port, cvar_port = werte[\"Portfolio A+B\"]\n", " print(\"-\" * 88)\n", " print(f\"VaR: Portfolio {var_port:7.2f} vs. Summe {var_summe:7.2f} -> \"\n", " f\"{'VERLETZT die Subadditivitaet!' if var_port > var_summe else 'subadditiv'}\")\n", " print(f\"CVaR: Portfolio {cvar_port:7.2f} vs. Summe {cvar_summe:7.2f} -> \"\n", " f\"{'subadditiv (kohaerent)' if cvar_port <= cvar_summe + 1e-6 else 'verletzt'}\")\n", "\n", " print(\"\\nDeutung: Einzeln betrachtet meldet der VaR fuer jede Anleihe einen\")\n", " print(\"GEWINN von 2 EUR - denn mit 96 % Wahrscheinlichkeit passiert nichts,\")\n", " print(\"und 4 % liegen unterhalb der 5-%-Schwelle. Im Portfolio steigt die\")\n", " print(\"Wahrscheinlichkeit mindestens eines Ausfalls auf 7,8 % und damit UEBER\")\n", " print(\"die Schwelle - der VaR springt auf 98. Er behauptet also, Streuung\")\n", " print(\"habe das Risiko erhoeht. Das ist oekonomisch unsinnig und der Grund,\")\n", " print(\"warum die Bankenaufsicht mit Basel III auf den Expected Shortfall\")\n", " print(\"umgestellt hat.\")\n", " print(\"=\" * 88)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Implementierung: CVaR-Portfolio mit Reibung\n", "\n", "`CVaR_Portfolio.py`\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!/usr/bin/env python3\n", "\n", "# CVaR_Portfolio.py\n", "\"\"\"\n", "Kapitel CVaR: CVaR-Optimierung mit L1-Transaktionskosten via CVXPY.\n", "\n", "Eigenschaften:\n", " * Spaltenreihenfolge erzwungen\n", " * Einheiten konsistent (alles taeglich, Annualisierung nur in der Ausgabe)\n", " * Szenario-Nebenbedingungen VEKTORISIERT statt in einer Python-Schleife\n", " (eine Matrixbedingung statt S einzelner Constraints - deutlich schneller)\n", " * Vergleich CVaR- gegen Varianz-Optimierung\n", " * Nachrechnung von VaR/CVaR aus den realisierten Szenarien\n", "\"\"\"\n", "\n", "import time\n", "\n", "import cvxpy as cp\n", "import numpy as np\n", "import pandas as pd\n", "import yfinance as yf\n", "\n", "TICKER = [\"AAPL\", \"MSFT\", \"NVDA\", \"AMZN\", \"JNJ\", \"PFE\", \"JPM\", \"GS\", \"XOM\", \"CVX\"]\n", "ALPHA = 0.95 # Konfidenzniveau: schlechteste 5 % der Tage\n", "MAX_GEWICHT = 0.25\n", "GEBUEHRENSATZ = 0.002 # 0,2 % je Einheit Turnover (Spread + Brokerage)\n", "RISIKOAVERSION = 1.5 # bezogen auf TAEGLICHE Groessen\n", "HANDELSTAGE = 252\n", "\n", "\n", "def lade_renditen():\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() # erzwingt eigene Spaltenreihenfolge\n", " else:\n", " kurse = roh[[\"Close\"]].dropna()\n", " kurse.columns = TICKER\n", " assert list(kurse.columns) == TICKER, \"Spaltenreihenfolge weicht ab!\"\n", " return kurse.pct_change().dropna()\n", "\n", "\n", "def optimiere_cvar(R, w_alt, vektorisiert=True):\n", " \"\"\"\n", " Maximiere: taegliche Rendite - lambda * CVaR - Transaktionskosten\n", " Alle Groessen TAEGLICH.\n", " \"\"\"\n", " S, N = R.shape\n", " mu_taeglich = R.mean(axis=0)\n", "\n", " w = cp.Variable(N, nonneg=True)\n", " gamma = cp.Variable() # wird im Optimum zum VaR\n", " u = cp.Variable(S, nonneg=True) # Ueberschuss ueber die Schwelle\n", "\n", " cvar = gamma + (1.0 / (S * (1.0 - ALPHA))) * cp.sum(u)\n", " turnover = cp.norm1(w - w_alt)\n", " kosten = GEBUEHRENSATZ * turnover\n", "\n", " ziel = cp.Maximize(mu_taeglich @ w - RISIKOAVERSION * cvar - kosten)\n", "\n", " bedingungen = [cp.sum(w) == 1, w <= MAX_GEWICHT]\n", " if vektorisiert:\n", " # EINE Matrixbedingung statt S einzelner - deutlich schneller\n", " bedingungen.append(u >= -(R @ w) - gamma)\n", " else:\n", " for s in range(S): # Alternative: S einzelne Constraints (langsamer)\n", " bedingungen.append(u[s] >= -R[s] @ w - gamma)\n", "\n", " problem = cp.Problem(ziel, bedingungen)\n", " problem.solve()\n", " if problem.status not in (\"optimal\", \"optimal_inaccurate\"):\n", " raise SystemExit(f\"CVaR-Problem nicht loesbar: {problem.status}\")\n", " return w.value, float(gamma.value), float(cvar.value), float(turnover.value)\n", "\n", "\n", "def optimiere_varianz(R, w_alt):\n", " \"\"\"Klassisches Mean-Variance zum Vergleich - ebenfalls taeglich gerechnet.\n", "\n", " ACHTUNG, DCP-Falle: Die Standardabweichung ist hier NICHT als\n", " cp.sqrt(cp.quad_form(w, sigma)) formulierbar. cp.sqrt ist konkav und\n", " verlangt ein konkaves Argument; quad_form ist konvex - CVXPY lehnt den\n", " Ausdruck mit einem DCPError ab, und zwar voellig zu Recht (Anhang\n", " Fehlerdiagnose). cp.psd_wrap() hilft dagegen nicht: Es behebt eine\n", " NUMERISCHE Beanstandung an sigma, keine Regelverletzung im Aufbau.\n", "\n", " Der Standardweg ist die Cholesky-Zerlegung sigma = L L^T. Damit gilt\n", " w' sigma w = ||L^T w||^2, also ist die Standardabweichung die 2-Norm\n", " eines AFFINEN Ausdrucks - konvex und damit regelkonform.\n", " \"\"\"\n", " N = R.shape[1]\n", " mu_taeglich = R.mean(axis=0)\n", " sigma = np.cov(R, rowvar=False, ddof=1)\n", " # Der winzige Diagonalzuschlag faengt den Fall ab, dass sigma numerisch\n", " # nur halbdefinit ist (mehr Titel als Handelstage, doppelte Spalten).\n", " L = np.linalg.cholesky(sigma + 1e-12 * np.eye(N))\n", "\n", " w = cp.Variable(N, nonneg=True)\n", " ziel = cp.Maximize(mu_taeglich @ w\n", " - RISIKOAVERSION * cp.norm2(L.T @ w)\n", " - GEBUEHRENSATZ * cp.norm1(w - w_alt))\n", " problem = cp.Problem(ziel, [cp.sum(w) == 1, w <= MAX_GEWICHT])\n", " problem.solve()\n", " if problem.status not in (\"optimal\", \"optimal_inaccurate\"):\n", " raise SystemExit(f\"Varianz-Problem nicht loesbar: {problem.status}\")\n", " return w.value\n", "\n", "\n", "def realisierte_kennzahlen(w, R):\n", " \"\"\"VaR und CVaR direkt aus den Szenarien - unabhaengige Gegenprobe.\"\"\"\n", " portfoliorenditen = R @ w\n", " verluste = -portfoliorenditen\n", " var = float(np.quantile(verluste, ALPHA))\n", " cvar = float(verluste[verluste >= var].mean())\n", " return var, cvar, float(portfoliorenditen.mean()), float(portfoliorenditen.std(ddof=1))\n", "\n", "\n", "if __name__ == \"__main__\":\n", " renditen = lade_renditen()\n", " R = renditen.values\n", " S, N = R.shape\n", " w_alt = np.ones(N) / N # Ausgangslage: Gleichgewichtung\n", "\n", " t0 = time.perf_counter()\n", " w_cvar, var_modell, cvar_modell, turnover = optimiere_cvar(R, w_alt, vektorisiert=True)\n", " dauer_vektor = time.perf_counter() - t0\n", "\n", " w_var = optimiere_varianz(R, w_alt)\n", "\n", " print(\"=\" * 90)\n", " print(\" CVaR-PORTFOLIO-OPTIMIERUNG MIT TRANSAKTIONSKOSTEN\")\n", " print(\"=\" * 90)\n", " print(f\"Datenbasis: {S} Handelstage, {N} Titel | Konfidenzniveau \"\n", " f\"{ALPHA*100:.0f} % | Loesungszeit {dauer_vektor:.2f} s\\n\")\n", "\n", " # --- Gegenprobe: Modellwerte gegen realisierte Szenariowerte ---------\n", " var_real, cvar_real, mu_real, sd_real = realisierte_kennzahlen(w_cvar, R)\n", " print(\"--- Gegenprobe: stimmen Modell und Szenarien ueberein? ---\")\n", " print(f\" VaR aus dem Modell (gamma): {var_modell*100:7.4f} % | \"\n", " f\"aus den Szenarien: {var_real*100:7.4f} %\")\n", " print(f\" CVaR aus dem Modell: {cvar_modell*100:7.4f} % | \"\n", " f\"aus den Szenarien: {cvar_real*100:7.4f} %\")\n", " assert abs(cvar_modell - cvar_real) < 1e-4, \"CVaR stimmt nicht mit den Szenarien!\"\n", " print(\" -> Der Rockafellar-Uryasev-Trick liefert exakt den empirischen CVaR.\")\n", "\n", " # --- Kennzahlen beider Portfolios ------------------------------------\n", " print(f\"\\n{'Portfolio':<24} {'Rendite p.a.':>13} {'Vola p.a.':>11} \"\n", " f\"{'VaR 95% (Tag)':>15} {'CVaR 95% (Tag)':>16} {'Turnover':>10}\")\n", " print(\"-\" * 90)\n", " for name, w in [(\"CVaR-optimiert\", w_cvar), (\"Varianz-optimiert\", w_var),\n", " (\"Gleichgewichtung\", w_alt)]:\n", " v, c, m, s = realisierte_kennzahlen(w, R)\n", " to = float(np.abs(w - w_alt).sum())\n", " print(f\"{name:<24} {m*HANDELSTAGE*100:>12.2f} % \"\n", " f\"{s*np.sqrt(HANDELSTAGE)*100:>10.2f} % {v*100:>14.3f} % \"\n", " f\"{c*100:>15.3f} % {to*100:>9.1f} %\")\n", "\n", " print(\"\\nHinweis zur Annualisierung: Renditen werden mit 252 skaliert,\")\n", " print(\"Volatilitaeten mit sqrt(252). Fuer VaR/CVaR ist eine solche Skalierung\")\n", " print(\"nur unter starken Annahmen (Unabhaengigkeit, kein Drift) zulaessig -\")\n", " print(\"sie werden hier deshalb bewusst als TAGESwerte ausgewiesen.\")\n", "\n", " # --- Allokationstabelle ------------------------------------------------\n", " print(\"\\n--- Allokation ---\")\n", " tabelle = pd.DataFrame({\n", " \"Ticker\": TICKER,\n", " \"vorher\": [f\"{v*100:5.1f} %\" for v in w_alt],\n", " \"CVaR-opt.\": [f\"{v*100:5.1f} %\" for v in w_cvar],\n", " \"Handel\": [f\"{(w_cvar[i]-w_alt[i])*100:+6.1f} %\" for i in range(N)],\n", " \"Varianz-opt.\": [f\"{v*100:5.1f} %\" for v in w_var],\n", " })\n", " print(tabelle.to_string(index=False))\n", " print(f\"\\nTurnover {turnover*100:.1f} % -> Transaktionskosten \"\n", " f\"{GEBUEHRENSATZ*turnover*100:.3f} % des Portfoliowerts\")\n", "\n", " # --- Laufzeitvergleich vektorisiert vs. Schleife ----------------------\n", " if S <= 600: # bei sehr vielen Szenarien zu langsam\n", " t0 = time.perf_counter()\n", " optimiere_cvar(R, w_alt, vektorisiert=False)\n", " dauer_schleife = time.perf_counter() - t0\n", " print(f\"\\n--- Laufzeit: {S} Nebenbedingungen aufbauen ---\")\n", " print(f\" vektorisiert (u >= -(R @ w) - gamma): {dauer_vektor:6.2f} s\")\n", " print(f\" Schleife ueber Szenarien (V01): {dauer_schleife:6.2f} s \"\n", " f\"({dauer_schleife/dauer_vektor:.1f}x langsamer)\")\n", " print(\"=\" * 90)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "name": "python", "version": "3.11" } }, "nbformat": 4, "nbformat_minor": 5 }