Gravity must be understood as a curvature of spacetime rather than space itself because the 1915 general theory of relativity, Einstein's new theory of gravity, is an extension of the 1905 special theory of relativity and the special theory of relativity introduces an inseparable connection between the space and the time and forces us to talk about them in a unified – talk about spacetime.
Space and time have to mix according to special relativity because the theory starts from two postulates, including the absolute constancy of the speed of light in the vacuum, and if space and time were separated, such a constancy would be incompatible with the other postulate, the identical form of the physical laws as seen by an arbitrary inertial observer. It makes no sense to discuss a better, post-Newtonian theory of gravity without taking special relativity into account; the general theory of relativity with its insights about the spacetime curvature is a result of the reconciliation of Newton's gravity and special relativity.
In fact, when one studies how Newton's approximate (inverse square) laws of gravity emerge from general relativity, it turns out that the "curvature of time", and not so much "curvature of space", as a function of space plays the decisive role in determining the gravitational fields at each point. Technically speaking, the rate of time at a given point is determined by $g_{00}$ which is approximately a linear function of the gravitational potential $\Phi$ known from Newton's theory.
Conceptual ideas may precede the mathematical formulation of some principles but one usually can't get too far if he avoids mathematics. Well over 99% of important insights in modern physics depend on mathematical equations and structures that may be at most translated to "awkward and confusing" words.
I don't think the way you take the covariant derivatives is correct.
I will show that, in fact, there is no time derivative of the lapse in the trace of extrinsic curvature.
We'll be using the four-dimensional point of view, in which fields living on $\mathbb{R} \times \Sigma$ are treated as defined on the whole spacetime manifold $\mathcal{M}$.
The extensions of tensor fields on $\Sigma$ to $\mathcal{M}$ we call spatial, as their contractions with the normal vector field $\mathbf{n}$ vanish, and the action of the orthogonal projector $\gamma^{\mu}_{\;\;\nu}$ on these fields is an identity. For example, $\beta$ is defined to be a spatial vector field that quantifies the difference between normal observers and coordinate time observers:
$$ \partial_{t} = \alpha n + \beta, \\[0.5em]
g(\beta, n) = 0.$$
From the latter, knowing that $n_{\mu}=(-\alpha, 0,0,0)$, we can derive:
$$ 0= \beta^{\mu}n_{\mu} = \beta^{t}n_{t} + \beta^{i}n_{i}= \beta^{t}n_{t}.$$
And therefore the temporal component of the shift vector is identically zero $$ \beta^{t}=0.$$
CONTENT OF THE EDIT:
Contract $K_{\mu\nu}=-\nabla_{\nu}n_{\mu} -a_{\mu}n_{\nu}$ with the inverse 4-metric:
$$ K=K_{\mu\nu}g^{\mu\nu} = -g^{\mu\nu}\nabla_{\nu}n_{\mu} = \nabla_{\mu}n^{\mu} = -\partial_{\mu}n^{\mu} - \Gamma^{\mu}_{\mu\lambda}n^{\lambda}.$$
Here $a_{\mu}n_{\nu}g^{\mu\nu}=0$, since $a_{\mu}=n^{\lambda}\nabla_{\lambda}n^{\mu}.$
Now let us introduce the notation: $g$, $\gamma$ are the determinants of the 4 and 3-metrics, respectively, satisfying the relation $\sqrt{-g}=\alpha \sqrt{\gamma}$ (ref 1, eq. 4.55) and we will use the fact that:
$$\Gamma^{\mu}_{\mu\lambda}=\frac{1}{2}g^{\mu\nu}g_{\mu\nu,\lambda}=\frac{1}{\sqrt{-g}}\frac{\partial \sqrt{-g}}{\partial x^{\lambda}}=\frac{1}{\alpha\sqrt{\gamma}}\frac{\partial (\alpha\sqrt{\gamma})}{\partial x^{\lambda}}, $$
as well as:
$$ n^{\mu}=(\frac{1}{\alpha}, \frac{-\beta^{i}}{\alpha})$$
Expanding the expression of interest:
$$-\partial_{\mu}n^{\mu} - \Gamma^{\mu}_{\mu\lambda}n^{\lambda} = -\partial_{t}n^{t} - \partial_{i}n^{i} - \frac{1}{\alpha\sqrt{\gamma}}\partial_{\lambda}(\alpha\sqrt{\gamma}) n^{\lambda} =\\
\frac{1}{\alpha^{2}}\partial_{t}\alpha - \partial_{i}n^{i} -\frac{n^{\lambda}}{\sqrt{\gamma}}\partial_{\lambda}\sqrt{\gamma} -\frac{1}{\alpha}n^{\lambda}\partial_{\lambda}\alpha = \\
\frac{1}{\alpha^{2}}\partial_{t}\alpha - \partial_{i}n^{i} - \frac{n^{\lambda}}{\sqrt{\gamma}}\partial_{\lambda}\sqrt{\gamma} - \frac{1}{\alpha}n^{t}\partial_{t}\alpha - \frac{1}{\alpha}n^{i}\partial_{i}\alpha=\\
\frac{1}{\alpha^{2}}\partial_{t}\alpha - \partial_{i}n^{i} - \frac{n^{\lambda}}{\sqrt{\gamma}}\partial_{\lambda}\sqrt{\gamma} - \frac{1}{\alpha^{2}}\partial_{t}\alpha - \frac{1}{\alpha}n^{i}\partial_{i}\alpha$$
As you see, the time derivatives of the lapse cancel out.
Going further:
$$
-\partial_{i}n^{i} - \frac{n^{\lambda}}{\sqrt{\gamma}}\partial_{\lambda}\sqrt{\gamma} - \frac{1}{\alpha}n^{i}\partial_{i}\alpha =\\
\partial_{i}(\frac{\beta^{i}}{\alpha}) - \frac{n^{t}}{\sqrt{\gamma}}\partial_{t}\sqrt{\gamma} - \frac{n^{i}}{\sqrt{\gamma}}\partial_{i}\sqrt{\gamma} + \frac{1}{\alpha^{2}}\beta^{i}\partial_{i}\alpha=\\
\frac{1}{\alpha}\partial_{i}\beta^{i} - \frac{1}{\alpha^{2}}\beta^{i}\partial_{i}\alpha - \frac{1}{\alpha\sqrt{\gamma}}\partial_{t}\sqrt{\gamma} + \frac{\beta^{i}}{\alpha\sqrt{\gamma}}\partial_{i}\sqrt{\gamma} + \frac{1}{\alpha^{2}}\beta^{i}\partial_{i}\alpha=\\
\frac{-1}{\alpha}\Big{[} \frac{1}{\sqrt{\gamma}}\partial_{t}\sqrt{\gamma} - \frac{1}{\sqrt{\gamma}}\partial_{i}(\sqrt{\gamma}\beta^{i})\Big{]}=\\
\frac{-1}{\alpha}\Big{[} \frac{1}{\sqrt{\gamma}}\partial_{t}\sqrt{\gamma} - D_{i}\beta^{i} \Big{]}$$
Where the last expression involving the shift vector is just an alternative formula for the divergence of a vector field.
We can get to the expression we have just derived by using an an alternative (but equivalent) formula for the extrinsic curvature tensor.
$$ K_{ij} = \frac{-1}{2\alpha}(\mathcal{L}_{m}\gamma)_{ij}=\frac{-1}{2\alpha}\big{[}\partial_{t}\gamma_{ij} - (\mathcal{L}_{\beta}\gamma)_{ij} \big{]}.$$
Upon contraction and using $(\mathcal{L}_{\beta}\gamma)_{ij}= D_{i}\beta_{j} + D_{j}\beta_{j}$:
$$ K=K_{ij}\gamma^{ij} = \frac{-1}{\alpha}\Big{(} \frac{1}{\sqrt{\gamma}} \frac{\partial \sqrt{\gamma}}{\partial t} - D_{i}\beta^{i}\Big{)} .$$
An easier derivation could be made from the start, based on the fact that $\mathbf{K}$ is a spatial tensor, and that $g^{\mu\nu}= \gamma^{\mu\nu} + n^{\mu}n^{\nu}$.
Contracting as before:
$$ K=K_{\mu\nu}g^{\mu\nu} = K_{\mu\nu}\gamma^{\mu\nu}=K_{ij}\gamma^{ij} = -\gamma^{ij}\nabla_{j}n_{i}$$
and once again, no time derivative of the lapse is present in the end.
In the above we use the fact that $\gamma^{t\mu}=0$ is zero, as it is a spatial tensor with an upper temporal index $\gamma^{\mu\nu}n_{\nu}=0$. Additionally, we use the fact that $n_{i}=0$ (but not the lowered-index version covariant derivative of its covariant derivative, $(\nabla_{j}n)_{i}$, mind you!
This $\gamma$ is, by the way, the extension of the 3-metric $\gamma_{ij}$ to a spatial tensor field in 4-dimensional spacetime.
As for the derivations of different slicings:
to obtain maximal slicing prescription for the lapse, you should contract the evolution equation for $K_{\mu\nu}$ and obtain an evolution equation for its trace; then, you enforce for all times $t$ the vanishing of $K$, i.e. $K=0$ and $\partial_{t}K=0$, thereby obtaining an elliptic equation for the lapse (ref 1, eq. 4.63, 4.64, 6.90, 9.15); the maximal slicing prescription requires that at each time step, the elliptic equation for the lapse on the hypersurface is solved for before we evolve further
to obtain the harmonic slicing evolution equation, we demand that the time coordinate of our foliation is a harmonic function: $\Box_{g}t=0$. (ref 1, section 9.2.3); in there, the time derivative of lapse appears in a natural way and the derivation answers your question regarding the minus sign, I think; you can generalize this prescription and obtain a family of slicings based on different choices of the $f(\alpha)$ function
The arbitrariness in either of these choices is now evident:
for maximal slicing, you choose that your foliation has $K=0$ everywhere, and you use a no-longer-an-evolution-equation to fix your lapse; inverting the reasoning, you pick such a lapse that will guarantee $K=0$ for each hypersurface
for harmonic slicing, the choice that the time coordinate of the slicing is a harmonic function is totally arbitrary
For details I refer you to the 3+1 Formalism and Bases of Numerical Relativity, specifically Chapter 9 for details on the lapse conditions.
DESCRIPTION OF THE EDIT: I owe you an apology - the answer was not correct. In fact, both parts of the full expression (the one with the Christoffel symbol as well) contribute to the absence of the temporal derivative of the lapse. From one part we get $\sim \partial_{t}\alpha$ and from the other $\sim - \partial_{t}\alpha$, so these cancel out.
In particular, I was wrong stating that $g^{t \mu} = \frac{\beta^{\mu}}{\alpha^{2}}$. This concerns the $g^{ti}$, but not $g^{tt}$. In the quoted part, you can even see that $g^{tt}=\frac{-1}{\alpha^{2}}$.
The fact that $\beta^{t}=0$ is of course true, but does not come into play in this calculation in the end. Also, one must be careful. The orthogonality of a tensor field does not mean that its temporal component vanishes. When the one of the indices of a spatial tensor field component is $t$, the component will be zero by virtue of this contraction: $T_{\;\;\;\mu\nu}^{\alpha\beta}n_{\alpha}=0$ (since $n_{i}=0$). However, for a covariant index, this is not longer the case, easily visible even in the case of $\beta$.
$$ 0 = \beta_{\mu}n^{\mu} = \beta_{t}n^{t} + \beta_{i}n^{i}$$
which leads to:
$$\beta_{t}=g_{t\mu}\beta^{\mu}=g_{ti}\beta^{i} = \beta_{i}\beta^{i}$$
Best Answer
Schwarzschild geometry in Schwarzschild coordinates $(t,r,\theta,\phi)$ is time-symmetric \begin{equation} ds^2=-\left(1-\frac{2GM}{c^2 r}\right)c^2dt^2+\left(1-\frac{2GM}{c^2 r}\right)^{-1}dr^2+r^2\left(d \theta^2 +\sin^2 \theta d \phi^2\right) \;. \end{equation}
Novikov coordinate system is defined by a set of geodesic clocks. The coordinate clocks are freely falling from some maximal radius $r_m$ towards $r=0$, where $r_m$ is different for each clock. All clocks start falling at the same Schwarzschild time $t_0$ and they are synchronized in such a manner that each clock shows $0$ at $r_m$. Novikov coordinate is defined to stay constant along the trajectory of each clock, while for time coordinate proper time is taken.
From now on the angular part metric will be omitted, since is stays the same. We also take $r_s=2M$ and $G=c=1$: \begin{equation}\label{eq:sch-met2} ds^2=-\left(1-\frac{r_s}{r}\right)dt^2+\left(1-\frac{r_s}{r}\right)^{-1}dr^2 \;. \end{equation}
Geodesics in Schwarzschild gometry
To get the equation of geodesics in Schwarzschild geometry we have to solve equations of motion of a free particle: \begin{equation}\label{eq:lagrangian} \mathcal{L}=\frac{1}{2}mg_{\mu\nu}\dot{x}^\mu\dot{x}^\nu \;, \end{equation} \begin{equation}\label{eq:dot} \dot{x}^\mu=\frac{dx^\mu}{d\tau}=u^\mu \;. \end{equation} \begin{equation}\label{eq:lagrangian2} \mathcal{L}=-\frac{m}{2}\left(1-\frac{r_s}{r}\right)\dot{t}^2+\left(1-\frac{r_s}{r}\right)^{-1}\dot{r}^2 \;, \end{equation} \begin{equation}\label{eq:EL} \frac{d}{d\tau}\frac{\partial\mathcal{L}}{\partial \dot{x}^\mu}-\frac{\partial\mathcal{L}}{\partial x^\mu}=0 \;, \end{equation} For $\mu=0$ we get a constant of motion \begin{equation}\label{eq:ConstOfMotion} \frac{\partial}{\partial\tau}\left[\left(1-\frac{r_s}{r}\right)\dot{t}\right]=0 \qquad \Rightarrow \qquad \left(1-\frac{r_s}{r}\right)\dot{t}=a \;, \end{equation}
For timelike geodesics: $ds^2=-d\tau^2$ the radial geodesic equation becomes \begin{equation}\label{eq:orbit} \left(\frac{d\tau}{dr}\right)^2=\frac{1}{a^2-\left(1-\frac{r_s}{r}\right)} \;. \end{equation} Maximal radius is ($dr/d\tau=0$) \begin{equation}\label{eq:maximal} r_m=\frac{r_s}{1-a^2} \;. \end{equation} We use $\frac{dt}{dr}=\frac{dt}{d\tau}\frac{d\tau}{dr}$ and obtain the following relations: \begin{eqnarray} \frac{d\tau}{dr} &=& \frac{\varepsilon}{\sqrt{\frac{r_s}{r}-\frac{r_s}{r_m}}} \;,\label{eq:orbit1} \ \frac{dt}{dr} &=& \frac{\varepsilon\sqrt{1-\frac{r_s}{r_m}}}{\left(1-\frac{r_s}{r}\right)\sqrt{\frac{r_s}{r}-\frac{r_s}{r_m}}} \;, \label{eq:orbit2} \end{eqnarray} where $\varepsilon$ is $+1$ or $-1$. For falling particles we choose $\varepsilon=-1$.
Novikov time coordinate
We first transform from $(r,t)$ to $(r,\tau)$. From last two equations we obtain for $d\tau(dt,dr)$ \begin{equation}\label{eq:bit} d\tau=\left(1-\frac{r_s}{r_m}\right)^{1/2}dt+\frac{\left(\frac{r_s}{r}-\frac{r_s}{r_m}\right)^{1/2}}{1-\frac{r_s}{r}}dr \;. \end{equation} where we assumed $t$ are $r$ known.
This can be integratet from $r$ to $r_m$, where we take into account that all clocks reach their maximum radius at $\tau_{0i}=0$. It follows \begin{equation}\label{eq:integral} \tau=\left(1-\frac{r_s}{r_m}\right)^{1/2}(t-t_0)+\int_{r_m}^{r}\frac{\left(\frac{r_s}{y}-\frac{r_s}{r_m}\right)^{1/2}}{1-\frac{r_s}{y}}dy \;. \end{equation}
maximal radius $r_m$ is here a function of $r$ and \tau$. Their implicit relationship is
\begin{equation}\label{eq:implicit1} \tau=-f(r,r_m)\;, \end{equation} where \begin{equation} f(r,r_m) = \int_{r_m}^{r}\frac{dy}{\sqrt{\frac{r_s}{y}-\frac{r_s}{r_m}}} \label{eq:integral3} = -\left[\frac{rr_m}{r_s}(r_m-r)\right]^{1/2}-\frac{r_m^{3/2}}{\sqrt{r_s}}\arccos\left[\left(\frac{r}{r_m}\right)^{1/2}\right] \;.\label{eq:f} \end{equation}
We can now eliminate coordinate $t$ from the line element \begin{equation}\label{eq:sch-met3} ds^2=-d\tau^2+\frac{1}{1-\frac{r_s}{r_m}}\left[- dr-\left(\frac{r_s}{r}-\frac{r_s}{r_m}\right)^{1/2}d\tau\right]^2 \;. \end{equation}
Novikov radial coordinate
For radial coordinate we take the maximal Schwarzschild radius $r_m$, which remains constant along the worldline of a geodesic clock. \begin{equation}\label{eq:relation2} - dr-\left(\frac{r_s}{r}-\frac{r_s}{r_m}\right)^{1/2}d\tau=\left(\frac{r_s}{r}-\frac{r_s}{r_m}\right)^{1/2}\frac{\partial f}{\partial r_m}dr_m \;. \end{equation} With this we can eliminate the other Schwrazschild coordinate $r$: \begin{equation}\label{eq:sch-met4} ds^2=-d\tau^2+\frac{\left[g(r,r_m)\right]^2}{1-\frac{r_s}{r_m}}dr_m^2 \;. \end{equation} Here we $g(r,r_m)$ is the following \begin{eqnarray} g(r,r_m)&=&-\left(\frac{r_s}{r}-\frac{r_s}{r_m}\right)^{1/2}\frac{\partial f}{\partial r_m} \label{eq:g}\ &=&1+\frac{1}{2}\left(1-\frac{r}{r_m}\right)-\frac{3}{4}\left(\frac{r_m}{r}-1\right)^{1/2}\left[\sin^{-1}\left(\frac{2r}{r_m}-1\right)-\frac{\pi}{2}\right] \;. \nonumber \end{eqnarray} $r$ is not a radial coordinate anymore, but a metric function of coordinates $r_m$ and $\tau$, which is given implicitly by equation ().
Novikov metric
By introducting $r_m$ the metric became diagonal as in Schwarzschild coordinates. It also stays diagonal by introducig a new radial coordinate, that is only functionally related to the old one. Novikov-s choice is $r^*$ with the following monotonic relation to $r_m$: \begin{equation}\label{eq:r*} r^*=\left(\frac{r_m}{r_s}-1\right)^{1/2}\;. \end{equation} The metric now becomes \begin{equation}\label{eq:novikov-met} ds^2=-d\tau^2+4r_s^2\left({r^*}^2+1\right)\left[g(r,r^*)\right]^2d{r^*}^2 \;. \end{equation} We can show that the following also holds \begin{equation}\label{eq:relation} 4Mg(r,r^*)=\frac{1}{r^*}\frac{\partial r}{\partial r^*}\;. \end{equation} With this the metric gets the form that is standard in literature [MTW,p. 826]: \begin{equation}\label{eq:novikov} ds^2=-d\tau^2+\left(\frac{{r^*}^2+1}{{r^*}^2}\right)\left(\frac{\partial r}{\partial r^*}\right)^2d{r^*}^2+r^2\left(d \theta^2 +\sin^2 \theta d \phi^2\right) \;, \end{equation}
where we also included the angular part.
Relations among coordinates
We now give the relations between Schwarzschild coordinates $(t,r)$ and Novikov coordinates $(\tau,r^*)$. The first one, $r=(\tau,r^*)$, is obtained from equations () and () \begin{equation}\label{eq:CoordRela1} \tau=r_s\left({r^*}^2+1\right)\left[\frac{r}{r_s}-\frac{(r/r_s)^2}{{r^*}^2+1}\right]^{1/2}+r_s\left({r^*}^2+1\right)^{3/2}\arccos\left[\left(\frac{r/r_s}{{r^*}^2+1}\right)^{1/2}\right] \;. \end{equation} The second one, $t=(\tau,r^*)$, is obtained by integration from ()
\begin{equation} t=r_s\ln\left|\frac{r^*+\left(\frac{r_s}{r}\left({r^}^2+1\right)-1\right)^{1/2}}{r^-\left(\frac{r_s}{r}\left({r^*}^2+1\right)-1\right)^{1/2}}\right|+r_sr^\left[\left({r^}^2+3\right)\arctan\left(\frac{r_s}{r}\left({r^}^2+1\right)-1\right)^{1/2}+\left({r^}^2+1\right)\frac{\sqrt{\frac{r_s}{r}\left({r^*}^2+1\right)-1}}{\frac{r_s}{r}\left({r^*}^2+1\right)}\right] \;. \end{equation}