# Dynamical analysis of $f\left(Q\right)$-cosmology

Andronikos Paliathanasis<EMAIL_ADDRESS>Institute of Systems Science,
Durban University of Technology, Durban 4000, South Africa  Departamento de
Matemáticas, Universidad Católica del Norte, Avda. Angamos 0610, Casilla 1280
Antofagasta, Chile

###### Abstract

We study the evolution of the physical variables in $f\left(Q\right)$-gravity
for two families of symmetric and flat connections in a spatially flat
Friedmann–Lemaître–Robertson–Walker geometry where the equation of motion for
the nonmetricity scalar is not trivially identity. From the analysis of
dynamics we found that the de Sitter universe is always an attractor while the
cosmological models admit scaling solutions which can describe the early
acceleration phase of our universe or the matter and the radiation epochs.

Cosmology; $f\left(Q\right)$-gravity; symmetric teleparallel; dynamical
analysis.

##  I Introduction

Symmetric teleparallel theory Nester:1998mp  is a gravitational theory
equivalent to General Relativity where the Levi-Civita connection is replaced
by a torsion-free and flat connection and the fundamental scalar for the
gravitational Lagrangian is the non-metricity scalar $Q$. Indeed, General
Relativity, symmetric teleparallel theory, and teleparallel equivalence of
General Relativity form a family of theories known as trinity of gravity tr1 ,
in which the gravitational Lagrangian is a linear function of the Ricci scalar
$R$, the non-metricity scalar $Q$ and the torsion scalar $T$ respectively. The
equivalency between these three theories is lost when non-linear components of
the above scalars are introduced in the gravitational Lagrangian.

The $f\left(X\right)-$theories of gravity where $X$ is a geometric scalar,
have been proposed by cosmologists mod1  in order to explain the observational
phenomena od1 ; od2 ; od3 , like the late-time acceleration phase of the
universe attributed to the dark matter. In $f\left(X\right)-$theories of
gravity, the dark energy has a geometric origin sdc1 . Indeed, the new degrees
of freedom which are provided in the field equations by the nonlinear
components of the $f$ function define an effective geometrodynamic fluid which
can drive the dynamics in order to explain the acceleration of the universe.
When the three scalar of the trinity of gravity are used on the
$f\left(X\right)$-theories we end with the $f\left(R\right)~{}$Buda , the
teleparallel $f\left(T\right)~{}$Ferraro , and the symmetric teleparallel
$f\left(Q\right)$ f6  gravitational theories. $f\left(R\right)$ and
$f\left(T\right)$ gravitational models have been widely studied in the
literature for cosmological studies rt1 ; rt2 ; rt3 ; rt4 ; rt5 ; rt6 , and
for the analysis of astrophysical objects rt7 ; rt8 ; rt9 ; rt10 ; for more
details we refer the reader to the review articles rt11 ; rt12 .

In this work we focus on the cosmological application of
$f\left(Q\right)$-gravity in a isotropic and homogeneous spatially flat
Friedmann–Lemaître–Robertson–Walker background. In symmetric teleparallel
theory of gravity it is possible to separate the inertial effects from gravity
because through an appropriate coordinate transformation we can define the
coincident gauge where all the components of the flat connect vanish. Only the
recent years $f\left(Q\right)$-gravity has drawn the attention of
cosmologists. In ff1 ; ff2 ; ff3 ; ff4 ; ff4a  $f\left(Q\right)$-gravity has
been considered as a geometric dark energy model; furthermore in ff5  the
authors show that $f\left(Q\right)$-gravity could challenge the $\Lambda$CDM
model. An inflationary scenario in $f\left(Q\right)$-gravity is presented in
ff6 . Anisotropic cosmological models were investigated in an1 ; an2  while
inhomogeneous spacetimes were studied in Heis1 ; ww4 ; ww5 . For other recent
studies of $f\left(Q\right)$-gravity we refer the reader to nn1 ; nn2 ; nn3 ;
nn4 ; nn5 ; nn6  and references therein.

We perform a detailed analysis for the phase-space of the field equations in
$f\left(Q\right)$-gravity for different families of flat connections. These
different families of connections provide the same field equations in
symmetric teleparallel equivalence of general relativity, however for
nonlinear functions $f\left(Q\right)$ there are introduced dynamical degrees
of freedom from the connection. In order to understand the effects of these
dynamical terms in the cosmological parameters we determine the asymptotic
solutions. This approach has been widely applied in gravitational theories
with many important results df1 ; df2 ; df3 ; df4 ; df5 ; df6 ; df7 ; df8 ;
df9 ; df10 . From this analysis we can construct constraints for the free
parameters of the given theory and for the initial value problem rt2 . The
structure of the paper is as follows.

In Section II we present the gravitational model of our consideration which is
that of symmetric teleparallel $f\left(Q\right)$-gravity. For the background
space we assume isotropic and homogeneous spatially flat
Friedmann–Lemaître–Robertson–Walker (FLRW) spacetime. We present the three
different families of symmetric and flat connections and the corresponding
field equations. For two families of connections, an arbitrary scalar is
introduced in the field equations. As a result the corresponding field
equations admit a different dynamical evolution. We focus our study on these
two connections. The main results of this work are presented in Section III
where we perform a detailed analysis of the phase-space for the field
equations. We apply a set of dimensionless variables and we determine the
admitted asymptotic solutions by the cosmological field equations. Finally, in
Section IV we summarize our results and we draw our conclusions.

##  II $f\left(Q\right)$-gravity

In symmetric teleparallel $f\left(Q\right)-$gravity the gravitational Action
Integral is Zhao

$S=\frac{1}{2}\int d^{4}x\sqrt{-g}f(Q)+\int
d^{4}x\sqrt{-g}\mathcal{L}_{M}+\lambda_{\kappa}^{\;\lambda\mu\nu}R_{\;\lambda\mu\nu}^{\kappa}+\tau_{\lambda}^{\;\mu\nu}\mathrm{T}_{\;\mu\nu}^{\lambda},$
(1)

where $g=\mathrm{det}(g_{\mu\nu})$, $g_{\mu\nu}$ is the metric tensor,
quantities $R_{\;\lambda\mu\nu}^{\kappa}$, $\mathrm{T}_{\;\mu\nu}^{\lambda}$
are the Riemann and torsion tensors respectively,
$\lambda_{\kappa}^{\;\lambda\mu\nu}$, $\tau_{\lambda}^{\;\mu\nu}$ are two
Lagrange multipliers which impose that the connection
$\Gamma_{\mu\nu}^{\kappa}$ is flat and symmetric, that is,
$R_{\;\lambda\mu\nu}^{\kappa}=0$ and $\mathrm{T}_{\;\mu\nu}^{\lambda}=0$.

$f\left(Q\right)$ is an arbitrary function of the non-metricity scalar $Q$ for
the connection $\Gamma_{\mu\nu}^{\kappa}$ in which

$Q=Q_{\lambda\mu\nu}P^{\lambda\mu\nu},$ (2)

with

$Q_{\lambda\mu\nu}=\frac{\partial g_{\mu\nu}}{\partial
x^{\lambda}}-\Gamma_{\;\lambda\mu}^{\sigma}g_{\sigma\nu}-\Gamma_{\;\lambda\nu}^{\sigma}g_{\mu\sigma},$
(3)

and $P_{\;\mu\nu}^{\lambda}$ to be the non-metricity conjugate tensor defined
as

$P_{\;\mu\nu}^{\lambda}=-\frac{1}{4}Q_{\;\mu\nu}^{\lambda}+\frac{1}{2}Q_{(\mu\phantom{\lambda}\nu)}^{\phantom{(\mu}\lambda\phantom{\nu)}}+\frac{1}{4}\left(Q^{\lambda}-\bar{Q}^{\lambda}\right)g_{\mu\nu}-\frac{1}{4}\delta_{\;(\mu}^{\lambda}Q_{\nu)},$
(4)

where $Q_{\lambda}=Q_{\lambda~{}~{}~{}\mu}^{~{}~{}~{}\mu}$,
$\bar{Q}_{\lambda}=Q_{~{}~{}\lambda\mu}^{\mu}$.

In the case of vacuum, variation with respect to the metric tensor of the
gravitational Action Integral provides the field equations

$\frac{2}{\sqrt{-g}}\nabla_{\lambda}\left(\sqrt{-g}f^{\prime}(Q)P_{\;\mu\nu}^{\lambda}\right)-\frac{1}{2}f(Q)g_{\mu\nu}+f^{\prime}(Q)\left(P_{\mu\rho\sigma}Q_{\nu}^{\;\rho\sigma}-2Q_{\rho\sigma\mu}P_{\phantom{\rho\sigma}\nu}^{\rho\sigma}\right)=0,$
(5)

where now prime denotes a derivative with respect to the non-metricity scalar,
i.e. $f^{\prime}(Q)=\frac{df\left(Q\right)}{dQ}$. Equivalently, the field
equations (5) can be written with the use of the Einstein-tensor $G_{\mu\nu}$
as follow

$f^{\prime}(Q)G_{\mu\nu}+\frac{1}{2}g_{\mu\nu}\left(f^{\prime}(Q)Q-f(Q)\right)+2f^{\prime\prime}(Q)\left(\nabla_{\lambda}Q\right)P_{\;\mu\nu}^{\lambda}=0,$
(6)

where now $G_{\mu\nu}=\tilde{R}_{\mu\nu}-\frac{1}{2}g_{\mu\nu}\tilde{R},$ with
$\tilde{R}_{\mu\nu}$ and $\tilde{R}$ are the Riemannian Ricci tensor and
scalar respectively which are constructed by the Levi-Civita connection.
Indeed, when $f\left(Q\right)$ is a linear function, the field equations (6)
reduce to that of Genereal Relativity.

Over and above , from the variation of the Action Integral with respect to the
connection we derive the equation of motion for the non-metricity scalar

$\nabla_{\mu}\nabla_{\nu}\left(\sqrt{-g}f^{\prime}(Q)P_{\phantom{\mu\nu}\sigma}^{\mu\nu}\right)=0.$
(7)

###  II.1 FLRW background

In this study for the background geometry we consider the spatially flat FLRW
geometry with scale factor $a\left(t\right)$ and line element

$ds^{2}=-dt^{2}+a(t)^{2}\left[dr^{2}+r^{2}\left(d\theta^{2}+\sin^{2}\theta
d\phi^{2}\right)\right],$ (8)

where now the Hubble function is defined as $H=\frac{\dot{a}}{a}$, in which a
dot means total derivative with respect to the time-parameter $t$, i.e.
$\dot{a}=\frac{da}{dt}$.

The connections which are symmetric, flat and inherit the six isometries of
the FLRW as collineations have been determined before in Hohmann ; Heis2 .
Indeed, there exist three family of connections with common components

$\displaystyle\Gamma_{\theta\theta}^{r}$
$\displaystyle=-r~{},~{}\Gamma_{\phi\phi}^{r}=-r\sin^{2}\theta~{},$
$\displaystyle\Gamma_{r\theta}^{\theta}$ $\displaystyle=\Gamma_{\theta
r}^{\theta}=\Gamma_{r\phi}^{\phi}=\Gamma_{\phi r}^{\phi}=\frac{1}{r}$
$\displaystyle\Gamma_{\phi\phi}^{\theta}$
$\displaystyle=-\sin\theta\cos\theta~{},~{}\Gamma_{\theta\phi}^{\phi}=\Gamma_{\phi\theta}^{\phi}=\cot\theta.$

The first connection, namely $\Gamma_{1}$ has the additional non-zero
component

$\Gamma_{\;tt}^{t}=\gamma(t).$

Furthermore, the second-connection, $\Gamma_{2},$ has the additional non-zero
components

$\Gamma_{\;tt}^{t}=\frac{\dot{\gamma}(t)}{\gamma(t)}+\gamma(t),\quad\Gamma_{\;tr}^{r}=\Gamma_{\;rt}^{r}=\Gamma_{\;t\theta}^{\theta}=\Gamma_{\;\theta
t}^{\theta}=\Gamma_{\;t\phi}^{\phi}=\Gamma_{\;\phi t}^{\phi}=\gamma(t).$

Finally, the third connection $\Gamma_{3}$ has the extra non-zero components

$\Gamma_{\;tt}^{t}=-\frac{\dot{\gamma}(t)}{\gamma(t)},\quad\Gamma_{\;rr}^{t}=\gamma(t),\quad\Gamma_{\;\theta\theta}^{t}=\gamma(t)r^{2},\quad\Gamma_{\;\phi\phi}^{t}=\gamma(t)r^{2}\sin^{2}\theta.$

where function $\gamma$ is an arbitrary functions of the time-parameter $t$.

We remark that connections $\Gamma_{1}$, $\Gamma_{2}$ and $\Gamma_{3}$ are the
same only for zero function $\gamma\left(t\right)$. Moreover, for each
connection the field equations (6) are different.

Indeed, for the first connection $\Gamma_{1}$ we derive the non-metricity
scalar $Q=-6H^{2}$ and the field equations are

$\displaystyle
3H^{2}f^{\prime}(Q)+\frac{1}{2}\left(f(Q)-Qf^{\prime}(Q)\right)$
$\displaystyle=0~{},$ (9) $\displaystyle
2\frac{d}{dt}\left(f^{\prime}(Q)H\right)$ $\displaystyle=0,$ (10)

where equation (7) is satisfied as identity. From (10) it follows that
$Qf\,^{\prime}\left(Q\right)=const$, which means that $H$ is a constant
parameter, i.e. $H=H_{0}$. Thus, in the case of vacuum for connection
$\Gamma_{1}$, the only possible solution is the de Sitter universe, $H_{0}\neq
0$, or the Minkowski spacetime, $H_{0}=0$. For the power-law function
$f\left(Q\right)=Q^{\alpha}$, equation (9) gives that that the unique solution
is the Minkowski scacetime where $H=0$.

For the second connection $\Gamma_{2}$ we calculate the non-metricity scalar
$Q=-6H^{2}+9\gamma H+3\dot{\gamma}$, while equations (6) and (7) read

$\displaystyle
3H^{2}f^{\prime}(Q)+\frac{1}{2}\left(f(Q)-Qf^{\prime}(Q)\right)+\frac{3\gamma\dot{Q}f^{\prime\prime}(Q)}{2}$
$\displaystyle=0,$ (11) $\displaystyle
2\frac{d}{dt}\left(f^{\prime}(Q)H\right)-3\gamma\dot{Q}f^{\prime\prime}(Q)$
$\displaystyle=0,$ (12)
$\displaystyle\dot{Q}^{2}f^{\prime\prime\prime}(Q)+\left[\ddot{Q}+3H\dot{Q}\right]f^{\prime\prime}(Q)$
$\displaystyle=0.$ (13)

From the latter equations it is clear that function $\gamma\left(t\right)$
plays an important role on the dynamical evolution of the field equations.

On the other-hand, for the third connection $\Gamma_{3}$ the non-metricity
scalar is $Q=-6H^{2}+\frac{3\gamma}{a^{2}}H+\frac{3\dot{\gamma}}{a^{2}},$ and
the field equations are

$\displaystyle
3H^{2}f^{\prime}(Q)+\frac{1}{2}\left(f(Q)-Qf^{\prime}(Q)\right)-\frac{3\gamma\dot{Q}f^{\prime\prime}(Q)}{2a^{2}}$
$\displaystyle=0,$ (14) $\displaystyle
2\frac{d}{dt}\left(f^{\prime}(Q)H\right)+\frac{\gamma\dot{Q}f^{\prime\prime}(Q)}{a^{2}}$
$\displaystyle=0,$ (15)
$\displaystyle\dot{Q}^{2}f^{\prime\prime\prime}(Q)+\left[\ddot{Q}+\dot{Q}\left(H+\frac{2\dot{\gamma}}{\gamma}\right)\right]f^{\prime\prime}(Q)$
$\displaystyle=0.$ (16)

where again function $\gamma\left(t\right)$ effects the dynamics of the
physical parameters.

##  III Dynamical analysis

We proceed our analysis with the analysis of phase-space for the cosmological
field equations for the $f\left(Q\right)$-gravity. The dynamical analysis for
the first connection with a matter source performed in dyn1  where it was
found that the universe admits as asymptotic solutions the de Sitter universe,
and the matter dominated eras. A similar analysis appeared later in dyn2 .
Indeed, when there is not any external matter source for the first connection
the only possible asymptotic solution is that of the de Sitter universe. While
this work was in preparation the authors in dyn3  performed an analysis for
the phase-space of the field equations for the second connection $\Gamma_{2}$
with matter source. For this case new stationary points exist where the scalar
$\gamma\left(t\right)$ contributes in the cosmological fluid.

Indeed, the existence of self-similar solutions was the subject of study in
self  and it was found that for connections $\Gamma_{2}$ and $\Gamma_{3}$
self-similar solutions exist when a power-law component $Q^{\alpha}$ dominates
in the $f\left(Q\right)$ function. Thus for simplicity of our analysis we
consider the power-law function $f\left(Q\right)$ and we investigate the
dynamics of the field equations in order to reconstruct the evolution of the
physical parameters.

We consider the power-law $f\left(Q\right)=Q^{\alpha}$ function with
$\alpha\neq 0,1$ and we define new dimensionless variables in the context of
the $H$-normalization approach. We rewrite the field equations as a system of
first-order differential equations and we determine the stationary points.
Each stationary point describes an asymptotic solution for the background
space. Furthermore we investigate the stability properties of the stationary
points in order to reconstruct the cosmological history provided by the
theory.

###  III.1 Connection $\Gamma_{2}$

In order to study the dynamical evolution for the field equations of the
second connection for the power-law theory $f\left(Q\right)=Q^{\alpha},$ we
introduce the new variables

$x=\frac{\alpha-1}{6\alpha}\frac{Q}{H^{2}}~{},~{}y=\frac{\left(1-\alpha\right)}{2}\frac{\gamma}{H^{2}}\frac{\dot{Q}}{Q}~{},~{}z=\frac{\gamma}{H}~{},~{}\tau=\ln
a.$

Hence, the field equations (11), (12) and (13) are expressed as follow

$\displaystyle\frac{dx}{d\tau}$
$\displaystyle=\frac{2}{\alpha-1}\frac{xy}{z}\left(1-2\alpha+3\left(\alpha-1\right)z\right),$
(17) $\displaystyle\frac{dy}{d\tau}$
$\displaystyle=\frac{2}{\alpha-1}\frac{y}{z}\left(\alpha
x+\left(\alpha-1\right)\left(y-1\right)\left(3z-1\right)\right),$ (18)
$\displaystyle\frac{dz}{d\tau}$
$\displaystyle=\frac{2\alpha}{\alpha-1}x+\left(y-1\right)\left(3z-2\right),$
(19)

with algebraic constraint

$1-x-y=0.$ (20)

Moreover, the equation of state parameter
$w_{eff}=-1-\frac{2}{3}\frac{\dot{H}}{H^{2}}$ for the effective cosmological
fluid reads

$w_{eff}\left(x,y,z\right)=-1+y\left(2-\frac{4}{3z}\right).$ (21)

It is important to mention the variables $x,~{}y$ and $z$ are constraint and
they can take values in all the space of real numbers. Indeed, they can take
values at the infinity. Furthermore, with the use of the algebraic equation
(20) we can reduce the dimension by one the dimension of the dynamical system.
Thus, in the two-dimensional plane $\left(x,z\right)$ the dynamical system
reads

$\displaystyle\frac{dx}{d\tau}$
$\displaystyle=\frac{2}{\alpha-1}\frac{x\left(1-x\right)}{z}\left(1-2\alpha+3\left(\alpha-1\right)z\right),$
(22) $\displaystyle\frac{dz}{d\tau}$
$\displaystyle=x\left(\frac{2\alpha}{\alpha-1}-\left(3z-2\right)\right).$ (23)

We note that $z$ is constraint as $z>0$ or $z\,<0$.

The stationary points $A=\left(x\left(A\right),z\left(A\right)\right)$ of the
latter dynamical system at the finite regime are

$A_{1}=\left(1,\frac{2\left(2\alpha-1\right)}{3\left(\alpha-1\right)}\right)~{},~{}A_{2}=\left(0,z_{2}\right).$

For point $A_{1}$ we calculate $w_{eff}\left(A_{1}\right)=-1$, which means
that the asymptotic solution is that of the de Sitter universe. On the other-
hand $A_{2}$ describes a family of points with
$w_{eff}\left(A_{2}\right)=1-\frac{4}{3z_{2}}$, which means that the
asymptotic solutions at points $A_{2}$ describe accelerated universes for
$0<z_{2}<1$, while for $z_{2}=\frac{2}{3}$ the de Sitter universe is
recovered.

As far as the stability of the stationary points is concerned for the
linearized system (22), (23) near the stationary point $A_{1}$ we derive the
two eigenvalues $e_{1}\left(A_{1}\right)=-3$ and $e_{2}\left(A_{2}\right)=-3$
from where we infer that the de Sitter universe is always a future attractor.
On the other hand, for the family of points $A_{2}$ we determine the
eigenvalues $e_{1}\left(A_{2}\right)=0$ and
$e_{2}\left(A_{2}\right)=\frac{2\left(1-2\alpha-3z\left(1+\alpha\right)\right)}{\left(a-1\right)z}$.
When $e_{2}\left(A_{2}\right)<0,$ that is, for
$\left\\{\alpha<\frac{1}{2},0<z_{2}<\frac{2\alpha-1}{3\left(\alpha-1\right)}\right\\},~{}\left\\{\frac{1}{2}<\alpha<1\text{,
}\frac{2\alpha-1}{3\left(\alpha-1\right)}<z_{2}<0\right\\}$and
$\left\\{\alpha>1,0<z_{2}<\frac{2\alpha-1}{3\left(\alpha-1\right)}\right\\}$,
there may exist a stable submanifold. The latter can be calculated with the
application of the center manifold theorem. However, such a discussion is a
mathematical tool and does not contribute in the physical discussion of the
theory. Instead in Fig. 1 we present the phase-space portrait for the
dynamical system (22), (23) in the finite regime. For the diagrams it is clear
that stationary points $A_{2}$ admit always a stable submanifold, however in
general the family of points are saddle points.

Figure 1: Phase-space portrait for the dynamical system (22), (23) for
different values of parameter $\alpha$. It is clear that the family of
stationary points admit a stable submanifold. Moreover, point
$\left(x,z\right)=\left(1,0\right)$ is a transition point where parameter $z$
changes sign.

####  III.1.1 Analysis at Infinity

In order to investigate the existence of stationary points at the infinity we
define the Poincare variables

$x=\frac{X}{\sqrt{1-X^{2}-Z^{2}}}~{},~{}z=\frac{Z}{\sqrt{1-X^{2}-Z^{2}}}\text{,~{}}d\tau=\sqrt{1-X^{2}-Z^{2}}d\sigma\text{.}$

Therefore, the system (22), (23) becomes

$\frac{dX}{d\sigma}=F_{1}\left(X,Z;\alpha\right)~{};~{}\frac{dZ}{d\sigma}=F_{2}\left(X,Z;\alpha\right)\text{.
}$ (24)

Infinity is reached when $1-X^{2}-Z^{2}=0$. Hence the physical accepted
stationary points at infinity are

$A_{1}^{\infty}=\left(0,1\right)\text{ and }A_{2}^{\infty}=\left(0,-1\right).$

Points $A_{1}^{\infty}$ and $A_{2}^{\infty}$ describe scaling solutions with
$w_{eff}\left(A_{1,2}^{\infty}\right)=1$, that is, they describe universes
dominated by a stiff fluid. The two eigenvalues of the linearized system close
to the stationary points are determined to be zero. In Fig. 2 we present the
phase-space portrait for the dynamical system in the Poincare variables from
where it is clear that the two stationary points at infinity regime describe
unstable solutions.

Figure 2: Phase-space portrait for the dynamical system (22), (23) for
different values of parameter $\alpha$ at the Poincare variables. We observe
that the de Sitter point $A_{1}$ is the unique future attractor.

In Fig. 3 we present the qualitative evolution for the effective equation of
state parameter $w_{eff}\left(\sigma\right)$ for different values of the free
parameter $\alpha$ and for initial condition such that the matter dominated
era to be described by the effective fluid. We observe that the power-law
model $f\left(Q\right)=Q^{\alpha}$ for the second connection can be used as
unified model for the dark sector of the universe. On the other hand, in Fig.
4 we present the qualitative evolution of $w_{eff}\left(\sigma\right)$ for
different values of the free parameters with initial conditions which describe
an unstable accelerated universe. From the qualitative evolution of
$w_{eff}\left(\sigma\right)$ there exist solution trajectories which describe
universes with phases acceleration, matter era and future acceleration.
Consequently, this cosmological model can be used as a toy model to unify the
inflation with the late-time acceleration and the matter era.

Figure 3: Qualitative evolution for the effective equation of state parameter
$w_{eff}\left(\sigma\right)~{}$of connection $\Gamma_{2}$ for different values
of the free parameter $\alpha$ and for initial condition such that the matter
dominated era to be described by the effective fluid. Figure 4: Qualitative
evolution for the effective equation of state parameter
$w_{eff}\left(\sigma\right)$ of connection $\Gamma_{2}~{}$for different values
of the free parameter $\alpha$ and for initial condition such that to describe
exit from inflation.

#### III.1.2

###  III.2 Connection $\Gamma_{3}$

For the third connection and the field equations (14), (15) and (16) we apply
the same procedure. In this case we define the dimensionless variables

$\hat{x}=\frac{\alpha-1}{6\alpha}\frac{Q}{H^{2}}~{},~{}\hat{y}=\frac{\left(1-\alpha\right)}{2}\frac{a^{2}\gamma}{H^{2}}\frac{\dot{Q}}{Q}~{},~{}\hat{z}=\frac{a^{2}\gamma}{H}~{},~{}\tau=\ln
a.$

In the new variables the field equations are

$\displaystyle\frac{d\hat{x}}{d\tau}$
$\displaystyle=\frac{2}{\alpha-1}\frac{\hat{x}\hat{y}}{\hat{z}}\left(1-2\alpha+\left(\alpha-1\right)\hat{z}\right),$
(25) $\displaystyle\frac{d\hat{y}}{d\tau}$
$\displaystyle=-\frac{2}{\alpha-1}\frac{\hat{y}}{\hat{z}}\left(\alpha\hat{x}+\left(\alpha-1\right)\left(1+\hat{y}\right)\left(1+\hat{z}\right)\right),$
(26) $\displaystyle\frac{d\hat{z}}{d\tau}$
$\displaystyle=2\left(1+\frac{\alpha}{\alpha-1}\hat{x}\right)-3\hat{z}-\hat{y}\left(2+\hat{z}\right),$
(27)

with algebraic equation

$1-\hat{x}+\hat{y}=0,$ (28)

and $w_{eff}=-1-\frac{2}{3}\frac{\hat{y}}{\hat{z}}\left(2+\hat{z}\right),$
where $\hat{z}>0$ or $\hat{z}<0$.

With the use of the latter constraint equation we end with the two-dimensional
dynamical system

$\displaystyle\frac{d\hat{x}}{d\tau}$
$\displaystyle=\frac{2}{\alpha-1}\frac{\hat{x}\left(\hat{x}-1\right)}{\hat{z}}\left(1-2\alpha+\left(\alpha-1\right)\hat{z}\right),$
(29) $\displaystyle\frac{d\hat{z}}{d\tau}$
$\displaystyle=2\left(1+\frac{\alpha}{\alpha-1}\hat{x}\right)-3\hat{z}+\left(1-\hat{x}\right)\left(2+\hat{z}\right).$
(30)

The stationary points
$B=\left(\hat{x}\left(B\right),\hat{z}\left(B\right)\right)~{}$for the
dynamical system (29), (30) are

$B_{1}=\left(1,\frac{2\left(2\alpha-1\right)}{3\left(\alpha-1\right)}\right)~{},~{}B_{2}=\left(0,2\right)\text{
and
}B_{3}=\left(\frac{2\left(3-4\alpha\right)}{1+2\alpha},\frac{1-2\alpha}{\alpha-1}\right).$

For the stationary point $B_{1}$ we calculate $w_{eff}\left(B_{1}\right)=-1$
and the eigenvalues for the linearized system $e_{1}\left(B_{1}\right)=-3$ and
$e_{2}\left(B_{2}\right)=-5$. Thus, the asymptotic solution at point $B_{1}$
describes a de Sitter universe and point $B_{1}$ is a future attractor. For
point $B_{2}$ it follows $w_{eff}\left(B_{2}\right)=\frac{1}{3}$,
$e_{1}\left(B_{2}\right)=-2$ and
$e_{1}\left(B_{2}\right)=\frac{4\alpha-3}{\alpha-1}$; thus $B_{2}$ describes a
universe dominated by radiation and the point is an attractor for
$\frac{3}{4}<\alpha<1$, otherwise is a saddle point. Finally, for the point
$B_{3}$ we calculate $w_{eff}\left(B_{3}\right)=\frac{7-6\alpha}{3+6\alpha}$,
$e_{\pm}\left(B_{3}\right)=1-\frac{5}{1+2\alpha}\pm\frac{2\sqrt{\alpha\left(41\alpha-14\right)-11}}{1+2\alpha}$.
Therefore point $B_{3}$ describe an accelerated universe for $\alpha>2$, or
$\alpha<-\frac{1}{2}$ and the point is an attractor for
$-\frac{1}{2}<\alpha<\frac{3}{4}$.

In Fig. 5 we present the phase-space portrait for the two-dimensional system
(29), (30) for various values of the free parameter $\alpha$.

Figure 5: Phase-space portrait for the dynamical system (29), (30) for
different values of parameter $\alpha$. The de Sitter point $B_{1}$ is always
an attractor.

####  III.2.1 Analysis at Infinity

We define the new Poincare variables

$\hat{x}=\frac{\hat{X}}{\sqrt{1-\hat{X}^{2}-\hat{Z}^{2}}}~{},~{}\hat{z}=\frac{Z}{\sqrt{1-\hat{X}^{2}-\hat{Z}^{2}}}\text{,~{}}d\tau=\sqrt{1-\hat{X}^{2}-\hat{Z}^{2}}d\sigma\text{,}$

thus the field equations (29), (30) are written in the form

$\frac{d\hat{X}}{d\sigma}=G_{1}\left(\hat{X},\hat{Z};\alpha\right)~{};~{}\frac{d\hat{Z}}{d\sigma}=G_{2}\left(\hat{X},\hat{Z};\alpha\right)\text{.
}$ (31)

The stationary points of the latter dynamical system at the infinity, i.e. on
the surface $1-\hat{X}^{2}-\hat{Z}^{2}=0$, are

$B_{1}^{\infty}=\left(0,1\right)\text{ and }B_{2}^{\infty}=\left(0,-1\right).$

The two points describe scaling solutions with
$w_{eff}\left(B_{1}^{\infty}\right)=w_{eff}\left(B_{2}^{\infty}\right)=-\frac{1}{3}.$
The eigenvalues of the linearized system around the stationary points are
zero. From the phase-space portraits of Fig. 6 we observe that the stationary
points $B_{1}^{\infty}$, $B_{2}^{\infty}$ are sources.

Figure 6: Phase-space portrait for the dynamical system (22), (23) for
different values of parameter $\alpha$ at the Poincare variables. We observe
that the stationary points at the infinity regime are sources.

In Fig. 7 we plot the qualitative evolution of the equation of state parameter
$w_{eff}\left(\sigma\right)$ for initial conditions near infinity. We remark
that the $w_{eff}$ reaches the saddle point $B_{3}$ which describes a
radiation-like spacetime and then the de Sitter point $B_{1}$ is the future
attractor. Indeed, the evolution of the physical parameters for the third
connection can describe a radiation epoch as we exist from inflation and and
then to effective fluid has the behaviour of the cosmological constant.

Figure 7: Qualitative evolution for the effective equation of state parameter
$w_{eff}\left(\sigma\right)$ of connection $\Gamma_{3}~{}$for different values
of the free parameter $\alpha$ and for initial condition such that to describe
exit from inflation.

The results of this Section for the dynamical analysis of the cosmological
field equations for the two connections $\Gamma_{2}$ and $\Gamma_{3}$ are
summarized in Table 1.

Table 1: Stationary points and physical properties of the asymptotic
solutions. Connection $\mathbf{\Gamma}_{2}$
---
Finite regime | Point | $\left(\mathbf{x,z}\right)$ | $\mathbf{w}_{eff}$ | Acceleration? | Attractor?
| $A_{1}$ | $\left(1,\frac{2\left(2\alpha-1\right)}{3\left(\alpha-1\right)}\right)$ | $-1$ |  Yes | Yes
| $A_{2}$ | $\left(0,z_{2}\right)$ | $1-\frac{4}{3z_{2}}$ | $0<z_{2}<1$ | No
Infinity regime | Point | $\left(\mathbf{X,Z}\right)$ | $\mathbf{w}_{eff}$ | Acceleration? | Attractor?
| $A_{1}^{\infty}$ | $\left(0,1\right)$ | $1$ | No | No
| $A_{2}^{\infty}$ | $\left(0,-1\right)$ | $1$ | No | No
Connection $\mathbf{\Gamma}_{3}$
Finite regime | Point | $\left(\mathbf{\hat{x},\hat{z}}\right)$ | $\mathbf{w}_{eff}$ | Acceleration? | Attractor?
| $B_{1}$ | $\left(1,\frac{2\left(2\alpha-1\right)}{3\left(\alpha-1\right)}\right)$ | $-1$ | Yes | Yes
| $B_{2}$ | $\left(0,2\right)$ | $\frac{1}{3}$ | No | $\frac{3}{4}<\alpha<1$
| $B_{3}$ | $\left(\frac{2\left(3-4\alpha\right)}{1+2\alpha},\frac{1-2\alpha}{\alpha-1}\right)$ | $\frac{7-6\alpha}{3+6\alpha}$ |  $\alpha>2$ or $\alpha<-\frac{1}{2}$ | $-\frac{1}{2}<\alpha<\frac{3}{4}$
Infinity regime | Point | $\left(\mathbf{\hat{X},\hat{Z}}\right)$ | $\mathbf{w}_{eff}$ | Acceleration? | Attractor?
| $B_{1}^{\infty}$ | $\left(0,1\right)$ | $-\frac{1}{3}$ | No | No
| $B_{2}^{\infty}$ | $\left(0,-1\right)$ | $-\frac{1}{3}$ | No | No

##  IV Conclusions

In this study we investigated the asymptotic behaviour for the cosmological
parameters in symmetric teleparallel $f\left(Q\right)$-gravity for two
torsion-free and flat connections which provide non-trivial dynamical degrees
of the freedom in the field equations. We focused in the simplest case of
power-law $f\left(Q\right)=Q^{\alpha}$ theory and that there is not any matter
external matter source. While for the usual applied connection in the
literature such consideration leads to the Minkowski spacetime that is not
true for the two connections of the FLRW geometry that we consider in this
work.

We defined new variables in the $H$-normalization approach and we wrote the
field equations in the equivalent form of two-dimensional dynamical system of
first-order differential equations. We determined the stationary points at the
finity and infinity regimes. The analysis at the infinity performed with the
introduction of Poincare variables. For both connections we found that the de
Sitter universe is always a future attractor for the theory. However, new
scaling solutions were determined at the finite and infinity regimes.

For the connection $\Gamma_{2}$ at the finite regime there exist a family
unstable scaling solutions which can describe acceleration, while at infinity
there exist stiff fluid solutions. On the other hand, for connection
$\Gamma_{2}$ at the finite regime there exists a stationary point which
describes an unstable radiation solution and a scaling solution which can
describe acceleration or matter era. At infinity the stationary points
describe solutions dominated by curvature-like fluid source.

In order to understand the evolution of the general solution we plotted phase-
space portraits to constrain the initial conditions from the evolution of the
trajectories in the two-dimensional space of the dynamical system. We solved
numerically the trajectories and from the qualitative evolution of the
effective equation of state parameter we found that the theory can explain the
evolution of the universe as we exit from the inflationary epoch and for
various sets of initial conditions it can be used as a toy model for the
description of the matter era and of the present acceleration.

In future work we plan to investigate further the physical properties of these
solutions of $f\left(Q\right)$-gravity for these two connections.

###### Acknowledgements.

This work was financially supported in part by the National Research
Foundation of South Africa (Grant Numbers 131604). The author thanks the
support of Vicerrectoría de Investigación y Desarrollo Tecnológico (Vridt) at
Universidad Católica del Norte through Núcleo de Investigación Geometría
Diferencial y Aplicaciones, Resolución Vridt No - 098/2022.

## References

  * (1) M. Hohmann, Phys. Rev. D 104, 124077 (2021)
  * (2) J. Beltrán-Jiménez, L. Heisenberg and T.S. Koivisto, Universe 5, 173 (2019)
  * (3) S. Nojiri, S.D. Odintsov and V.K. Oikonomou, Phys. Rept. 692, 1 (2017)
  * (4) A.G. Riess et al., Astron. J. 116, 1009 (1998)
  * (5) M. Tegmark et al., Astrophys. J. 606, 702 (2004)
  * (6) E. Komatsu et al., Astrophys. J. Suppl. Ser. 180, 330 (2009)
  * (7) S. Capozziello, C.A. Mantica and L.G. Molinari, EPL 137, 19001 (2022)
  * (8) H.A. Buchdahl, Mon. Not. Roy. Astron. Soc. 150, 1 (1970)
  * (9) R. Ferraro and F. Fiorini, Phys. Rev. D 75, 084031 (2007)
  * (10) J. B. Jimenez, L. Heisenberg and T. Koivisto, Phys. Rev. D 98, 044048 (2018)
  * (11) S. M. Carrol, V. Duvvuri, M. Trodden and M. S. Turner, Phys. Rev. D., 70, 043528 (2004);
  * (12) L. Amendola, D. Polarski and S. Tsujikawa, Phys. Rev. Lett. 98, 131302 (2007)
  * (13) A. A. Starobinsky, Phys. Lett. B 91, 99 (1980)
  * (14) T. Clifton and J.D. Barrow, Class. Quant. Grav. 23, 2951 (2006)
  * (15) R.C. Nunes, A. Bonilla, S. Pan and E.N. Saridakis, EPJC 77, 230 (2017)
  * (16) S. Nesseris, S. Basilakos, E.N. Saridakis and L. Perivolaropoulos, Phys. Rev. D 88, 103010 (2013)
  * (17) A. de la Cruz-Dombriz, A. Dobado and A. L. Maroto, Phys. Rev. D 80, 124011 (2009)
  * (18) E. Santos, Astroph. Sp. Sci. 341, 411 (2012)
  * (19) M. E. Rodrigues, M. J. S. Houndjo, J. Tossa, D. Momeni and R. Myrzakulov, JCAP 11, 024 (2013)
  * (20) J. Aftergood and A. DeBenedictis, Phys. Rev. D 90, 124006 (2014)
  * (21) T. P. Sotiriou and V. Faraoni, Rev. Mod. Phys. 82 451, (2010)
  * (22) S. Bahamonde, K. F. Dialektopoulos, C. Escamilla-Rivera, G. Farrugia, V. Gakis, M. Hohmann, J. L. Said, J. Mifsud and E. Di Valentino, Rep. Prog. Phys. 86, 026901 (2023)
  * (23) L. Atayde and N. Frusciante, Phys. Rev. D 104, 064052 (2021)
  * (24) R. Solanki, A. De and P.K. Sahoo, Phys. Dark Universe 36, 100996 (2022)
  * (25) A. Lymperis, JCAP 11, 018 (2022)
  * (26) S.H. Shekh, Phys. Dark Univ. 33, 100850 (2021)
  * (27) S.A. Narawade and B. Mishra, Phantom Cosmological Model with Observational Constraints in f(Q) Gravity, to appear in Annalen de Physik, (2023) DOI: 10.1002/andp.202200626
  * (28) F.K. Anagnastopoulos, S. Basilakos and E.N. Saridakis, Phys. Lett. B 822, 136634 (2021)
  * (29) S. Capozziello and M. Shokri, Phys. Dark Univ. 37, 101113 (2022)
  * (30) A. De, S. Mandal, J.T. Beh, T.-H. Loo and P.K. Sahoo, EPJC 82, 72 (2022)
  * (31) W. Wang, H. Chen and T. Katsuragawa, Phys. Rev. D 105, 024060 (2022)
  * (32) F. D’ Ambrosio, S. D. B. Fell, L. Heisenberg and S. Kuhn, Phys. Rev. D 105, 024042 (2022)
  * (33) R.-H. Lin and X.-H. Zhai, Phys. Rev. D 103, 124001 (2021)
  * (34) S. Mandal, G. Mustafa, Z. Hassan and P.K.Sahoo, Phys. Dark Universe 35, 100934 (2022)
  * (35) K. Hu, T. Katsuragawa and T. Qiu, ADM formulation and Hamiltonian analysis of f(Q) gravity, [arXiv:2204.12826] (2022)
  * (36) N. Dimakis, A. Paliathanasis and T. Christodoulakis, Class. Quantum Grav. 38, 225003 (2021)
  * (37) A. De, S. Mandal, J.T. Beh, T.-H. Loo and P.K. Sahoo, Eur. Phys. J. C 82, 72 (2022)
  * (38) F. Esposito, S. Carloni and S. Vignolo, Class. Quantum Grav. 39, 235014 (2022)
  * (39) N. Dimakis, A. Paliathanasis and T. Christodoulakis, Class. Quantum Grav. 38, 225003 (2021)
  * (40) M. Koussour, S. Dahmani, M. Bennai and T. Quali, Eur. Phys. J. Plus 138, 179 (2023)
  * (41) W. Khyllep, J. Dutta, E.N. Saridakis and K. Yesmakhanova, Phys. Rev. D 107, 044022 (2023)
  * (42) M. Hohmann, Phys. Rev. D 104, 124077 (2021)
  * (43) K. Hu, T. Katsuragawa and T. Qiu, Phys. Rev. D 106, 044025 (2022)
  * (44) G. Subramaniam, A. De, T.-H. Loo and Y.K. Goh, How different connections in flat FLRW geometry impact energy conditions in f(Q) theory?, (2023) [arXiv:2304.02300]
  * (45) A.A. Coley, Phys. Rev. D 62, 023517 (2000)
  * (46) R. Lazkoz and G. Leon, Phys. Lett. B 638, 303 (2006)
  * (47) J. De-Santiago, J.L. Cervantes-Cota and D. Wands, Phys. Rev. D 87, 023502 (2013)
  * (48) A. Cid, G. Leon and Y. Leyva, JCAP 02, 027 (2016)
  * (49) S. Chatzidakis, A. Giacomini, P.G.L. Leach, G. Leon, A. Paliathanasis and S. Pan, JHEAp 36, 141 (2022)
  * (50) H.R. Kausar, I. Noureen and M.U. Shahzad, Eur. Phys. J. Plus 130, 204 (2015)
  * (51) C.R Fadragas and G. Leon, Class. Quantum Grav. 31, 195011 (2014)
  * (52) T. Nakamura, T. Ikeda, R. Saito, N. Tanahashi and C.-M. Yoo, Phys. Rev. D 103, 024009 (2021)
  * (53) S.Kr. Biswas, S. Chakraborty, J. Dutta and S. Chakraborty, Phys. Rev. D 95, 103009 (2017)
  * (54) M. Kerachian, G. Acquaviva and G. Lukes-Gerakopoulos, Phys. Rev. D 101, 043535 (2020)
  * (55) D. Zhao, Covariant formulation of f(Q) theory, Eur. Phys. J. C 82, 303 (2022)
  * (56) M. Hohmann, General covariant symmetric teleparallel cosmology, Phys. Rev. D 104 124077 (2021)
  * (57) F. D’ Ambrosio, L. Heisenberg and S. Kuhn, Revisiting cosmologies in teleparallelism, Class. Quantum Grav. 39 025013 (2022)
  * (58) W. Khyllep, A. Paliathanasis and J. Dutta, Phys. Rev. D 103, 103521 (2021)
  * (59) C. Bohmer, E. Jensko and R. Lazkoz, Universe 9, 166 (2023)
  * (60) H. Shabani, A. De and T.-H. Loo, Phase-space analysis of a novel cosmological model in f(Q) theory, (2023) [arXiv:2304.02949]
  * (61) N. Dimakis, M. Roumeliotis, A. Paliathanasis, P.S. Apostolopoulos and T. Christodoulakis, Phys. Rev. D 106, 123516 (2022)


# Prediction Error Estimation in Random Forests

Ian Krupkin
Department of Mathematics & Statistics, Pomona College
and
Johanna Hardin
Department of Mathematics & Statistics, Pomona College

###### Abstract

In this paper, error estimates of classification Random Forests are
quantitatively assessed. Based on the initial theoretical framework built by
Bates et al., (2023), the true error rate and expected error rate are
theoretically and empirically investigated in the context of a variety of
error estimation methods common to Random Forests. We show that in the
classification case, Random Forests’ estimates of prediction error is closer
on average to the true error rate instead of the average prediction error.
This is opposite the findings of Bates et al., (2023) which were given for
logistic regression. We further show that this result holds across different
error estimation strategies such as cross-validation, bagging, and data
splitting.

Keywords: true error rate, expected error rate, cross validation,
bootstrapping, bagging, data splitting

##  1 Introduction

As evidenced by the 2017 State of Data Science and Machine Learning report by
Kaggle, almost half of data scientists use Random Forests at work (Kaggle,,
2017). Random Forests (Breiman,, 2001) have become a popular classification
tool in a variety of fields, especially because of their excellent performance
in very complex data settings. The fact that out-of-bag (OOB) errors are
theoretically and computationally simple improvements over a train-test split,
lead to their ubiquity. When deploying a predictive model, it is important to
understand its prediction accuracy on future test points, so both good point
estimates and an understanding of the variability of the estimates are
essential. When Random Forests are implemented, the OOB error is a widely-used
approach for point and interval estimate tasks, but in spite of OOB’s seeming
simplicity, its properties remain opaque. In the past, the OOB error has been
affirmed to be an unbiased estimate of the true error rate (Zhang et al.,,
2010; Goldstein et al.,, 2011). Nonetheless, it has been shown that for two-
class classification problems the OOB error can overestimate the true
prediction error (Bylander,, 2002; Mitchell,, 2011). It was later argued that
the use of stratified subsampling with sampling fractions that are
proportional to response class sizes of the training data yielded almost
unbiased error rates (Janitza and Hornung,, 2018). The present work is
primarily concerned with OOB errors, but also addresses other common methods
such as data splitting and cross validation, as well as their combination with
OOB errors.

Despite the apparent straightforwardness of data splitting, cross validation,
and bagging, the formal properties of these modeling techniques are subtle.
When calculating an error rate, the question of “what are we estimating?”
rightfully reappears often. For the linear model, fit by ordinary least
squares, it was proven that the most popular estimates of prediction error
(including data splitting, bagging, cross-validation, etc.) are closer in
expected value to the average prediction error of models fit on other unseen
training sets drawn from the same population than the prediction error of the
model at hand (Bates et al.,, 2023). In other words, although the error of the
model fit on the training data may seem like a reasonable estimand, it is not
the closest error target. To our knowledge, in the case of Random Forests,
ours is the first study investigating the prediction error as a metric for the
accuracy of the model on training, cross validated, test, and future data.

In an innovative paper, Bates et al., (2023) provide a framework for
understanding and measuring prediction error in the case of least squares
linear models. They prove that cross validation does not estimate the
prediction error for the model at hand, fit to the training data. Instead, the
prediction error “estimates the average prediction error of models fit on
other unseen training sets drawn from the same population.” We illustrate and
extend their work to the setting of Random Forests where empirical evidence
for Random Forests shows the opposite of what Bates et al., (2023) proved for
linear models.

The main contribution of our work is two-fold: (i) the switch in direction of
the proximity of Random Forests’ error estimates to the error targets
presented by Bates et al., (2023) is investigated in comparison to the
proximity of logistic regression to its error targets, through studies with
different numbers of observations and predictors and (ii) the performance of a
variety of error estimation strategies is explored.

The paper is structured as follows: In Section 2, we set up notation and
introduce the two different error targets: true error rate ($\text{Err}_{XY}$)
and expected true error (Err). Subsequently, Section 3, introduces simulation-
based studies. The descriptions include an outline of the simulated data, the
considered settings, and several model building workflows that will be
investigated. In Section 4, we present the results of the studies. The results
are discussed in Section 5 alongside recommendations.

All code used for the simulations is available at
https://github.com/iankrupkin/Prediction-Error-Estimation-in-Random-Forests.

##  2 What Prediction Error are we Estimating?

Before turning to our main method in the next section, we introduce our
notation and review topics related to error targets. We consider the
supervised learning setting where we have $p$ features and $n$ observations,
denoted by ${\bf X}_{i}=(X_{i,1},\dots,X_{i,p})$, the $i^{th}$ row of ${\bf
X}\in\mathbb{R}^{n}\times\mathbb{R}^{p}$, and response ${\bf
Y}=(Y_{1},\dots,Y_{n})\in\mathbb{R}^{n}$. We assume that the data points
$({\bf X}_{i},Y_{i})$ for $i=1,\dots,n$ are independent and identically
distributed from some underlying distribution $P$ on $\mathbb{R}^{p+1}$. We
wish to understand the performance of our fitted model when generalized to
unseen data points, which can be formalized by a loss function on an observed
response $y$ compared to a predicted response $\hat{y}$:

$\ell(\hat{y},y)\geq 0$

such that $\ell(y,y)=0$ for all $y$. The form of $\ell$ need not be specified
and could be squared error loss, misclassification error, cross-entropy, etc.
Now consider a model $f(\cdot)$ parameterized by $\theta$. Let $\hat{f}({\bf
x},\hat{\theta})$ be the function that predicts $y$ from ${\bf
x}\in\mathbb{R}^{p}$ using the model with parameters $\theta$, which take
values in the space $\Theta$. Let $\mathcal{A}$ be a model-fitting algorithm
that takes any number of data points and returns a parameter vector
$\hat{\theta}\in\Theta$. Hence, $\hat{\theta}=\mathcal{A}({\bf X},Y)$ is the
fitted value of the parameter based on the observed data ${\bf X}$ and $Y$.
Let $({\bf X}_{n+1},Y_{n+1})\sim P$ be another independent test point from the
same distribution. Using the training data, we are interested in finding the
function $\hat{f}({\bf x},\theta)$ that minimizes the loss $\ell(\hat{f}({\bf
X}_{n+1},\theta),Y_{n+1})$. Note that $\ell(\hat{f}({\bf
X}_{n+1},\theta),Y_{n+1})$ is a random and unknown object, and our target is
one of two quantities:

$\displaystyle\text{True Error Rate:\;}\text{Err}_{XY}$
$\displaystyle=E[\ell(\hat{f}({\bf X}_{n+1},\hat{\theta}),Y_{n+1})|({\bf
X},Y)]$ (1) Expected Error Rate: Err $\displaystyle=E[\ell(\hat{f}({\bf
X}_{n+1},\hat{\theta}),Y_{n+1})]=E[\text{Err}_{XY}]$ (2)

These are the two most natural quantities of interest to the analyst. Known as
the _true error rate_ , $\text{Err}_{XY}$ is the expected test error of the
model that was fit on our actual training set. Err is the expected average
error of the fitting algorithm run on the same-sized data sets drawn from the
underlying distribution $P$, and called the _expected error rate_. It is
important to note that the random variable Err is a constant with respect to
$({\bf X},Y)$, while $\text{Err}_{XY}$ is a function of $({\bf X},Y)$ (Bates
et al.,, 2023; Rajanala et al.,, 2022).

Depending on the situation, one may prefer to estimate $\text{Err}_{XY}$, or
Err may be the quantity of interest. The former quantity is of the most
interest to a practitioner deploying a specific model, whereas the latter may
be of interest to a researcher comparing different fitting algorithms. To
illustrate this difference consider the following examples.

Suppose Statistician A is trying to estimate the average height of penguins,
found in the wild, based on a sample of 100 researched penguins. Statistician
A will use the original sample to build a model to estimate the average height
of the next sample of penguins. When presenting the model, they will be
interested in the true error rate ($\text{Err}_{XY}$) of their specific model
because they will want to know how the specific model they built will perform
on the next data set.

However, if Statistician B is trying to accomplish the same task of estimating
the average height of penguins, based on a finite sample, but is unsure of the
structure of the model to utilize, they will be interested in a slightly
different error metric. Statistician B will run a variety of fitting
algorithms to build numerous models based on the sample available and will
need to compare the models. They will want to estimate the expected error rate
(Err) to know the average error of the fitting algorithm run on same-sized
data sets drawn from the underlying distribution. Statistician B is less
interested in the performance of the model which was built using the sample at
hand (which is Statistician A’s target), but rather the performance of the
process of arriving at the model.

The estimation of the theoretical quantities, $\text{Err}_{XY}$ and Err,
deepens the understanding of the difference between the two. As mentioned
above, the _true error rate_ , $\text{Err}_{XY}$, is the test error of the
model that was fit on our actual training set. Hence, the estimation of this
quantity is the error produced by the model on a new theoretically infinitely
large test set. As illustrated in Figure 1, the sample is used to create the
model, and then the model is used to predict on a large test set from the
underlying population. The missclassification rate on the test set will be the
_true error rate_.

Subsequently, Err is the average error of the fitting algorithm run on the
same-sized data sets drawn from the underlying distribution $P$, and called
the _expected error rate_. When calculating the _expected error rate_ , the
average of the _true error rate_ , but using a new model every time, is taken.
As seen in Figure 1, the entire model fitting process is repeated to obtain
each $\text{Err}_{XY}$ from one sample.

In other words, the difference between the estimation of $\text{Err}_{XY}$ and
Err is that the former uses one model, while the latter averages over many
models. As seen from Equation 1, $\text{Err}_{XY}$ is conditional on the data,
while Err is unconditional. Note that Err averages over everything that is
random, including the randomness in the training set that produced the model.

Figure 1: Empirical estimation of $\text{Err}_{XY}$ and Err.

While it may initially appear that the quantity $\text{Err}_{XY}$ is easier to
estimate, since it concerns the model at hand, it has been observed that, in
some settings, the cross-validation estimate of error is weakly correlated
with $\text{Err}_{XY}$ (given a particular population at hand) (Yousef,,
2019). The disassociation issue is mainly attributed to data re-usage.

Let $({\bf X}_{\infty},Y_{\infty})$ represent a data set of unlimited size
enabling the best possible model $f(\cdot)$ to be chosen. Theoretically,
$\text{Err}_{XY}$ can be decomposed into four parts, as see in Faraway,
(2014). Let $({\bf X}_{\infty},Y_{\infty})$ represent a black box with the
ability to generate data sets of unlimited size and quantity but no sense of
the model generating mechanism. The true model $f$ is given in contrast to the
model found from the original data set, $\hat{f}_{n}$. The parameter is either
given by the infinite data, $\theta_{\infty}$, or estimated using the original
data, $\hat{\theta}_{n}$, or a new (independent) data set,
$\hat{\theta}_{n+1}.$

$\displaystyle\text{Err}_{XY}$ $\displaystyle=E[\ell(f({\bf
X}_{\infty},\theta_{\infty}),Y_{\infty})|({\bf X},Y)]$ best performance (3)
$\displaystyle+E[\ell(\hat{f}_{n}({\bf
X}_{\infty},\theta_{\infty}),Y_{\infty})-\ell(f({\bf
X}_{\infty},\theta_{\infty}),Y_{\infty})|({\bf X},Y)]$ model selection cost
(4) $\displaystyle+E[\ell(\hat{f}_{n}({\bf
X}_{n+1},\hat{\theta}_{n+1}),Y_{n+1})-\ell(\hat{f}_{n}({\bf
X}_{\infty},\theta_{\infty}),Y_{\infty})|({\bf X},Y)]$ parameter estimation
cost (5) $\displaystyle+E[\ell(\hat{f}_{n}({\bf
X},\hat{\theta}_{n}),Y)-\ell(\hat{f}_{n}({\bf
X}_{n+1},\hat{\theta}_{n+1}),Y_{n+1})|({\bf X},Y)]$ data re-use cost (6)

Term (3) represents the best performance of the prediction error; it is the
expected loss on the correct model using all the possible data. Term (4)
represents the difference between the loss for the true model on infinite data
and the loss for the selected model on infinite data. Term (5) represents the
differences in loss using the data model but estimating the parameter with a
single (independent) data set versus infinite data. Term (6) represents the
difference in loss between the data model and the parameter estimates using
the original data versus a single (independent) data set.

The most interesting component is the final term (6), caused by data re-usage,
which has a non-zero expectation when the same data points are used for both
model selection and parameter estimation. If one uses a validation set
approach, term (6) will have an expectation of zero because each observation
is only used once. However, as estimated empirically, when using a full data
approach, term (6) can be large and easily cancel out any advantages the full
data has in model selection and parameter estimation (Faraway,, 2014). Thus,
the full data strategy will have lower model selection and parameter
estimation costs than the validation set strategy due to the higher number of
observations used to complete the model selection and parameter estimation
processes, but they can be outweighed by the data re-use cost.

The difference between the full data and validation set strategies, seen in
comparing terms (4) & (5), is bounded and well understood as an effect of
sample size (Faraway,, 2014). Despite suffering in model selection and
parameter estimation costs, the validation set strategy will have a lower data
re-use cost than the full data strategy, and we know the data re-use cost term
(6) could be very large. Therefore, we would like to investigate the
situations when the data re-use cost outweighs the model selection and
parameter estimation costs, in the Random Forest context. This comes in the
form of analyzing various model fitting approaches and resulting estimates of
error. Specifically the use of OOB errors compared to validation set and cross
validation strategies is investigated.

##  3 Methods

In this section, some simulation-based studies are described. Simulated data
are used to study the behavior of modeling strategies in simple settings, in
which all predictor variables are uncorrelated. The results provide insight
into the mechanisms which lead to different targets in error estimates.

###  3.1 Data Generation and Settings

The bias of error estimates in different data settings with numeric predictor
variables is systematically investigated by means of simulation studies in
balanced binary two-class response variable data. The settings considered are:

  * •

Different number of predictors, $p\in\\{10,100\\}$.

  * •

Different number of observations such that $n<p,n>p,n\gg p$.

As done when modeling real data, several Random Forests with different $mtry$
values are constructed for each setting ($mtry$ is the parameter that
determines the number of randomly chosen variables to be considered for each
split on a tree). The values for $mtry$ range from $mtry=1$ all the way up to
$mtry=p$. Note that for $mtry=1$ there is no selection of an optimal predictor
variable for a split, while for $mtry=p$ the Random Forest method coincides
with the bagging procedure which selects the best predictor variable from the
entire set of predictors.

The number of trees chosen is a trade-off between accuracy and computational
speed. More trees are necessary when using a large number of predictor
variables. The OOB error stabilizes at around 250 trees in convergence studies
(Goldstein et al.,, 2010), and they concluded that 1000 trees might be
sufficiently large for their genome-wide data set of more than 300,000
predictor variables. Also in high-dimensional settings, Random Forests with
500 trees and 1000 trees yield very similar OOB errors (Genuer et al.,, 2008).
In accordance with these findings the number of trees is set to 500 in all
studies of this paper. Each setting was repeated 1000 times to obtain stable
results.

Only numeric predictor variables are considered in the studies. Both
predictors associated with the response and predictors not associated with the
response were considered, with all predictors still independently distributed.
The predictors not associated with the response follow a standard normal
distribution. The distribution of predictors with association is different for
each response class. The predictor values for observations from class 1 are
always drawn from a standard normal distribution. The predictor values for
observations from class 2 are drawn from a normal distribution with variance 1
and a mean different from zero. Figure 2 gives an overview of the distribution
of predictors in the response classes. Let us consider the setting with $p=10$
as an example. The first two predictors $X_{1}$ and $X_{2}$ are associated
with the response, while the other predictors $X_{3},\dots,X_{10}$ are noise.
Hence, $X_{3},\dots,X_{10}$ follow a standard normal distribution, while the
distributions of $X_{1}$ and $X_{2}$ depend on the class to which the
observations belong. If the observation comes from class 1, the distribution
of $X_{1}$ and $X_{2}$ is $N(0,1)$, and $X_{1}$ and $X_{2}$ are distributed
$N(0.75,1)$ for class 2. Randomly drawing the mean separately for $X_{1}$ and
$X_{2}$ and for each repetition of the study insures that predictors with
different effect strengths are considered.

$\begin{array}[]{|l|l|l|l|}\hline\cr\text{Number
of}&\text{Predictors}&\text{Class 1}&\text{Class 2}\\\
\text{Predictors}&&N(\mu_{1},1)&N(\mu_{2},1)\\\ \hline\cr
p=10&X_{1}&\mu_{1}=0&\mu_{2}\sim N(0.75,1)\\\
\cline{2-4}\cr&X_{2}&\mu_{1}=0&\mu_{2}\sim N(0.75,1)\\\
\cline{2-4}\cr&X_{3},\dots,X_{10}&\mu_{1}=0&\mu_{2}=0\\\ \cline{1-4}\cr
p=100&X_{1}&\mu_{1}=0&\mu_{2}\sim N(0.75,1)\\\
\cline{2-4}\cr&\vdots&\vdots&\vdots\\\
\cline{2-4}\cr&X_{10}&\mu_{1}=0&\mu_{2}\sim N(0.75,1)\\\
\cline{2-4}\cr&X_{11},\dots,X_{100}&\mu_{1}=0&\mu_{2}=0\\\
\hline\cr\end{array}$ Figure 2: Distribution of predictors in class 1 and
class 2 of the simulated data set up as in Janitza and Hornung, (2018).

It is important to note that the settings are simplistic because all
predictors are uncorrelated. Although assuming no correlations between any of
the predictors is not realistic, such settings are important to investigate in
order to understand the mechanisms which lead to different targets in error
estimation.

###  3.2 Strategies for Error Estimation

The modeling process consists of parameter estimation followed by error
estimation. An important point of consideration when completing the two
estimation steps is the choice of which subset of observations will be used in
each operation. Often, data for parameter estimation and data for error
estimation are collected at the same time, thus resulting in a single sample
that needs to be apportioned to both parameter and error estimation. Finding
the optimal model complexity requires an external test set (Hastie et al.,,
2001). In an ideal world, to avoid “data snooping”, one needs one data set for
model building, one for parameter estimation, and then after a model is
accepted, another data set for error estimation. However, rarely are three
independent data sets available, so one may need to do the best one can with
the data available. Hence, when modeling it is important to outline the
strategy that will be used to construct the model and then estimate its error.

Various strategies were chosen to separately target the parameter estimation
and error estimation steps in the modeling process, and thus each strategy
consists of two parts. In the descriptions that follow, in-fold represents the
observations used to fit the model using cross validation; out-of-fold
represents the observations that were held out of the model fitting using
cross validation. In-bag are the observations that were used to fit the model
using Random Forests; out-of-bag (OOB) represents the observations that were
held out of the model fitting when using Random Forests. A more complete
description of cross validation and Random Forests can be found at James et
al., (2013).

For each of the considered strategies, there were three aspects to consider
(see Figure 3).

  1. 1.

Model Fitting: Data1 is used to fit a series of models to consider. For
example, with Random Forests, there might be 5 models fit
($\hat{f}_{1},\hat{f}_{2},\hat{f}_{3},\hat{f}_{4},\hat{f}_{5}$), each with a
different value of $mtry.$

  2. 2.

Parameter Estimation: Data2 is used to select the value of the parameter. In
the above case, the value of $mtry$ is selected by assessing Data2 on each of
$\hat{f}_{1},\hat{f}_{2},\hat{f}_{3},\hat{f}_{4},\hat{f}_{5}$.

  3. 3.

Error Estimation: Once the model and parameter have been selected, Data3 is
used to estimate the prediction error of the selected model.

Figure 3: For each strategy, the data is divided differently into three data
sets: Data1, Data2, and Data3. See Figure 4 for details on how each strategy
divides up the original data set.

The following strategies were considered:

  * •

Logistic Regression CV Error (LGCV): The logistic model is built on the in-
fold data set and the error of the model is estimated via the out-of-fold data
using 4-fold cross-validation. Logistic regression models are run only for
$n\geq p.$

  * •

Full Data Set CV Error (FDCV): Parameters are set prior to model building with
$mtry=\sqrt{p}$. The Random Forest is built on the in-fold data set and the
error of the model is estimated via out-of-fold data using 4-fold cross-
validation.

  * •

Full Data Set OOB Error (FDO): Parameter and error estimation is done on the
same data set. $mtry$ is chosen by using the OOB error rate. Hence, the Random
Forest with the lowest OOB error rate is chosen and the OOB error is returned
as the error estimate.

  * •

Split Data Set OOB Error (SDO): The sample is divided into training and
testing sets. Parameter estimation is done on the training set, using the OOB
error rate to select $mtry$. The error of the Random Forest, built on the
entire training set, is then estimated by predicting on the testing set with
an accuracy measure as the error estimate.

  * •

Split Data Set CV Error (SDCV): The sample is divided into training and
testing sets. Parameter estimation is done on the training set, using 4-fold
cross-validated error estimates to select $mtry$. The error of the Random
Forest, built on the entire training set, is then estimated by predicting on
the testing set with an accuracy measure as the error estimate.

  * •

Split Data Set Test Error (SDT): The sample is divided into three independent
training, validation, and testing sets. Parameter estimation of $mtry$ is done
on the training and validation sets. The Random Forest is fit on the training
set, and the fitted forest is used to predict the responses for the
observations in the validation set. The model with the highest accuracy is
chosen. The error of this Random Forest, built solely on the training set, is
estimated by predicting on the testing set and obtaining an accuracy measure,
and thus uses observations that are not part of the set of observations that
are considered for constructing the Random Forest.

$\begin{array}[]{|c|c|c|c|}\hline\cr\text{Strategy}&\text{Model
Fitting}&\text{Parameter Estimation}&\text{Error Estimation}\\\
&(\text{Data}_{1})&(\text{Data}_{2})&(\text{Data}_{3})\\\
\hline\cr\text{LGCV}&\text{In-Fold observations}&\text{N/A}&\text{Out-of-Fold
observations}\\\ &\text{from full data set}&&\text{from full data set}\\\
\hline\cr\text{FDCV}&\text{In-Fold observations}&\text{N/A}&\text{Out-of-Fold
observations}\\\ &\text{from full data set}&&\text{from full data set}\\\
\hline\cr\text{FDO}&\text{In-Bag observations}&\text{OOB
observations}&\text{OOB observations}\\\ &\text{from full data set}&\text{from
full data set}&\text{from full data set}\\\ \hline\cr\text{SDO}&\text{In-Bag
observations}&\text{OOB observations}&\text{Test data set}\\\ &\text{from
training data set}&\text{from training data set}&\\\
\hline\cr\text{SDCV}&\text{In-Fold observations}&\text{Out-Fold
observations}&\text{Test data set}\\\ &\text{from training data
set}&\text{from training data set}&\\\ \hline\cr\text{SDT}&\text{Training data
set}&\text{Validation data set}&\text{Test data set}\\\ \hline\cr\end{array}$

Figure 4: Data used in each step of the error estimation strategies

##  4 Results

###  4.1 Distance to Target Errors

Recall, when creating statistical models, one wants estimates of error that
are both close to the truth and low. But that begs the question: close to what
truth, $\text{Err}_{XY}$ or Err? In an effort to compare the results of our
simulations to those of Bates et al., (2023), Figure 5 shows the distance of
$\widehat{\text{Err}}^{\text{(LGCV)}}$ from Err compared to its distance from
$\text{Err}_{XY}$. Similar to Bates et al.,, we can see that in the logistic
regression model,

$|\widehat{\text{Err}}^{\text{(LGCV)}}-\text{Err}_{XY}|>|\widehat{\text{Err}}^{\text{(LGCV)}}-\text{Err}|$

The difference lessens as $n\rightarrow\infty$. Regardless, repeated
simulations consistently show that $\widehat{\text{Err}}^{\text{(LGCV)}}$ is,
on average, closer to Err than $\text{Err}_{XY}$.

Figure 5: Simulation results comparing the error of cross-validated estimates
of a logistic regression model when estimating Err to its error when
estimating $\text{Err}_{XY}$: the mean absolute deviation between
$\widehat{\text{Err}}^{\text{(LGCV)}}$ and Err or $\text{Err}_{XY}$. The dark
grey vertical line in each panel is where $n=p$. Notice that across all n,
$\widehat{\text{Err}}^{\text{(LGCV)}}$ is closer to Err than
$\text{Err}_{XY}$. Note that logistic regression models are run only for
$n\geq p$. Figure 6: Simulation results comparing the error of CV estimates of
a Random Forest model when estimating Err to its error when estimating
$\text{Err}_{XY}$: the mean absolute deviation between
$\widehat{\text{Err}}^{\text{(FDCV)}}$ and Err or $\text{Err}_{XY}$. The dark
grey vertical line in each panel is where $n=p$. Notice that across n,
$\widehat{\text{Err}}^{\text{(FDCV)}}$ is closer to $\text{Err}_{XY}$ than
Err.

In Figure 6 we see that this relationship has flipped for Random Forests.
$\widehat{\text{Err}}^{\text{(FDCV)}}$ is closer to $\text{Err}_{XY}$ than
Err. As a reminder, the difference between
$\widehat{\text{Err}}^{\text{(LGCV)}}$ and
$\widehat{\text{Err}}^{\text{(FDCV)}}$ is that the former is an error estimate
for a logistic model while the latter a Random Forest. Both are cross-
validated estimates on the in-fold data set with no parameter tuning. Thus,
the relationship highlighted by Bates et al., (2023) seems to be specific to
generalized linear models as they investigated linear and logistic regression
models.

The flip in relationship may be attributed to the difference in the way each
model utilizes the data. In logistic regression, the coefficients are
estimated via maximum likelihood estimation, thus possibly leading to over-
fitting and biased estimates of error due to the model optimizing for the
specific data set. On the other hand, bagging and other resampling techniques
can be used to reduce the variance in model predictions. In Random Forests,
the bias of the full forest is equivalent to the bias of a single decision
tree (which itself has low bias and high variance) (Hastie et al.,, 2001).
However, by creating many trees and then averaging them, the variance of the
final forest can be greatly reduced over that of a single tree. In practice,
the only limitation on the size of the forest is computing time as an infinite
number of trees could be trained without ever increasing bias and with a
continual (if asymptotically declining) decrease in the variance. As a result,
the logistic regression model may be less informative on the “next” sample,
than a Random Forest. Hence, $\widehat{\text{Err}}$ is closer to
$\text{Err}_{XY}$ than Err for Random Forests because its resampling methods
build the model on data more akin to wild data.

We further explore the difference between $\widehat{\text{Err}}$ for logistic
regression and Random Forests in an experiment with $n=50$ observations and
$p=10$ features; see Figure 7. In the right plot, there seems to exist a
pairing between $\widehat{\text{Err}}^{\text{(FDCV)}}$ and $\text{Err}_{XY}$,
where high estimates of $\widehat{\text{Err}}^{\text{(FDCV)}}$ are paired with
high estimates of $\text{Err}_{XY}$ and vice-versa (i.e., very few of the
linking lines cross). In a logistic regression model, there does not seem to
exist this pairing as seen in the left plot (i.e., most of the linking lines
cross). In other words, $\widehat{\text{Err}}^{\text{(LGCV)}}$ is neither
closer nor farther to Err than to $\text{Err}_{XY}$, and
$\widehat{\text{Err}}^{\text{FDCV}}$ tends to be closer to $\text{Err}_{XY}$
than to Err.

Figure 7: Random sample of 20 iterations of $\widehat{\text{Err}}$ linked with
the corresponding $\text{Err}_{XY}$ for logistic regression (left plot) as
compared to Random Forests (right plot). Each of 20 data sets has $n=50$ and
$p=10$. The lines connecting the boxplots represent the $\widehat{\text{Err}}$
and $\text{Err}_{XY}$ values for a given data set and analysis method. The
crossing lines on the left plot (logistic regression) indicate that that
$\widehat{\text{Err}}$ and $\text{Err}_{XY}$ are not correlated; the non-
crossing lines on the right plot (Random Forest) indicate that
$\widehat{\text{Err}}$ and $\text{Err}_{XY}$ are correlated. Solid black
horizontal line is Err.

Returning to the remaining strategies that all use Random Forests, Figures 8
\- 11 show that across the error estimation strategies, $\widehat{\text{Err}}$
is closer to $\text{Err}_{XY}$ than Err. Despite this relationship, the
differences in mean absolute deviations, from $\text{Err}_{XY}$ and Err, tend
to be quite small.

Figure 8: Simulation results comparing the error of OOB estimates of a Random
Forest model when estimating Err to its error when estimating
$\text{Err}_{XY}$: the mean absolute deviation between
$\widehat{\text{Err}}^{\text{(FDO)}}$ and Err or $\text{Err}_{XY}$. The dark
grey vertical line in each panel is where $n=p$. Figure 9: Simulation results
comparing the error of OOB estimates of a Random Forest model when estimating
Err to its error when estimating $\text{Err}_{XY}$: the mean absolute
deviation between $\widehat{\text{Err}}^{\text{(SDO)}}$ and Err or
$\text{Err}_{XY}$. The dark grey vertical line in each panel is where $n=p$.
Figure 10: Simulation results comparing the error of CV estimates of a Random
Forest model when estimating Err to its error when estimating
$\text{Err}_{XY}$: the mean absolute deviation between
$\widehat{\text{Err}}^{\text{(SDCV)}}$ and Err or $\text{Err}_{XY}$. The dark
grey vertical line in each panel is where $n=p$. Figure 11: Simulation results
comparing the error of validation set estimates of a Random Forest model when
estimating Err to its error when estimating $\text{Err}_{XY}$: the mean
absolute deviation between $\widehat{\text{Err}}^{\text{(SDT)}}$ and Err or
$\text{Err}_{XY}$. The dark grey vertical line in each panel is where $n=p$.

###  4.2 Distance Across Error Estimation Strategies

Section 4.1 detailed our investigation of the distance of the estimates of
error to $\text{Err}_{XY}$ and Err. We showed that in Random Forests, the
error estimates are closer to $\text{Err}_{XY}$ than to Err. Here we assess a
follow-up question: how close is $\widehat{\text{Err}}$ to $\text{Err}_{XY}$?
Figure 12 compares the strategies according to the expected value of
$|\widehat{\text{Err}}-\text{Err}_{XY}|$.

In the case of $p=10$ features, the strategies that utilize the in-fold or in-
bag data set to train the model (LGCV, FDCV, and FDO) seem to outperform the
split data approaches. The strategies without independent validation sets seem
to not over-fit compared the error estimates with train/test splits that
suffer a drop in performance. It is important to mention that FDO is the only
strategy of the three (LGCV, FDCV, and FDO) that tunes parameters.

In the case of $p=100$ features, FDCV seems to be the best candidate when
$n<p$ as well as $n>p$. As was the case when $p=10$ features, FDO performs
better as $n$ increases, and its performance is equivalent to FDCV when $n\geq
p$. In contrast to $p=10$ features, FDCV and FDO perform better than LGCV when
$n>p$. Once again, it is important to note that, out of the strategies that
build the model on the entire data set, FDO tunes the model’s parameters,
compared to FDCV which does not.

Figure 12: Mean absolute deviation between $\widehat{\text{Err}}$ and
$\text{Err}_{XY}$ across error estimation strategies. The dark grey vertical
line in each panel is where $n=p$, note that logistic regression models were
run only for $n\geq p$.

##  5 Discussion

This investigation had two main components. First, we discussed the difference
in error targets presented by Bates et al., (2023). In their work, they found
that in the special case of the generalized linear model using unregularized
OLS for model-fitting, common estimates of prediction error — cross-
validation, bootstrap, data splitting, and covariance penalties — should be
viewed as estimates of the average prediction error, averaged across other
hypothetical data sets from the same distribution. Our primary result is that
in the classification case, Random Forests’ estimates of prediction error can
be taken as an estimate of the true error rate instead of as an estimate of
the average prediction error, which is the opposite of the result of Bates et
al., (2023) whose work was on logistic regression. In simulations the result
held across error estimation strategies such as cross-validation, bagging, and
data splitting (See Figures 6 \- 11). The result was present regardless of
$n,p$. Nonetheless, we wish to be clear that the estimates of prediction error
were a good approximation of both the true error rate and expected error rate
in the data splitting cases.

A fundamental open question is to understand the size of the gap of estimates
of prediction error with the true error rate and expected error rate. The
present work focused on to which target the estimate is closer. Moreover, it
is necessary to understand under what conditions the gap is large, making it
necessary to modify the method of error estimation depending on the target.
Roughly speaking, we expect the gap to be small when $n/p$ is large. In our
simulations, the difference between estimates of prediction error with the
true error rate and expected error rate was always smaller than 0.01. As $n$
increases, however, the difference decreases. Other future directions include
the investigation of the relationship among $\widehat{\text{Err}}$,
$\text{Err}_{XY}$, and Err in correlated data and imbalanced data.

Second, we discussed the performance of a variety of error estimation
strategies. The models built on the entire sample were closer to the true
error rate compared to those built on a training set with error estimates
obtained from a testing set. Therefore, the data strategies that do not use a
holdout set seem more appealing choice for model building, regardless if
parameter tuning is to be performed or not. Empirically, the strategies that
use resampling techniques as opposed to a holdout set are favorable.

Whereas Bates et al., (2023) show $\widehat{\text{Err}}$ to be closer to Err
than $\text{Err}_{XY}$ in generalized linear models, we show
$\widehat{\text{Err}}$ to be closer to $\text{Err}_{XY}$ than Err in Random
Forests. Additionally, resampling techniques seem to outperform data splitting
models in Random Forests.

## References

  * Bates et al., (2023) Bates, S., Hastie, T., and Tibshirani, R. (2023).  Cross-validation: what does it estimate and how well does it do it?
  * Breiman, (2001) Breiman, L. (2001).  Random forests.  Machine Learning, 45:5–32.
  * Bylander, (2002) Bylander, T. (2002).  Estimating generalization error on two-class datasets using out-of-bag estimates.  Machine Learning, 48(1-3):287–297.  Copyright - Kluwer Academic Publishers 2002.
  * Faraway, (2014) Faraway, J. J. (2014).  Does data splitting improve prediction?  Statistics and Computing, 26(1–2):49–60.
  * Genuer et al., (2008) Genuer, R., Poggi, J.-M., and Tuleau, C. (2008).  Random forests: some methodological insights.
  * Goldstein et al., (2010) Goldstein, B., Hubbard, A., Cutler, A., and Barcellos, L. (2010).  An application of random forests to a genome-wide association dataset: Methodological considerations & new findings.  BMC genetics, 11:49.
  * Goldstein et al., (2011) Goldstein, B. A., Polley, E. C., and Briggs, F. B. S. (2011).  Random forests for genetic association studies.  Statistical Applications in Genetics and Molecular Biology, 10(1).
  * Hastie et al., (2001) Hastie, T., Tibshirani, R., and Friedman, J. (2001).  The Elements of Statistical Learning.  Springer Series in Statistics. Springer New York Inc., New York, NY, USA.
  * James et al., (2013) James, G., Witten, D., Hastie, T., and Tibshirani, R. (2013).  An Introduction to Statistical Learning: with Applications in R.  Springer.
  * Janitza and Hornung, (2018) Janitza, S. and Hornung, R. (2018).  On the overestimation of random forest’s out-of-bag error.  PLOS ONE, 13(8):1–31.
  * Kaggle, (2017) Kaggle (2017).  The state of data science & machine learning.
  * Mitchell, (2011) Mitchell, M. (2011).  Bias of the random forest out-of-bag (oob) error for certain input parameters.  Open Journal of Statistics, 01:205–211.
  * Rajanala et al., (2022) Rajanala, S., Bates, S., Hastie, T., and Tibshirani, R. (2022).  Confidence intervals for the generalisation error of random forests.
  * Yousef, (2019) Yousef, W. A. (2019).  A leisurely look at versions and variants of the cross validation estimator.
  * Zhang et al., (2010) Zhang, G.-Y., Zhang, C.-X., and Zhang, J.-S. (2010).  Out-of-bag estimation of the optimal hyperparameter in subbag ensemble method.  Communications in Statistics - Simulation and Computation, 39(10):1877–1892.