Models of non-relativistic quantum gravity: the good, the bad and the healthy

Horava's proposal for non-relativistic quantum gravity introduces a preferred time foliation of space-time which violates the local Lorentz invariance. The foliation is encoded in a dynamical scalar field which we call `khronon'. The dynamics of the khronon field is sensitive to the symmetries and other details of the particular implementations of the proposal. In this paper we examine several consistency issues present in three non-relativistic gravity theories: Horava's projectable theory, the healthy non-projectable extension, and a new extension related to ghost condensation. We find that the only model which is free from instabilities and strong coupling is the non-projectable one. We elaborate on the phenomenology of the latter model including a discussion of the couplings of the khronon to matter. In particular, we obtain the parameters of the post-Newtonian expansion in this model and show that they are compatible with current observations.


Introduction
Recently P. Hořava argued that it may be possible to construct a consistent renormalizable theory of quantum gravity within the framework of perturbative quantum field theory (QFT) [1] (see also [2]). The proposal exploits the improved ultraviolet (UV) behavior of non-relativistic QFTs possessing an UV fixed point with anisotropic scaling of space and time. In this type of theories the UV behavior of the field propagators is improved thanks to terms with higher spatial derivatives. At the same time the number of time derivatives in the Lagrangian remains equal to two allowing to bypass the problems with ghosts arising in Lorentz invariant higher-derivative theories of gravity [3]. Evidently, by breaking the symmetry between space and time one sacrifices Lorentz invariance (LI). The latter is no longer a fundamental property of the theory and may only emerge at low energies as an approximate symmetry. For matter fields in flat space this does not pose an immediate problem, other than the stringent observational constraint that Lorentz invariance needs to emerge to an extremely high accuracy, see [4,5,6,7,8,9] and references therein for the discussion of experimental bounds on Lorentz violating extensions of the Standard Model.
However, abandoning relativistic invariance has dramatic effects for gravity. The reason is that in general relativity (GR) Lorentz symmetry is part of the local gauge group whose role is to remove the unphysical degrees of freedom contained in the metric, leaving only two massless helicity-2 modes for the graviton. Thus one expects that abandoning LI will lead to the appearance of new degrees of freedom in the theory. As already pointed out in [1] this is indeed the case: in addition to the helicity-2 modes the Hořava model propagates a helicity-0 excitation 1 . The physical meaning of the new mode can be understood as follows. From the geometrical point of view the introduction of the preferred time coordinate amounts to equipping the space-time manifold with a foliation by space-like surfaces. In the gravitational theory the foliation inevitably becomes dynamical together with the geometry of the manifold. The helicity-0 mode is nothing but the excitation of this foliation structure. In this sense, the extra mode describes fluctuations of the global time, so we coin for it the name 'khronon' 2 .
Importantly, the new mode does not have a mass gap and thus cannot be consistently decoupled at low energies [11,12]. This implies that, contrary to the original expectation [1], the theory cannot flow to GR at low energies. Instead, one may entertain the possibility that at low energies the theory reduces to a (Lorentz-violating) model of modified gravity, with the modifications being small enough not to contradict the experimental data. The studies of modified gravity models, both in Lorentz-invariant [14,15,16,17] and Lorentzviolating [18,19,20] contexts, have shown that the properties of the extra degrees of freedom can make them fail as phenomenologically acceptable effective field theories (EFT). This is precisely what happens in the original realization [1] of Hořava's proposal where the behavior of the khronon turns out to be pathological [21,22,12,23,24].
It is worth clarifying whether these pathologies are completely generic and invalidate the approach to quantum gravity proposed in [1], or they represent merely a failure of a specific realization of the general framework. The purpose of the present paper is to address this question.
To this end we consider three different models of non-relativistic gravity. First we will reexamine the 'projectable' version of the original proposal [1] outlining its problems associated with the scalar sector. Then we will consider the two extensions of the original 1 The helicity-0 mode is absent in a special case when the model obeys anisotropic Weyl invariance. However, this symmetry is not compatible with observations. For another proposal to eliminate the helicity-0 mode see [10]. 2 From Greek χρoνoς -time. This should not be confused with the term 'chronon' appearing in the sense of "fundamental interval of time" in some theories [13].
proposal suggested in [12] to remedy these problems. The first one is based on a smaller symmetry group and can be viewed as a (power-counting) renormalizable version of ghost condensation [18]. Such a model allows to alleviate the problems of the low-energy EFT. However, we find that this leads to other troubles. The reason is the unavoidable presence in this formulation of a second helicity-0 mode at high energies that leads either to fast instabilities or to the break-down of the perturbative description. Hence, we conclude that this attempt does not offer a promising candidate theory of quantum gravity.
The second extension is the one presented in [25]. In contrast with the two previous models, in this case the scalar sector is free from pathologies. A preliminary study of this healthy extension was reported in [25,26]. Here we analyze it in more detail and confirm that the model is compatible with phenomenological constraints for suitable choices of parameters. This suggests that it may serve as a starting point for constructing a renormalizable theory of quantum gravity.
All the models we will consider contain flat space-time as a consistent background and are naively power-counting renormalizable. However, the behavior of the extra mode(s) is drastically different in the three cases, the differences stemming from the symmetries or other details of the particular realization. Hence, before going to our analysis let us spend a few words concerning the main features of the different possible implementations of non-relativistic quantum gravity. One can distinguish the following: • The choice of the anisotropic scaling to implement power-counting renormalizability.
Namely, one postulates the scaling transformations with a given critical exponent z, together with the scaling weights of the different fields. Then one classifies the operators in the theory according to their dimensions with respect to this scaling. The theory is (power-counting) renormalizable if the action contains only a finite number of terms of zero (marginal operators) and negative (relevant operators) scaling dimensions 3 . The case of a relativistic QFT corresponds to z = 1. In this paper we stick to the choice z = 3. As discussed in [1], this is the minimal value of z that allows to construct a power-counting renormalizable Lagrangian for gravity in (3 + 1) dimensions. Larger values of z are also possible and lead to super-renormalizable models.
• The subgroup of the four-dimensional diffeomorphisms (Diff) under which the theory is invariant. The distinction between space and time enforces a preferred frame, and so a preferred time coordinate. As stated above, this corresponds to endowing the space-time manifold with an additional structure: a preferred foliation by space-like surfaces. In particular this means that the arbitrary reparameterizations of time t →t(t, x) are not an invariance of the theory. Instead, the following unbroken symmetries have been considered in the literature: x →x(t, x) and t →t(t). (1.2) We will refer to these transformations as 'foliation-preserving Diffs', or FDiffs for short. This is the largest possible unbroken gauge group. It is the one originally considered in [1]. (ii) x →x(t, x) and t →t = t + const.
We shall refer to this as the 'restricted-foliation-preserving Diffs' (RFDiffs). This symmetry arises in a number of effective field theories of modified gravity, such as the shift-symmetric k−essence [27] or the ghost condensation [18]. The invariance under time translations implies existence of a conserved energy.
This is the unbroken group in potential-driven inflation and in non-shift-symmetric k−inflation around time-dependent spatially homogeneous solutions. It serves as the basis of the effective field theory of inflation [28,29]. The action in this case contains explicit time dependence and there is no energy conservation.
This option is realized in Einstein-aether theory [30] (see [31] for recent review) and gauged ghost condensation [32]. It leads to the appearance of propagating helicity-1 degrees of freedom (in addition to the helicity-2 and helicity-0 modes).
Here ξ(t) is and arbitrary time-dependent three-vector. This is the symmetry group of Lorentz violating massive gravity [19,33].
In principle it may be possible to construct power-counting renormalizable theories of gravity with any of the above unbroken symmetries. A natural expectation is that the larger the unbroken gauge group, the more constrained the model is, and the fewer degrees of freedom it contains. Investigating all the possibilities is beyond the scope of this article. We shall limit the analysis to the FDiff-invariant and RFDiff-invariant theories. In these cases there are only extra scalar modes, but their number grows when relaxing the symmetry. We will see that the covariant form of these theories involves only a single scalar field (a khronon) in addition to the metric. 4 It seems therefore appropriate to refer to these as "khrono-metric" theories.
• Finally, one may impose additional restrictions on the action that do not follow from a symmetry. Examples of such restrictions are the 'projectability' and the 'detailed balance' conditions of the original paper [1]. In the present work we do not impose the detailed balance condition, and we consider models both with and without the projectability property.
The paper is organized as follows. In Sec. 2 we introduce the basic notations and tools for our analysis. In Sec. 3 we consider the 'projectable' version of the original Hořava model and discuss problems associated with the extra scalar mode. In Sec. 4 we elaborate on an attempt to fix these problems by relaxing the symmetry from FDiffs down to RFDiffs. We show that the model obtained in this way does not provide a consistent candidate for quantum gravity. The relation between this model and the ghost condensation [18] is discussed. In Sec. 5 we return to the case of full FDiff-invariance, now without the projectability condition. In this case the original action of Hořava's model must be supplemented by additional terms [25] allowed by the symmetry and power-counting renormalizability. In Sec. 5.1 we demonstrate that these terms make the scalar graviton stable and weakly coupled. In the rest of Sec. 5 we analyze the phenomenological bounds on the healthy model. In doing this we exploit the analogy [26,34] between the low-energy limit of the model and the Einstein-aether theory [30,31]. In Sec. 6 we summarize our results and discuss future directions. Some details of the analysis are deferred to the Appendices.
Readers interested in the phenomenology of the healthy model may skip Secs. 3, 4 and go directly from Sec. 2 to Sec. 5.

Three theories under scrutiny
In this section we introduce the basic ingredients of the Hořava-type theories that we are going to consider. The field content includes the spatial metric γ ij , the shift N i and the lapse N entering into the (3+1) (ADM) decomposition of the 4-dimensional metric, These fields transform in the standard way under the 4-dimensional coordinate transformations.
We first consider the FDiff-invariant case. Let us also impose the additional requirement that the lapse is 'projectable', i.e. that it does not depend on space coordinates, (2.1) Note that this restriction is compatible with the transformation rules for the lapse under FDiffs, One writes down the following action with two time derivatives [1], Here M P is the Planck mass; λ is a dimensionless constant; K ij is the extrinsic curvature tensor for the surfaces of constant time, 5 K is its trace, K ≡ K ij γ ij ; and the 'potential' V I depends on the spatial metric γ ij via the 3-dimensional Ricci tensor R ij and its spatial covariant derivatives. One notices that the first two terms in (2.2), which comprise the kinetic part of the action, are invariant under the scaling (1.1) with z = 3 provided γ ij , N i and N scale as In other words, the kinetic terms are marginal with respect to the anisotropic scaling with z = 3. A power-counting renormalizable theory is obtained by including in the action all possible marginal and relevant terms. This corresponds to considering the most general potential V I containing local operators with the scaling dimensions up to 6. There is a finite number of these terms, whose classification was performed in [11]. For our purposes it suffices to write schematically, 6 where ξ, A n , B n are dimensionless coupling constants and dots stand for all the possible inequivalent operators of the given dimension. Note that we have introduced here the mass M * for the scale suppressing higher-order operators; this scale may or may not coincide with M P . In what follows we set the parameter ξ to 1 which can always be achieved by a constant rescaling of the time coordinate. In the absence of matter this does not affect the physical content of the model. We will return to the general case ξ = 1 in Sec. 5.4. Omitting higher-derivative terms in the potential and setting λ = 1 one formally recovers the action of GR [1]. 7 However, as we discuss below the limit λ → 1, A n , B n → 0 is discontinuous: instead of flowing to GR the model (2.2), (2.4) becomes pathological in this limit. This model will be analyzed in Sec. 3.

Model II: non-projectable FDiff gravity
Let us now relax the projectability condition (2.1), so that the lapse N is allowed to depend both on time and space. In this case, the object 5 Throughout the paper we use lower-case Latin letters i, j, . . . for 3-dimensional indices. They are raised and lowered using the spatial metrics γij, γ ij . The covariant derivatives carrying these indices are understood accordingly. 6 Recall that in three dimensions the Riemann tensor is completely determined in terms of the Ricci tensor. 7 Still, the resulting theory is different from GR because of the projectability restriction (2.1).
transforms covariantly under FDiffs and has scaling dimension 1. To obtain a (powercounting) renormalizable theory one should allow the potential to depend on a i . Again, this dependence must be limited to local operators of dimension up to 6. Note that a i can be consistently excluded from the kinetic part of the Lagrangian. Indeed, recall that each time derivative raises the dimension of an operator by 3. Thus operators of dimensions less or equal 6 can contain together with a i at most one time derivative. It is straightforward to check that, up to integration by parts in the action, there are three possible combinations (all having dimension 5): These operators are odd under T and CPT transformations. While in general it might be interesting to study the effect of such operators, in the present paper we forbid them by assuming the CPT invariance. Thus we arrive to the following action Note that due to the positive scaling dimension of a i the potential again contains only finite number of terms. The projectable model can be recovered from (2.6) by taking the limit α → ∞ which enforces the spatial gradient of the lapse to vanish, ∂ i N = 0. The potential (2.6) was first proposed in [25]. The model with the action (2.5), (2.6) will be studied in Sec. 5.

Model III: RFDiff gravity
The third model we are going to consider corresponds to relaxing the symmetry from FDiffs to RFDiffs. This has a dramatic effect on the possible structure of interactions in the theory. The reason is that the lapse N is now a scalar under the symmetry group 8 . At the same time its scaling dimension is 0. Therefore, all dimensionless couplings in the Lagrangian may acquire arbitrary dependence on it without spoiling the power-counting renormalizability. Besides, the reduced symmetry allows and power-counting renormalizability requires to include into the action a kinetic term for N . Thus the most general action reads, 8 Note that the reduced gauge invariance is compatible with restricting the lapse N to a fixed value, say N = 1. Then N drops out of the action. The only difference of the resulting theory from the projectable model is the absence of the integral Hamiltonian constraint, so locally the two theories are equivalent (cf. the discussion in the beginning of Sec. 3.2). Thus all results about local dynamics of the projectable model apply to this case.
Note that the form of the kinetic term for N is fixed by the invariance under RFDiffs (1.3). Besides, relaxing the symmetry down to RFDiffs allows to include the standard potential terms V (N ) for the lapse. The presence of arbitrary functions of N in the Lagrangian makes this theory perhaps less attractive, since it involves an infinite number of coupling constants. In spite of this, we find it instructive to explore what happens in this type of extension. We study this model in Sec. 4.
Let us end this section by introducing a few notations that will be common in the analysis of all the above models. We are going to study the dynamics of small metric perturbations in these models above flat background. Since the most worrisome modes are the helicity-0 excitations we concentrate on scalar perturbations of the metric which we parameterize as follows: where ∆ is the spatial Laplacian. Finally, the following mass scales will appear in the analysis: (2.10)

Stückelberg formalism and the khronon field
In our study we will make extensive use of the Stückelberg formalism. This formalism allows to single out explicitly the extra degrees of freedom appearing due to breaking of Diff invariance and study their dynamics in a transparent way. In the present context, it amounts to rewriting the action of the theory in a generally covariant form at the expense of introducing a compensator field which transforms non-homogeneously under the broken part of the 4-dimensional Diffs. It is worth stressing that this procedure does not introduce new degrees of freedom into the theory, but just makes explicit those already present. See [15,19,33] for previous applications of the Stückelberg formalism in modified gravity theories.
In the present paper we follow the approach of [12]. One encodes the foliation structure in a scalar field ϕ(x). Namely, the foliation surfaces are identified with the surfaces of constant ϕ, ϕ = const . (2.11) Clearly, to define a regular foliation the field ϕ must possess non-zero time-like gradient. In this sense ϕ defines an absolute time, hence we call this field 'khronon'. This identification allows a covariant definition of the foliation structure. The actions (2.2), (2.5), (2.7) correspond to the frame where the coordinate time coincides with ϕ, We will refer to this choice of coordinates as the 'unitary gauge'. The action for ϕ in a generic frame is obtained by realizing that the objects appearing in (2.2), (2.5), (2.7) are the standard geometrical quantities that characterize the embedding of the hypersurfaces (2.11) into space-time. The central role in the construction of these quantities is played by the unit normal vector 9 (2.14) Note that u µ is invariant under the reparameterizations of ϕ, where f is an arbitrary monotonic function. This reflects the invariance of the foliation structure under reparameterizations of ϕ. The time-dependent VEV of the khronon field ϕ breaks the product of the reparameterizations (2.15) and general covariance down to the diagonal subgroup which in the unitary gauge coincides with the FDiffs (1.2). One concludes that the covariant form of a FDiff-symmetric theory must be invariant under (2.15). On the other hand, in the case of RFDiffs the symmetry of the khronon action is reduced to the shift symmetry ϕ →φ = ϕ + const , (2.16) allowing general dependence of the action on the derivatives of ϕ.
Other geometrical quantities associated to the foliation are constructed out of u µ and its derivatives. We have the following expressions for the spatial projector: the extrinsic curvature: the acceleration of the congruence of curves normal to the foliation: 2.18) and the intrinsic Riemann tensor: where in the last equation (4) R α βγδ is the 4-dimensional Riemann tensor. Now it is straightforward to obtain the covariant form of the actions (2.2), (2.5), (2.7) by identifying the quantities appearing in the ADM decomposition with the appropriate combinations of ∂ µ ϕ, u µ , P µν , K µν , etc. Namely, writing the covariant objects in the unitary gauge we obtain: In this way one ends up with an action describing GR plus a derivatively coupled khronon field ϕ.
Postponing the detailed analysis to the following sections, let us discuss here an important issue about the counting of degrees of freedom in the Stückelberg picture. Note that the higher spatial derivatives of Hořava theories in the unitary gauge translate in the covariant picture into higher covariant derivatives of the khronon field. Consider, for example, the covariant form of the action (2.5) where K ≡ K µ µ , and K µν and a µ depend on the khronon field via the expressions (2.17), (2.18), (2.13); dots stand for the higher-order terms 10 . It is easy to see that the second and the third terms in (2.19) contain four derivatives 11 of ϕ. Thus in a general coordinate frame the equations of motion for ϕ are fourth-order in time. Naively, this implies presence of ghosts and hence inconsistency of the theory. However, as discussed in [12] this conclusion would be incorrect due to the following reason. Though the action (2.19) has a covariant form, the physical content of the model does depend on the choice of the coordinate frame. The counting of degrees of freedom of the model must be performed in the preferred frame which is set by the preferred foliation. In this frame the number of time derivatives in the equations following from (2.19) is reduced to two, and the ghost modes are absent.
While the above statement is obvious in the unitary gauge where perturbations of the khronon field are set to zero, it is rather non-trivial for general gauges where ϕ is allowed to fluctuate. Let us formulate this statement more precisely. Consider perturbations of the field ϕ around a backgroundφ, ϕ =φ + χ .
We do not make any assumptions about the backgroundφ, in particular, it does not need to obey equations of motion. For the sake of the argument we treat the metric as external assuming that it remains unperturbed; again, the background valueḡ µν (x) can be arbitrary.
Let us now choose a coordinate frame such that the time coordinate coincides with the background value of the khronon field, t =φ . (2.20) Then in this frame the linearized equations of motion for χ following from (2.19) are second order in time. This proposition is proved in Appendix A. As explained there, the proof uses in an essential way the invariance of the action under the transformations (2.15). As a consequence it also applies to the projectable model (2.2) which, being FDiff-invariant, obeys this symmetry. On the other hand, the above statement does not hold for the model (2.7) whose covariant form is not invariant under (2.15). We will see in Sec. 4 that this makes the equation of motion for the khronon in the model (2.7) fourth-order in time; as one can anticipate, this will lead to certain pathologies.

Gradient instability
In this section we will study the projectable case (2.2), where the lapse N is assumed to obey the condition (2.1). Let us analyze the dynamics of small perturbations around flat background. At this stage it is convenient to work in the unitary gauge. We concentrate on the scalar perturbations of the metric and use the decomposition (2.9), where due to the projectability condition φ depends only on time. Substituting these expression into the action (2.2) we obtain the following quadratic Lagrangian: where the constants f 1 and g 1 are related to the coefficients in the potential (2.4); the precise form of this relation is not important for us. Note that the perturbation of the lapse φ drops out of the quadratic Lagrangian because of the projectability condition: when φ does not depend on the space coordinates all terms containing it are total derivatives. Integrating out the non-dynamical fields B and E, one finds Let us discuss this result. The positivity of the kinetic term for the field ψ, needed for the positivity of kinetic energy, imposes the restriction: On the other hand, the dispersion relation following from (3.2) has the form We see that for λ satisfying (3.3) the field ψ is tachyonic at small values of spatial momenta. This instability can be cut off by the second term in the bracket at p ∼ M * . One estimates the maximal rate of instability as To obtain an acceptable phenomenology, we must require that the instability does not develop during the life-time of the Universe, 12 where H 0 is the present Hubble expansion rate. Recall that the physical meaning of M * is that of the scale suppressing higher-derivative operators in the gravitational action. In other words, it is the genuine scale of quantum gravity in the model at hand. A theoretically motivated value for M * would be a microscopic scale of order or somewhat below the Planck mass. However, to be as general as possible, we choose not to impose any a priori assumptions about the value of M * . Then we are left with the experimental lower bounds on M * . The direct bound is rather mild and comes from the tests of Newton's law at the distances ∼ 10µm [37] and reads, Even taking M * at the lower end of the allowed range, the stability condition (3.6) implies that the parameter λ must be extremely close to 1, This is unappealing from the theoretical viewpoint. In particular, it is unclear if (3.8) can be preserved by radiative corrections. Most importantly, though, we are going to see that the condition (3.8) introduces a more important problem: the scalar mode becomes strongly coupled at unacceptably low energies.

Strong coupling
In this section we extend the study of the projectable model beyond the linear level and consider self-interaction of the scalar mode. As we will show, these considerations rule out the model as a weakly coupled alternative to GR. We will use the Stückelberg formalism described in Sec. 2.2. To construct the appropriate khronon action one notices [12] that the projectable model (2.2) with the condition (2.1) is locally equivalent to the theory with the same action where the lapse N is set to a constant; for example one can choose 13 N = 1. Indeed, the only difference between the two theories is in the equation obtained from the variation of the action with respect to the lapse. But this equation in the projectable model has the form of a Hamiltonian constraint integrated over the whole space and does not affect the local physics, cf. [38]. In Ref. [12] the condition N = 1 was implemented by introducing a Lagrange multiplier, and it was shown that the classical low-energy dynamics of the khronon field are equivalent to that of a pressureless fluid [38]. In the present paper we adopt a different approach more suitable for the analysis of the quantum properties of the system. Namely, we implement the constraint by adding to the action (2.2) a potential term and considering the limit M → ∞. In this limit the potential becomes very steep and enforces N = 1. We will keep the coefficient M finite in the intermediate calculations and take the limit only at the last step. Now we are ready to write down the covariant form of the action for the projectable model. Using the prescription of Sec. 2.2 we obtain, (3.9) where ϕ is the khronon field, M λ and X are defined respectively in (2.10) and (2.14). In this expression we have omitted the terms with higher spatial derivatives from the potential (2.4). Note that (3.9) coincides with the action of the ghost condensate model [18]. In particular, the last term in (3.9) is a covariant realization of the term with four spatial derivatives of the ghost condensate from Ref. [18]. However, whereas in [18] M and M λ are assumed to be of the same order, we are interested in the case M λ ≪ M .
To get a first feeling about the khronon dynamics it is instructive to neglect the backreaction of the khronon on the metric. Formally, this is achieved by considering the limit The perturbations of the metric get frozen in this limit, so the khronon dynamics effectively decouples from that of gravity. Following the terminology adopted in massive gravity we refer to this regime as 'decoupling limit'. Note that the phenomenological constraint (3.8) implies M λ ≪ M P so that one of the decoupling limit assumptions is automatically satisfied. On the other hand the second assumption M ≪ M P is in conflict with the fact that we want to take M eventually to infinity in order to recover the projectable model. This means that the proper analysis of the projectable model will require going beyond the previous decoupling regime. However, working in the decoupling limit is still instructive. First, it will allow to outline the steps of the analysis in a simplified setting. Second, we will use the decoupling limit results in the next section where we extend the model to a version of the ghost condensation. With this in mind, let us proceed for the moment with the assumptions (3.10). Setting the metric perturbations to zero, h µν = 0, and writing ϕ = t + χ we obtain: 11) where we have written only the most relevant (cubic) interaction terms. This action describes a non-relativistic scalar with derivative self-coupling. The dispersion relation for χ reads, We can compare this dispersion relation with that obtained previously for the scalar mode of the projectable model, Eq. (3.4). In deriving (3.11) we neglected the higher-derivative terms. This explains the absence in (3.12) of terms with higher powers of momentum suppressed by M * . Moreover, in the projectable limit M → ∞ (3.12) reduces to ω 2 = 0 which coincides with the |λ−1| → 0 limit of (3.4). The latter corresponds to the decoupling limit (3.10).
From the presence of derivative interactions in (3.11) one expects the field χ to become strongly coupled above certain energy/momentum scale. The easiest way to identify this scale is to perform power-counting. Due to the non-relativistic form of the dispersion relation for the χ-field the power counting rules in our case are different from the standard ones. To identify the proper rules we follow the approach of [18,12] (see also [39]). One notices that the quadratic part of the action is invariant under the scaling transformations The interaction terms in (3.11) scale as b 1/2 and b 5/2 implying that the strength of the interaction grows at short distances. To simplify the further analysis let us choose the units in such a way that the quadratic part of the action takes canonical form, i.e. it contains only order-one numerical coefficients. This is achieved by the rescaling which yields, Now the momentum and frequency cutoffs (strong coupling scales) of the theory are identified as the appropriate powers of the scale suppressing the interaction terms, where the subscript "dec" reminds that these expressions are obtained in the decoupling limit. Hereω is the frequency corresponding to the rescaled timet; going back to the physical frequency one obtains Note that the cutoffs in momentum and energy are different reflecting the non-relativistic nature of the theory; they are related by the dispersion relation (3.12).
In the projectable limit M → ∞ both cutoffs (3.13), (3.14) vanish which would imply that the model is strongly coupled at all scales. However, as already pointed above, this reasoning has a caveat: the expressions (3.13), (3.14) are obtained under the assumptions M P → ∞ (which in particular implies M ≪ M P ). We are going to see that a finite value of M P raises the cutoff from zero; but it is still too low to be phenomenologically acceptable, cf. [23].
To obtain the correct value of the cutoff we have to take into account the fluctuations of the metric. Let us first see how this affects the analysis at the quadratic level. Taking the metric perturbations in the form (2.9) and expanding the action (3.9) to quadratic order we obtain, Note that the first line here is nothing but the quadratic part of the standard Einstein-Hilbert action. If one fixes in (3.15) the gauge χ = 0 and takes the limit M → ∞ forcing φ to be constant one recovers the unitary gauge Lagrangian 14 (3.1). For our present purposes it is more convenient to choose instead the gauge where khronon perturbations do not vanish. To integrate out the non-propagating degrees of freedom recall that we are working in the regime M λ ≪ M P . From the dispersion relation (3.4) we know that in this case the frequency of scalar perturbations is much smaller than the spatial momentum, Thus we can neglect the termψ 2 compared to ψ∆ψ in (3.15). Then the equation of motion of ψ implies ψ = φ; substituting this into the action we obtain, S Varying with respect to φ one finds and hence

S
(2) This gives the dispersion relation which reduces to (3.12) in the decoupling limit M P → ∞. On the other hand, in the projectable limit M → ∞ it correctly reproduces the first term in the exact dispersion relation (3.4).
Let us turn to the interactions. One has to compare contributions coming from various terms in the action (3.9). We start by considering the second term in (3.9). Due to the inequality (3.16) the leading interaction is given by the cubic term with the smallest number of time derivatives. A simple calculation yields where in the second line we substituted φ from (3.17). Note that in the decoupling limit this interaction term reproduces the third term in the action (3.11). Next, let us estimate the contributions from the first and third terms in (3.9) and show that they are suppressed compared (3.20). This can be done without explicit calculation on purely dimensional grounds. Indeed, (4) R does not explicitly depend on the khronon and has mass dimension two. Thus the possible leading contribution (the one without time derivatives when written in terms of the metric) has the schematic form where in the last relation we used Eq. (3.17) and took the limit M → ∞. This contains more time derivatives than (3.20) and hence is suppressed. Finally, the schematic form of the leading contribution from the third term in (3.9) is M 2 λχ (∆χ) 2 . This is clearly suppressed compared to (3.20) by the ratio (M λ /M P ) 2 .
In the projectable limit M → ∞ both the quadratic action (3.18) and the interaction term (3.20) simplify. Combining them together we obtain It is now straightforward to find the strong coupling scales of the model. Upon bringing the quadratic part of the action into canonical form by the rescaling one finds that the interaction term is suppressed by the parameter which thus sets the momentum cutoff of the theory. The energy cutoff contains additional factor M λ M −1 P and reads Note that these cutoffs go down to zero in the naive GR limit λ → 1 (M λ → 0). This agrees with the results of [23].
Using the definition of M λ , Eq. (2.10), we find that the stability bound (3.8) implies This corresponds to the distance of order 10 13 cm, within which the theory is strongly coupled and the perturbative analysis breaks down. Note that the strong coupling cannot be resolved by the higher-derivative terms: such terms would lead to sizable modifications of the Newton's law that are forbidden at these distances. We conclude that the projectable model fails to provide adequate description of gravity within distances ∼ 10 13 cm.

Extensionà la ghost condensation
The analysis of the previous section suggests a possible way to address the problems of the projectable model. The idea is to consider a modification of the theory with finite value of the parameter M in the khronon action (3.9) [12]. This approach implies reducing the symmetry of the theory from FDiffs down to RFDiffs as the khronon action is no longer invariant under the general reparameterizations (2.15) but only under the shifts (2.16). We presently explore this option and show that in spite of the possibility to improve the behavior of the khronon at low energies the theory still encounters serious problems in the UV completion. As already pointed out, at finite values of M the low-energy khronon action (3.9) coincides with the action of the ghost condensate model [18]. The latter is a consistent effective theory describing low-energy modification of GR. From the EFT point of view it is natural to assume that the scales M and M λ appearing in (3.9) are of the same order. Then both the momentum and frequency cutoffs obtained in the decoupling limit (3.13), (3.14) are finite and of order M . From the dispersion relation (3.19) one observes that the instability of the scalar graviton at low values of momenta is still present in this model but is truncated by the quartic term at This corresponds to the maximal instability rate This rate is slower than the present rate of Hubble expansion for which would make the instability harmless. 15 The problems appear, however, when one tries to UV complete the model to a renormalizable theory of gravity. First, a necessary requirement to obtain a weakly coupled UV completion is that the strong coupling observed within low-energy EFT must be resolved by the higher-derivative operators, cf. [26]. This implies that the scale M * suppressing higher derivatives should be lower than M , M * < M . Thus we obtain Note that this is many orders of magnitude smaller than M P so that the model contains two very different scales. For the gravitational sector alone, this seems to be only an esthetic problem: the values (4.2) are compatible with the direct gravitational bound (3.7). Besides, we expect this hierarchy to be stable under radiative corrections. The scale which effectively truncates the power-law divergencies in the model is the lower scale M * ; therefore, the corrections to both M P and M * are small 16 . However, the real tension arises when we take into account coupling of the theory to matter. The size of the higher-derivative Lorentz violating terms in the action of the matter sector is experimentally constrained from below. The leading effect of these terms is the modification of dispersion relations of the matter fields at high energies due to contributions with higher powers of spatial momentum, see to enter into radiative corrections for M * thus giving rise to a fine-tuning problem 17 . A more detailed study of this issue is beyond the scope of the present article.
The second problem stems from the reduction of the symmetry from FDiffs down to RFDiffs. To get a renormalizable theory one must consider the most general RFDiffinvariant action (2.7) which contains, among other terms, the term with time derivatives of the lapse, 18 where M λ 1 is defined in (2.10). The easiest way to see that this term leads to pathologies is to adopt the Stückelberg picture. The covariant form of (4.3) reads Let us work in the decoupling limit; the conditions for the validity of this approximation will be specified below. One writes ϕ = t + χ and expands (4.4) up to quadratic order in χ: 15 It is argued in [40] that nonlinear dynamics of the ghost condensate suppresses the exponential growth of perturbations predicted by the linear theory, making the model phenomenologically viable for M up to 100 GeV. 16 See also the discussion of a similar issue in the next section. 17 The fine tuning may be absent if there is a symmetry protecting M * from this kind of corrections. 18 Even if this term is absent in the bare action, one expects it to be generated by quantum corrections.
This contains four time derivatives of the khronon perturbations which implies the presence of a second helicity-0 mode besides the excitation studied so far. As the new mode appears due to higher time derivatives one is tempted to qualify it as a ghost. In fact, we are going to show that this mode also exhibits a gradient instability, so it is more appropriate to call it tachyonic ghost. In Appendix B we discuss the difference between the tachyonic ghost and an ordinary ghost with stable dispersion relation 19 . In particular, we demonstrate there that the presence of the tachyonic ghost is compatible with the fact that in the unitary gauge the action for the model, Eq. (2.7), is only second order in time derivatives.
To proceed we need the quadratic action for the khronon perturbations. This is obtained in the usual way by first covariantizing the terms entering (2.7), (2.8) and then taking the decoupling limit. This yields, In deriving this expression we have omitted the contributions of terms with higher spatial derivatives; thus the domain of validity of (4.5) is restricted to spatial momenta smaller than M * . The last term in (4.5) represents the contribution of the term αa i a i in the potential (2.8) (recall the definition of M α in (2.10)). The action (4.5) describes two modes with dispersion relations At small values of momenta we find, In the first expression one recognizes the dispersion relation (3.12) for the previously encountered gapless mode. The second mode (4.7) is unstable. It possesses a frequency gap of order M 2 /M λ 1 . For the choice M λ 1 ∼ M , natural from the low-energy effective theory point of view, the frequency of this mode lies above the scale M . Thus this mode would be simply discarded as unphysical in the EFT considerations with cutoff M , such as in [18]. However, our case is different: we are looking for a UV-complete model and thus have to take into account all possible excitations of the system. According to the above expressions the instability rate of the new mode grows with momentum, the fastest instability occurring at |p| ∼ M * . [At larger momenta terms with higher spatial derivatives in the full action (2.7) can, in principle, stabilize the mode.] A lower estimate for the instability rate is obtained by keeping only the last term under the square root in Eq. (4.6). This yields, (4.8) 19 Of course, both type of ghosts signal pathologies of the theory.
where we have assumed Requiring the instability to be slower than the Hubble rate H 0 and taking into account the experimental bound (3.7) on M * one obtains, M λ 1 10 60 eV ≈ 10 32 M P . (4.10) Having such a huge value for M λ 1 is rather unsatisfactory. First, it means that one introduces into the theory, besides M * and M P , one more hierarchical scale M λ 1 ≫ M P . We do not know whether this hierarchy is technically natural or not. More importantly, this term also gives rise to interactions which are enhanced by the large value of M λ 1 . We will demonstrate shortly that this reintroduces the strong coupling problem. Before addressing the interactions let us make a step back and discuss the conditions for the validity of the decoupling limit. The simplest requirement would be that all the parameters in the khronon action (4.5) are smaller than the Planck mass. This is automatically satisfied by the parameters M λ , M , M α that must be small for the stability of the gapless (ghost-condensate) mode. However, the value (4.10) of M λ 1 clearly violates this condition. Thus additional considerations are needed to establish whether the decoupling limit holds or not for M λ 1 as large as (4.10). A rigorous method to do this is to verify that the terms in the Lagrangian describing mixing between the khronon field and the metric components are small compared to the other terms. Instead of following this route we will take a shortcut and compare the decoupling limit dispersion relations (4.6) for the scalar modes with the exact expressions; their coincidence will serve as a criterion for decoupling. The exact dispersion relations are obtained in Appendix B within the unitary gauge, see Eq. (B.2). In the case M λ 1 ≫ M P they simplify to 20 , On the other hand, the decoupling limit result (4.6) at large M λ 1 becomes The two expressions coincide and thus the decoupling holds, provided (4.13) 20 Note that these expressions immediately imply a lower estimate for the instability rate of the new mode which is complementary to (4.8), This estimate does not depend on M λ 1 . Evaluating it at |p| ∼ M * and requiring the instability to be slower than the Universe expansion rate we obtain the upper bound This is marginally compatible with the requirement (4.9) and the experimental bound (3.7).
Note that this allows for M λ 1 to be much larger than M P . Taking M λ ∼ 0.1 eV as in the previous estimates (see also the footnote 20) one finds that the decoupling limit is applicable up to M λ 1 of order (4.10).
The previous considerations show that to estimate the scale of strong coupling we can concentrate on the self-interaction of the khronon neglecting its mixing with the metric. Clearly, for large M λ 1 the leading interactions come from the term (4.4). At the cubic level we obtain the following contribution, This must be combined with the quadratic action (4.5). The form of the dispersion relation (4.12) shows that at large M λ 1 the latter action is dominated by the first and the third terms. Retaining only these two terms we obtain the khronon action in the limit of interest The quadratic part is brought to canonical form by the rescaling From the scale suppressing the interaction term in the resulting action one reads off the momentum and frequency cutoffs, Substituting here M λ from (4.11) and M λ 1 from (4.10) one obtains 21 Λ p 10 −16 eV. Such a low cutoff is phenomenologically unacceptable. This shows that the initial hope, namely that the model could be UV complete and weakly coupled, is not met in reality 22 . The general conclusion of this section is that reducing the symmetry of the theory from FDiff to RFDiffs cannot possibly improve the scalar sector. The reason is that the smaller symmetry allows for the new operators such as (4.3) which bring in additional difficulties in the form of tachyonic ghosts. In retrospect, we learn that the symmetry under the FDiffs (or the khronon reparameterizations (2.15) in covariant language) plays the quite important role of preventing this kind of pathologies. In the next section, we return to the other model compatible with the FDiffs -the general non-projectable case -where this problem is automatically turned away. 23 21 Note that this is close to the result in the projectable case, Eq. (3.21). 22 One may wonder if it is possible to cure the model by choosing some of the parameters in the action (4.5) negative. A simple reasoning demonstrates that this is not the case. The analysis in the unitary gauge (Appendix B) shows that negative values of M 2 λ 1 , M 2 λ imply negative energies in the UV. Next, making M 2 α negative does not essentially improve the behavior of the tachyonic ghost. Finally, the choice M 4 < 0, while removing the instability from the mode (4.7), destabilizes the old gapless mode. This brings us back to the situation studied in the previous section. 23 Clearly, the symmetry (2.15) must be non-anomalous in order that this statement remains valid in the quantum theory. This question will not be addressed here.

Stability and absence of strong coupling
We presently consider the non-projectable FDiff-invariant model with the action (2.5), (2.6). We are going to show that, remarkably, the problems associated to the additional modes which plague the two models considered previously are absent in this case [25,26].
Let us start by studying the spectrum of linear perturbations around flat background. As usual, we concentrate on the scalar sector and use the decomposition (2.9). In the unitary gauge the quadratic Lagrangian reads where f n , g n are related to the coefficients in front of the higher derivative operators in the potential (2.6); the precise form of this relation is not important for our purposes.
Integrating out the non-dynamical fields B, E and φ we obtain, where the polynomials P , Q have the form, The Lagrangian (5.2) describes a single mode which is free of pathologies provided that two conditions are satisfied. First, the positivity of the kinetic term can be achieved by choosing 24 λ > 1. Second, from the dispersion relation of the propagating mode ψ, one reads off the condition to avoid exponential instabilities, This puts certain restrictions on the coefficients α, f n , g n , that are presented in Appendix C.
In particular, at low energies the dispersion relation takes the form, Thus stability requires 25 Note that the dispersion relation (5.7) describes a gapless mode propagating with constant velocity which is generically different from one (the velocity of the helicity-2 modes, i.e. gravitons). This signals that in the model at hand Lorentz symmetry is broken down to arbitrary low energies. In the UV the dispersion relation (5.5) takes the form ω 2 ∝ p 6 which obeys the anisotropic scaling with z = 3. This is compatible with the power-counting arguments in favor of renormalizability.
It is important to realize that the healthy behavior of the scalar mode can be achieved simultaneously with the stability in the sector of the helicity-2 perturbations. Consider operators in the action (2.5), (2.6) which contribute at the quadratic level. Upon integrating by parts and using Bianchi identities one obtains a list of 10 inequivalent combinations, The dispersion relation of the helicity-2 modes depends only on the coefficients in front of the operators in the first column of this list. After fixing these to ensure stability of the helicity-2 modes, we still have the freedom to choose the coefficients of the remaining operators in the list. The number of free parameters matches with the number of coefficients α, f n , g n in the scalar Lagrangian implying that we have freedom to adjust the latter coefficients to satisfy (5.6).
To get more insight into the dynamics of the model (in particular, at the nonlinear level), we make use of the Stückelberg formalism. The covariant form of the model action was given before, see (2.19). It is convenient to rewrite it as where u µ is related to the khronon field ϕ by (2.13) and dots stand for the terms with higher derivatives which are not important at low energies. The action (5.9) is closely related [26,34] to the action of the Einstein-aether theory [30] (see [31] for a review). The latter is an effective theory describing breaking of Lorentz-invariance by a time-like vector field with unit norm. The difference between our case and the Einstein-aether theory is that in our model the unit vector u µ is by construction hypersurface orthogonal, i.e. it is completely characterized by the scalar khronon field ϕ. The similarity between (5.9) and the Einstein-aether theory will be exploited in the next subsections where we study phenomenological consequences of the model. Let us use (5.9) to study interactions of the khronon perturbations at low energies. It is convenient to introduce the scales M λ and M α as defined in (2.10). These scales characterize the khronon action. We will assume them to be much smaller than M P , so that the metric perturbations are frozen out (in other words, we will be working in the decoupling limit). This assumption is justified by the phenomenological bounds (5.35) that will be obtained in Sec. 5.4, and which constrain the dimensionless parameters |λ − 1|, α to be much smaller than one 26 . Writing down ϕ = t + χ and expanding (5.9) up to cubic order in χ we obtain, Let us analyze this expression. The action (5.10) describes a propagating mode with dispersion relation This coincides with the exact dispersion relation (5.7) in the limit |λ − 1|, α ≪ 1. The form of the action (5.10) is uniquely fixed by the reparameterization symmetry (2.15) of the khronon and by the Lorentz symmetry, both non-linearly realized on the khronon perturbation χ. Up to quadratic order the reparameterization transformations read where ǫ(t) is an arbitrary function of time. Under boosts the field χ transforms as where ε i is a 3-dimensional vector characterizing the boost. It is straightforward to check that (5.10) is the only action invariant under these symmetries up to cubic order. Finally, one observes that (5.10) contains an interaction term with three time derivatives; it produces a contribution with third time derivative in the equation of motion for χ. One may be worried that this leads to appearance of a new unwanted degree of freedom. However, from the analysis in the unitary gauge, we know that this degree of freedom is spurious and it should be possible to eliminate it with an appropriate choice of variables (which corresponds to fixing the foliation consistently). To see how this is done explicitly, one considers the change of variable χ =χ +χχ .
This substitution can be interpreted as the change of the time foliation, t → t − χ. Indeed, up to cubic orderχ(t, x) = χ(t − χ, x). In terms of the new variable the action takes the form We observe that a term with three time derivatives is still present but its structure has changed: its contribution to the equation of motion contains now only second time derivative. One concludes that no new degrees of freedom appear. It is worth emphasizing that the existence of the change of variables with the above properties follows from the general statement proved in Appendix A that the equations of motion for the khronon field are second order in time for an appropriate choice of the coordinate system. The low-energy action (5.12) contains derivative couplings which become stronger as the energy / momentum grows. Let us estimate the scale where the strength of the interactions would become of order one. Performing the rescaling which casts the quadratic part of the action into canonical form we obtain, From the scales suppressing the interaction terms one reads out the momentum and frequency cutoffs of the low-energy description: If (5.9) were the full action of the theory, it would become inconsistent at energies / momenta above these scales [41]. Note that Λ ω and Λ p are related by the dispersion relation (5.11) and are different if the khronon velocity differs from 1. We see that the scales (5.14), (5.15) are low whenever there is a large discrepancy between M α and M λ . On the other hand, in the case M α ∼ M λ the cutoffs essentially coincide with M α , We concentrate on this latter case in what follows. Of course, the action (5.9), and hence (5.13), represents only the low-energy part of the full action (2.5). So the existence of a finite cutoff for the low-energy theory (5.9) does not imply any inconsistency. At energies/momenta larger than M * one has to take into account the higher-derivative terms in the full action, and any conclusions drawn from (5.9) become invalid. By construction, the role of the higher-derivative terms in (2.5) is to modify the power-counting rules at high energies in such a way that all the interactions become marginal under the anisotropic scaling. One concludes that strong coupling is avoided provided the scale M * of higher-derivative operators is lower than 27 Λ p,ω [26]. Let us emphasize that, contrary to the claim in [41], the choice of M * parametrically below M P does not introduce a fine-tuning in the model. In fact, having M * somewhat below M P is technically natural. From the point of view of the low-energy theory, the reason is that the cutoff for the power-law growing of the couplings is set by the lower scale M * . Thus neither M P nor M * receive large quantum corrections.

Phenomenological considerations
The absence of notorious pathologies in the model (2.5) makes it worth to have a closer look on its phenomenological consequences. We are interested in the phenomenology of the model at energies much lower than M * . The natural language for this analysis is provided by the Stückelberg formalism where at low energies the theory reduces to GR plus the khronon field. Then all non-trivial effects of the model are clearly interpreted as due to the presence of the khronon. Before proceeding we have to specify the coupling of the khronon to the fields of the Standard Model, to which we collectively refer as "matter". The FDiff symmetry requires that matter couples to the khronon via geometrical objects, such as u µ , a µ , K µν , etc. Possible interactions fall into two classes having qualitatively different phenomenological consequences. The first class consists of couplings which contain derivatives of the vector u µ ; examples of this type of couplings are: a µψ γ µ ψ , K µνψ γ µ ∂ ν ψ , (5.18) where ψ is some fermionic matter field. Importantly, the combinations of the khronon field entering into the operators (5.18) have vanishing VEV in the flat background. Therefore, these operators start linear in the khronon perturbation χ and describe its derivative interaction with matter. It is easy to see that in terms of canonically normalized fields these couplings have dimensions larger than four 28 and are suppressed by the high-energy scale M * . Given that M * is large, one expects the effect of these operators to be highly suppressed at the energies accessible to the present-day experiments. It is still possible that couplings of the type (5.18) may be probed in some precision measurements. This issue is, however, beyond the scope of the present article. The second, more 'dangerous', class of operators describe coupling of the matter fields directly to the vector u µ ; examples are: The crucial property of these operators is that they give rise to Lorentz-violating effects within the Standard Model as they couple matter fields to the VEV of u µ . Such effects and constraints on them have been extensively studied in the literature, see [6,7,8,9] for reviews. The effects of the operators (5.19) fall into two categories: those that become pronounced only at high energies, and those that persist at all energy scales. Effects of the first category are associated to Lorentz-violating operators of dimensions larger than four such as the last operator in (5.19). The strength of these effects is characterized by the scale suppressing the higher-order operators which we will collectively denote as M The experimental data constrain the scale M (mat) * from below. A rather robust bound comes from astrophysical observations and exploits the fact that the higher-order operators lead to modification of the dispersion relations of the matter fields, in particular, photons, at high momenta 29 where η is a dimensionless coefficient. This would produce a frequency dependent delay in the arrival times of γ-rays emitted by a distant source. The absence of such a time-lag in the signals coming from active galactic nuclei [44] and γ-ray bursts [45]  In the simple case when M * (mat) and M * are of the same order this translates into the lower bound on M * . Note that this is several orders of magnitude below the upper limit (5.17) imposed by the absence of strong coupling. It has been argued that considerably stronger constraints on the scale M (mat) * come from the physics of ultra-high-energy cosmic rays (UHECR) [46,47,48,49]. These constraints are, however, less robust than (5.21) as they are sensitive to a number of assumptions about the sources of cosmic rays and the UHECR chemical composition (e.g. the constraints essentially disappear if UHECR of the highest energies consist of heavy nuclei [49]). Thus we do not discuss these constraints in the present paper.
Finally, violation of Lorentz symmetry in the matter sector by operators of dimensions 3 and four 4, such as the first and the second operators in (5.19), would lead to sizable effects even at low energies 30 . On the other hand, the experimental constraints on these effects are 29 Terms with odd powers of momentum can be forbidden by imposing the CPT invariance. In particular, this allows to avoid bounds from the absence of vacuum birefringence discussed in [7,8]. It is also worth mentioning that modification of dispersion relation of the matter fields at high momenta is strongly suppressed in supersymmetric Lorentz-violating extensions of the Standard Model [42,43]. 30 It is worth stressing that without additional assumptions the khronon -matter couplings generically would be non-universal, i.e. species-dependent, implying, in particular, violation of the weak equivalence principle. extremely tight [9]. One arrives to the conclusion that breaking of Lorentz symmetry in the matter sector at the level of dimension 3 and 4 operators must be highly suppressed. For some of the couplings this may be achieved by imposing discrete symmetries. For example, the first operator listed in (5.19) can be forbidden by requiring CPT invariance. However, a stronger mechanism is needed to suppress all dimension 3 and 4 operators. We will mention one possibility in the Discussion section. For the moment let us just assume that such a mechanism exists. Then to the leading approximation the coupling of the khronon field to matter at the lowest-derivative level must be encoded in a universal effective metric where β is a dimensionless constant. Clearly, this kind of coupling preserves Lorentz invariance of the matter sector. We concentrate on this type of coupling in what follows.

Velocity-dependent forces and instantaneous interaction
Before engaging into a systematic analysis of the observational constraints on the universally coupled khronon let us make some preliminary estimates of the type of effects induced by the coupling (5.22). When β is small one expands the interaction terms and obtains to leading order in β where T µν is the matter energy-momentum tensor. Expanding the khronon field as usual, ϕ = t + χ, one obtains the action for small perturbations including the source term From this expression we read off the khronon exchange amplitude between two sources with energy-momentum tensors T µν (1) and T µν (2) , This amplitude does not encapsulate all the new interactions due to the presence of the khronon field. In the systematic treatment Eq. (5.23) must be supplemented by the amplitudes involving khronon-graviton mixing and by the contributions coming from the modification of the metric propagator due to the Lorentz symmetry breaking. This will be done in the next section within the formalism of the parameterized post-Newtonian (PPN) expansion. Here we note that the direct khronon exchange (5.23) dominates over other effects in the case α ∼ |λ − 1| ≪ β . The simplicity of the amplitude (5.23) makes the analysis in this case much more transparent than in the full PPN treatment. Thus we first concentrate on the case (5.24). Let us consider explicitly the khronon-induced interaction for two point-like nonrelativistic sources. In this case T 0i (a) = m a v i a δ(x − x a ), where x a , v a , a = 1, 2 are the coordinates and velocities of the two particles. In the non-relativistic regime one can neglect the time derivatives in the khronon propagator. After integration over spatial coordinates one obtains from (5.23) the following contribution into the two-particle Lagrangian, where r = x 1 −x 2 andr = r/r. This interaction is quite peculiar: it depends simultaneously on the distance between the particles and on their velocities with respect to the preferred frame. Note that due the velocity dependence this term contributes non-trivially into the expressions for the conserved energy and momentum of the system. An interesting special case of Eq. (5.25) corresponds to the situation when one of the particles (say, particle 2) is much heavier than the other so that its velocity is approximately conserved (like, e.g., in the Sun -Earth system). One can show that in this case (5.25) reduces to up to a total time-derivative. We observe that the dependence on the test particle velocity has dropped out. Equation (5.26) has a simple interpretation as a direction-dependent contribution into the gravitational potential of the heavy source. In the next section we will discuss the phenomenological constraints on this type of contributions within the PPN framework.
An important remark is in order. The form of the khronon exchange amplitude (5.23) makes it manifest that the model involves instantaneous interaction. This is due to the factor ∆ −1 in the khronon propagator. This factor may cancel out in some special cases when the sources contain sufficient number of spatial derivatives, but not in general. Appearance of instantaneous interactions is a common feature of modified gravity models (more generally, gauge theories) with broken Lorentz invariance [19,50,51,52]. In a theory that does not aspire to a Lorentz invariant UV completion this does not pose an obvious obstruction. In particular, it certainly does not imply any problems with causality which is defined with respect to the preferred time slicing. Instantaneous interactions do introduce a kind of non-locality since far-away sources may affect the immediate future in any given local domain. However, the strength of the non-local effects seems to decay with distance and is suppressed both by post-Newtonian factors of v 2 and the (small) model parameters. Hence, it is not obvious if the presence of the instantaneous interaction leads to any significant constraints on the model. A more detailed investigation of related issues (the possibility to have any direct observational limits on instantaneous interactions, the implications for BH physics, etc.) is left for future.

Universal coupling and Post-Newtonian Parameterization
In the universally coupled case the effects of the khronon field are naturally interpreted as modification of the (universal) gravitational interaction between matter particles. The size of allowed effects is constrained by the existing tests of GR [37]. So we now turn to the bounds that the model has to satisfy in order to pass these tests.
The khronon sector is completely described by three parameters: (λ − 1), α in the action (5.9) and β in the effective metric (5.22). It is convenient to make a field redefinition by trading the "Einstein-frame" metric g µν in favor of the effective metricg µν to which matter couples minimally. The result of this redefinition is readily obtained if we assume the parameters (λ − 1), α, β to be much smaller than one; we will see that this assumption is justified by the phenomenological bounds which indeed require the above parameters to be small. Then to the leading order we write Using the identity and we have omitted tildes over the new variables (we will always work with the redefined metric in the rest of this section so this will not lead to confusion). Note that the appearance of the new parameter β in the khronon action (5.27) can be traced back to the presence of the parameter ξ in front of the 3-dimensional scalar curvature in the unitary-gauge potential, see (2.6), (2.4). In the pure gravity theory this parameter does not have physical meaning as it can be eliminated, say, by a rescaling of the time coordinate, cf. the discussion after Eq. (2.4). This is no longer true in the presence of matter: as we are going to see physical observables depend on the value of β. As already mentioned before, the action (5.27) has the same form as the action of the Einstein-aether theory [30] which has been extensively studied as a phenomenological model for violation of Lorentz invariance [31]. The difference of our model from Einsteinaether is that in our case the aether vector u µ is by construction hypersurface-orthogonal. As a consequence, in our case aether propagates a single longitudinal degree of freedom (khronon), while in general there are additional transverse modes. This implies that the results about Einstein-aether theory that are insensitive to the presence of transverse modes are also valid for our model 31 . 31 The relation with Einstein-aether can also be used to derive the result of metric redefinition (5.22) in the general case when the parameters (λ − 1), α, β are not small. In Einstein-aether the substitution (5.22) leads to the change of the coefficients in the aether Lagrangian; the corresponding formulas have been worked out in [53]. The most general Einstein-aether action contains an additional term ∇µuν ∇ µ u ν compared to (5.27). In our case when the aether vector uµ is hypersurface-orthogonal the four terms in the aether action are not independent [34] implying that the term ∇µuν∇ µ u ν can be eliminated. In this way The first set of constraints obtained in this way is related to the fact that the velocities of the modes propagated by the action (5.27) are in general different from 1, the maximal velocity of matter. Using the formulas for the Einstein-aether theory [31] we obtain to first non-trivial order in parameters the velocity of helicity-2 modes (graviton): and that of the helicity-0 mode (khronon): Note that (5.29b) coincides with the velocity entering the dispersion relation (5.11). If the velocities of graviton or khronon are smaller than 1 relativistic matter particles will quickly loose their energy via vacuum Cherenkov radiation [54]. This is strongly contrained by the existence of high-energy cosmic rays. Thus we conclude that the graviton and khronon velocities must be larger or equal to 1 which yields the bounds Another constraint is obtained from the comparison of the gravitational constants appearing in the Newton law and the Friedmann equation governing the cosmological expansion. Again, the transverse aether modes do not play any role in these considerations, so we can directly apply the results from the Einstein-aether theory to our case. The Newton constant, which is defined as the coefficient in the Newton law for the gravitational force between two static masses, is related to the parameters appearing in the action (5.27) as follows [31]: On the other hand, the cosmological expansion in the Einstein-aether theory is described by the standard Friedmann equation, where H is the Hubble rate, ρ -the energy density of the Universe, but with a different proportionality coefficient 32 [31], . (5.32) one arrives at the action (5.27) with which coincides with (5.28) for β ≪ 1. 32 The expressions (5.30), (5.32) for the case β = 0 were derived directly from the action (2.5) in [25].
The discrepancy between G N and G cosm is constrained by Big Bang nucleosynthesis [55], Barring accidental cancellations, this yields an order-of-magnitude bound on the parameters of the model, α, β, λ ′ 0.1 .
Stringent limits on any alternative theory of gravity come from the observational constraints on the parameters of the parameterized post-Newtonian (PPN) formalism. Remarkably, in the Einstein-aether theory all PPN parameters except two are the same as in GR [56]. The non-trivial parameters are called α P P N 1 , α P P N 2 and describe preferred frame effects related to breaking of Lorentz symmetry. We now argue that in the khrono-metric theory (5.27) all PPN parameters except α P P N 1 , α P P N 2 are the same as in Einstein-aether. Indeed, besides α P P N 1 , α P P N 2 and the parameters which vanish automatically in any theory described by a Lagrangian, there are three more PPN parameters: β P P N , γ P P N and ξ P P N .
The key point is that these three are determined from spherically symmetric solutions which are identical in the khrono-metric and Einstein-aether theories, see Appendix D. Thus in the khrono-metric model at hand these parameters have the same values as in GR, β P P N = γ P P N = 1 , ξ P P N = 0 .
One cannot use the results for Einstein-aether in the case of the parameters α P P N 1 , α P P N 2 : they describe effects related to the motion of the source with respect to the preferred frame and are contaminated in Einstein-aether by the contributions of the transverse aether modes. These parameters can be defined as the coefficients in the linearized metric produced by a point source of mass m in its rest frame [57]: where r is the distance from the source and v i is the velocity of the source with respect to the preferred frame. Note that the contribution proportional to α P P N 2 has the form of the direction-dependent gravitational potential encountered in Sec. 5.3. The current Solar system limits on these parameters are [37]: Derivation of α P P N 1 , α P P N 2 for the model (5.27) in the case α, β, λ ′ ≪ 1 is given in Appendix E; the result is Note that both parameters vanish if α − 2β = 0. Another interesting case is β = 0, λ ′ = α; this corresponds to the situation when the velocities of all modes in the theory are equal to 1, see Eqs. (5.29). In this case the parameter α P P N 2 which is most tightly constrained vanishes. Barring these special cases and assuming α, β, λ ′ to be of the same order we obtain the bound α , β , λ ′ 10 −7 ÷ 10 −6 . (5.35) To the best of our knowledge, this is the strongest constraint on the parameters of the model that can be obtained at present. As discussed in Sec. 5.1, the bound (5.35) combined with the requirement that the theory is weakly coupled implies the upper limit (5.17) on the scale suppressing higherorder operators in the gravitational action.

Summary and discussion
In this paper we have investigated the self-consistency issues related to the scalar graviton modes in Hořava's approach to quantum gravity. We have considered three models of non-relativistic gravity differing by the symmetry group and the requirement or not of the projectability condition. In our study we extensively used the Stückelberg formalism which makes the scalar modes explicit by encoding them into the khronon field: a scalar field with time-dependent VEV. This facilitates a lot the analysis of the scalar mode dynamics.
In two of the considered models (the projectable version of the original proposal [1] and a possible extension based on a smaller symmetry group), the scalar modes were found to exhibit pathological behavior, such as instability and strong coupling. These pathologies undermine the perturbative analysis. In particular, they invalidate the naive power-counting argument for renormalizability of these models.
Several qualitative lessons can be extracted from these studies. While it is relatively easy to make the scalar modes well-behaved in the UV, this is much harder to achieve in the IR. The scalar gravitons tend to develop gradient instabilities which can be suppressed only by pushing the model parameters to extreme values. This, in turn, introduces strong coupling. Thus the primary goal in constructing a consistent non-relativistic gravity model is to stabilize the scalar modes.
The analysis of Sec. 4 teaches us that reducing the symmetry to a smaller group than the FDiffs does not give any advantage in achieving this goal. Though this approach allows to improve the behavior of the scalar graviton of the projectable Hořava's model, it introduces yet other scalar modes which bring the pathologies back. Even relaxing the symmetry to the RFDiffs, that have as many local generators as the FDiffs, turns out to have the quite dramatic consequence of allowing for new operators that lead to tachyonic ghosts and therefore fast instabilities. Even if not present at tree level these operators would be generated by quantum effects. From this we conclude that the timereparameterization symmetry contained in the FDiffs plays the quite important role of preventing the appearance of tachyonic ghosts in the model.
Remarkably, the scalar mode is stable and also free from other pathologies in the third model which we analyzed. The latter is a natural extension of the original non-projectable Hořava model obtained by including into the action all terms compatible with invariance under FDiffs and renormalizability by power counting. For appropriate choice of parameters the unique propagating scalar mode possesses stable dispersion relation in the entire range of spatial momenta. Moreover, the dispersion relation has nice properties both at low and high momenta. In the first case it is linear, ω 2 ∝ p 2 , implying that the mode remains stable in any sufficiently smooth background 33 . On the other hand, the asymptotic form ω 2 ∝ p 6 of the dispersion at high momenta is compatible with the anisotropic scaling postulated in the UV. Therefore, the presence of this mode does not spoil the powercounting arguments which strongly suggest that the model is renormalizable. This implies, in particular, that strong coupling is avoided in the model by construction provided the coupling constants are chosen small enough.
Encouraged by these results we studied some phenomenological aspects of the healthy model. Making use of the Stückelberg formalism we have demonstrated that at low energies the theory reduces to GR interacting with an additional scalar field -the khronon. The time-dependent VEV of the khronon field breaks Lorentz invariance down to arbitrarily low energies. We observed that the structure of the low-energy limit of the model is similar to that of the Einstein-aether theory [30,31], even though the two theories are not exactly equivalent. The difference is due to transverse modes present in the Einstein-aether theory and absent in the healthy model of this paper. The transverse modes do not affect the form of homogeneous isotropic cosmological solutions and spherically symmetric solutions. This allowed us to directly apply to our case the bounds on Einstein-aether theory coming from the expansion history of the Universe, as well as to conclude that all but two PPN parameters in the healthy model coincide with those of GR. We have calculated the values of the remaining two PPN parameters, α P P N 1 , α P P N 2 , which characterize preferred frame effects. Current observational bounds on α P P N 1 , α P P N 2 significantly constrain the parameter space of the healthy model (but do not rule it out). Combined with the requirement that the theory is weakly coupled these bounds translate into the upper bound on the scale of quantum gravity M * 10 15 ÷ 10 16 GeV. Within the simplifying assumption that M * coincides with the scale of Lorentz symmetry breaking in the matter sector, astroparticle data constrain it from below: M * 10 10 ÷ 10 11 GeV. [We stress though that this lower bound must be taken cautiously: its validity depends on the details of the matter sector.] Thus the net result of our phenomenological study is that the healthy model can be tested experimentally using existing techniques, which we believe makes this model quite attractive.
The above results suggest that the healthy model of this paper can serve as a starting point for construction of a viable renormalizable theory of gravity. Admittedly, a lot of open issues remain. Let us discuss some of them.
The renormalizability of the model is yet to be demonstrated beyond the naive power counting. In principle, this amounts to an explicit analysis of the loop corrections to the action of the form (2.5), (2.6). In practice though, this may turn out to be a formidable task, given the large (of order 100) number of allowed terms in the Lagrangian and the complications related to the gauge symmetry of the action under FDiffs. A particularly subtle issue is a proper treatment of the time-reparameterization invariance. For consistency of the model this symmetry must be free of anomalies. Otherwise the symmetry would be reduced to that of the model of Sec. 4, with the corresponding re-appearance of pathologies. An interesting development along these lines is the renormalization of the energy-momentum tensor of test fields in curved backgrounds in Hořava-type theories [58]. It is shown that, in contrast to the relativistic case, this does not require counterterms with more than two time derivatives.
Even if the model proves to be renormalizable, there will be a question if it is UVcomplete or not. In other words, if it possesses a weakly coupled UV fixed point. An answer to this question requires the study of the renormalization group (RG) flow of the theory. Let us mention in this connection interesting recent results [59] about quantum electrodynamics (QED) in five space-time dimensions. It is demonstrated that 5d QED can be UV completed within Lorentz-violating framework by adding operators with higher spatial derivatives. The resulting theory is weakly coupled at all energies, possesses a weakly coupled UV fixed point with anisotropic scaling exponent z = 2 and flows to the usual 5d QED in IR (though, in general, with different velocities of photons and electrons).
It will be interesting to assess the quantum properties of the healthy model beyond perturbation theory. A possible approach to this difficult problem would be to use the canonical formalism. It is worth noting that inclusion of the terms with derivatives of the lapse N into the potential (2.6) significantly improves the canonical structure of the theory compared to the original non-projectable version of the proposal. Indeed, due to these terms the lapse is no longer a Lagrange multiplier and the analysis [22,24] unveiling the pathological structure of the constraints in the original Hořava model does not apply. Instead, in the model of this paper the Hamiltonian constraint obtained as the variation of the action with respect to the lapse, H ≡ δS δN = 0, has non-vanishing Poisson bracket with the equation π N = 0, where π N is the canonically conjugate momentum for N . In other words, these two constraints form a second class pair and must be used to eliminate the variables N , π N from the phase space [25,60]. Presence of second class constraints is an interesting feature of the model and its implications both for classical and quantum dynamics of the theory deserve a detailed study.
The non-perturbative dynamics of the model can be also addressed semiclassically by developing perturbative expansion in the backgrounds of classical non-linear solutions. Clearly, the first step here is to find explicitly such solutions. Let us point out in this context the existence of cosmological solutions which are free from the initial singularities 34 [61,62,63]. It will be interesting to understand if these solutions are stable. Another important issue is the structure of black hole solutions in the healthy model. Given that the high frequency modes in this model propagate with arbitrarily high velocities (moreover, as discussed in Sec. 5.3, the theory involves instantaneous interactions), it is a priori unclear if the notion of black hole as a region of space surrounded by an event horizon makes sense. On the other hand, there certainly must exist some low-energy notion of horizon defined using low-frequency modes which have finite velocity. It will be interesting to understand the physical meaning of black hole entropy associated to this low-energy horizon. The situation is even more intriguing because generically, due to the violation of Lorentz invariance, different low-energy species will propagate with different velocities thus giving rise to a number of nested horizons. As discussed in [64,65] one would expect this to lead to break down of black hole thermodynamics, which on general grounds implies violation of unitarity of the underlying theory. Thus it is vital for the proposal to understand if and how this paradox is resolved. It is worth mentioning that some spherically symmetric solutions of the healthy model have been reported in [66].
From the phenomenological perspective the major challenge for the non-relativistic gravity framework is a mechanism for emergence of Lorentz invariance in the matter sector at low energies. As discussed above this must happen with very high accuracy to satisfy existing experimental bounds. This poses a severe fine-tuning problem which is aggravated by the fact that Lorentz-violating parameters, such as differences of velocities of various matter species, run with the energy scale [67,68]. A promising mechanism to avoid fine-tuning is related to supersymmetry. Indeed, it has been argued [42,43] that within Lorentz-violating extensions of the Minimal Supersymmetric Standard Model (MSSM) it is impossible to write any Lorentz-violating operator of dimension less than 5. Dimension 5 Lorentz-breaking operators can be further forbidden by imposing discrete symmetries, e.g. CPT; in this case Lorentz breaking starts at dimension 6. In other words, given supersymmetry, Lorentz invariance emerges as an accidental low-energy symmetry. Eventual supersymmetry breaking gives rise to Lorentz-violating effects, but these are suppressed by the small ratio of the soft supersymmetry breaking masses m sof t to the scale M * of the higher-order Lorentz-violating operators. For example, the coefficients of dimension-4 Lorentz-violating operators generated from the operators of dimension 6 upon supersymmetry breaking are of order (m sof t /M * ) 2 . These are comfortably within the experimental bounds for m sof t ∼ 1 TeV, M * ∼ 10 15 GeV. Interestingly, supersymmetry also suppresses the contributions with higher powers of momentum in the dispersion relations of matter particles [42,43] which weakens significantly the lower limits on M * coming from astrophysical observations. Needless to say, realization of this scenario for the emergence of Lorentz invariance requires supersymmetrizing the non-relativistic gravitational action (2.5), (2.6). The first step in this direction would be to supersymmetrize the low-energy limit of this action, Eq. (5.27). This essentially amounts to constructing a supersymmetric action for the khronon sector as the gravitational part coincides with GR and can be made supersymmetric in the standard way.
The present work makes only first steps in the study of observational consequences of the healthy model. Many topics which have not been touched in this paper deserve a thorough investigation. These include implications of the model for emission and propagation of gravitational waves, dynamics of binary pulsars, spectrum of CMB perturbations and structure formation. Interesting results about the dynamics of the cosmological perturba-tions in the healthy model were reported recently in [69,70].
We plan to return to some of the above questions in future.
Note added. When this paper was in preparation the article [71] appeared which also uses the Stückelberg formalism to study the properties of the healthy model. Where there is an overlap, our results agree. However, Ref. [71] pushes the analysis using the decoupling limit to the regime of high energies (higher than M * ) and argues that the terms with the higher spatial derivatives do not resolve strong coupling of the scalar mode. This seems to contradict the conclusions of the present paper. In fact, pushing the decoupling limit to high energies, as done in [71], corresponds to going beyond its range of applicability. Indeed, the limit considered in [71] corresponds to throwing away all the higher-derivative mixings between the khronon field and the metric. On the other hand, as admitted in the published version of the [71], these terms are crucial to ensure the correct UV behavior ω 2 ∝ p 6 of the dispersion relation of the scalar mode, required by the anisotropic scaling with z = 3 and satisfied in the full theory, see Eq. (5.5). Instead, for the decoupling limit considered in [71] the UV asymptotic of the dispersion relation for the khronon field is ω 2 = const. Thus the loop integrals over spatial momenta instead of being suppressed become even more divergent than in the relativistic case. Clearly, this only represents a failure of the limit considered in [71] and does not imply the inconsistency of the theory.

A. Proof of the no ghost theorem
In this appendix we prove the proposition formulated at the end of Sec. 2.2. For simplicity, let us first consider the case when the Lagrangian depends only on the normal vector u µ and its first derivatives. These derivatives enter into the Lagrangian through the extrinsic curvature K µν and the acceleration vector a µ , This case covers all terms in the general Hořava-type Lagrangian (2.5) except those involving spatial derivatives of the 3-dimensional tensor R ij or of the acceleration a i .
We start by expanding u µ up to quadratic order, Here bar refers to the background values of the fields. The crucial observation is that in the background (2.20) the first variation (A.1) of the normal vector does not contain time derivatives of χ, while the second variation (A.2) contains only first time derivative. The next step is to consider variations of the extrinsic curvature and acceleration One observes that in the background (2.20) δ (1) K µν does not contain time derivatives of χ; δ (1) a µ , δ (2) K µν contain one time derivative, and δ (2) a µ -two time derivatives. The quadratic Lagrangian consists of two types of terms: the terms containing the first variations of the fields squared, and the terms linear in the second variations. For example, writing down explicitly the terms with the variation of the acceleration we obtain Clearly, both these terms contain at most two time derivatives of χ. In the case whenφ does not satisfy equations of motion one has to consider also the linear variation of the Lagrangian L (1) . Using the same reasoning as before one concludes that L (1) contains at most first time derivative of χ. This completes the proof. It is straightforward to generalize the above proof to include the dependence of the Lagrangian on higher spatial derivatives. In covariant language a spatial derivative translates into the operator P • ∇ • P , which is purely spatial (i.e. does not introduce further time derivatives). For example, the object ∇ i a j takes the covariant form P λ µ ∇ λ (P ρ ν a ρ ). Using the same reasoning as before one can show that the first (second) variation of this type of objects contain at most one (two) time derivatives of χ.
Note that the above proof uses in an essential way the invariance of the action under the transformations (2.15): this forces the Lagrangian to depend on the khronon field only via u µ and its derivatives. Recall that this invariance stems from the invariance of the original action (2.5) under FDiffs. As a consequence the proof also goes through for the projectable model (2.2) which obeys this symmetry. Indeed, as pointed out in Sec. 2.1, the projectable model is recovered from the theory (2.5), (2.6) by taking the limit α → ∞. Thus the order of time derivatives in the covariant equations for these two theories coincide.
On the other hand, the proof does not apply to models with reduced symmetry, such as the theory (2.7), because in these cases the Lagrangian contains additional dependence on the ϕ gradients. We will see in Sec. 4 that this makes the equation of motion for the khronon fourth-order in time; as one can anticipate, this will lead to certain pathologies of the theory.

B. Two faces of the tachyonic ghost
The purpose of this appendix is to clarify the following puzzle. In Sec. 4 we have shown using the Stückelberg formalism that the RFDiff model includes in its spectrum a tachyonic ghost. On the other hand, in the unitary gauge the action (2.7) of the model contains only two time derivatives and according to the standard lore one would not expect any ghosts in this picture. This seems to contradict the results of the Stückelberg analysis. We are going to show that the contradiction is removed when one properly formulates the physical questions to assess the effects of the tachyonic ghost.
We start by writing the unitary gauge Lagrangian for the sector of scalar perturbations around Minkowski background. Substituting the decomposition (2.9) in the action (2.7) and integrating out the non-dynamical fields B and E one obtains instability also has a negative kinetic term (tachyonic ghost). The standard definition of the ghost as "a field whose Hamiltonian is not positive definite" is not very useful: the Hamiltonian is not sign-definite in both cases. Moreover, the two cases can be related by a canonical transformation. To illustrate this point consider a toy Lagrangian representing a single Fourier mode of a tachyonic ghost, The corresponding Hamiltonian reads where π is the canonically conjugate momentum for η. The canonical transformatioñ π = pη ,η = π/p casts this into the Hamiltonian of a simple tachyon, One may still try to distinguish if a mode is ghost or not by the sign of the residue at the pole in the one-particle exchange amplitude considered as function of ω 2 ; this must be negative for the mode to qualify as a ghost. However, in the case of the tachyonic ghost the sign depends on the type of the coupling to the source used to define the amplitude. Taking again the toy model (B.3) as an example consider two couplings: source = ηΣ 1 and L (2) source =ηΣ 2 .
In the first case the one-particle exchange amplitude reads A Σ 1 ∝ Σ 1 −1 ω 2 + p 2 Σ 1 and the residue is negative. However, in the second case −ω 2 ω 2 + p 2 Σ 2 and the residue becomes positive: −ω 2 = p 2 > 0. Note that this ambiguity is related to the instability of the mode; it is absent for the case of a ghost with a stable dispersion relation when the pole lies at positive ω 2 . Thus we conclude that for a mode with gradient instability there is no unambiguous way to tell if it is a ghost or not. Rather, this notion makes sense only for a given coupling of the mode to the source.
It is instructive to trace explicitly the agreement of one-khronon exchange amplitudes calculated in the Stückelberg and unitary gauge pictures. To this aim we need to specify the source. This is easier to do on the Stückelberg side where we couple the khronon to a scalar field Σ. The khronon field must enter with derivatives so we write, In terms of the khronon perturbations this takes the form S source = d 4 x (−χ + ∆χ) Σ .
From this expression and the quadratic action (4.5) one reads off the khronon exchange amplitude One observes that this amplitude is a sum of two contributions, with positive and negative residues at the poles. In this sense one of the khronon modes is indeed a ghost.
Let us now see how the same amplitude is recovered in the unitary gauge. Fixing ϕ = t and substituting the decomposition (2.9) into the source term (B.4) we obtain The next step is to integrate out the non-dynamical fields B, E. Importantly, this produces a contribution into the Lagrangian which is quadratic in Σ. Omitting the higher-derivative terms the resulting Lagrangian reads where for simplicity we have taken the limit M P → ∞ (this corresponds to the decoupling limit in the Stückelberg picture) and introduced the rescaled variablẽ ψ = (2M 2 P /M 2 λ ) ψ .
This Lagrangian leads to the following propagators: Finally, the khronon exchange amplitude reads, Combining everything together one obtains the result (B.5). Thus we find that the unitary gauge calculation leads (in a quite non-trivial way) to the same result as in the Stückelberg formalism, as, of course, it should be. Technically, in the unitary gauge the ghost pole appears due to the structure of the coupling to the source.

D. Spherically symmetric solutions in Einstein-aether and khrono-metric theories
In this appendix we demonstrate that spherically symmetric solutions of the khrono-metric model (5.27) are identical to those of the Einstein-aether theory. Let us consider the equation of motion for the khronon field coming from varying the action (5.27) with respect to the field χ ∇ µ J µ = 0 , (D.1) where At the same time the equation of motion for the aether is obtained by varying (5.27) with respect to the field u ν and reads In deriving this equation one has to take into account the constraint u µ u µ = 1: it leads to the appearance of the projector P µ ν on the l.h.s. Finally, the energy-momentum tensor appearing in the Einstein equations for the khrono-metric theory can be written as where the second term comes from the explicit dependence of the vector u µ in the khronon theory on the metric, see Eqs.(2.13), (2.14). This coincides with the energy-momentum tensor of the aether. To obtain the second term in this case one again has to take into account the constraint u µ u µ = 1 [34]. Any spherically symmetric configuration of aether is automatically hypersurface-orthogonal implying that any spherically symmetric solution of the Einstein-aether theory is a solution for the khrono-metric theories [34]. The converse is less obvious as the equation of motion (D.1) of the khronon field contains an additional derivative compared to the aether equation (D.2) and thus, a priori, admits more solutions. However, for spherically symmetric configurations (D.1) implies (D.2). Indeed, the current J µ obeys the relation Hence in the unitary gauge (2.12) its component J 0 identically vanishes imlying that the corresponding charge is identically zero. On the other hand, the time derivative of Q is equal to the flux of the spatial component J i of the current through the 2-sphere at spatial infinity 35 . In spherically symmetric situation this implies that the current itself is zero which brings us to the equation (D.2). Combining this result with the equality of the energy-momentum tensors we conclude that spherically symmetric solutions in the khronon and Einsteinaether theories are indeed identical.
E. PPN parameters α P P N 1 , α P P N 2 for the healthy model In this appendix we derive the formulas for the PPN parameters α P P N 1 , α P P N 2 in the khronon theory (5.27). We assume the couplings α, β, λ ′ to be small and perform calculations to the leading order in these couplings. We consider the metric produced by a point source of mass m in its rest frame. This frame does not coincide with the frame defined by the preferred foliation. Hence the background valueφ of the khronon field in this frame differs from the coordinate time. Using the reparameterization symmetry (2.15) we fix where v i is the velocity of the source with respect to the preferred frame. This corresponds to the background value of the vector u μ The source perturbs the metric and the khronon. One writes, g µν = η µν + h µν , ϕ =φ + χ .
To the leading order h µν is given by the standard Newtonian expressions, Our goal is to find corrections to (E.1) in powers of v.
It is convenient to introduce the following notations Here we have imposed the harmonic gauge, and have used the fact that the metric is static. Note that the first order in the post-Newtonian approximation requires to find h 00 , h 0i and h ij components of the metric with the accuracy O(v 2 ), O(v) and O(1) respectively. This implies that we need to determine T µν χ and T χ to order O(v 2 ), T 0i χ -to order O(v), and T ij χ -to order O(1). Expanding up to terms O(v 2 ) we obtain: where we have used that derivatives act on static configurations. Substituting these expressions together with (E.1) into (E.3) we find where on the l.h.s we have neglected terms O(v 2 ∆ 2 χ) as they are of higher order in v.
From (E.9) we find meaning that χ is of order O(v). Let us estimate the orders of the variations (E.7). One finds by inspection that δSχ δΦ is of order O(1), δSχ δV λ -at most of order O(v), δSχ δH λρ -at most of order O(v 2 ). This allows to simplify the khronon energy-momentum tensor. To the required orders we have: One evaluates the variations appearing in these formulas, Inserting (E.11) into (E.10) and substituting the result into (E.8) we find the equations for the first order corrections to the metric, ∆h It is straightforward to solve these equations for the explicit form (E.2) of the function φ(r). Note that the second term on the r.h.s. of (E.12b) can be removed by a timeindependent gauge transformation. Combining the result with the Newtonian expressions (E.1) we obtain the metric (5.33) with G N = 1 8πM ′2 P 1 + α 2 (E. 13) and the PPN parameters α P P N 1 , α P P N 2 quoted in (5.34). The expression (E.13) coincides with (5.30) to linear order in α.