Aims. We investigate the contribution of shot-noise and sample variance to uncertainties in the cosmological parameter constraints inferred from cluster number counts, in the context of the Euclid survey. Methods. By analysing 1000 Euclid-like light cones, produced with the PINOCCHIO approximate method, we validated the analytical model of Hu & Kravtsov (2003, ApJ, 584, 702) for the covariance matrix, which takes into account both sources of statistical error. Then, we used such a covariance to define the likelihood function that is better equipped to extract cosmological information from cluster number counts at the level of precision that will be reached by the future Euclid photometric catalogs of galaxy clusters. We also studied the impact of the cosmology dependence of the covariance matrix on the parameter constraints. Results. The analytical covariance matrix reproduces the variance measured from simulations within the 10 percent; such a difference has no sizeable effect on the error of cosmological parameter constraints at this level of statistics. Also, we find that the Gaussian likelihood with full covariance is the only model that provides an unbiased inference of cosmological parameters without underestimating the errors, and that the cosmology-dependence of the covariance must be taken into account.