  We present a comprehensive analysis of the dispersion relations for the
doubly-virtual process $\gamma^*\gamma^*\to\pi\pi$. Starting from the
Bardeen-Tung-Tarrach amplitudes, we first derive the kernel functions that
define the system of Roy-Steiner equations for the partial-wave helicity
amplitudes. We then formulate the solution of these partial-wave dispersion
relations in terms of Omn\`es functions, with special attention paid to the
role of subtraction constants as critical for the application to hadronic
light-by-light scattering. In particular, we explain for the first time why for
some amplitudes the standard Muskhelishvili-Omn\`es solution applies, while for
others a modified approach based on their left-hand cut is required unless
subtractions are introduced. In the doubly-virtual case, the analytic structure
of the vector-resonance partial waves then gives rise to anomalous thresholds,
even for space-like virtualities. We develop a strategy to account for these
effects in the numerical solution, illustrated in terms of the $D$-waves in
$\gamma^*\gamma^*\to\pi\pi$, which allows us to predict the doubly-virtual
responses of the $f_2(1270)$ resonance. In general, our results form the basis
for the incorporation of two-meson intermediate states into hadronic
light-by-light scattering beyond the $S$-wave contribution.
