I have came up with this type of equations during solving Quantum field theory problems for finding n-point correlation functions. This equation most probably can be solved via perturbation series method but, Is there any general method to do this kind of problem in WKB with more than one perturbation parameter? $[P_4(z)\dfrac{\mathrm{d}^4}{\mathrm{d}z^4}+P_3(z)\dfrac{\mathrm{d}^3}{\mathrm{d}z^3}+P_2(z)\dfrac{\mathrm{d}^2}{\mathrm{d}z^2}+P_1(z)\dfrac{\mathrm{d}}{\mathrm{d}z}]G(z)+\epsilon_1\Delta^2(z)+\epsilon_2\Delta(z)=0$
Here $P_n(z)$ is the $n$-th degree polynomial and $\Delta(z)$ is some function of z (in this case bare propagator). Range of z is 0 to 1 and $\epsilon_i$'s are perturbation parameters.