Assessing whether disc evolution is driven by viscosity or MHD winds is crucial to determine the conditions in which planets form. The rich statistical information encoded in nearby star-forming regions provides an excellent opportunity to probe disc evolution. In this paper, we apply the Bayesian inference tool DEMOS, designed to exploit the statistical properties of disc populations and infer initial conditions and model parameters, to the Lupus disc population. We consider the observed distributions of stellar mass M_⋆, disc mass M (from continuum flux), mass accretion rate Ṁ and CO disc size R_(90%, CO). We constrain all parameters for both the viscous model, with internal photoevaporation, and the MHD winds model, with unprecedented precision. The inferred initial conditions suggest that positive correlations between initial disc properties and M_⋆ must be present. Both models yield a typical R₀ of 3-5 au for a disc around a 0.3 M_⊙ star, with positive stellar-mass dependence. We infer α (or αDW) of 3-4 × 10⁻⁴ for discs around 0.3 M_⊙ stars, with an approximately quadratic dependence on M_⋆. For the viscous-photoevaporative model, the typical wind mass-loss rate is Ṁw ∼ 10⁻⁹ M_⊙ yr⁻¹, scaling with M_⋆ with an exponent of ∼ 2.5. The inferred mean age is 1.6 Myr for the MHD winds model and 0.9 Myr for the viscous-photoevaporative model. Despite a slight preference for the latter, Bayesian evidence shows that both models reproduce the data comparably well, and neither can be ruled out. However, their best-fit synthetic populations predict substantially different evolutionary pathways. Applying DEMOS to multiple populations, especially evolved ones like Upper Scorpius, will allow to place stronger constraints on disc evolution mechanisms.