Metric variation of a scalar-curvature coupling ($E_{ab}=\Phi G_{ab}+{\nabla }_a{\nabla }_b\Phi -g_{ab}\square \Phi$)

With the curvature sign convention in which $\delta R_{bc}={\nabla }_c\delta {\Gamma }^a{}_{ba}-{\nabla }_a\delta {\Gamma }^a{}_{bc}$, the inverse-metric variation of $\int \sqrt{-g}\Phi R$ is $\int \sqrt{-g}{E}_{ab}\delta g^{ab}$ plus a boundary term, for a scalar multiplier held fixed in that variation. The Palatini identity expresses the curvature variation as divergences of $\delta \Gamma$. Two integrations by parts move the derivatives onto $\Phi$, yielding the displayed terms. For constant $\Phi$ they vanish. The opposite curvature convention instead has $E_{ab}=\Phi G_{ab}+g_{ab}\square \Phi -{\nabla }_a{\nabla }_b\Phi$.