  We develop a quantitative theory of stochastic homogenization for linear,
uniformly parabolic equations with coefficients depending on space and time.
Inspired by recent works in the elliptic setting, our analysis is focused on
certain subadditive quantities derived from a variational interpretation of
parabolic equations. These subadditive quantities are intimately connected to
spatial averages of the fluxes and gradients of solutions. We implement a
renormalization-type scheme to obtain an algebraic rate for their convergence,
which is essentially a quantification of the weak convergence of the gradients
and fluxes of solutions to their homogenized limits. As a consequence, we
obtain estimates of the homogenization error for the Cauchy-Dirichlet problem
which are optimal in stochastic integrability. We also develop a higher
regularity theory for solutions of the heterogeneous equation, including a
uniform $C^{0,1}$-type estimate and a Liouville theorem of every finite order.
