Monte Carlo Methods • Part 1 of 3
📝Draft

Monte Carlo Prediction

Estimating value functions by averaging complete returns

In Policy Evaluation, we computed VπV^\pi by sweeping over every state and applying Bellman backups—possible only because we knew the model. Monte Carlo prediction solves the same problem with none of that machinery. The whole algorithm fits in a sentence:

📖Monte Carlo Prediction

To estimate Vπ(s)V^\pi(s), run episodes under policy π\pi, record the return that followed each visit to ss, and average the returns. As the number of visits grows, the average converges to Vπ(s)V^\pi(s).

That works because of how VπV^\pi is defined: the expected return from ss under π\pi. An expectation can always be estimated by a sample mean. Everything in this section is about doing that carefully—and understanding what the estimate costs us in noise.

The Running Example

We’ll use the standard 4×4 GridWorld from Value Functions, so every estimate can be checked against exact values we already computed with DP:

S . . .
. X . .
. . . .
. . . G

S = start (0,0) · X = wall (1,1) · G = goal (3,3)

−1
per step (bumps too)
+10
on entering the goal
0.9
discount γ
0
V(terminal)

Rewards land on the transition: bumping a wall or edge means you stay put but still pay −1, and the final step into G nets 1+10=+9-1 + 10 = +9. The policy we’re evaluating is the uniform random policy—each of the four moves with probability 1/4. From the linear-solve in the Value Functions chapter, the exact answer at the start state is Vπ(start)=8.87V^\pi(\text{start}) = -8.87.

Under this policy the agent needs a while to stumble into the goal: across 10,000 simulated episodes, the mean episode length is 60 steps (median 46, and one unlucky episode took 592).

Averaging Returns, Concretely

Run one random-policy episode from the start and compute its discounted return G0G_0. Then another, and another. Here are the first five returns an actual simulation produced (seed 0), and the running average after each:

−9.94
avg −9.94
−10.00
avg −9.97
−5.17
avg −8.37
−9.72
avg −8.71
−9.99
avg −8.96

Individual returns are all over the place—one lucky episode came home at −5.17, others languished near −10. But after just five episodes the average sits at −8.96, already close to the true −8.87. That’s the law of large numbers at work: averages of independent samples of a random quantity converge to its expectation. Each return is a noisy, unbiased sample of Vπ(start)V^\pi(\text{start}); noise cancels, signal accumulates.

Mathematical Details

Formally, the state value is defined as an expectation over episodes generated by π\pi:

Vπ(s)=Eπ[GtSt=s],Gt=Rt+1+γRt+2+γ2Rt+3+V^\pi(s) = \mathbb{E}_\pi[G_t \mid S_t = s], \qquad G_t = R_{t+1} + \gamma R_{t+2} + \gamma^2 R_{t+3} + \cdots

Given nn observed returns G(1),,G(n)G^{(1)}, \ldots, G^{(n)} following visits to ss, the Monte Carlo estimate is the sample mean:

V^(s)=1ni=1nG(i)\hat{V}(s) = \frac{1}{n} \sum_{i=1}^{n} G^{(i)}

Because each G(i)G^{(i)} is drawn from the distribution whose mean is Vπ(s)V^\pi(s), the estimator is unbiased, and by the law of large numbers V^(s)Vπ(s)\hat{V}(s) \to V^\pi(s) as nn \to \infty. Its standard error shrinks like σ/n\sigma / \sqrt{n}, where σ\sigma is the standard deviation of the return—so quadrupling the data halves the noise.

Note what’s absent: no transition probabilities, no sum over successor states, no reference to other states’ values. Each state’s estimate stands on its own samples.

First-Visit vs Every-Visit

One subtlety: in GridWorld (unlike blackjack), a random walk can pass through the same state several times within one episode. Which of those visits get to contribute a return?

First-visit MC

Per episode, only the return following the first occurrence of ss counts.

One independent sample per episode → a plain i.i.d. average → unbiased. The cleaner theory; most convergence results are stated for this variant.

Every-visit MC

The return following every occurrence of ss counts.

More samples, but returns from the same episode overlap and are correlated → slightly biased in finite samples, yet still consistent: it converges to Vπ(s)V^\pi(s) as episodes accumulate.

In practice the two behave similarly, and every-visit’s extra samples sometimes make it more accurate despite the bias—you’ll see that in the numbers below. First-visit is the standard textbook default, so that’s what we implement.

The Incremental Update

Storing every return and re-averaging would work, but there’s no need. You met the fix in the bandits chapter: the incremental mean.

Keep a count N(s)N(s) of returns seen for state ss. When a new return GG arrives, nudge the estimate toward it:

new estimate = old estimate + step size × (target − old estimate)

with step size 1/N(s)1/N(s). This is exactly the bandit update with the reward RR replaced by the return GG. Bandits averaged immediate rewards per arm; MC averages full returns per state. Same machinery, bigger target.

Mathematical Details

N(s)N(s)+1N(s) \leftarrow N(s) + 1

V(s)V(s)+1N(s)[GV(s)]V(s) \leftarrow V(s) + \frac{1}{N(s)} \big[ G - V(s) \big]

With step size 1/N(s)1/N(s) this reproduces the sample mean exactly. Replacing it with a constant α\alpha gives constant-α MC:

V(s)V(s)+α[GV(s)]V(s) \leftarrow V(s) + \alpha \big[ G - V(s) \big]

which forgets old data exponentially. That’s what you want when the target moves—as it will in the next section, where the policy improves while we evaluate it. Note the shape of this update: target minus estimate, scaled by a step size. Every learning rule in the rest of this book—TD, Q-learning, DQN—has this shape; only the target changes.

Implementation

</>Implementation

Complete, runnable first-visit MC prediction for the random policy on our GridWorld:

import numpy as np
from collections import defaultdict

ROWS, COLS = 4, 4
WALL, GOAL, START = (1, 1), (3, 3), (0, 0)
MOVES = {'up': (-1, 0), 'down': (1, 0), 'left': (0, -1), 'right': (0, 1)}
GAMMA = 0.9


def step(state, action):
    """One environment step. Returns (next_state, reward, done)."""
    dr, dc = MOVES[action]
    nxt = (state[0] + dr, state[1] + dc)
    if not (0 <= nxt[0] < ROWS and 0 <= nxt[1] < COLS) or nxt == WALL:
        nxt = state  # bounce off wall or edge, still pay the step cost
    reward = -1.0 + (10.0 if nxt == GOAL else 0.0)
    return nxt, reward, nxt == GOAL


def generate_episode(policy, rng):
    """Roll out one complete episode from the start state."""
    state, trajectory, done = START, [], False
    while not done:
        action = policy(state, rng)
        next_state, reward, done = step(state, action)
        trajectory.append((state, reward))
        state = next_state
    return trajectory


def random_policy(state, rng):
    return list(MOVES)[rng.integers(4)]


def first_visit_mc(num_episodes, gamma=GAMMA, seed=0):
    """First-visit Monte Carlo prediction for V under the random policy."""
    rng = np.random.default_rng(seed)
    V = defaultdict(float)   # value estimates
    N = defaultdict(int)     # first-visit counts
    for _ in range(num_episodes):
        trajectory = generate_episode(random_policy, rng)
        # Walk backward to compute the return G_t after each state
        G = 0.0
        returns = []
        for state, reward in reversed(trajectory):
            G = reward + gamma * G
            returns.append((state, G))
        returns.reverse()  # time order: (S_0, G_0), (S_1, G_1), ...
        # First-visit: only the return from each state's FIRST occurrence counts
        seen = set()
        for state, G in returns:
            if state not in seen:
                seen.add(state)
                N[state] += 1
                V[state] += (G - V[state]) / N[state]  # incremental mean
    return V


V = first_visit_mc(num_episodes=10_000)
print("First-visit MC estimates after 10,000 episodes:")
for r in range(ROWS):
    print('  '.join(
        '  WALL' if (r, c) == WALL else '  GOAL' if (r, c) == GOAL
        else f'{V[(r, c)]:6.2f}' for c in range(COLS)))

Output:

First-visit MC estimates after 10,000 episodes:
 -8.87   -8.62   -7.78   -7.22
 -8.62    WALL   -6.46   -5.50
 -7.72   -6.35   -4.20   -0.78
 -7.24   -5.61   -0.87    GOAL

The backward pass is the one implementation idea worth internalizing: computing Gt=Rt+1+γGt+1G_t = R_{t+1} + \gamma G_{t+1} from the end of the trajectory makes all returns available in a single sweep, instead of an O(T2)O(T^2) forward recomputation.

For comparison, the exact values from the DP linear solve:

 -8.87  -8.61  -7.75  -7.24
 -8.61   WALL  -6.38  -5.51
 -7.75  -6.38  -4.28  -0.91
 -7.24  -5.51  -0.91   GOAL

Every cell matches to within a few hundredths—learned purely from wandering, with no knowledge of the grid’s dynamics.

How Fast Does It Converge?

Running both MC variants against the exact values, here’s the RMS error over all 14 non-terminal states (single run, seed 0):

EpisodesFirst-visit RMS errorEvery-visit RMS error
1000.5520.381
1,0000.2260.079
10,0000.0560.040

Two things to notice. Error falls roughly like 1/n1/\sqrt{n}, as the math predicted. And every-visit—the “biased” variant—is actually more accurate here, because its extra correlated samples still carry information. Bias and error are not the same thing.

Unbiased, but Noisy

MC’s defining tradeoff: each return is an honest, unbiased sample of the value—but a single return bundles together every random event in the rest of the episode. In our grid, that’s 60 random steps on average. One episode tells you very little; you pay for the missing model with variance.

Here’s the first-visit estimate of Vπ(start)V^\pi(\text{start}) (true value −8.87), with the spread across 30 independent runs:

EpisodesEstimate (seed 0)Std across 30 seeds
10−9.230.52
100−8.690.16
1,000−8.740.07
10,000−8.870.02

No run is systematically high or low—that’s the unbiasedness. But at 10 episodes the run-to-run spread (±0.52) is larger than the gaps between many states’ true values. High variance is the tax MC pays for using complete, model-free returns; the full accounting of that tradeoff against TD’s low-variance-but-biased targets is in TD vs Monte Carlo.

The Classic Example: Blackjack

📌Blackjack (Sutton & Barto, Chapter 5)

The traditional showcase for MC prediction—and the reason the method feels so natural—is blackjack. A state is (your card sum, the dealer’s showing card, whether you hold a usable ace): about 200 states in all. The policy under evaluation is simple, e.g. “stick on 20 or 21, otherwise hit.” Rewards are +1, −1, or 0 at the end of the hand.

Writing down the transition model is miserable: the probability of moving from sum 14 to sum 19 depends on the composition of the remaining deck, and the reward requires marginalizing over the dealer’s entire drawing procedure. But simulating a hand is trivial—deal cards, follow the rule, see who wins. Sutton and Barto estimate the full value function by simulating half a million hands, producing the well-known value surfaces over player sum and dealer card.

That asymmetry—easy to sample, hard to model—is the signature of problems where Monte Carlo methods shine.

What’s Next

We can now measure how good a policy is, using nothing but experience. Measurement isn’t the goal, though—improvement is. The next section closes the loop: use MC estimates to make the policy greedier, and handle the exploration problem that immediately appears. That’s Monte Carlo Control.