The COVID-19 outbreak was a major event that greatly impacted the economy and the health systems around the world. Understanding the behavior of the virus and being able to perform long-term and short-term future predictions of the daily new cases is a working field for machine learning methods and mathematical models. This paper compares Verhulst’s, Gompertz´s, and SIR models from the point of view of their efficiency to describe the behavior of COVID-19 in Spain. These mathematical models are used to predict the future of the pandemic by first solving the corresponding inverse problems to identify the model parameters in each wave separately, using as observed data the daily cases in the past. The posterior distributions of the model parameters are then inferred via the Metropolis–Hastings algorithm, comparing the robustness of each prediction model and making different representations to visualize the results obtained concerning the posterior distribution of the model parameters and their predictions. The knowledge acquired is used to perform predictions about the evolution of both the daily number of infected cases and the total number of cases during each wave. As a main conclusion, predictive models are incomplete without a corresponding uncertainty analysis of the corresponding inverse problem. The invariance of the output (posterior prediction) with respect to the forward predictive model that is used shows that the methodology shown in this paper can be used to adopt decisions in real practice (public health).