We propose a multiscale method based on a finite element heterogeneous multiscale method (in space) and the implicit Euler integrator (in time) to solve nonlinear monotone parabolic problems with multiple scales due to spatial heterogeneities varying rapidly at a microscopic scale. The multiscale method approximates the solution at the scale of interest at computational cost independent of the small scale by performing numerical upscaling (coupling of macro and micro finite element methods). Optimal a priori error estimates in the L^2(H^1) and C^0(L^2) norm are derived taking into account the error due to time discretization as well as macro and micro spatial discretizations. Further, we present numerical simulations to illustrate the theoretical error estimates and the applicability of the multiscale method to practical problems.