We propose a twofold adjustment for Event Studies considering spatiotemporal data in a multivariate time series framework where the data are characterized by spatial and temporal dependence. The first adjustment consists of modeling the spatiotemporal dynamics of the data by implementing several geostatistical models capable of handling both spatial and temporal components, as well as estimating the relationship between the response variable and a set of exogenous factors. With the second adjustment, we propose to use cross-sectional-adjusted test statistics directly accounting for spatial cross-correlation. The proposed methods are applied to the case of NO22 concentrations observed in Northern Italy during the first wave of the COVID-19 pandemic. The key findings are as follows. First, all the considered geostatistical models estimate larger reductions in the major metropolitan and congested areas, while smaller reductions are estimated in rural plains and in the mountains. Second, the models are nearly equivalent in terms of fitting and are capable of identifying the true event window. Third, by using spatiotemporal models we ensure the residuals are uncorrelated across space and time, thus allowing Event Studies test statistics to provide reliable and realistic estimates. Fourth, as expected, all test statistics show significant reductions in NO22 concentrations starting from the first few days of lockdown. Supplementary materials accompanying this paper appear online.