A factorization of the matrix elements of the Dyall Hamiltonian in N-electron valence state perturbation theory allowing their evaluation with a computational effort comparable to the one needed for the construction of the third-order reduced density matrix at the most is presented. Thus, the computational bottleneck arising from explicit evaluation of the fourth-order density matrix is avoided. It is also shown that the residual terms arising in the case of an approximate complete active space configuration interaction solution and containing even the fifth-order density matrix for two excitation classes can be evaluated with little additional effort by choosing again a favorable factorization of the corresponding matrix elements. An analogous argument is also provided for avoiding the fourth-order density matrix in complete active space second-order perturbation theory. Practical calculations indicate that such an approach leads to a considerable gain in computational efficiency without any compromise in numerical accuracy or stability.