Context. Determining photometric redshifts (photo-zs) of extragalactic sources to a high accuracy is paramount to measure distances in wide-field cosmological experiments. With only photometric information at hand, photo-zs are prone to systematic uncertainties in the intervening extinction and the unknown underlying spectral-energy distribution of different astrophysical sources, leading to degeneracies in the modern machine learning algorithm that impacts the level of accuracy for photo-z estimates.

Aims. Here, we aim to resolve these model degeneracies and obtain a clear separation between intrinsic physical properties of astrophysical sources and extrinsic systematics. Furthermore, we aim to have meaningful estimates of the full photo-z probability distribution, and their uncertainties.

Methods. We performed a probabilistic photo-z determination using mixture density networks (MDN). The training data set is composed of optical (griz photometric bands) point-spread-function and model magnitudes and extinction measurements from the SDSS-DR15 and WISE mid-infrared (3.4 mu m and 4.6 mu m) model magnitudes. We used infinite Gaussian mixture models to classify the objects in our data set as stars, galaxies, or quasars, and to determine the number of MDN components to achieve optimal performance.

Results. The fraction of objects that are correctly split into the main classes of stars, galaxies, and quasars is 94%. Furthermore, our method improves the bias of photometric redshift estimation (i.e., the mean Delta z=(z(p)-z(s))/(1+z(s))) by one order of magnitude compared to the SDSS photo-z, and it decreases the fraction of 3 sigma outliers (i.e., 3xrms(Delta z) < Delta z). The relative, root-mean-square systematic uncertainty in our resulting photo-zs is down to 1.7% for benchmark samples of low-redshift galaxies (z(s)

Conclusions. We have demonstrated the feasibility of machine-learning-based methods that produce full probability distributions for photo-z estimates with a performance that is competitive with state-of-the art techniques. Our method can be applied to wide-field surveys where extinction can vary significantly across the sky and with sparse spectroscopic calibration samples. The code is publicly available.