L07 — Correlazioni a finestre mobili e quanto salgono nei crolli

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}")