In abelian case, the photon propagator under Pauli- Villars regularization goes as $$\frac{1}{k^2}\rightarrow\frac{1}{k^2}- \frac{1}{k^2-{\Lambda^2}}$$ which is basically $1/(k^2- \Lambda^{-2}k^4)$. The presence of the extra term in the propagator can be obtained by supplementing the Lagrangian with an extra term $$\frac{1}{2\Lambda^2}\Box A_\mu\Box A^\mu$$ which using Lorenz gauge can be converted into $$\frac{1}{2\Lambda^2}\partial^\nu F_{\nu\mu}\partial_\sigma F^{\sigma\mu}$$ as $$\partial_\sigma F^{\sigma\mu}=\partial_\sigma(\partial^\sigma A^\mu-\partial^\mu A^\sigma)$$ and using Lorenz gauge second term can be dropped. Now the expression $\frac{1}{2\Lambda^2}\partial^\nu F_{\nu\mu}\partial_\sigma F^{\sigma\mu}$ is U(1) invariant as $F^{\nu\mu}$ is U(1) invariant so Pauli -Villars regularization preserves gauge invariance in abelian case.
In Yang -Mills case, a term like $\frac{1}{2\Lambda^2}\Box A_\mu\Box A^\mu$ violates gauge invariance due to presence of simple derivative instead of covariant derivatives. It is not possible to convert this term by any gauge transformation to a form which has an explicit gauge $covariance$. Gauge invariance can only be obtained by modifying the Lagrangian appropriately for Yang-Mills case which does not follow from Pauli- Villars procedure.
The requirement of covariant derivative in non-Abelian case for gauge invariance comes form the fact that in non-Abelian case the Yang-Mills field strength transform homogeneously under gauge transformation as $F\rightarrow UFU^{-1}$ as opposed to the Abelian case where it is gauge invariant. Covariant derivative has the same transformation law as the field under gauge transformation as $D\rightarrow UDU^{-1}$, so the combination $D_\mu F_{\mu\nu}$ will transform homogeneously which will not be the case if we use a partial derivative instead of covariant derivative. Under a trace term in the Lagrangian which includes a couple of terms like this one can use cyclic property of the trace to show that the term (like the one mentioned below) is unaffected by gauge transformation.
A requirement of gauge invariant term in the Lagrangian will involve higher covariant derivative term like $tr(\frac{1}{\Lambda^2}D_\alpha F_{\mu\nu}D_\alpha F_{\mu\nu})$ .
This modified form of Pauli- Villars regularization was used by A. A. Slavnov in Nucl. Phys. B, 31, 301 (1971).