Abstract We present a comprehensive analysis of the dispersion relations for the doubly-virtual process γ ∗ γ ∗ → ππ. 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ès 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ès 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 γ ∗ γ ∗ → ππ, 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.