We provide the full asymptotic description of the quasi-normal modes (resonances) in any strip of fixed width for Dirac fields in slowly rotating Kerr-Newman-de Sitter black holes. The resonances split in a way similar to the Zeeman effect. The method is based on the extension to Dirac operators of techniques applied by Dyatlov in [Quasi-normal modes and exponential energy decay for the Kerr-de Sitter black hole, Commun. Math. Phys. 306(1) (2011) 119-163; Asymptotic distribution of quasi-normal modes for Kerr-de Sitter black holes. Ann. Henri Poincare 13(5) (2012) 1101-1166] to the (uncharged) Kerr-de Sitter black holes. We show that the mass of the Dirac field does not have an effect on the two leading terms in the expansions of resonances. We give an expansion of the solution of the evolution equation for the Dirac fields in the outer region of the slowly rotating Kerr-Newman-de Sitter black hole which implies the exponential decay of the local energy. Moreover, using the r-normal hyperbolicity of the trapped set and applying the techniques from [Asymptotics of linear waves and resonances with applications to black holes, Commun. Math. Phys. 335 (2015) 1445-1485: Resonance projectors and asymptotics for r-normally hyperbolic trapped sets. J. Amer. Math. Soc. 28 (2015) 311-381]. we give location of the resonance free band and the Weyl-type formula for the resonances in the band near the real axis.