We present Zephyr, a novel method that integrates cutting-edge normalizing flow techniques into a mixture density estimation framework, enabling effective utilization of the heterogeneous training data for photometric redshift inference. Compared to previous methods, Zephyr demonstrates enhanced robustness for both point estimation and distribution reconstruction by leveraging normalizing flows for density estimation and incorporating careful uncertainty quantification. Moreover, Zephyr offers unique interpretability to disentangle contributions from multi-source training data, which can facilitate future weak lensing analysis by providing an additional quality assessment. As probabilistic generative deep learning techniques gain increasing prominence in astronomy, Zephyr may serve as an inspiration for handling miscellaneous dataset issues, achieving good interpretability, and robustly accounting for uncertainties in heterogeneous training data.