Backpropagation Through Time

A Recurrent Neural Network computes a hidden state $h_t$ from the current input $x_t$ and the previous hidden state $h_{t-1}$, using the same parameters $W_{hh}, W_{xh}, b_h$ at every time step:

Because $h_t$ depends on $h_{t-1}$, which depends on $h_{t-2}$, and so on back to $h_0$, the computation graph is not a simple chain of independent layers — it is a single recurrent cell applied $T$ times, with parameters tied across every application. Unrolling means writing out this graph explicitly for time steps $1, \dots, T$ as if it were a $T$-layer feedforward network, where "layer $t$" is the RNN cell instantiated with the actual $x_t$ and connected to layer $t-1$'s output. Once unrolled, the graph is an ordinary computational graph with a clear topological order, and reverse-mode automatic differentiation (backpropagation) applies unchanged. The one subtlety — and it is the entire story of BPTT — is parameter tying. In a normal deep feedforward network, layer $\ell$ has its own weight matrix $W^{(\ell)}$, and $\partial L/\partial W^{(\ell)}$ only involves that layer. In the unrolled RNN, every layer uses the identical $W_{hh}$. Backpropagation still computes a local gradient contribution at each time step, but because $W_{hh}$ is one tensor, not $T$ separate tensors, its total gradient is the sum of the local contributions from every time step:

where the notation $|_{(t)}$ denotes "the gradient contribution computed using $W_{hh}$ as it appears at time step $t$'s copy in the unrolled graph." This sum-over-time is the defining structural feature of BPTT, distinguishing it from vanilla backpropagation in a non-recurrent network. The other defining feature is what happens to the hidden-state gradient (often called $\delta_t = \partial L/\partial h_t$) as it flows backward through time. Because $h_t = \phi(W_{hh}h_{t-1} + \dots)$, the chain rule gives

so propagating error from time step $T$ back to time step $1$ multiplies by this Jacobian $T-1$ times. Repeated multiplication by the same (or similar) matrix is exactly the setting in which linear algebra tells us eigenvalues dominate: if the spectral radius of $\operatorname{diag}(\phi')W_{hh}$ is consistently below 1, the product shrinks geometrically (vanishing gradients); if consistently above 1, it grows geometrically (exploding gradients). This is developed rigorously in the Complexity & Correctness section.

BPTT proceeds in three conceptual phases: unroll, forward pass with caching, backward pass with gradient accumulation. 1. Unroll. Fix a sequence length $T$ (or a truncation window, see below). Instantiate the recurrent cell $T$ times, one per time step, all sharing the same parameters $W_{xh}, W_{hh}, W_{hy}, b_h, b_y$. 2. Forward pass. Run the RNN forward from $h_0$ (usually zero) to $h_T$, computing and caching every intermediate quantity needed for the backward pass: each pre-activation $z_t = W_{hh}h_{t-1} + W_{xh}x_t + b_h$, each hidden state $h_t$, each output $\hat y_t$, and each per-step loss $L_t$ (e.g. cross-entropy against target $y_t$). The total loss is typically $L = \sum_{t=1}^T L_t$. 3. Backward pass. Walk time backward from $t = T$ down to $t = 1$. At each step, combine the gradient flowing in from the local output loss with the gradient flowing in from the future hidden state (i.e., from $t+1$), producing $\delta_t = \partial L / \partial h_t$. Use $\delta_t$ to accumulate gradients into the shared parameter matrices, then propagate $\delta_t$ further backward into $\delta_{t-1}$. The crucial bookkeeping variable is the carried hidden-state gradient, which plays the same role that "gradient flowing into a layer's input" plays in ordinary backprop, except here it also folds in the direct loss gradient at that time step. ```text PSEUDOCODE: Backpropagation Through Time (single sequence, length T) Inputs: x[1..T] # input sequence y[1..T] # target sequence W_xh, W_hh, W_hy # shared weight matrices b_h, b_y # shared biases h0 # initial hidden state (often zero vector)

# ---------- FORWARD PASS (cache everything needed for backward) ---------- h[0] = h0 for t = 1 to T: z[t] = W_xh @ x[t] + W_hh @ h[t-1] + b_h h[t] = phi(z[t]) # e.g. tanh o[t] = W_hy @ h[t] + b_y yhat[t] = psi(o[t]) # e.g. softmax L[t] = loss(yhat[t], y[t]) L_total = sum(L[t] for t in 1..T)

# ---------- BACKWARD PASS ----------

# Gradient accumulators for shared parameters, all initialized to zero dW_xh = 0; dW_hh = 0; dW_hy = 0 db_h = 0; db_y = 0 dh_next = 0 # gradient flowing back from step t+1 into h[t] for t = T downto 1: do_t = d(loss)/d(o[t]) # local output gradient, e.g. yhat[t]-y[t] dW_hy += do_t @ h[t].T db_y += do_t dh_t = W_hy.T @ do_t + dh_next # loss gradient + future-time gradient dz_t = dh_t * phi_prime(z[t]) # elementwise, backprop through phi dW_xh += dz_t @ x[t].T dW_hh += dz_t @ h[t-1].T db_h += dz_t dh_next = W_hh.T @ dz_t # this becomes "dh_next" for step t-1 return dW_xh, dW_hh, dW_hy, db_h, db_y ``` Two lines deserve emphasis because they encode the two structural facts from the previous section. The line `dW_hh += dz_t @ h[t-1].T` (rather than `=`) implements gradient summation over shared parameters. The line `dh_t = W_hy.T @ do_t + dh_next` implements the merging of local loss gradient and carried-forward-in-time gradient — this addition is exactly the multivariate chain rule applied to a node ($h_t$) that has two downstream consumers: the output $\hat y_t$ and the next hidden state $h_{t+1}$. After computing gradients, an ordinary optimizer step (SGD, Adam, etc.) updates $W_{xh}, W_{hh}, W_{hy}, b_h, b_y$ using the accumulated gradients, exactly as in feedforward backprop.

Truncated BPTT. For long or streaming sequences, caching all $T$ hidden states and backpropagating through all of them is often infeasible. Truncated BPTT (TBPTT) runs the forward pass continuously but only backpropagates through the last $k$ time steps (a fixed window), detaching (stopping gradient flow at) the boundary. This trades exactness — gradients no longer reflect dependencies older than $k$ steps — for bounded memory $O(k)$ and bounded per-step compute, which is essential for training on very long or infinite streams. ```text PSEUDOCODE: Truncated BPTT with window k, chunk size k (common variant: k1 = k2 = k) h = h0 for chunk_start = 1, 1+k, 1+2k, ... up to T: chunk = x[chunk_start .. chunk_start+k-1] h = h.detach() # cut gradient flow into earlier chunks run forward pass over chunk, starting from h, caching intermediates compute chunk loss backward pass over just this chunk (as in full BPTT pseudocode above) accumulate/apply parameter gradients h = final hidden state of this chunk # carried forward, but gradient-detached ```

Time complexity. The forward pass performs $O(1)$ matrix–vector work per time step (dominated by $W_{hh}h_{t-1}$, $W_{xh}x_t$, and $W_{hy}h_t$, each costing $O(d^2)$ or $O(d \cdot |x|)$ for hidden size $d$), so the forward pass over a length-$T$ sequence is $O(T d^2)$. The backward pass performs the same amount of matrix work per time step (each line in the backward loop is one matrix–vector product), so it is also $O(T d^2)$. Total: $\Theta(T d^2)$ time, i.e. linear in sequence length, same asymptotic order as the forward pass — this is the general fact that reverse-mode autodiff costs a small constant factor (roughly 2–3×) times the forward pass, regardless of graph shape.

Space complexity. The forward pass must cache every $z_t, h_t, o_t$ (and whatever the elementwise nonlinearity needs to compute its local derivative) because the backward pass visits time steps in reverse order and needs these values. This gives $\Theta(T d)$ memory for full BPTT — memory that grows linearly with sequence length, which is the practical reason long sequences are expensive to train on and the reason truncated BPTT exists (it caps memory at $\Theta(k d)$ regardless of total stream length).

Correctness argument. BPTT's correctness reduces entirely to two facts, both standard in reverse-mode differentiation: 1. Unrolling is a faithful reformulation. The unrolled graph, with $T$ copies of the cell connected in a chain and tied weights, computes exactly the same function of $(x_1,\dots,x_T, h_0, W)$ as the original recurrence — it is the same computation written out explicitly rather than expressed with a self-reference. Therefore any correct differentiation procedure on the unrolled graph yields the correct gradient of the original recurrent function. 2. The multivariate chain rule handles shared parameters and fan-out correctly. When a single scalar (or tensor) parameter feeds into multiple downstream computations — here, $W_{hh}$ feeding into $z_1, z_2, \dots, z_T$, and $h_t$ feeding into both $\hat y_t$ and $h_{t+1}$ — the total derivative is the sum of the derivatives along each path, by the multivariable chain rule:

The pseudocode's `+=` accumulation and the `dh_t = ... + dh_next` line are literal implementations of these two identities. Since reverse-mode autodiff is known to compute exact gradients on any directed acyclic computational graph (a standard result inherited from ordinary backpropagation), and the unrolled RNN graph is exactly such a DAG, BPTT computes the exact gradient of the sequence loss with respect to every shared parameter.

Why gradients vanish or explode. Unrolling the recursive Jacobian relation gives, for the gradient flowing from time $T$ back to time $t$,

This is a product of $T - t$ matrices. If $W_{hh}$ were a scalar $w$ and $\phi'$ were bounded by a constant $c$, this product behaves like $(cw)^{T-t}$: for $|cw| < 1$ it shrinks exponentially in $T - t$ (vanishing gradient — early time steps get essentially zero learning signal from distant losses), and for $|cw| > 1$ it grows exponentially (exploding gradient — parameter updates become enormous and unstable, or overflow to NaN). For matrices, the rigorous statement uses the spectral radius $\rho$ of $\operatorname{diag}(\phi'(z))W_{hh}$ evaluated along the trajectory: if $\rho < 1$ consistently, the largest singular values of the product decay geometrically in $T-t$; if $\rho > 1$, they grow geometrically. With $\tanh$ or sigmoid activations, $\phi'(z) \in (0, 1]$ (and is often much less than 1 away from the origin), which biases the product toward vanishing — this is why plain ("vanilla") RNNs empirically struggle to learn dependencies spanning more than a few dozen time steps, a phenomenon documented well before gated architectures existed.

Tradeoffs, honestly stated.

• Full BPTT is exact but $O(T)$ in both time and memory — infeasible for very long sequences (e.g. streaming audio, book-length text) without truncation.
• Truncated BPTT bounds memory and compute but introduces bias: gradients cannot represent dependencies longer than the truncation window, so the network is structurally prevented from learning them, no matter how much data it sees.
• Gradient clipping (rescaling the gradient norm if it exceeds a threshold) is a cheap, effective mitigation for exploding gradients, but does nothing for vanishing gradients — it only bounds the magnitude of instability, not the loss of long-range signal. This asymmetry is precisely why architectural fixes were needed.

Consider a toy scalar RNN (hidden size 1, so all "matrices" are scalars) with $\phi = \tanh$, no biases, mapped through a linear readout $\hat y_t = w_{hy}h_t$, trained with squared-error loss $L_t = \tfrac12(\hat y_t - y_t)^2$, over $T = 3$ steps. Let

Inputs $x = (1.0,\ 1.0,\ 1.0)$, targets $y = (0.5,\ 0.5,\ 0.5)$.

Forward pass.

Backward pass. Use $\phi'(z_t) = 1 - h_t^2$: $\phi'(z_1)=1-0.2135=0.7865$, $\phi'(z_2)=1-0.6137=0.3863$, $\phi'(z_3)=1-0.7989=0.2011$. Start with $dh_{\text{next}} = 0$ and step $t = 3 \to 1$:

• t = 3: $do_3 = \hat y_3 - y_3 = 0.3938$. $\; dw_{hy} \mathrel{+}= do_3 \cdot h_3 = 0.3938(0.8938)=0.3520$.

$dh_3 = w_{hy}\cdot do_3 + dh_{\text{next}} = 0.3938 + 0 = 0.3938$. $dz_3 = dh_3\cdot \phi'(z_3) = 0.3938(0.2011)=0.0792$. $dw_{xh}\mathrel{+}=dz_3\cdot x_3 = 0.0792$; $\;dw_{hh}\mathrel{+}=dz_3\cdot h_2 = 0.0792(0.7834)=0.0620$. $dh_{\text{next}} \leftarrow w_{hh}\cdot dz_3 = 1.2(0.0792)=0.0951$ (passed to $t=2$).

• t = 2: $do_2 = 0.2834$. $\;dw_{hy}\mathrel{+}=0.2834(0.7834)=0.2220$ → running total $0.3520+0.2220=0.5740$.

$dh_2 = w_{hy}\cdot do_2 + dh_{\text{next}} = 0.2834 + 0.0951 = 0.3785$. $dz_2 = 0.3785(0.3863)=0.1462$. $dw_{xh}\mathrel{+}=0.1462$ → running total $0.0792+0.1462=0.2254$; $\;dw_{hh}\mathrel{+}=0.1462(0.4621)=0.0676$ → running total $0.0620+0.0676=0.1296$. $dh_{\text{next}} \leftarrow 1.2(0.1462)=0.1754$ (passed to $t=1$).

• t = 1: $do_1 = -0.0379$. $\;dw_{hy}\mathrel{+}=-0.0379(0.4621)=-0.0175$ → running total $0.5740-0.0175=0.5565$.

$dh_1 = -0.0379 + 0.1754 = 0.1375$. $dz_1 = 0.1375(0.7865)=0.1082$. $dw_{xh}\mathrel{+}=0.1082$ → final $dw_{xh}=0.2254+0.1082=0.3336$; $\;dw_{hh}\mathrel{+}=0.1082(h_0)=0.1082(0)=0$ → final $dw_{hh}=0.1296+0=0.1296$.

Final accumulated gradients: $\partial L/\partial w_{xh} \approx 0.334$, $\partial L/\partial w_{hh}\approx 0.130$, $\partial L/\partial w_{hy}\approx 0.557$. Notice concretely how the recurrent-weight gradient $dw_{hh}$ is a sum of three terms ($0.0620 + 0.0676 + 0$), one per time step, exactly as the theory predicted — and notice how $dh_{\text{next}}$ carries information from $t=3$ back to $t=1$, shrinking somewhat at each hop because $\phi'(z_t) < 1$ and getting rescaled by $w_{hh}=1.2$; with a smaller $w_{hh}$ or deeper unrolling this shrinkage compounds into the vanishing-gradient regime described above.

The unrolling idea and the general reverse-mode differentiation machinery behind BPTT trace back to the broader development of backpropagation for feedforward networks, but its specific application to recurrent networks — treating time as depth — was articulated explicitly in the mid-1980s and formalized soon after. Paul Werbos described the "backpropagation through time" formulation in the context of dynamical systems, and the technique was subsequently popularized within the connectionist community alongside the broader adoption of backpropagation itself.

BPTT is the training engine behind every sequence model built on recurrence, and its failure modes directly motivated the architectures that dominate sequence modeling today. Language models built on LSTM and GRU cells — used historically in machine translation, speech recognition, and text generation before the rise of attention-based models — are trained with BPTT (or its truncated variant) at their core; the gating mechanisms in those architectures exist specifically to create more stable gradient paths through time, directly addressing the vanishing-gradient analysis derived above. Speech recognition systems that process audio streams frame-by-frame use truncated BPTT to keep training tractable over very long recordings. Time-series forecasting systems in finance and demand prediction, sequence tagging in NLP pipelines, and any online/streaming RNN deployment (e.g. on-device keyboard prediction) rely on the exact mechanics — accumulation over shared weights, caching for the backward pass, and truncation windows — described here. Even in the transformer era, understanding BPTT remains foundational: it is the conceptual bridge that explains why attention mechanisms and gated recurrent units were considered necessary advances, and it remains the training algorithm for any modern architecture that retains explicit recurrent state (e.g. certain streaming or state-space sequence models).

Ask the tutor: How exactly do LSTM gates modify the recurrent Jacobian to prevent gradients from vanishing the way they do in a vanilla RNN?