The Einstein-Hilbert actions is: $$ \begin{align} S=\int_M R(g)\text{dvol}_g=\int_M\sqrt{-g}Rd^4x \end{align} $$ I'm looking to vary the action in a coordinate free way without appealing to the Euler-Lagrange equations. For example, in Yang-Mills theory, we fix a principal bundle $P$ over a Pseudo-Riemannian manifold $M$, with compact structure group $G$. Then, with an ad-invariant inner product on $\mathfrak{g}$, we get an induced bundle metric on $\text{Ad}(P)$, the vector bundle associated to $P$ with fibre $V=\mathfrak{g}$, and the representation the adjoint representation. Then we write the action: $$ \begin{align} S[A]=\int_M\langle F_M^A,F_M^A \rangle_{\text{Ad}(P)}\text{dvol}_g \end{align} $$ Where $F_M^A\in \Omega^2(M,\text{Ad}(P))$, is the curvature two form corresponding to $A$. Since the set of connections is an affine space over $\Omega^1(M,\text{Ad}(P))$, we are able to define critical points as those that satisfy:
$$ \begin{align} \frac{d}{dt}\Big|_{t=0}S[A+t\alpha_m]=0 \end{align} $$ for all $\alpha_M\in \Omega(M,\text{Ad}(P))$. So how do we do something similar for the Einstein Hilbert action? I am mainly confused because I don't think the set of pseudo riemannian metrics on some manifold is an affine space in anyway. Is there any less physicsy text that treats GR this way? i.e. not appealing to some calculation in coordinates?