76 - L'equazione di Bellman derivata per intero, e perché la ricompensa futura va scontata

Se c'è una sola formula da capire bene in tutto il reinforcement learning, è l'equazione di Bellman. Ogni algoritmo che vedremo dopo (Q-Learning, SARSA, Deep Q-Learning) è un modo di risolverla, o di avvicinarsi a risolverla, quando non si conosce l'ambiente. Il suo contenuto è un'idea quasi banale che ha conseguenze enormi: il valore di una situazione è la ricompensa che prendo subito più il valore della situazione in cui mi ritrovo dopo. Qui la ricaviamo per intero, partendo dalla definizione di ritorno e senza saltare passaggi; spieghiamo perché serve lo sconto γ, con tre argomenti diversi e con numeri; la usiamo per calcolare, senza simulare nulla, la tabella dei valori dell'articolo sul processo decisionale di Markov; e ricaviamo la versione per la politica ottima, che è quella che davvero interessa.

Il vocabolario che serve

Riusiamo i termini dell'articolo sull'MDP (stato, azione, politica, ritorno, γ, P, R) e ne aggiungiamo alcuni. Il valore atteso E[X] è la media che si otterrebbe ripetendo l'esperimento infinite volte, e E[X | Y] è la media di X sapendo il valore di Y. La funzione valore dello stato v_π(s) = E_π[G_t | S_t = s] risponde a «quanto è buona questa situazione per me, se gioco secondo π?»; la funzione valore dell'azione q_π(s, a) = E_π[G_t | S_t = s, A_t = a] dice quanto vale fare a in s e poi seguire π.

L'equazione di Bellman lega il valore di uno stato ai valori degli stati successivi; prende il nome da Richard Bellman, al quale Sutton e Barto attribuiscono il termine «programmazione dinamica». Un backup è l'operazione di calcolare il valore di uno stato a partire da quelli dei successori, e un diagramma di backup la disegna; uno sweep è una passata di aggiornamento su tutti gli stati. Un punto fisso è un valore che l'aggiornamento lascia invariato: il valore corretto è il punto fisso dell'equazione di Bellman. Infine v*(s) è il valore massimo ottenibile da s fra tutte le politiche, e una politica avida (greedy) rispetto a un valore sceglie, in ogni stato, l'azione che massimizza «ricompensa più valore del successivo».

L'analogia: il tempo di viaggio rimasto

Sei in auto e il navigatore ti dice «arrivo fra 42 minuti». Come fa a saperlo? Non simula il viaggio intero dall'inizio: ragiona così. «Da qui al prossimo incrocio ci sono 3 minuti. All'incrocio, io so già (dalla stessa logica) che mancano 39 minuti. Quindi da qui mancano 3 + 39 = 42.» Il tempo rimasto da qui è il costo del primo pezzetto più il tempo rimasto da dove arrivo. Questo è tutto Bellman: il valore di un posto è il costo (la ricompensa) del prossimo passo più il valore del posto successivo. Se aggiungi incertezza («dall'incrocio posso andare a sinistra con il 70% o a destra con il 30%») il tempo rimasto è una media ponderata sulle possibilità. E c'è un dettaglio che richiederà lo sconto: un'ora fra dieci anni vale per te meno di un'ora adesso.

Prima di tutto: perché la ricompensa futura va scontata?

Nell'articolo precedente abbiamo definito il ritorno scontato G_t = R_{t+1} + γR_{t+2} + γ²R_{t+3} + … con 0 ≤ γ ≤ 1 senza dire perché. Ci sono tre ragioni indipendenti, e conviene conoscerle tutte perché ciascuna dice una cosa diversa su come scegliere γ.

Ragione 1: la somma deve avere un valore finito (la serie geometrica)

Se un problema non finisce mai (o può durare molto) e le ricompense sono, diciamo, sempre 1, la somma 1 + 1 + 1 + … diverge: non si possono confrontare due politiche con un ritorno infinito. Con lo sconto la somma è una serie geometrica, che si calcola in due righe. Poniamo S = 1 + γ + γ² + γ³ + …. Moltiplicando per γ: γS = γ + γ² + γ³ + …, cioè γS = S − 1. Quindi S − γS = 1, da cui

S = 1 / (1 − γ)         (per 0 ≤ γ < 1)

In generale, se ogni ricompensa è al più R_max in valore assoluto, il ritorno non può superare R_max / (1 − γ): è finito. Con γ = 0.9 la somma di tanti «1» vale 10; con γ = 0.99, 100. Il numero 1/(1 − γ) è, in un senso preciso, l'«orizzonte» dell'agente: quanti passi, grosso modo, conta nel decidere.

Ragione 2: lo sconto è la probabilità che la storia continui

Immagina che a ogni passo l'episodio possa finire all'improvviso (l'auto si guasta, il giocatore abbandona, il robot si scarica): continua con probabilità γ e finisce con probabilità 1 − γ. La ricompensa del passo k la incassi solo se sei sopravvissuto ai k passi precedenti, e la probabilità è γ^k. Allora il valore atteso della somma non scontata delle ricompense è

E[ Σ_k R_{k+1} · 1{sono ancora vivo al passo k} ] = Σ_k γ^k · E[R_{k+1}]

che è proprio il ritorno scontato. Quindi γ si può leggere come «la probabilità, per passo, che il mondo non finisca». Lo verifichiamo con una simulazione nel codice qui sotto: ricompensa 1 per passo, episodio che continua con probabilità γ.

Ragione 3: il futuro è più incerto e (spesso) preferiamo l'immediato

Un euro oggi vale più di un euro fra un anno, perché fra un anno potresti non averlo, o potresti averlo investito. Allo stesso modo, nell'RL, una ricompensa lontana dipende da tante cose che possono andare diversamente. Più γ è vicino a 0, più l'agente è «miope»; più è vicino a 1, più è «lungimirante». Il valore giusto dipende dal problema: per un robot che deve afferrare un oggetto bastano pochi passi di orizzonte, per giocare a scacchi conta una partita intera.

# Why a discount? Three views of the same number.
import numpy as np
rng = np.random.default_rng(0)
for gamma in (0.5, 0.9, 0.99):
    closed = 1 / (1 - gamma)                                            # sum of gamma^k for k = 0, 1, 2, ...  (geometric series)
    horizon = int(np.ceil(np.log(0.01) / np.log(gamma)))                # steps after which gamma^k < 1%
    # view 3: reward 1 per step, but at every step the episode SURVIVES with probability gamma. Undiscounted total reward:
    lengths = []
    for _ in range(20000):
        n = 1
        while rng.random() < gamma:                                    # survive another step with prob. gamma
            n += 1
        lengths.append(n)
    print(f"gamma={gamma}:  discounted sum of 1s = {closed:6.2f}   effective horizon ~{horizon:3d} steps   mean undiscounted reward if the episode ends with prob 1-gamma per step: {np.mean(lengths):6.2f}")

# Output (due esecuzioni consecutive, stessa uscita):
# gamma=0.5:  discounted sum of 1s =   2.00   effective horizon ~  7 steps   mean undiscounted reward if the episode ends with prob 1-gamma per step:   1.99
# gamma=0.9:  discounted sum of 1s =  10.00   effective horizon ~ 44 steps   mean undiscounted reward if the episode ends with prob 1-gamma per step:  10.04
# gamma=0.99:  discounted sum of 1s = 100.00   effective horizon ~459 steps   mean undiscounted reward if the episode ends with prob 1-gamma per step:  99.79

Per ogni γ il listato calcola tre cose. La prima è closed = 1/(1 − γ), la somma geometrica: 2, 10, 100. La seconda è horizon, il numero di passi dopo il quale il peso γ^k scende sotto l'1%, cioè log(0.01)/log(γ): 7, 44, 459 passi, una finestra di attenzione che si allarga molto rapidamente con γ. La terza è la simulazione della ragione 2: 20.000 episodi in cui, a ogni passo, il mondo continua con probabilità γ (rng.random() < gamma); la durata media, cioè la ricompensa totale non scontata (1 per passo), è 1.99, 10.04, 99.79 e coincide con 1/(1 − γ) entro il rumore di campionamento. Le tre ragioni descrivono lo stesso oggetto da tre punti di vista. Il seme fisso (default_rng(0)) rende l'uscita riproducibile; con un altro seme le due medie campionarie cambiano nell'ultima cifra.

La derivazione dell'equazione di Bellman

Ora il pezzo centrale. Partiamo dalle definizioni: il ritorno G_t = R_{t+1} + γ G_{t+1} (la ricorsione ricavata nell'articolo precedente) e il valore v_π(s) = E_π[G_t | S_t = s]. Procediamo per passi, dichiarando ogni volta che cosa si usa.

Passo 1: sostituire la ricorsione del ritorno.

v_π(s) = E_π[ G_t | S_t = s ]
       = E_π[ R_{t+1} + γ G_{t+1} | S_t = s ]

Passo 2: separare con la linearità del valore atteso (la media di una somma è la somma delle medie, e le costanti escono):

v_π(s) = E_π[ R_{t+1} | S_t = s ]  +  γ · E_π[ G_{t+1} | S_t = s ]

Passo 3: aprire il primo termine. La ricompensa immediata dipende dall'azione che scelgo (con probabilità π(a|s)) e dal prossimo stato in cui il mondo mi manda (con probabilità P(s′|s,a)). Per definizione di valore atteso come somma pesata sulle possibilità:

E_π[ R_{t+1} | S_t = s ] = Σ_a π(a|s) Σ_{s′} P(s′|s,a) · R(s,a,s′)

Passo 4: aprire il secondo termine con la legge della probabilità totale per i valori attesi. Per calcolare la media di G_{t+1} sapendo solo S_t = s, si condiziona sui possibili valori di (A_t, S_{t+1}) e si pesano per la loro probabilità:

E_π[ G_{t+1} | S_t = s ] = Σ_a π(a|s) Σ_{s′} P(s′|s,a) · E_π[ G_{t+1} | S_t = s, A_t = a, S_{t+1} = s′ ]

Passo 5: la proprietà di Markov. Il ritorno futuro G_{t+1} dipende dal passato solo attraverso lo stato presente, che qui è S_{t+1} = s′. Quindi, una volta noto s′, sapere anche che prima eravamo in s e abbiamo fatto a non aggiunge nulla:

E_π[ G_{t+1} | S_t = s, A_t = a, S_{t+1} = s′ ] = E_π[ G_{t+1} | S_{t+1} = s′ ] = v_π(s′)

L'ultima uguaglianza è la definizione di v_π applicata allo stato s′. Questo è il passo dove la proprietà di Markov entra davvero: senza di essa l'equazione non vale.

Passo 6: rimettere insieme. Mettiamo i risultati dei passi 3, 4 e 5 nel passo 2:

v_π(s) = Σ_a π(a|s) Σ_{s′} P(s′|s,a) · [ R(s,a,s′) + γ · v_π(s′) ]

Questa è l'equazione di Bellman per v_π. A parole: il valore di uno stato è la media, sulle mie scelte (π) e sulle risposte del mondo (P), di «ricompensa subito più γ per il valore dello stato in cui finisco».

Il diagramma di backup disegna esattamente questa somma. In alto lo stato s; sotto, i nodi neri sono le azioni, con le probabilità π(a|s) sui rami; sotto ancora, gli stati successivi con le probabilità P(s′|s,a); in fondo, ciò che si media, r + γ·v(s′). «Backup» perché il valore viaggia dal basso (i successori, già noti) verso l'alto (lo stato che vogliamo valutare).

La versione per l'azione e il legame fra le due

Con lo stesso ragionamento si ricava l'equazione per q_π. Il valore di fare a in s non ha più la media sulle azioni al primo livello (l'azione è già data), ma ce l'ha al livello successivo:

q_π(s,a) = Σ_{s′} P(s′|s,a) · [ R(s,a,s′) + γ · Σ_{a′} π(a′|s′) · q_π(s′,a′) ]
v_π(s)   = Σ_a π(a|s) · q_π(s,a)

La seconda riga dice che il valore dello stato è la media pesata dei valori delle azioni; la prima, che il valore di un'azione guarda un passo avanti e poi usa la stessa politica. Le due sono due viste della stessa struttura.

Usarla: valutazione iterativa di una politica

Se conosciamo P (la tabella dell'MDP), l'equazione di Bellman diventa una ricetta di calcolo: non sappiamo v_π ma sappiamo che deve soddisfare quell'uguaglianza; partiamo da una stima qualunque (zero ovunque) e applichiamo ripetutamente il lato destro come aggiornamento, finché non cambia più. Questo è il punto fisso. Il banco di prova è la griglia 4×4 dell'articolo sull'MDP, che qui riporto in forma compatta perché il listato sia eseguibile da solo: sedici stati, due terminali agli angoli (0 e 15), ricompensa −1 a ogni passo, P[s][a] come lista di quattro elementi (probabilità, stato_successivo, ricompensa, finito). Con slip=0.0 (il valore predefinito) ogni azione ha un solo esito.

import numpy as np

# Actions: 0=up, 1=right, 2=down, 3=left
MOVES = {0: (-1, 0), 1: (0, 1), 2: (1, 0), 3: (0, -1)}

class GridWorld:
    """A small MDP written out explicitly: P[s][a] = list of (probability, next_state, reward, done)."""
    def __init__(self, n=4, terminals=(0, 15), step_reward=-1.0, slip=0.0):
        self.n, self.terminals, self.n_states, self.n_actions = n, set(terminals), n * n, 4
        self.P = {s: {a: self._transitions(s, a, step_reward, slip) for a in range(4)} for s in range(self.n_states)}

    def _move(self, s, a):
        r, c = divmod(s, self.n)
        dr, dc = MOVES[a]
        r2, c2 = r + dr, c + dc
        if not (0 <= r2 < self.n and 0 <= c2 < self.n):       # bumping into a wall: stay where you are
            r2, c2 = r, c
        return r2 * self.n + c2

    def _transitions(self, s, a, step_reward, slip):
        if s in self.terminals:                                # terminal states: absorbing, no reward
            return [(1.0, s, 0.0, True)]
        outcomes = {}                                          # next_state -> probability
        for a2, p in [(a, 1.0 - slip), ((a + 1) % 4, slip / 2), ((a - 1) % 4, slip / 2)]:   # intended, or slip to a side
            if p > 0:
                s2 = self._move(s, a2)
                outcomes[s2] = outcomes.get(s2, 0.0) + p
        return [(p, s2, step_reward, s2 in self.terminals) for s2, p in outcomes.items()]

Questa è la stessa griglia, per esempio, che Sutton e Barto usano come Esempio 4.1 (compito episodico senza sconto, ricompensa −1 su tutte le transizioni, politica casuale equiprobabile), dove il valore finale è l'opposto del numero atteso di passi fino al termine. Con questa classe, la valutazione iterativa si scrive così, con γ = 1 (episodi finiti).

def policy_evaluation(env, pi, gamma=1.0, tol=1e-10, history=()):
    """Iterative policy evaluation: turn the Bellman equation into an update and sweep until nothing changes.
    pi[s][a] = probability of action a in state s."""
    V = np.zeros(env.n_states)
    snapshots = {}
    for sweep in range(1, 10_000):
        delta = 0.0
        for s in range(env.n_states):
            v_new = sum(pi[s][a] * sum(p * (r + gamma * V[s2] * (not done))        # no future value after a terminal step
                                       for p, s2, r, done in env.P[s][a])
                        for a in range(env.n_actions))
            delta = max(delta, abs(v_new - V[s]))
            V[s] = v_new                                                          # in-place update
        if sweep in history:
            snapshots[sweep] = V.copy()
        if delta < tol:
            return V, sweep, snapshots
    return V, sweep, snapshots

env = GridWorld()
random_pi = np.full((env.n_states, 4), 0.25)                   # equiprobable random policy
V, sweeps, snaps = policy_evaluation(env, random_pi, gamma=1.0, history=(1, 2, 3, 10))
for k, Vk in snaps.items():
    print(f"after {k:2d} sweeps:", np.round(Vk.reshape(4, 4), 1).tolist()[1])      # second row, as a taste
print(f"converged after {sweeps} sweeps. v_pi:")
print(np.round(V.reshape(4, 4), 2))

# Output (due esecuzioni consecutive, stessa uscita):
# after  1 sweeps: [-1.0, -1.5, -1.7, -1.8]
# after  2 sweeps: [-1.9, -2.8, -3.2, -3.4]
# after  3 sweeps: [-2.8, -4.0, -4.7, -4.9]
# after 10 sweeps: [-7.8, -10.4, -11.8, -11.9]
# converged after 272 sweeps. v_pi:
# [[  0. -14. -20. -22.]
#  [-14. -18. -20. -20.]
#  [-20. -20. -18. -14.]
#  [-22. -20. -14.   0.]]

policy_evaluation parte da V a zero. In ogni sweep scorre tutti gli stati e calcola la formula di Bellman per ciascuno: la doppia somma Σ_a π(a|s) Σ P·(r + γ·V[s2]); il fattore (not done) annulla il valore del successore se l'episodio è finito, perché dopo uno stato terminale non c'è futuro. L'aggiornamento è «sul posto» (V[s] = v_new): gli stati aggiornati prima nello stesso sweep sono già usati dagli altri. delta è la massima variazione dello sweep; quando scende sotto tol = 10⁻¹⁰, ci si ferma.

Dopo 1 sweep il valore di ogni stato è circa −1: si è «visto» solo la ricompensa del primo passo. Dopo 2 sweep, fra −2 e −3: l'informazione sul costo si è propagata di un passo. Dopo 10 sweep si è fra −8 e −12 e il valore sta ancora scendendo, perché le camminate casuali sono lunghe. Dopo 272 sweep converge, e la tabella finale (0, −14, −20, −22 …) è quella che il libro riporta per questo esempio nella Figura 4.1. Si vede l'idea del backup: il valore «si diffonde» a partire dai terminali (dove è 0) verso l'esterno, un passo per sweep.

Senza iterare: un sistema lineare

L'equazione di Bellman per una politica fissata è lineare in v: scritta in forma matriciale è v = r_π + γ P_π v, dove P_π[s, s′] è la probabilità di passare da s a s′ seguendo π e r_π[s] è la ricompensa immediata attesa. Si risolve con (I − γ P_π) v = r_π, un sistema di 16 equazioni in 16 incognite, con una chiamata a numpy.linalg.solve.

# The same numbers WITHOUT iterating: v = r_pi + gamma * P_pi v   ->   (I - gamma P_pi) v = r_pi   (a linear system)
def policy_matrices(env, pi, gamma=1.0):
    n = env.n_states
    P_pi = np.zeros((n, n)); r_pi = np.zeros(n)
    for s in range(n):
        for a in range(env.n_actions):
            for p, s2, r, done in env.P[s][a]:
                r_pi[s] += pi[s][a] * p * r                           # expected immediate reward under pi
                if not done:
                    P_pi[s, s2] += pi[s][a] * p                       # probability of moving s -> s2 (terminal: no continuation)
    return P_pi, r_pi

P_pi, r_pi = policy_matrices(env, random_pi)
V_exact = np.linalg.solve(np.eye(env.n_states) - 1.0 * P_pi, r_pi)
print("linear system vs iteration, max difference:", np.abs(V_exact - V).max())
residual = np.abs(r_pi + P_pi @ V_exact - V_exact).max()
print("Bellman residual  max| r + P v - v | :", residual)

# Output (due esecuzioni consecutive, stessa uscita):
# linear system vs iteration, max difference: 1.0726104449076956e-09
# Bellman residual  max| r + P v - v | : 3.552713678800501e-15

policy_matrices costruisce P_π e r_π accumulando, per ogni stato, le probabilità e le ricompense pesate per la politica; salta il termine di continuazione quando l'episodio finisce (if not done). Poi np.linalg.solve(np.eye(n) − γ·P_π, r_π) risolve il sistema. Il risultato coincide con l'iterazione (differenza massima dell'ordine di 10⁻⁹) e il residuo di Bellman, cioè quanto la soluzione viola l'equazione, è dell'ordine di 10⁻¹⁵: zero a precisione macchina. Il fatto che il sistema sia risolubile con γ = 1 dipende dal fatto che gli episodi terminano: la matrice I − P_π è invertibile perché dai terminali non si continua. L'equazione di Bellman ha dunque una soluzione che si calcola in due modi indipendenti, e i due coincidono.

# Compare with the Monte Carlo table of the previous lesson
V_mc = np.array([[0.0, -13.8, -19.6, -21.7],
                 [-14.0, -18.4, -19.8, -20.4],
                 [-19.9, -19.9, -18.4, -14.1],
                 [-22.1, -20.3, -14.2, 0.0]]).ravel()
print("largest gap between Monte Carlo and the exact value:", round(np.abs(V_mc - V_exact).max(), 2))

# Output (due esecuzioni consecutive, stessa uscita):
# largest gap between Monte Carlo and the exact value: 0.4

Il confronto con l'articolo precedente: la stima Monte Carlo, ottenuta mediando 4000 episodi per stato, differisce dal valore esatto al massimo di 0.4, su valori che arrivano a −22. Un errore standard di circa 0.3 per cella, misurato nell'articolo sull'MDP, rende plausibile che su 14 celle il massimo scarto sia 0.4: è rumore di campionamento. Il metodo Monte Carlo deve simulare migliaia di episodi per stato e restituisce una stima rumorosa; Bellman con il modello dà il valore esatto in una frazione di secondo (l'intero listato di questo articolo gira in circa due secondi su questa macchina, una prova). Il prezzo è che richiede di conoscere P.

Che cosa fa lo sconto al calcolo e alla risposta

Ripetiamo la valutazione della politica casuale per diversi valori di γ, contando gli sweep necessari e guardando due stati: l'1, adiacente al traguardo, e il 10, lontano.

# How does gamma change the work and the answer?  Policy evaluation of the random policy for different discount factors
print(f"{'gamma':>6s} {'sweeps':>7s}   value of state 1 (next to the goal)   value of state 10 (far from the goal)")
for gamma in (0.0, 0.5, 0.9, 0.99, 1.0):
    Vg, sweeps, _ = policy_evaluation(env, random_pi, gamma=gamma)
    print(f"{gamma:6.2f} {sweeps:7d}   {Vg[1]:30.2f}   {Vg[10]:30.2f}")

# Output (due esecuzioni consecutive, stessa uscita):
#  gamma  sweeps   value of state 1 (next to the goal)   value of state 10 (far from the goal)
#   0.00       2                            -1.00                            -1.00
#   0.50      23                            -1.69                            -1.92
#   0.90      95                            -5.28                            -6.61
#   0.99     231                           -11.95                           -15.32
#   1.00     272                           -14.00                           -18.00

Due effetti, entrambi misurati. Il costo del calcolo cresce con γ: da 2 sweep a γ = 0, 23 a 0.5, 95 a 0.9, 231 a 0.99, 272 a 1. Con γ più grande l'informazione deve propagarsi su un orizzonte più lungo, quindi servono più sweep a parità di tolleranza. La risposta cambia con γ: con γ = 0 l'agente è completamente miope e tutti gli stati valgono −1 (guarda solo la ricompensa del passo successivo, e quella è −1 ovunque), quindi non distingue «vicino al traguardo» da «lontano». Con γ più alto i valori si differenziano: a 0.99 lo stato 1 vale −11.95 e il 10 vale −15.32, a 1 valgono −14 e −18. Il valore di γ è parte della definizione del problema: ciò che l'agente considera «buono» dipende da quanto guarda lontano.

L'equazione di Bellman di ottimalità

La valutazione dice quanto vale una politica. Ciò che davvero vogliamo è la migliore. Sia v*(s) il valore massimo ottenibile da s fra tutte le politiche. In s l'agente ottimo sceglie l'azione che massimizza «ricompensa subito più valore (ottimo) di dove arrivo»; non c'è motivo di mediare su azioni peggiori. Sostituendo, nell'equazione di Bellman, la media sulle azioni con il massimo:

v*(s) = max_a Σ_{s′} P(s′|s,a) · [ R(s,a,s′) + γ · v*(s′) ]
q*(s,a) = Σ_{s′} P(s′|s,a) · [ R(s,a,s′) + γ · max_{a′} q*(s′,a′) ]

Nel caso di q* la struttura è proprio quella che si trasformerà in Q-Learning. Qui il massimo rende l'equazione non lineare: non si può risolvere con un sistema lineare, ma si può usare lo stesso trucco dell'aggiornamento ripetuto, chiamato value iteration.

def value_iteration(env, gamma=1.0, tol=1e-10):
    """Bellman OPTIMALITY equation as an update:  V(s) <- max_a sum_{s'} P (r + gamma V(s'))."""
    V = np.zeros(env.n_states)
    for sweep in range(1, 10_000):
        delta = 0.0
        for s in range(env.n_states):
            q = [sum(p * (r + gamma * V[s2] * (not done)) for p, s2, r, done in env.P[s][a]) for a in range(env.n_actions)]
            delta = max(delta, abs(max(q) - V[s]))
            V[s] = max(q)
        if delta < tol:
            return V, sweep
    return V, sweep

def greedy_policy(env, V, gamma=1.0):
    arrows = "^>v<"
    out = []
    for s in range(env.n_states):
        if s in env.terminals:
            out.append("G"); continue
        q = [sum(p * (r + gamma * V[s2] * (not done)) for p, s2, r, done in env.P[s][a]) for a in range(env.n_actions)]
        best = np.flatnonzero(np.isclose(q, max(q)))                 # all the actions that tie for best
        out.append("".join(arrows[a] for a in best))
    return out

V_star, sweeps = value_iteration(env)
print(f"value iteration converged in {sweeps} sweeps. v*:")
print(np.round(V_star.reshape(4, 4), 1))
pol = greedy_policy(env, V_star)
for r in range(4):
    print("  ".join(f"{pol[r * 4 + c]:>4s}" for c in range(4)))

# Output (due esecuzioni consecutive, stessa uscita):
# value iteration converged in 4 sweeps. v*:
# [[ 0. -1. -2. -3.]
#  [-1. -2. -3. -2.]
#  [-2. -3. -2. -1.]
#  [-3. -2. -1.  0.]]
#    G     <     <    v<
#    ^    ^<  ^>v<     v
#    ^  ^>v<    >v     v
#   ^>     >     >     G

value_iteration ha la stessa struttura della valutazione, con max(q) al posto della media pesata sulla politica. Converge in 4 sweep (contro 272): massimizzando si propaga subito l'informazione «la strada più breve» invece di aspettare la camminata casuale. Il risultato v* è, cambiato di segno, la distanza in passi dal terminale più vicino: −1 per le celle adiacenti a un angolo, −2, −3 per le più lontane. Le frecce mostrano la politica avida rispetto a v*: greedy_policy calcola, per ogni stato, il valore di ogni azione q e prende quelle che massimizzano (np.isclose per trovare i pareggi). Dove più azioni sono equivalenti vengono mostrate tutte: lo stato 6 ha quattro frecce perché ogni direzione porta in una cella da −2 e quindi dà lo stesso valore −3. Questa politica è ottima per questa griglia deterministica.

Perché non basta: i limiti di questo approccio

Il primo limite è che serve il modello: occorre conoscere P(s′|s,a) e R, cioè l'intera tabella, e nei problemi reali (un videogioco, un robot) non la si ha. Il secondo è che bisogna scorrere tutti gli stati: con 16 stati è istantaneo, con gli stati di una partita a scacchi o di un videogioco la tabella non entra in memoria. Il terzo è che questa è pianificazione e non apprendimento: non c'è esperienza, si ragiona sul modello. Queste tre limitazioni definiscono la strada: se non si conosce P, si può stimare l'aggiornamento di Bellman da campioni di esperienza, ed è il Q-Learning; se gli stati sono troppi per una tabella, si usa una rete neurale per approssimare q, ed è la Deep Q-Learning.

I compromessi

Su γ: alto significa orizzonte lungo e comportamento lungimirante, ma calcolo più lento (231 sweep a 0.99 contro 23 a 0.5) e apprendimento più sensibile al rumore; basso significa miope ma veloce. Fra valutazione iterativa e sistema lineare: il sistema lineare è esatto e rapido per decine o poche migliaia di stati (il costo di una risoluzione diretta cresce con il cubo del numero di stati), mentre l'iterazione scala meglio e si adatta a problemi grandi. E fra pianificare e imparare: con il modello si fa programmazione dinamica, senza si impara dai campioni.

Trappole reali

1. Dimenticare il termine «not done». Se il valore del successore viene aggiunto anche dopo uno stato terminale, l'agente «vede» futuro dove non c'è (e valori gonfiati).
2. Scambiare γ con un parametro di apprendimento. γ è una proprietà del problema (quanto guarda lontano l'obiettivo), non un iperparametro da ottimizzare sulla ricompensa ottenuta: cambiandolo cambi il problema.
3. Usare γ = 1 su problemi senza fine. Il ritorno diverge. Con episodi finiti e raggiungibili è lecito.
4. Confondere v e q. v(s) media sulle azioni secondo la politica; q(s, a) fissa la prima azione.
5. Verificare un'implementazione senza il residuo di Bellman. È il test più semplice: calcola max |r + γ P v − v|; se non è vicino a zero, l'aggiornamento è sbagliato.

Le prime due trappole producono valori numericamente plausibili e sbagliati: un futuro in più dopo il terminale gonfia le stime, un γ scelto per far salire la ricompensa cambia il problema che stai risolvendo. La quinta è il rimedio generale: il residuo di Bellman è un controllo che non dipende dall'algoritmo che ha prodotto v, quindi lo puoi applicare a qualunque implementazione, comprese quelle degli articoli successivi.

Esercizi con soluzione

1. Uno stato ha due azioni equiprobabili. La prima dà ricompensa 2 e porta in uno stato di valore 10, la seconda dà −1 e porta in uno di valore 20. Con γ = 0.5, quanto vale lo stato?
   Soluzione: 0.5·(2 + 0.5·10) + 0.5·(−1 + 0.5·20) = 0.5·7 + 0.5·9 = 8.
2. Con γ = 0.95, quanto vale la somma di ricompense costanti pari a 1?
   Soluzione: 1/(1 − 0.95) = 20.
3. In che passo della derivazione serve la proprietà di Markov?
   Soluzione: nel passo 5: dire che E[G_{t+1} | S_t=s, A_t=a, S_{t+1}=s′] = v_π(s′), cioè che il passato non conta una volta noto lo stato attuale.
4. Perché la value iteration (4 sweep) converge molto più in fretta della valutazione della politica casuale (272)?
   Soluzione: perché il massimo sceglie subito l'azione verso il traguardo; la politica casuale vaga e il suo valore dipende da camminate lunghe, che richiedono molti sweep per propagarsi.

I primi due esercizi si risolvono con la formula e un'addizione: servono a controllare di aver capito la media ponderata e la serie geometrica. Gli ultimi due sono di interpretazione e sono quelli che contano: se sai dire dove entra Markov e perché il massimo accelera la propagazione, hai capito perché l'equazione funziona, non solo come si usa.

Domande frequenti

Che differenza c'è fra v e q?

v_π(s) è il valore di uno stato seguendo la politica; q_π(s, a) è il valore di fare una specifica azione nello stato e poi seguire la politica. Il legame è v_π(s) = Σ_a π(a|s) q_π(s, a).

Serve davvero conoscere il modello per usare Bellman?

Per calcolare il valore come in questo articolo sì, perché serve P. Quando P non è noto si stima l'aggiornamento da campioni di esperienza, ed è ciò che fa il Q-Learning.

Perché la value iteration converge in 4 sweep e la valutazione in 272?

Perché la value iteration massimizza e propaga subito la strada più breve, mentre la politica casuale vaga e il suo valore dipende da camminate lunghe.

Che cosa sai fare adesso

Sai derivare l'equazione di Bellman, v_π(s) = Σ_a π(a|s) Σ_{s′} P(s′|s,a)[R + γ v_π(s′)], in sei passi dalla ricorsione del ritorno, e sai dire dove entrano la linearità del valore atteso, la legge della probabilità totale e la proprietà di Markov. Sai spiegare lo sconto con tre letture equivalenti e vedere che controlla sia la risposta (con γ = 0 tutti gli stati valgono −1) sia il costo del calcolo (da 2 a 272 sweep). Sulla griglia, valutazione iterativa e sistema lineare danno la stessa tabella; l'equazione di ottimalità, con un massimo al posto della media, si risolve in 4 sweep e dà v* = −distanza dal traguardo. Tutto questo richiede il modello P. Nel prossimo articolo vediamo che cosa succede quando la ricompensa arriva tardi e quanto è facile scriverla male (Ricompense ritardate, reward shaping e reward hacking); poi impariamo a fare a meno del modello.