The Einstein-Hilbert action reads:

$$ S=\int{\sqrt{-g}d^4x\:R} $$

Suppose $\Sigma$ is a 3-dimensional manifold, and suppose our universe is a 4-dimensional oriented manifold $M$ with a smooth function $t$ defined on it, such that each level set $\Sigma_{t}$ of $t$ is a hypersurface that is diffeomorphic to $\Sigma$. We define a smooth vector field $t^{a}$ on $M$ such that $t^{a} \nabla_{a} t = 1$. So far, we have not specified a metric on $M$. Now, suppose we equip $M$ with a metric $g_{ab}$, such that each level set $\Sigma_{t}$ is a spacelike hypersurface, and $t^{a}$ is a future-directed timelike vector field. In this case, $\left\{ \Sigma_{t} \right\}$ and $t^{a}$ define a ADM decomposition of the spacetime $(M, g_{ab})$. The vector field $t^{a}$ can then be expressed as

$$ t^{a} = N n^{a} + N^{a} $$

where $n^{a}$ is the future-directed unit vector normal to $\Sigma_{t}$, and $N$ and $N^{a}$ are known as the lapse function and shift vector field, respectively. Let $h_{ab}$ be the induced metric of $g_{ab}$ on $\Sigma_{t}$, defined by $h_{ab} = g_{ab} + n_{a} n_{b}$. Define

$$ h^{ab} \equiv g^{ab} + n^{a} n^{b},\tag{1} $$

then

$$ g^{ab} = h^{ab} - n^{a} n^{b} = h^{ab} - N^{-2} \left( t^{a} - N^{a} \right) \left( t^{b} - N^{b} \right).\tag{2} $$

Let us first explain how $h^{ab}$ can be determined from $h_{ab}$. Firstly, $h^{ab}$ in equation (1) satisfies

$$ \left( h^{ab} h_{bc} \right) w^{c} = w^{a} \quad \forall w^{c} \in W_{q}\tag{3} $$

however, there are many $h^{ab}$ satisfying this equation. Moreover, $h^{ab}$ should be a spatial tensor, meaning that $h^{ab} \nabla_{b} t = 0$. This condition uniquely specifies $h^{ab}$ among the many $h^{ab}$ that satisfy equation (2).With $h^{ab}$, we can then obtain $N^{a}$ from $N_{a}$ using $N^{a} = h^{ab} N_{b}$. Substituting this into equation (2) gives $g^{ab}$. Now we can rewrite Einstein-Hilbert Lagrangian density, $\mathscr{L} = (-g)^{1/2} R$, in terms of $N$, $N_{a}$, $h_{ab}$, and their temporal and spatial derivatives. Using $G_{ab} = R_{ab} - \frac{1}{2} R g_{ab}$, we have

$$ R = 2 \left( G_{ab} n^{a} n^{b} - R_{ab} n^{a} n^{b} \right) = \left( {}^{3} R - K_{ab} K^{ab} + K^{2} \right) - 2 R_{ab} n^{a} n^{b}\tag{4} $$

The second term on the right-hand side can be computed as follows:

$$ \begin{aligned} n^aR_{ab}n^b&=n^aR_{acb}^{\,\,c}n^b=-n^a\left( \nabla _a\nabla _c-\nabla _c\nabla _a \right) n^c\\ &=-\nabla _a\left( n^a\nabla _cn^c \right) +\left( \nabla _an^a \right) \nabla _cn^c+\nabla _c\left( n^a\nabla _an^c \right) -\left( \nabla _cn^a \right) \nabla _an^c\\ &=-\nabla _a\left( n^a\nabla _cn^c \right) +K^2+\nabla _c\left( n^a\nabla _an^c \right) -K_{ac}K^{ac}\\ \end{aligned} $$

where in the last step we used $K^{ac} = K^{ca}$. Substituting this into equation (4) and then into $\mathscr{L} = \sqrt{-g} R$, and using $\sqrt { - g } = N \sqrt { h }$, we obtain

$$ \mathscr{L} = \sqrt{h} N \left( {}^{3} R + K_{ab} K^{ab} - K^{2} \right) + 2 \sqrt{-g} \left[ \nabla_{a} \left( n^{a} \nabla_{c} n^{c} \right) - \nabla_{c} \left( n^{a} \nabla_{a} n^{c} \right) \right]. $$

When integrating $\mathscr{L}$, the last two terms on the right-hand side can be converted into boundary integrals using Gauss’s theorem. Boundary terms will appear frequently in subsequent derivations. To maintain focus on the main calculation, we will temporarily ignore all boundary terms and return to them later for a unified discussion. Ignoring the boundary terms, the above expression simplifies to

$$ \mathscr{L} = \sqrt{h} N \left( {}^{3} R + K_{ab} K^{ab} - K^{2} \right). $$