{'fraction_non_alphanumeric': 0.03444492469328723, 'fraction_numerical': 0.010705854972237923, 'mean_word_length': 4.42630344578748, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
On stability of Friedmann-Lema\^itre-Robertson-Walker solutions in doubled geometries

Motivated by the models of geometry with discrete spaces as additional dimensions we investigate the stability of cosmological solutions in models with two metrics of the Friedmann-Lema\^itre-Robertson-Walker type. We propose an effective gravity action that couples the two metrics in a similar manner as in the bimetric theory of gravity and analyse whether standard solutions with identical metrics are stable under small perturbations.

2 ical bimetric theory (none of the metrics can be thought of a background metric) the obtained potential is much similar to the bimetric case (though it is expressed as a rational function and not a polynomial in the eigenvalues of the metrics ratio). Moreover, the symmetric coupling to the matter and radiation makes it closer to the symmetric bimetric theory, where both metrics couple (in the same way) to matter and radiation.
The paper is organized as follows: we present the assumptions of our model (the structure of the two-sheeted geometry) and the methods of deriving the leading two terms of the spectral action using the pseudodifferential calculus and the Wodzicki residue. After computing the Euclidean action functional for flat as well for curved geometries we perform the Wick rotation and obtain a set of nonlinear differential equations for the four functions that describe the model. In the rest we focus on the stability of the symmetric solutions, which are the standard Friedmann-Lemaître-Robertson-Walker geometries for both sheets and analyse small perturbations for the typical cosmological solutions of dark-energy, matter and radiation-dominated universes. In the last section we briefly discuss the possible physical consequences and argue why the model is physically viable.

II. ALMOST COMMUTATIVE FRIEDMANN-LEMAÎTRE-ROBERTSON-WALKER MODELS
A. Almost-commutative geometries The Gelfand-Naimark equivalence between topological spaces and commutative C * -algebras was further enriched by A. Connes in order to include noncommutative algebras and also to describe more than only the topology. In his formulation of noncommutative geometry [13] the crucial role is played by a spectral triple which is a system (A, H, D) consisting of an unital * -algebra A, Hilbert space H and a Dirac type selfadjoint operator acting on H. Usually, more additional structure is assumed (e.g. the existence of grading-type operator γ, and an anti-unitary operator J, called real structure) and further compatibility conditions between all these elements. The canonical commutative example is (C ∞ (M ), L 2 (M, S), D M ), where M is a manifold equipped with a spin structure, L 2 (M, S) is the Hilbert space of square-integrable spinors, and D M = iγ µ (∂ µ + ω µ ) is the canonical Dirac operator expressed in the terms of the connection ω µ on the spinor bundle.
It turns out that crucial from the applications in particle physics point of view are triples with algebras being tensor products of the above one which some finite-dimensional matrix algebras A F . The Hilbert space is the tensor product of L 2 (M, S) with some finitedimensional Hilbert space H F on which A F is represented, and its dimension determines the number of fermionic degrees of freedom in the theory. Grading operators and real structures are also composed in an appropriate way in order to define analogous objects on the resulting triple. The Dirac operator, however, is not just the simple tensor product of D M and D F , but has the following form: The resulting triple forms the so-called almost-commutative geometry and have been the backbone of multiple models applied to the physics of elementary particles (see [14,15]). The starting point to consider physical models based on spectral triples is the spectral action.
Its bosonic part is given by where Λ is some cut-off parameter and f is some smooth approximation of the characteristic function of the interval [0, 1]. In the case of particle physics models it reproduces the bosonic part of the Lagrangian of such theories minimally coupled to gravity, together with the standard Hilbert-Einstein action for the metric.
B. The classical geometry We consider geometries described by the generalized Friedmann-Lemaître-Robertson-Walker metric, and a(t), b(t) are positive (sufficiently smooth) functions. The orthogonal coframe {θ a } for ds 2 is defined so that ds 2 = θ a θ a . It allows us to immediately compute the spin connection ω which is determined by dθ a = ω ab ∧ θ b . Then, the Dirac operator is, in a local coordinates, given by where γ a 's are gamma matrices chosen to be antihermitian and so that γ a γ b +γ b γ a = −2δ ab I. Instead of the original Dirac operator we can equivalently analyse the operator, which is conformally rescaled, D h = h −1 Dh, with the scale factor h(t) = a(t) −3/2 b(t) −1/2 . This assures that we can work with the Hilbert space of spinors, where the scalar product does not depend on a(t) and b(t).

C. The two-sheet almost commutative model
We consider a generalised almost-commutative geometry, which is described by a productlike spectral triple of the spectral triple over the manifold with the Friedmann-Lemaître-Robertson-Walker metric and the triple over two points. However, instead of the usual product Dirac operator, we take a more general one, where D 1 , D 2 are both of the form (II.5), yet with possibly different scaling functions a and b, and Φ being a priori a field (which can be later restricted to be constant).

4
The choice of the full Dirac operator with the γ in the off-diagonal part is motivated by the fact that in the case of D 1 = D 2 it yields a usual almost-commutative product geometry. Note that, in principle one can study generalized objects with arbitrary order-zero operators on the off-diagonal of D, so the only thing we require of γ is that it anticommutes with γ a matrices and is not necessarily the chirality grading operator of the Euclidean spin geometry of the manifold M . In order to have the full Dirac operator hermitian we must require that γ is hermitian and, consequently, we have to normalize γ 2 = 1. However, we shall relax this assumption and consider also models with γ 2 = −1. This allows for much more flexibility, in particular, for the models that are derived from higher-dimensional Kaluza-Klein type geometries and would lead to some more realistic effective physical situations. One of the interesting possibilities is that when passing to the Lorentzian signature for the manifold M we can as well choose the Lorentzian signature for the discrete degrees of freedom. This possibility has been discussed for finite geometries, albeit in a different context of the Standard Model in [16], in the natural language of Krein spaces. What is important for our consideration is that the only difference will be that the operator D will be only Krein self-adjoint, meaning that γ will be antiselfadjoint and γ 2 = −1. To accommodate for both possibilities in the discrete degrees of freedom we do not fix γ 2 and we allow that (after normalization) To simplify the presentation in the paper we introduce the following matrices: D. The spectral action For the geometry described by a given Dirac operator D the main object of interest is the Laplace-type operator D 2 , which is a second-order differential operator acting on the sections of the doubled spinor bundle. Its symbol σ D 2 (x, ξ) consists of three parts a 0 +a 1 +a 2 , each of a k (x, ξ) being homogeneous of degree k in ξ's. Then we compute the symbol of its inverse, where b k (x, ξ) is homogeneous of order −2 − k in ξ (we briefly review the mathematical details of how the computations of the symbols are performed in the Appendix A) and use it to compute the first two terms of the spectral action for the considered model. It can be expressed in terms of Wodzicki residua [17,18] as, Λ 4 Tr Tr Cl b 2 0 + cΛ 2 Tr Tr Cl b 2 , (II.9) where Tr Cl denotes the trace performed over the Clifford algebra and Tr is the trace over the matrices M 2 (C) that are used in the mild noncommutativity introduced in the model.

E. Flat geometries.
Although the topology of the flat case in physics is not exactly toroidal, from the point of view of local behaviour it is identical to such, which was already analysed for b = 1 in [12]. In this section we generalize those results to the case with arbitrary function b(t), so we consider here toroidal Friedmann-Lemaître-Robertson-Walker geometries described by the following metric in the coordinate system (t, x) = (t, x 1 , x 2 , x 3 ): (II.10) Hence an orthogonal frame for ds 2 is of the form, while the matrix of connection 1-forms is, As a result, the (single) Dirac operator takes the following form, and after the conformal rescaling so that the full Dirac operator acting on the doubled Hilbert space of spinors is, The resulting Laplace-type operator in this model is of the following form: (II. 16) The symbol σ(D 2 ) = a 0 + a 1 + a 2 is given by where we denoted by ξ 2 = ξ 2 1 + ξ 2 2 + ξ 2 3 . Now, computing the symbol of D −2 using the prescription presented in Appendix A, we obtain b 0 (D) and b 2 (D). Then, taking the trace over the Clifford algebra and the matrices M 2 (C), and integrating over the cosphere bundle |ξ| 2 = 1, we compute the Wodzicki residue that gives us the Euclidean spectral action of the considered model. The final result is, where the scalar curvature for the flat spatial geometry is,

F. The non-flat case
In this subsection we concentrate on the case with positive (k = 1) curvatures, with the negative (k = −1) case that can be treated in a similar manner. Although the effective Lagrangian and the equations of motion are local and hence the dynamical terms are expected to be unchanged, we derive them explicitly using appropriate coordinates. For the case of k = 1 we use the spherical coordinates (t, χ, θ, φ), so that the metric is then described by: (II.20) The orthogonal frame is given by Therefore the only nonvanishing components for the spin connection ω are [19,20]: Now, for the Dirac operator we get explicitly 24) where in this case After the conformal rescaling by using h(t) = a(t) −3/2 b(t) −1/2 we end up with the following Dirac operator Therefore, for the doubled model that we are considering, the Dirac operator is As a result we have (II.28) In order to compute its symbol σ D 2 = a 2 + a 1 + a 0 we first notice that the symbol of D 2 3 is given by: As a result, for the operator D 2 , we have Using the prescription presented in Appendix A, we first compute the symbols σ D −2 = b 0 + b 1 + b 2 + ..., then we proceed, in an exactly similar manner as in the case of the toroidal geometry, to compute the spectral action. The result is, where now R(a, b) denotes the scalar of curvature for spherical spatial geometries, Note that the action functional differs (II.30) from (II.18) only through the last term that arises from the scalar curvature of the spherical spatial geometry, where a relevant term depending on k = 1 is added. The above result can be generalized to the negative curvature case (we omit straightforward but tedious parametrization and computation of symbols). In fact, taking R(a, b) as R(a, b, k), depending on the space curvature k, we have a general action functional for all geometries in the doubled spacetime Friedmann-Lemaître-Robertson-Walker models.

G. The interactions of the metrics
Before we pass to the equations of motions an their stability, let us briefly compare the effective potential describing the interaction between the two metrics to the bimetric gravity models [1,4,21]. Certainly, apart from the fact that we have an action for two metrics, there is a much deeper symmetry between the two, since neither plays a role of a "background" metric. In fact, the usual solution in the case of vanishing α gives both metric totally independent of each other. Introducing the variables, which depend only on the entries of the matrix X a c = g ab 2 g 1 bc , we can express the interactions between the metrics as proportional to: where the function V g −1 2 g 1 is of the form: which can be efficiently expressed as a rational function of the symmetric polynomials in √ X.
We stress that the resulting model possesses features that are characteristic to bimetric gravity models: the potential V depends on the metrics only through √ X and satisfies (II.32). On the other hand, in the usual bimetric models such potential is a polynomial in eigenvalues of √ X rather than a rational function. It was proposed in [22] that the construction presented here might results in the derivation of bimetric theories out of the geometric data. The above result suggest that indeed this class of models resembles some characteristics of bimetric gravity models, but is a different one. We postpone for the future research the detailed analysis of these differences and their cosmological implications.

H. The equations of motion
The action functional (II.30) depends on the field B only via b 1 and b 2 but not their derivatives. As a result, b 1 and b 2 are not dynamical and its Euler-Lagrange equations give rise to the constraints of the model. Moreover, due to the reparametrization invariance we can fix one of these functions or relate them with each other.
Furthermore, the action functional was derived for the Euclidean model and to pass to physical situation we need to perform Wick rotation, as described in the [12]. In our case, this will affect only the square of the time derivative of the scaling factors a i (t), which will change signs. Consequently, the action and the equation of motion for the rest of this paper are in the Lorentzian signature of the metric (−, +, +, +). Let us remind that the discrete degrees of freedom of the geometry might be Riemannian or pseudo-Riemannian, which results in the appropriate choice of the sign κ.
After integration by parts and omitting the boundary terms that are full derivatives in t, we obtain the following action for the pure gravity Friedmann-Lemaître-Robertson-Walker doubled geometries for the Lorentzian signature and arbitrary spatial curvature k, where we have factored out the overall constant so that the dynamical term appears only with a numerical factor, denoted the effective cosmological constant by Λ e and introduced the effective coupling between the two metrics by α: Unlike the bare cut-off parameter Λ, here the effective cosmological constant can vanish or be negative for a particular model. We shall use the above convention with Λ e and α throughout the rest of the paper.

10
The four Euler-Lagrange equations take the following form, In the above equations we use the convention that (i, i ) = {(1, 2), (2, 1)}, and H b,j = ∂ta j a j b j are the generalized Hubble parameters. Before we analyse the inclusion of matter fields and possible solutions, let us observe that in the flat k = 0 case, inserting Λ e from first two equations in last two, one obtains with some rational function L(a 1 , a 2 , b 1 , b 2 ), so in particular, whenever a 1 b 2 = a 2 b 1 both H b,1 and H b,2 must be constant.

III. INTERACTION WITH MATTER FIELDS AND RADIATION
The equation of motion derived in the previous section describe the empty universe in the doubled model. Here, we can ask how they are modified by the presence of the matter fields. The crucial point is to see how the effective matter and radiation action depend on the components of the metrics described in terms of fields a 1 , a 2 , b 1 , b 2 . The main difficulty is the passage from the microscopic action for spinor and gauge fields to the effective averaged energy-momentum tensor in the Einstein equations.
The microscopic action for the spinor fields in the doubled universe will be the usual fermionic action ΨDΨ. Since both components of the spinor couple to the respective Dirac operators D 1 and D 2 on each of the single sheets separately, and the Ψ field is, by assumption, independent of the metric fields, we conclude that the resulting action will be split into separate actions that do not mix the metric components on each of the single universes. Similar argument can be used for the radiation energy-momentum tensor that originates from the gauge fields over the considered model. As the model has two U (1) symmetries there are two gauge fields that couple to the Higgs field. A linear combination of them will become a massive one, due to spontaneous symmetry breaking of the Higgs field, whereas another linear combination will correspond to the massless photons. Again, the effective Yang-Mills action for the photon field will not mix the metric components over the two sheets and therefore we shall have independent tensor-energy components for each equation.
These heuristic arguments suggest that the effective equations of motion are modified by the respective components of the overall energy-momentum tensor T 0 0 and T 1 1 , which depend separately on a 1 , b 1 and a 2 , b 2 , As in the conventional cosmology we consider the model of the perfect fluid, i.e. the stress-energy tensor is taken to be of the form where P is refered to as pressure, while ρ is called energy density. For the generalized Friedmann-Lemaître-Robertson-Walker metric, the vector u µ is 1 b(t) , 0, 0, 0 , so that u µ u µ = −1. As a result, T 0 0 = −ρ and T 1 1 = P . Furthermore, the continuity equation ∇ µ T µν = 0 reduces to the standard one: We assume that the thermodynamics of the matter content is characterized by the following equation of state: From the continuity equation we immediately infer that then where η is the proportionality constant, exactly as in the standard cosmology. The resulting Einstein equations for the double-sheeted universe are of the following form: We stress that the above model is a straightforward generalization of the classical one for the doubled theory, with the only difference that we allow two different scaling factors and the interaction between them derived from the spectral action. Indeed, for a 1 = a 2 ≡ a, b 1 = b 2 = 1, or α = 0, equations of motion reduces to the usual Friedmann equations yielding the well-known solutions.
In what follows we shall aim to analyse the possibility of small perturbations of the classical solutions of Friedmann-Lemaître-Robertson-Walker models, trying to answer the question whether the double-sheeted universe is stable. However, before we start the computations to see when the small perturbations of the classical solution are possible, let us observe that thanks to the reparametrisation invariance of the time variable in the equations (II.34), (II.35), we can decide to fix either b 1 or b 2 or relate them with each other. There are, in principle, many choices of the possible parametrizations and we choose a particular one, which is motivated by the existing symmetry with respect to the exchange between left and right modes in the geometric setup we consider. We shall set b 1 (t) + b 2 (t) = 2, therefore, effectively, one can introduce a new function, and derive the equations of motion for a 1 (t), a 2 (t) and b(t). Taken an appropriate linear combination of the derivatives of (II.34) and the equations (II.35) we shall obtain three nonlinear first order differential equations for these functions. Despite the fact that a full analysis of these equations is complicated and can be done possibly only numerically, we can obtain some significant results.
We shall finish this section by a remark that one cannot a'priori assume that both lapse functions are identically 1. Indeed, we shall see that such solutions (in the linearised regime) are not possible. Moreover, a far more general argument holds also for the full equations in the case of the empty universe. Then, there are no solutions with b(t) = 0 apart from a 1 (t) = a 2 (t). The argument is quite simple and relies on algebraic manipulation of the equations (II.34). Indeed, assuming b 1 (t) = 1 = b 2 (t) and subtracting the equations we obtain the following relation between a 1 and a 2 , which is true only if a 1 = a 2 (since both a 1 , a 2 are positive functions) or 1 a 1 + 1 a 2 = const. The latter condition can be solved, and when used in either of the first two equations (II.34) it leads to the constant solutions for a 1 and a 2 . Therefore, the functional based on the action (II.33) has extremal points only if the time scaling factor differs for the two metrics, so b(t) = 0. We leave aside the interpretation of this observation and its potential physical consequences to see whether the solutions that differ from the standard ones allow physically feasible models.

IV. PERTURBATIVE SOLUTIONS
In what follows we study infinitesimal perturbations of the classical solutions of the Einstein equations in different scenarios like empty universe with a cosmological constant, with 13 and without curvature, the matter dominated (i.e. for w = 0) flat universe with a vanishing cosmological constant Λ e = 0 and the radiation dominated (i.e. for w = 1 3 ) flat universe. Our working assumption is that we look for small perturbations around the symmetric, product, geometry of the form, and linearise the equations of motion, taking the first terms in .
In the zeroth order, we obtain (from all equations, as expected): whereas the first order yields the following set of linear equations for r 1 , r 2 and s, for the function a(t), which already satisfies the equation (IV.2): where we have introduced Λ e = 6λ 2 for simplicity, and denote the time derivative by a dot. Note that for a given background solution a(t) we have a homogeneous equation for the sum r 1 (t) + r 2 (t), which has a simple solution that, however satisfies reasonable initial conditions r 1 (t 0 ) = r 2 (t 0 ) = 0 if and only if it is constantly 0. Therefore, we may freely restrict ourselves to the case r 1 (t) = r(t) = −r 2 (t), and final set of perturbative equations,

A. The empty universe
In the case of an empty, or dark-energy dominated universe, we have the simple case of η = 0 and cosmological solutions depending only on the curvature k and the cosmological constant Λ e .

De Sitter universe (k = 0)
The solution of (IV.2) is, 5) and the equations of motion for r, s are: (IV.6) Solving this system of linear equations we obtain, where On the other hand, we see that the perturbative solutions cannot be extended to −∞ as, independently of the value of the parameters, they then become much bigger than the de Sitter solution. This puts the limits of applicability of the perturbative expansion which is entirely consistent with the dark-energy dominated universe solutions. As a last remark we note that even independently of the value of α perturbations, which are decaying exponentially, are possible for certain values of initial parameters. For example, if at t = 0 we set, then C 1 = C 3 = 0 and the perturbations will be exponentially damped for all range of parameters.

Geometries with positive and negative curvatures k = ±1
We start with the easier case of negative curvature, for which the solution of (IV.2) is: and in what follows we choose t 0 = 0 to simplify the notation. It is convenient to change the variables and write the equations (IV.4) in τ = sinh(λt). Then we obtain, (IV.11) The above set of equations can be solved explicitly, and the solution for s(τ ) is given by, where 2 F 1 is the hypergeometric function, G 2,0 2,2 is the generalized Meijer's function [23] and Since the solution is of the Big-Bang cosmology type we shall look for the small t (small τ ) behaviour of solutions. Both functions are defined in the region τ 2 < 1 and can be extended analytically to the other values of τ 2 , yet τ 2 = 1 is the point at which they are discontinuous or singular. Additionally the Meijer's function has a pole at 0 of order at least 2 unless the parameter ζ is quantized, The possibility of having both solutions regular at τ = 0 means that there exists a nonzero perturbation of the standard solution, which has both perturbations vanishing at the initial time s(0) = r(0) = 0. However, the fact that τ = 1 is a singular point of the Meijer's function restricts the possibility of extending the assumed linearised perturbation beyond certain time frame. The long-time behaviour of the solutions that are arbitrary (not necessarily vanishing) at t = 0 is similar to the flat case and governed by value of ζ, with asymptotically vanishing solutions for the same range of parameters α, λ as in the k = 0 situation.
Finally, for the positive curvature, k = 1, the pure dark energy solution is, 14) and the small perturbations at t 0 = 0 are, again changing the variable to τ = sinh(λt), (IV. 15) which, similarly as in the previous situation, has the solutions that are expressed in terms of the hypergeometric function 2 F 1 : From the fact that in this case we deduce that if we require r(0) = 0 then c 1 = 0. One can easily check that then also s(0) = 0, however, both solutions will grow with t. On the other hand, the exponentially decreasing solution requires c 2 = 0.

B. Matter dominated universe
In a completely similar manner we consider the limit in a matter-dominated universe, in which we put Λ e = 0 and w = 0, while η = 0. We start with the Einstein-de Sitter universe, k = 0. In this case the standard solution, gives the following equations for r(t) and s(t), s(t). (IV.17)

V. SUMMARY AND OUTLOOK
The simplest almost-commutative geometry of the two-sheeted universe, which is motivated by the Connes-Lott idea [11] is an interesting model to study its potential relevance not only for the particle physics but also for its implication to the large-scale structure of the Universe. We have shown that an abstract model, with a more general type of metric structure that is not necessarily a product structure allows a two-metric theory, which is very similar to the bimetric theory of gravity. Although we are aware that both the interaction structure as well as the interpretation of the model's origin are quite different there are striking similarities in the potential term of the action. It shall be noted that models originating from quantum deformations of spacetime have a similar feature of two metrics although their origins are different [25].
Leaving the full model that was developed for the particle interactions [8,9] aside and concentrating first on a simplified one, we have focused on a primary question of stability of classical Friedmann-Lemaître-Robertson-Walker solutions. To be more precise, our idea was to check whether for some range of parameters a small perturbation in the Dirac operators making the full one, and hence the metrics different from each other on the two sheets of the Universe, will diverge or collapse.
Our conclusion is that for the considered range of models, including flat and curved spatial geometries with dark-energy, radiation or matter dominance there exist a range of parameters so that the symmetric solution (product geometry) is dynamically stable. Our analysis confirms but hugely extends the earlier indications [12] by allowing both the scale factors as well as lapse functions to vary. The stability of the cosmological solutions suggests that the models with two metrics are admissible from the physical point of view and are an interesting modification of geometry that may be used in future models.
This has an important bearing on the physical consequences of the model. First of all, cosmological observables like redshift and observable Hubble constant will be related to the background standard Friedmann-Lemaître-Robertson-Walker solution. This follows from the fact that both light and matter will couple (as argued in section III) to both metrics and, taking into account that in most models the difference between metrics is decreasing as the Universe evolves, only the average (background) scale factor a(t) will determine the observable redshift. However, one can speculate that a possible sign of the fluctuating two metrics might be seen in physical effects that couple only to one metric (as might be the case of massless Majorana particles) or couple to metrics in a nonlinear way.
The constructed (simplified) model is predominantly based on the idea that allowed to explain the appearance of Higgs field and Higgs quartic symmetric-breaking potential from purely geometric considerations as a form of generalized gauge theory. Transferring this concept to the theory of metric and generalized general relativity appears to be a natural an well-motivated physical step. Unlike in the bimetric theory, here the interaction terms between the two metrics are completely determined by the structure of the theory yet are not computable in full generality. This prevents us from an analysis of the possibility of ghost-free sectors in the way it was done for bimetric theories [7,26] . Nevertheless, since the model has strong features similar to bimetric gravity (as we have stressed in II G), in particular, even though the effective interaction potential between the two metric is not a symmetric polynomial of g −1 1 g 2 but rather a rational function, where the nominator and denominator are of this form, we expect that a similar result will hold. Apart from the fundamental questions of physical consistency and interpretation of the degrees of freedom of the theory there are still several questions that remain open. First of all, in case of small deviations from the product geometry it is interesting whether they might have some observable physical consequences both in the pure gravity sector as well as in the sector of the matter and radiation. Though this might be considered as pure speculation, such fluctuations of the metrics, if existing in the radiation era, might be linked to some parity anisotropies [27] in the Cosmic Microwave Background radiation. Another possible sector of the theory to explore are solutions with singularities like black holes. All such ideas need to be explored carefully in future studies.