The following general master equation (from this paper 'Dynamical quantum correlations of Ising models on arbitrary lattice and their resilience to decoherence') describes the various types of decoherence relevant to trapped ions, Rydberg atoms and condensed matter systems \begin{equation} \dot{\rho} = i \mathscr{H}(\rho) - \mathscr{L}_{ud}(\rho)-\mathscr{L}_{du}(\rho)-\mathscr{L}_{el}(\rho)~~~~~~~~~~~~~(1) \end{equation} where $$\mathscr{H}(\rho) = [\mathcal{H},\rho],\\ \mathscr{L}_{ud}(\rho) = \frac{\Gamma_{ud}}{2}\sum_{j}(\hat{\sigma}^{+}_j\hat{\sigma}^{-}_{j}\rho + \rho \hat{\sigma}^{+}_{j}\hat{\sigma}^{-}_{j}-2\hat{\sigma}^{-}_{j}\rho \hat{\sigma}^{+}_{j})\\ \mathscr{L}_{du}(\rho) = \frac{\Gamma_{ud}}{2}\sum_{j}(\hat{\sigma}^{+}_j\hat{\sigma}^{-}_{j}\rho + \rho \hat{\sigma}^{+}_{j}\hat{\sigma}^{-}_{j}-2\hat{\sigma}^{-}_{j}\rho \hat{\sigma}^{+}_{j})\\ \mathscr{L}_{el}(\rho) = \frac{\Gamma_{el}}{8}\sum_{j}(2\rho-2\hat{\sigma_j^z}\rho\hat{\sigma}_j^z) $$ The first term involving a commutator describes coherent evolution due to the Ising interaction, and the various terms having subscripts ‘ud’, ‘du’ and ‘el’ correspond respectively to spontaneous relaxation, spontaneous excitation and dephasing. Equation (1) has the formal solution $\rho(t) = \mathscr{U}(t)\rho(0)$, with $$\mathscr{U}(t) = \text{exp}[-t(i\mathscr{H + \mathscr{L}_{ud}+\mathscr{L}_{du} + \mathscr{L}_{el})}]~~~~~~~~~~~~~~~~~(2)$$ An immediate simplification follows from the observation that $$[\mathscr{L}_{el}, \mathscr{H}] = [\mathscr{L}_{el},\mathscr{L}_{du}] = [\mathscr{L}_{el}, \mathscr{L}_{ud}] = 0~~~~~~~~~~~~~~~~~~(3)$$ That the last two commutators vanish is less obvious than the firsts, but physically it has a very clear meaning: spontaneous relaxation/excitation on a site $j$ causes the $j$'th spin to have a well defined value of $\sigma^{z}_j$, and thus to be unentangled with the rest of the system. Since the dephasing jump operator $\hat{\sigma}_{j}^{z}$ changes the relative phase between the states $| \sigma_{j}^{z} = \pm 1 \rangle$, whether spontaneous relaxation/excitation occurs before or after a dephasing event only affects the sign of the overall wave function, which is irrelevant.
As a result, we can write $$\mathscr{U}(t) = e^{-t\mathscr{L}_{el}}e^{-t(i\mathscr{H} + \mathscr{L}_{ud} + \mathscr{L}_{du})}~~~~~~~~~~~~~~(4)$$
Question: I might be missing something simple, but can anyone see why equation (2) reduces to equation (4) as result of the commutativity in equation (3)?