From (5) the effective energy density and pressure of dark energy can be expressed as 8 πGρ$_{DE}$ = 1 2 ( F ′ R − F ) − 3 H $^{2}$F ′ − 3 H ˙ F ′ = − 3 H ˙ R F ′′ + 3( H 2 + ˙ H ) F ′ − 1 2 F , (8) 8 πG ( ρ$_{DE}$ + P$_{DE}$ ) = 2 ˙ H F ′ − H ˙ F ′ + ¨ F ′ , (9) respectively, where R = 12 H $^{2}$+ 6 ˙ H . We define the DE equation of state parameter w$_{DE}$ by the ratio w$_{DE}$ ≡ P$_{DE}$/ρ$_{DE}$ . With the appropriate initial condition after cosmic inflation mentioned above, F takes an asymptotically constant value F = − λR$_{s}$ at high redshift (apart from the R $^{2}$/ 6 M 2 term which we neglect here). In this regime, evolution of the Universe is the same as that obtained from the Einstein action with a cosmological constant Λ( ∞ ) = λR$_{s}$/ 2. The scale factor therefore evolves as a = a$_{i}$ ( 16 πGρ$_{i}$ λR$_{s}$ ) 1 3 sinh 2 3 (√ 3 λR$_{s}$ 8 t ) ∼ = a$_{i}$ ( t t$_{i}$ ) 2 3 , (10) where the suffix i denotes quantities at an initial time t = t$_{i}$ . The time dependence of ρ$_{DE}$ is mainly governed by the first term in the right-most expression of (8) initially. Since ˙ R < 0 and F ′′ > 0 for stability, this means that the effective energy density of dark energy increases with time in this regime. Therefore, DE exhibits the phantom behaviour, w$_{DE}$ < − 1, during the matter-dominated stage with z > 1, which lasts only temporarily because the late-time asymptotic de Sitter stage has an effective cosmological constant smaller than Λ( ∞ ). So, ρ$_{DE}$ stops growing after the end of the matter- dominated stage and begins to decrease. Indeed, as shown in Ref. [14], the late-time asymptotic de Sitter solution has a curvature R ≡ R$_{1}$ ≡ x$_{1}$R$_{s}$ where x$_{1}$ is the maximal solution of the equation, λ = x (1 + x $^{2}$) n +1 2 [(1 + x $^{2}$) n +1 − 1 − ( n + 1) x $^{2}$] . (11) It satisfies the inequality x$_{1}$ < 2 λ , so that Λ( R$_{1}$ ) = R$_{1}$/ 4 < Λ( ∞ ). These inequalities are saturated in the limit n ≫ 1 for fixed x$_{1}$ , or x$_{1}$ ≫ 1 for fixed n . In these cases cosmic evolution is indistinguishable from the standard ΛCDM model. Thus, this model naturally realizes crossing of the phantom boundary w$_{DE}$ = − 1 in a recent epoch. Note that phantom behaviour of DE is generic in its models based on the scalar-tensor gravity [23] which includes the f ( R ) theory. Here we see that it is realized in all simplest stable f ( R ) models of present DE. 5