We study data-driven prediction of coarse-grained dynamics in multiscale PDE systems. Adopting a closure-free operator-learning viewpoint, we apply a linear coarse-graining map and learn a surrogate evolution operator for the resolved field directly from filtered high-fidelity trajectories. Motivated by the Mori--Zwanzig formalism, we propose a spatiotemporal neural operator mapping a resolved history slab on to a resolved future slab. Spatial mixing uses Fourier convolution, while temporal mixing uses a causal kernel operator with position-attention weights on time lags. This temporal kernel supports training and inference with histories on uniform time grids of varying resolution and enables evaluation at arbitrary query times within the prediction window. To improve rollout robustness and suppress nonconservative artifacts, we embed a flux-form inductive bias by parameterizing the windowed update in explicit divergence form. We also provide a data-driven guideline for selecting the memory length via the decorrelation time of a closure-injection diagnostic computed from filtered trajectories. We validate on coarse-grained the viscous Burgers' equation, the Kuramoto--Sivashinsky, and two-dimensional turbulent flows, obtaining stable autoregressive rollouts with improved long-horizon accuracy and statistical fidelity, and robustness to changes in the temporal discretization of the input history.