L07 — Rolling-window correlations and how far they rise in crashes

import datetime as dt

import matplotlib.pyplot as plt
import numpy as np
import polars as pl

from cvbook.dati import carica
from cvbook.metriche import drawdown, rendimenti

NOMI = ["btcusdt", "ethusdt", "solusdt"]  # ← PROVA / TRY: aggiungi "ftsemib" o "eni"
                                          # (aggiungili anche a avvio.prepara([...]))
DA = dt.date(2020, 9, 1)  # PROVA / TRY: sposta la data d'inizio

serie = [
    carica(n).filter(pl.col("data") >= DA).select(["data", "chiusura"]).rename({"chiusura": n})
    for n in NOMI
]
tabella = serie[0]
for s in serie[1:]:
    tabella = tabella.join(s, on="data")
tabella = tabella.sort("data")

date = tabella["data"].to_list()[1:]
M = np.column_stack([rendimenti(tabella[n].to_numpy()) for n in NOMI])
print(f"{len(NOMI)} serie allineate, {len(M)} giorni comuni, dal {date[0]} al {date[-1]}")

coppie = np.triu_indices(len(NOMI), 1)


def correlazione_media(blocco: np.ndarray) -> float:
    return float(np.corrcoef(blocco.T)[coppie].mean())


primo = M[:, 0]
peggiori = np.argsort(primo)[: len(primo) // 20]  # il 5% dei giorni peggiori
                                                   # PROVA / TRY: // 5 per il 20% (vedi esercizio 3)

print(f"correlazione media su tutto il periodo:        {correlazione_media(M):.3f}")
print(f"correlazione media sul 5% dei giorni peggiori: {correlazione_media(M[peggiori]):.3f}")
print("\nSembra che nei crolli la correlazione SCENDA. E' falso, ed e' il "
      "risultato del modo in cui abbiamo selezionato i giorni.")

FINESTRA = 60  # PROVA / TRY: 20 · 60 · 200 (vedi esercizio 2)

corr = np.array([correlazione_media(M[i - FINESTRA:i]) for i in range(FINESTRA, len(M))])
date_c = date[FINESTRA:]
dd = drawdown(np.concatenate([[1.0], np.cumprod(1 + primo)]))[FINESTRA + 1:]
brutti = dd < -0.30

with avvio.figura("schermo"):
    fig, ax = plt.subplots(figsize=(10, 4))
    ax.fill_between(date_c, 0, 1, where=brutti, step="mid", alpha=0.25,
                    label="mercato oltre il 30% sotto il massimo")
    ax.plot(date_c, corr, linewidth=1.3, label="correlazione media a 60 giorni")
    ax.axhline(corr.mean(), linestyle=":", linewidth=1.2)
    ax.set_ylim(0, 1)
    ax.set_ylabel("Correlazione media")
    ax.legend(loc="lower left")
    fig.autofmt_xdate()
    plt.show()

print(f"minimo {corr.min():.2f}   massimo {corr.max():.2f}   media {corr.mean():.2f}")
print(f"media nei periodi difficili: {corr[brutti].mean():.2f}")
print(f"media nel resto del tempo:   {corr[~brutti].mean():.2f}")

def oscillazione_residua(n: int, rho: float) -> float:
    return float(np.sqrt(1 / n + (1 - 1 / n) * rho))


# Due stimatori della stessa parola. `correlazione_media(M)` è la correlazione
# **statica** su tutto il periodo: 0,675. La media delle finestre mobili della
# cella precedente è 0,73, ed è quella che il capitolo usa in ogni suo conto.
# Non sono la stessa cosa e non vanno confuse — è l'errore che il capitolo
# stesso ha dovuto correggere. Qui si usa quella del capitolo.
#
# Two estimators of the same word: the static correlation over the whole period
# (0.675) and the average of the rolling windows (0.73). The chapter's numbers
# all use the second one, so this table uses it too.
rho_misurata = float(corr.mean())
rho_statica = correlazione_media(M)
print(f"correlazione media delle finestre mobili (quella del capitolo): {rho_misurata:.2f}")
print(f"correlazione statica sull'intero periodo:                       {rho_statica:.2f}\n")
etichetta_misurata = f"a corr. {rho_misurata:.2f}".replace(".", ",")
print(f"{'asset':>6s} {'indipendenti':>14s} {'a corr. 0,3':>14s} "
      f"{etichetta_misurata:>14s} {'a corr. 0,9':>14s}")
for n in (2, 3, 4, 5, 10, 15, 30):
    valori = "".join(f"{oscillazione_residua(n, rho) * 100:13.1f}%"
                     for rho in (0.0, 0.3, rho_misurata, 0.9))
    print(f"{n:6d} {valori}")

# `len(NOMI)` e non un 3 scritto a mano: la cella sopra invita ad aggiungere
# serie, e con un 3 fisso il numero effettivo restava quello di tre asset
# qualunque cosa il lettore ci mettesse dentro.
quanti = len(NOMI)
n_eff = 1 / (1 / quanti + (1 - 1 / quanti) * rho_misurata)
print(f"\nnumero di asset DAVVERO indipendenti equivalenti ai tuoi {quanti}: {n_eff:.2f}")