We present a Bayesian analysis to constrain the equation of state of dense nucleonic matter by exploiting the available data from symmetric nuclear matter at saturation and from observations of compact X-ray sources and from the gravitational wave event GW170817. For the first time, such analysis is performed by using a class of models, the relativistic mean field models, which allow to consistently construct an equation of state in a wide range of densities, isospin asymmetries and temperatures. The selected class of models contains five nuclear physics empirical parameters at saturation for which we construct the joint posterior distributions. By exploring different types of priors, we find that the equations of state with the largest evidence are the ones featuring a strong reduction of the effective mass of the nucleons in dense matter which can be interpreted as an indication of a phase transition to a chiral symmetry restored phase. Those equations of state in turn predict $R_{1.4} sim 12$ km. Finally, we present a preliminary investigation on the effect of including $Lambda$ hyperons showing that they appear in stars more massive than about $1.6 M_{odot}$ and lead to radii larger than about $R_{1.4} sim 14$ km. Within the model here explored, the formation of such particles provide a poor agreement with the constraints from GW170817.