We introduce a variational multiscale closure modeling strategy for the numerical stabilization of proper orthogonal decomposition reduced-order models of convection-dominated equations. As a first step, the new model is analyzed and tested for convection dominated convection-diffusion equations. The numerical analysis of the finite element discretization of the model is presented. Numerical tests show the increased numerical accuracy over a standard reduced-order model and illustrate the theoretical convergence rates.