In our contribution we propose a variational Bayesian projection-based Reduced Order Modeling (ROM) based on a specific type of autoencoders called Variational AutoEncoders (VAEs). The importance and the use of this ROM are demonstrated for fluid dynamics applications based on the resolution of the incompressible and unsteady Navier-Stokes equations. The online stage can be repeated for the required number of unsteady solution samples, thanks to the projection of the residual over the posterior distribution of the random temporal coefficients given by the VAE latent space. We show that the proposed ROM is very practical to use thanks to the sharpness of the mean confidence interval around the high-fidelity solution, the length of this interval being defined by a multiple of the standard deviation of the solution samples.