# Scaled Brownian motion with random anomalous diffusion exponent

Hubert Woszczek<EMAIL_ADDRESS>Aleksei Chechkin  Agnieszka
Wyłomańska

###### Abstract

The scaled Brownian motion (SBM) is regarded as one of the paradigmatic random
processes, featuring the anomalous diffusion property characterized by the
diffusion exponent. It is a Gaussian, self-similar process with independent
increments, which has found applications across various fields, from
turbulence and stochastic hydrology to biophysics. In our paper, inspired by
recent single particle tracking biological experiments, we introduce a process
termed the scaled Brownian motion with random exponent (SBMRE), which
preserves SBM characteristics at the level of individual trajectories, albeit
with randomly varying anomalous diffusion exponents across the trajectories.
We discuss the main probabilistic properties of SBMRE, including its
probability density function (pdf), and the q-th absolute moment.
Additionally, we present the expected value of the time-averaged mean squared
displacement (TAMSD) and the ergodicity breaking parameter. Furthermore, we
analyze the pdf of the first hitting time in a semi-infinite domain, the
martingale property of SBMRE, and its stochastic exponential. As special
cases, we consider two distributions of the anomalous diffusion exponent,
namely the two-point and beta distributions, and discuss the asymptotics of
the presented characteristics in such cases. Theoretical results for SBMRE are
validated through numerical simulations and compared with the corresponding
characteristics for SBM.

###### keywords:

scaled Brownian motion , anomalous diffusion exponent , random diffusion
exponent , ergodicity breaking , stochastic exponent

###### MSC:

60G22 , 60G65

[inst1]organization=Faculty of Pure and Applied Mathematics, Hugo Steinhaus
Cenetr, Wroclaw University of Science and Technology,addressline=Wyspiańskiego
27, city=Wrocław, postcode=50-370, country=Poland

[inst2]organization=Institute for Physics & Astronomy, University of
Potsdam,addressline=, city=Potsdam-Golm, postcode=14476 , country=Germany

[inst3]organization=Akhiezer Institute for Theoretical Physics, National
Science Center ‘Kharkov Institute of Physics and
Technology’,,addressline=Akademicheskaya st.1, city=Kharkov, postcode=61108,
country=Ukraine

##  1 Introduction

The stochastic processes of anomalous diffusion behavior are ubiquitous in
many areas of interest. They are characterized by the non-linear, mostly
power-law, second moment, namely $\mathbb{E}\left[X^{2}(t)\right]\sim t^{\mu}$
[1]. The $\mu$ parameter is called the anomalous diffusion exponent. Depending
on the $\mu$ value, one can distinguish between sub-diffusive ($\mu<1$) and
super-diffusive ($\mu>1$) behavior. When $\mu=1$, the process exhibits linear
time dependence that is the hallmark of the normal (ordinary) Brownian
diffusion [2].

Mathematically, ordinary Brownian motion, or Wiener process, is the Gaussian
self-similar process with stationary independent increments [3]. There are two
generic Gaussian processes exhibiting anomalous diffusion behavior. A classic
example is the fractional Brownian motion (FBM), introduced by Kolmogorov [4]
and further examined by Mandelbrot and van Ness [5], particularly within the
context of economic time series analysis. FBM is characterized by the Hurst
exponent $H\in(0,1)$, which dictates its anomalous diffusion behavior as
$\mu=2H$. It is the only Gaussian self-similar process with stationary power-
law correlated increments, closely associated with what is known as long-range
phenomena [6]. FBM has found numerous interesting applications across various
areas, including hydrology [7, 8], telecommunications and signal processing
[9, 10, 11], image analysis [12, 13], economics [14, 15, 16], biological
systems such as single-particle tracking experiments [17, 18, 19, 20, 21], and
many others.

Another Gaussian process demonstrating anomalous diffusion behavior is the
scaled Brownian motion (SBM) that was introduced implicitly in the context of
turbulent diffusion in [22]. It was then rediscovered in [23] as a simple
time-changed subordination of the Brownian motion by a deterministic power
function. This process has similar characteristics as FBM, namely it is
Gaussian and self-similar, however, in contrast to FBM its increments are
independent (which implies Markovianity) and non-stationary.

Many authors have examined the SBM from a theoretical perspective. For
instance, in [24, 25], discussions revolve around the asymptotics of the
supremum of SBM, while [26] delves into the Lamperti transformation for SBM.
In [27] SBM is considered as the mean-field model for continuous-time random
walk. In [28] highly nonstationary behavior of the SBM in confining potential
is demonstrated, while the authors of [28, 27, 29] analyze ergodicity breaking
in SBM. In [30] the power spectral density of SBM is explored. The literature
also addresses various modifications of the SBM. For instance, the papers [31]
and [32] explore ultraslow and underdamped versions of SBM, respectively,
while [33, 34] discuss the SBM with resetting. The authors of [35] propose a
modification of SBM to render it stationary. Further studies include SBM with
random diffusion coefficient [36], and scaled geometric Brownian motion [37].
Recently, the authors of [38] explore FBM and SBM on a sphere, investigating
the effects of long-time correlations on navigation strategies. The SBM has
also found intriguing applications, particularly in physical sciences and
biology, as evidenced by works such as [28, 39, 40, 31, 41]. Interested
readers can also explore [42, 43, 44] for different approaches proposed for
testing SBM.

With the recent developments in single particle tracking technique it became
clear that the ”standard” anomalous diffusion models, like e.g., FBM and SBM,
may be inadequate for describing certain complex systems. For example, modern
experiments suggest that the motion of biological cells displays anomalous
diffusion behavior at the individual trajectory level. However, there is
variability in the anomalous diffusion exponent from one trajectory to
another, see e.g. [45, 46, 47, 48]. In such a scenario the natural
generalization of the standard anomalous diffusion model is to take into
account the randomness of the parameters characterizing the random process,
such as diffusivity and/or anomalous diffusion exponent. The idea of such
doubly stochastic behavior is the essence of the superstatistics approach [49,
50] and recently has got its further development in the concept of ”diffusing
diffusivity” suggested in [51] and then explored in a set of papers, see e.g.,
the review [52] and references therein.

The multifractional generalization of FBM, which allows the Hurst exponent to
be a stationary random process was discussed in mathematical literature [53,
54]. Its simplified version, called FBMRE, which considers the Hurst exponent
as a variable randomly changing from trajectory to trajectory is analysed in
[55], where interesting effects like accelerating diffusion and persistence
transitions in course of time is found. We also refer the reader to recent
publications that discuss FBM-based models with varying scenarios for random
parameters [56, 57, 58, 59].

In this paper, we investigate scaled Brownian motion with a random exponent
(SBMRE), which was introduced in [60]. We extend the mathematical description
of SBMRE for general distribution of the anomalous diffusion exponent. We
present the probability density function (pdf) and analyse q-th absolute
moment of SBMRE, expected value of the time-averaged mean squared displacement
(TAMSD), ergodicity breaking parameter, and the first hitting time pdf.
Finally, we also discuss martingale property of SBMRE and its stochastic
exponential. As a special case we examine two generic distributions of the
anomalous diffusion exponent, namely two-point and beta distributions. In such
cases we present the abovementioned characteristics of SBMRE and analyze their
asymptotics. Our theoretical results are supported by numerical analysis
demonstrating the specific behavior of the considered stochastic process and
its differences in contrast to the classical SBM.

The rest of the paper is organized as follows. In Section 2 we recall all the
definitions needed. In Section 3 we remind the readers important properties of
SBM while in Section 4 we present the same characteristics for SBMRE for
general distribution of the anomalous diffusive exponent. In Section 5 we
present the results for two exemplary distributions of the anomalous diffusive
exponent. In Section 6 we confirm the theoretical results with numerical
simulations. Last section concludes the paper.

##  2 Preliminaries

In this section, we introduce all definitions and notations used in the
following sections. We assume that all processes take values in $\mathbb{R}$.
If we do not express it explicitly, we assume, that corresponding probability
space is
$\left(\Omega,\mathcal{F},\mathbb{P},\left\\{\mathcal{F}_{t}\right\\}\right)$.

######  Definition 2.1.

Let $\left\\{X\left(t\right)\right\\}_{t\geq 0}$ be the stochastic process .
Then, the time averaged mean square displacement (TAMSD) of the process
$\left\\{X\left(t\right)\right\\}_{t\geq 0}$ is defined as follows [19]:

$\delta\left(\tau\right)=\frac{1}{T-\tau}\int_{0}^{T-\tau}\left(X\left(t+\tau\right)-X\left(t\right)\right)^{2}dt,$
(1)

where $\tau\in\left[0,T\right)$ is the time lag of measured time series,
represented by the stochastic process and $T\in\left(0,\infty\right)$ is a
time horizon (the trajectory length).

Next, we introduce an important concept, called an ergodicity breaking
parameter, which helps to quantify the ergodicity of a given process.

######  Definition 2.2.

Let $\left\\{X\left(t\right)\right\\}_{t\geq 0}$ be the stochastic process.
Then, the ergodicity breaking parameter of the process
$\left\\{X\left(t\right)\right\\}_{t\geq 0}$ is defined as follows [61, 62]:

$EB\left(\tau\right)=\frac{\mathbb{E}\left[\delta^{2}\left(\tau\right)\right]-\mathbb{E}\left[\delta\left(\tau\right)\right]^{2}}{\mathbb{E}\left[\delta\left(\tau\right)\right]^{2}}=\frac{\mathcal{N}\left(\tau\right)}{\mathcal{D}\left(\tau\right)}.$
(2)

Next, we define the hitting time in the barrier $b>0$ of a given stochastic
process.

######  Definition 2.3.

Let $\left\\{X\left(t\right)\right\\}_{t\geq 0}$ be the stochastic process.
Then the hitting time in the barrier $b>0$ of the process
$\left\\{X\left(t\right)\right\\}_{t\geq 0}$ is a random variable defined as
follows:

$\tau_{b}=\inf\left\\{t\geq 0:X\left(t\right)=b\right\\}.$ (3)

The next two definitions are the quadratic variation and the stochastic
exponent. These are important tools when we focus on the martingale properties
of the process.

######  Definition 2.4.

Let $\left\\{X\left(t\right)\right\\}_{t\geq 0}$ be the semimartingale (for
more details, see [63]). Then the quadratic variation of the process
$\left\\{X\left(t\right)\right\\}_{t\geq 0}$ is defined by [63]:

$\left[X\right]_{t}=\lim_{\Delta t_{k}\to 0}\sum_{t_{k}\leq
t}\left(X\left(t_{k+1}\right)-X\left(t_{k}\right)\right)^{2},$ (4)

where $0=t_{1}<t_{2}<\ldots<t_{n}=t$, $\Delta t_{k}=t_{k+1}-t_{k}$ and limit
is taken in probability.

######  Definition 2.5.

Let $\left\\{X\left(t\right)\right\\}_{t\geq 0}$ be the semimatringale. Then
the stochastic exponential of the process
$\left\\{X\left(t\right)\right\\}_{t\geq 0}$ is given by [63]:

$Y\left(t\right)=\exp\left\\{\lambda
X\left(t\right)-\frac{\lambda^{2}}{2}\left[X\right]_{t}\right\\},\;\lambda>0.$
(5)

The last two definitions are the big $\mathcal{O}$ and the little $o$
notation. They are important for providing asymptotic formulas. In both
definitions, we consider functions $f\left(\cdot\right)$ and
$g\left(\cdot\right)$ whose domain and co-domain are $\mathbb{R}$.

######  Definition 2.6.

We say, that $f\left(x\right)\in\mathcal{O}\left(g\left(x\right)\right)$ as
$x\to\infty$, if there exist $c>0$ and $x_{0}\in\mathbb{R}$, such that for all
$x\geq x_{0}$ we have [64]:

$\left|f\left(x\right)\right|\leq c\cdot g\left(x\right).$ (6)

######  Definition 2.7.

We say, that $f\left(x\right)\in o\left(g\left(x\right)\right)$ as
$x\to\infty$, if for all $c>0$ there exists $x_{0}\in\mathbb{R}$, such that
for all $x\geq x_{0}$ we have [64]:

$\left|f\left(x\right)\right|\leq c\cdot g\left(x\right).$ (7)

##  3 Scaled Brownian motion

In this section we recall the definition of the scaled Brownian motion and
describe its main properties.

######  Definition 3.1.

Let $\left\\{B\left(t\right)\right\\}_{t\geq 0}$ be the standard Brownian
motion. Then, the scaled Brownian motion (SBM) is defined as follows:

$B_{\alpha}\left(t\right)=B\left(t^{\alpha}\right),$ (8)

where $\alpha>0$ is an anomalous diffusion exponent. It is a Gaussian self-
similar process with independent increments (Markovianity), zero mean, and
autocovariance
$Cov\left(B_{\alpha}\left(s\right),B_{\alpha}\left(t\right)\right)=\min\left\\{t^{\alpha},s^{\alpha}\right\\}$.
SBM is subdiffusive process for $\alpha<1$ and superdiffusive process for
$\alpha>1$. For $\alpha=1$ the SBM becomes an ordinary Brownian motion. The
probability density function (pdf) is given by

$p\left(x,t\right)=\frac{1}{\sqrt{2\pi
t^{\alpha}}}\exp\left\\{-\frac{x^{2}}{2t^{\alpha}}\right\\},\quad
x\in\mathbb{R}.$ (9)

The SBM can also be represented by the following stochastic differential
equation (SDE) [23]:

$dB_{\alpha}\left(t\right)=\sqrt{\alpha}t^{\frac{\alpha-1}{2}}dB\left(t\right),\quad
B_{\alpha}\left(0\right)=0.$ (10)

The solution of (10) reads

$B_{\alpha}\left(t\right)=\int_{0}^{t}\sqrt{\alpha}s^{\frac{\alpha-1}{2}}dB\left(s\right).$
(11)

In the next part of this section we present the main properties of SBM.

######  Fact 3.1.

Let $\left\\{B_{\alpha}\left(t\right)\right\\}_{t\geq 0}$ be the SBM. Then,
the expected value of the TAMSD for SBM is given by [28]:

$\mathbb{E}\left[\delta_{\alpha}\left(\tau\right)\right]=\frac{T^{\alpha+1}-\tau^{\alpha+1}-\left(T-\tau\right)^{\alpha+1}}{\left(\alpha+1\right)\left(T-\tau\right)}.$
(12)

######  Fact 3.2.

Let $\left\\{B_{\alpha}\left(t\right)\right\\}_{t\geq 0}$ be the SBM. Then,
the ergodicity breaking parameter for SBM is given by [29]

$EB_{\alpha}\left(\tau\right)=\frac{\mathcal{N}_{\alpha}\left(\tau\right)}{\mathcal{D}_{\alpha}\left(\tau\right)},$
(13)

where

$\displaystyle\mathcal{N}_{\alpha}\left(\tau\right)=\frac{4\tau^{2\alpha+2}}{\left(T-\tau\right)^{2}}\left[\frac{\left(T/\tau-1\right)^{2\alpha+1}}{2\alpha+1}+\frac{\left(3\alpha+1\right)\left(T/\tau-1\right)^{2\alpha+2}}{2\left(\alpha+1\right)^{2}\left(2\alpha+1\right)}-\frac{2\left(T/\tau\right)^{\alpha+1}\left(T/\tau-1\right)^{\alpha+1}}{\left(\alpha+1\right)^{2}}\right.+$
(14)
$\displaystyle+\left.\frac{\left(T/\tau\right)^{2\alpha+2}}{2\left(\alpha+1\right)\left(2\alpha+1\right)}-\frac{\left(2\alpha^{2}+\alpha+1\right)}{2\left(\alpha+1\right)^{2}\left(2\alpha+1\right)}+\frac{2}{\alpha+1}\int_{0}^{T/\tau-1}x^{\alpha+1}\left(x+1\right)^{\alpha}dx\right]$

and $\mathcal{D}_{\alpha}\left(\tau\right)$ is the square of (12).

######  Fact 3.3.

Let $\left\\{B_{\alpha}\left(t\right)\right\\}_{t\geq 0}$ be the SBM. Let
$\tau_{b}=\inf\left\\{t\geq 0:B_{\alpha}\left(t\right)=b\right\\}$ be the time
of first hitting a point $b>0$. Then, the pdf of $\tau_{b}$ is given by [33]

$f_{\tau_{b}}\left(t\right)=\frac{\alpha
b}{\sqrt{2\pi}}e^{-\frac{b^{2}}{2t^{\alpha}}}t^{-1-\frac{\alpha}{2}},\quad
t\geq 0.$ (15)

######  Theorem 3.1.

Let $\left\\{B_{\alpha}\left(t\right)\right\\}_{t\geq 0}$ be the SBM. Then,
SBM is a martingale with respect to its natural filtration. Its quadratic
variation $\left[B_{\alpha}\right]_{t}=t^{\alpha}$. The stochastic exponential
of SBM is also a martingale with respect to its natural filtration.

###### Proof.

We start with the proof that SBM is a martingale. The first moment of absolute
value of SBM is, obviously, finite. For $s<t$ we have:

$\mathbb{E}\left[B_{\alpha}\left(t\right)\;|\;\mathcal{F}_{s}\right]=\mathbb{E}\left[B_{\alpha}\left(t\right)-B_{\alpha}\left(s\right)+B_{\alpha}\left(s\right)\;|\;\mathcal{F}_{s}\right]=\mathbb{E}\left[B_{\alpha}\left(t\right)-B_{\alpha}\left(s\right)\right]+B_{\alpha}\left(s\right)=B_{\alpha}\left(s\right).$
(16)

Next, using the SDE representation (10), we calculate the quadratic variation
of the SBM. Namely, we have:

$\left[B_{\alpha}\right]_{t}=\int_{0}^{t}\left(\sqrt{\alpha}s^{\frac{\alpha-1}{2}}\right)^{2}ds=\int_{0}^{t}\alpha
s^{\alpha-1}ds=t^{\alpha}.$ (17)

Now, we will show that the stochastic exponential of SBM is a martingale.
Indeed, we have:

$\mathbb{E}\left[e^{B_{\alpha}\left(t\right)}\;|\;\mathcal{F}_{s}\right]=\mathbb{E}\left[e^{B_{\alpha}\left(t\right)-B_{\alpha}\left(s\right)}e^{B_{\alpha}\left(s\right)}\;|\;\mathcal{F}_{s}\right]=e^{B_{\alpha}\left(s\right)}\mathbb{E}\left[e^{B_{\alpha}\left(t\right)-B_{\alpha}\left(s\right)}\right]=e^{\frac{t^{\alpha}-s^{\alpha}}{2}}e^{B_{\alpha}\left(s\right)}.$
(18)

The last equality is a straightforward consequence of the Gaussianity.
Finally, we have the following:

$\mathbb{E}\left[e^{B_{\alpha}\left(t\right)-\frac{t^{\alpha}}{2}}\;|\;\mathcal{F}_{s}\right]=e^{B_{\alpha}\left(s\right)-\frac{s^{\alpha}}{2}}.$
(19)

As a consequence, it also gives us the integrability of stochastic
exponential. ∎

##  4 Scaled Brownian motion with random anomalous diffusion exponent

In this section, we define scaled Brownian motion with random anomalous
diffusion exponent and investigate its main properties.

######  Definition 4.1.

Let $\mathcal{A}$ be some positive random variable with values from the
interval $\left(0,K\right)$ with pdf $f_{\mathcal{A}}\left(\cdot\right)$ and a
moment generating function $M_{\mathcal{A}}\left(\cdot\right)$:

$M_{\mathcal{A}}(s)=\mathbb{E}\left[e^{s\mathcal{A}}\right]=\int_{0}^{K}e^{sa}f_{\mathcal{A}}\left(a\right)da.$
(20)

$K>0$ is finite constant.

Then, the scaled Brownian motion with random anomalous diffusion exponent
(SBMRE) is defined as follows:

$B_{\mathcal{A}}\left(t\right)=B\left(t^{\mathcal{A}}\right).$ (21)

We assume that $\left\\{B\left(t\right)\right\\}_{t\geq 0}$ and $\mathcal{A}$
are independent. SBMRE can also be represented as the following SDE:

$dB_{\mathcal{A}}\left(t\right)=\sqrt{\mathcal{A}}t^{\frac{\mathcal{A}-1}{2}}dB\left(t\right),\quad
B_{\mathcal{A}}\left(0\right)=0.$ (22)

The solution of (22) is given by:

$B_{\mathcal{A}}\left(t\right)=\int_{0}^{t}\sqrt{\mathcal{A}}s^{\frac{\mathcal{A}-1}{2}}dB\left(s\right).$
(23)

The representation via SDE (22) is a direct consequence of random time change
for Brownian motion (for the details see Chapter 8 in [63]).

######  Fact 4.1.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the SBMRE.
Then, its pdf is given by the following formula:

$p_{B_{\mathcal{A}}}\left(x,t\right)=\int_{0}^{K}\frac{1}{\sqrt{2\pi
t^{a}}}exp\left\\{\frac{-x^{2}}{2t^{a}}\right\\}f_{\mathcal{A}}\left(a\right)da.$
(24)

###### Proof.

It is a direct consequence of the law of total probability. ∎

######  Fact 4.2.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the SBMRE.
Then, the moment of order q of the absolute value of
$\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ is given by:

$\mathbb{E}\left[\left|B_{\mathcal{A}}\left(t\right)\right|^{q}\right]=c_{q}M_{\mathcal{A}}\left(\frac{1}{2}q\log
t\right),$ (25)

where $c_{q}=\frac{2^{q/2}\Gamma\left(\frac{q+1}{2}\right)}{\sqrt{\pi}}$.

###### Proof.

$\displaystyle\mathbb{E}\left[\left|B_{\mathcal{A}}\left(t\right)\right|^{q}\right]=\mathbb{E}\left[\mathbb{E}\left[\left|B_{\mathcal{A}}\left(t\right)\right|^{q}\;|\;\mathcal{A}\right]\right]=\mathbb{E}\left[\mathbb{E}\left[\left|B\left(t^{\mathcal{A}}\right)\right|^{q}\;|\;\mathcal{A}\right]\right]=\mathbb{E}\left[\mathbb{E}\left[t^{\frac{1}{2}\mathcal{A}q}\left|B\left(1\right)\right|^{q}\;|\;\mathcal{A}\right]\right]=\mathbb{E}\left[t^{\frac{1}{2}\mathcal{A}q}Z^{q}\right]=$
(26) $\displaystyle=c_{q}\mathbb{E}\left[e^{\frac{1}{2}q\mathcal{A}\log
t}\right]=c_{q}M_{\mathcal{A}}\left(\frac{1}{2}q\log t\right),$

where $Z=\left|B\left(1\right)\right|^{q}$ and $\Gamma\left(\cdot\right)$ is
the gamma function. ∎

######  Corollary 4.1.

The second moment of $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq
0}$ is given by:

$\mathbb{E}\left[B_{\mathcal{A}}^{2}\left(t\right)\right]=M_{\mathcal{A}}\left(\log
t\right).$ (27)

######  Fact 4.3.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the SBMRE.
Then, the autocovariance function of SBMRE is as follows:

$Cov\left(B_{\mathcal{A}}\left(s\right),B_{\mathcal{A}}\left(t\right)\right)=M_{\mathcal{A}}\left(\log\left(\min\left\\{s,t\right\\}\right)\right).$
(28)

###### Proof.

$\displaystyle
Cov\left(B_{\mathcal{A}}\left(s\right),B_{\mathcal{A}}\left(t\right)\right)=\mathbb{E}\left[B_{\mathcal{A}}\left(s\right)B_{\mathcal{A}}\left(t\right)\right]=\mathbb{E}\left[\mathbb{E}\left[B_{\mathcal{A}}\left(s\right)B_{\mathcal{A}}\left(t\right)\mathcal{\;}|\;A\right]\right]=\mathbb{E}\left[\min\left\\{s^{\mathcal{A}},t^{\mathcal{A}}\right\\}\right]=$
(29)
$\displaystyle=\mathbb{E}\left[\min\left\\{s,t\right\\}^{\mathcal{A}}\right]=\mathbb{E}\left[\exp\left\\{\mathcal{A}\log\left(\min\left\\{s,t\right\\}\right\\}\right)\right]=M_{\mathcal{A}}\left(\log\left(\min\left\\{s,t\right\\}\right)\right).$

∎

######  Fact 4.4.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the SBMRE.
Then, the expected value of TAMSD for SBMRE is given by:

$\mathbb{E}\left[\delta_{\mathcal{A}}\left(\tau\right)\right]=\frac{1}{T-\tau}\int_{0}^{T-\tau}\left(M_{\mathcal{A}}\left(\log\left(t+\tau\right)\right)-M_{\mathcal{A}}\left(\log
t\right)\right)dt.$ (30)

###### Proof.

$\displaystyle\mathbb{E}\left[\delta_{\mathcal{A}}\left(\tau\right)\right]=\mathbb{E}\left[\frac{1}{T-\tau}\int_{0}^{T-\tau}\left(B_{\mathcal{A}}\left(t+\tau\right)-B_{\mathcal{A}}\left(t\right)\right)^{2}dt\right]=\mathbb{E}\left[\frac{1}{T-\tau}\int_{0}^{T-\tau}(B^{2}_{\mathcal{A}}\left(t+\tau\right)-2B_{\mathcal{A}}\left(t+\tau\right)B_{\mathcal{A}}\left(t\right)+B^{2}_{\mathcal{A}}\left(t\right)dt\right]=$
(31)
$\displaystyle=\frac{1}{T-\tau}\int_{0}^{T-\tau}\left(M_{\mathcal{A}}\left(\log\left(t+\tau\right)\right)-M_{\mathcal{A}}\left(\log
t\right)\right)dt.$

The change in the order of integration is made according to the Fubinis
theorem [65]. ∎

######  Fact 4.5.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the SBMRE.
Then, the ergodicity breaking parameter of SBMRE is given by:

$EB_{\mathcal{A}}\left(\tau\right)=\frac{4\int_{0}^{\tau}\int_{0}^{T-\tau-
t_{2}}\left(M_{\mathcal{A}}\left(\log\left(t_{1}+\tau\right)\right)-M_{\mathcal{A}}\left(\log\left(t_{1}+t_{2}\right)\right)\right)^{2}dt_{1}dt_{2}}{\int_{0}^{T-\tau}\int_{0}^{T-\tau}\left(M_{\mathcal{A}}\left(\log\left(t_{1}+\tau\right)\right)-M_{\mathcal{A}}\left(\log
t_{1}\right)\right)\left(M_{\mathcal{A}}\left(\log\left(t_{2}+\tau\right)\right)-M_{\mathcal{A}}\left(\log
t_{2}\right)\right)dt_{1}dt_{2}}.$ (32)

###### Proof.

Using $(17)$ from [29] we obtain

$\mathcal{N}_{\mathcal{A}}\left(\tau\right)=\frac{4}{\left(T-\tau\right)^{2}}\int_{0}^{\tau}\int_{0}^{T-\tau-
t_{2}}\left(M_{\mathcal{A}}\left(\log\left(t_{1}+\tau\right)\right)-M_{\mathcal{A}}\left(\log\left(t_{1}+t_{2}\right)\right)\right)^{2}dt_{1}dt_{2}.$
(33)

By taking a square of (1) and changing the square of the integral into a
double integral, we obtain:

$\mathcal{D}_{\mathcal{A}}\left(\tau\right)=\frac{1}{(T-\tau)^{2}}\int_{0}^{T-\tau}\int_{0}^{T-\tau}\left(M_{\mathcal{A}}\left(\log\left(t_{1}+\tau\right)\right)-M_{\mathcal{A}}\left(\log
t_{1}\right)\right)\left(M_{\mathcal{A}}\left(\log\left(t_{2}+\tau\right)\right)-M_{\mathcal{A}}\left(\log
t_{2}\right)\right)dt_{1}dt_{2}.$ (34)

∎

Then, using (2) we receive (32).

######  Fact 4.6.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the SBMRE.
Then, the distribution of the first hitting time in barrier b for SBMRE
$\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ has the following
pdf:

$f_{\tau_{b}}\left(t\right)=\int_{0}^{K}\frac{ab}{\sqrt{2\pi}}e^{-\frac{b^{2}}{2t^{a}}}t^{-1-\frac{a}{2}}f_{\mathcal{A}}\left(a\right)da.$
(35)

###### Proof.

We use (15) and the law of total probability. ∎

Next we discuss the martingale properties or SBMRE. In the proof, we follow
the same steps as in the proof of Theorem $2.1$ in [66]. First, we define the
filtration $\left\\{\mathcal{F}_{\tau}\right\\}_{\tau\geq 0}$ given by the
following formula:

$\mathcal{F}_{\tau}=\bigcap_{u>\tau}\left\\{\sigma\left(B\left(y\right)\;:\;0\leq
y\leq u\right)\vee\sigma\left(y^{\mathcal{A}}:\;y\geq 0\right)\right\\}.$ (36)

$\mathcal{F}_{2}\vee\mathcal{F}_{2}$ denotes the $\sigma-$ algebra generated
by the union of the $\sigma-$ algebras $\mathcal{F}_{1},\;\mathcal{F}_{2}$.

######  Theorem 4.1.

Let $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ be the
SBMRE.Suppose that the density of the random variable $\mathcal{A}$ is
continuous. Then, $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$
is a martingale with respect to filtration
$\left\\{\mathcal{F}_{t^{\mathcal{A}}}\right\\}_{t\geq 0}$, where
$\left\\{\mathcal{F}_{\tau}\right\\}_{\tau\geq 0}$ is given by (36). The
quadratic variation of $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq
0}$ is given by $\left[B_{\mathcal{A}}\right]_{t}=t^{\mathcal{A}}$.
Furthermore, the stochastic exponential of
$\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ is also
$\left\\{\mathcal{F}_{t^{\mathcal{A}}}\right\\}_{t\geq 0}-$martingale.

###### Proof.

We start with observations that $\left\\{\mathcal{F}_{\tau}\right\\}$ is
right-continuous, $B\left(\tau\right)$ is
$\left\\{\mathcal{F}_{\tau}\right\\}-$martingale, and for every fixed
$t_{0}>0$ the random variable $t_{0}^{\mathcal{A}}$ is the stopping time with
respect to $\left\\{\mathcal{F}_{t^{\mathcal{A}}}\right\\}$. These facts
ensure that
$\mathcal{G}_{t}=\left\\{\mathcal{F}_{t_{0}^{\mathcal{A}}}\right\\}$ is well
defined. Next, we introduce a sequence of
$\left\\{\mathcal{F}_{\tau}\right\\}$-stopping times defined by

$T_{n}=\inf\left\\{\tau>0:\;\left|B\left(\tau\right)\right|=n\right\\}.$ (37)

Using Doob’s optional sampling theorem [67] we obtain

$\mathbb{E}\left[B\left(T_{n}\wedge
t^{\mathcal{A}}\right)\;|\;\mathcal{G}_{s}\right]=B\left(T_{n}\wedge
s^{\mathcal{A}}\right)\quad s<t.$ (38)

Let us note that as $n\to\infty$ the left-hand side and the right-hand side of
(38) converge to
$\mathbb{E}\left[B\left(t^{\mathcal{A}}\right)\;|\;\mathcal{G}_{s}\right]$ and
$B\left(s^{\mathcal{A}}\right)$, respectively. Thus,
$\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ is
$\left\\{\mathcal{G}_{t}\right\\}-$ martingale.

From the SDE representation of SBMRE (22) we have that the quadratic variation
of $\left\\{B_{\mathcal{A}}\left(t\right)\right\\}_{t\geq 0}$ is given by:

$\left[B_{\mathcal{A}}\right]_{t}=\int_{0}^{t}\left(\sqrt{\mathcal{A}}s^{\frac{\mathcal{A}-1}{2}}\right)^{2}ds=\int_{0}^{t}\mathcal{A}s^{\mathcal{A}-1}ds=t^{\mathcal{A}}.$
(39)

Now we show that the stochastic exponential
$\left\\{Y\left(t\right)\right\\}_{t\geq 0}$ of SBMRE is also a martingale.
Using the Proposition $3.4$ from [67] we conclude that
$\left\\{Y\left(t\right)\right\\}_{t\geq 0}$ is a local martingale. Hence, it
is sufficient to show that $\sup_{0\leq u\leq t}Y\left(u\right)$ is integrable
to conclude that $\left\\{Y\left(t\right)\right\\}_{t\geq 0}$ is a martingale.
We have the following:

$\sup_{0\leq u\leq t}Y\left(u\right)\leq\sup_{0\leq u\leq
t}\exp\left\\{\lambda B\left(u^{\mathcal{A}}\right)\right\\}.$ (40)

In addition, $\exp\left\\{\lambda B\left(u^{\mathcal{A}}\right)\right\\}$ is a
positive submartingale. Thus, using Doob’s maximal inequality, we get:

$\mathbb{E}\left[\left(\sup_{0\leq u\leq t}\exp\left\\{\lambda
B\left(u^{\mathcal{A}}\right)\right\\}\right)^{2}\right]\leq
4\mathbb{E}\left[\exp\left\\{2\lambda
B\left(u^{\mathcal{A}}\right)\right\\}\right].$ (41)

Now, using the properties of Gaussian random variable and the conditioning on
$\sigma\left(y^{\mathcal{A}}:\;y\geq 0\right)$, we have:

$\displaystyle\mathbb{E}\left[\exp\left\\{2\lambda
B\left(t^{\mathcal{A}}\right)\right\\}\right]\leq\mathbb{E}\left[\exp\left\\{2\lambda^{2}t^{\mathcal{A}}\right\\}\right]=\sum_{n=0}^{\infty}\frac{\left(2\lambda^{2}\right)^{n}\mathbb{E}\left[t^{n\mathcal{A}}\right]}{n!}=\sum_{n=0}^{\infty}\frac{\left(2\lambda^{2}\right)^{n}\mathcal{M}_{\mathcal{A}}\left(n\log\left(t\right)\right)}{n!}=$
(42)
$\displaystyle=\sum_{n=0}^{\infty}\frac{\left(2\lambda^{2}\right)^{n}\int_{0}^{K}e^{an\log\left(t\right)}f_{\mathcal{A}}\left(a\right)da}{n!}=\int_{0}^{K}f_{\mathcal{A}}\left(a\right)\sum_{n=0}^{\infty}\frac{\left(2\lambda^{2}\right)^{n}e^{an\log\left(t\right)}}{n!}da=\int_{0}^{K}e^{2t^{a}\lambda^{2}}f_{\mathcal{A}}\left(a\right)da\leq
e^{2t^{c}\lambda^{2}}f_{\mathcal{A}}\left(c\right),$

where $c\in\left(0,K\right)$. The last inequality is a straightforward
consequence of the mean value theorem and the change in order of integration
and summation is performed according to the dominated convergence theorem
[65]. Finally, $\mathbb{E}\left[\sup_{0\leq u\leq
t}Y\left(u\right)\right]<\infty$, thus
$\left\\{Y\left(t\right)\right\\}_{t\geq 0}$ is a martingale. ∎

##  5 Selected distributions of anomalous diffusion exponent

In this section, we apply the general formulas calculated in the previous
section to characterize the SBMRE with two-point and beta distributions of
random anomalous diffusion exponent.

######  Example 5.1.

First, we consider the case of the simplest distribution of the random
variable $\mathcal{A}$, that is, the two-point distribution concentrated at
two points $A_{1},A_{1}\in\left(0,K\right)$, $A_{1}<A_{2}$ with pdf given by:

$f_{\mathcal{A}}(a)=p\delta\left(a-A_{1}\right)+\left(1-p\right)\delta\left(a-A_{2}\right),$
(43)

where the weight parameter $p\in\left(0,1\right)$ and
$\delta\left(\cdot\right)$ is the Dirac delta function. The moment generating
function of the two-point distribution is given by:

$M_{\mathcal{A}}(s)=pe^{A_{1}s}+\left(1-p\right)e^{A_{2}s}.$ (44)

Pdf of SBMRE in this case is given by:

$p_{B_{\mathcal{A}}}\left(x,t\right)=\frac{p}{\sqrt{2\pi
t^{A_{1}}}}\exp\left\\{-\frac{x^{2}}{2t^{A_{1}}}\right\\}+\frac{1-p}{\sqrt{2\pi
t^{A_{2}}}}\exp\left\\{-\frac{x^{2}}{2t^{A_{2}}}\right\\},\quad
x\in\mathbb{R}.$ (45)

It is easy to see that (45) is a mixture of pdfs of two Gaussian distributed
random variables with zero mean and variances equal to $t^{A_{1}}$ and
$t^{A_{2}}$, respectively.

In the considered case the second moment of SBMRE has the following form:

$\mathbb{E}\left[B_{\mathcal{A}}^{2}\left(t\right)\right]=pt^{A_{1}}+\left(1-p\right)t^{A_{2}}.$
(46)

Notify that for $t\ll 1$ $t^{A_{1}}\in o\left(t^{A_{2}}\right)$ and for $t\gg
1$ $t^{A_{2}}\in o\left(t^{A_{1}}\right)$. Thus we have the following
asymptotics:

$\mathbb{E}\left[B_{\mathcal{A}}^{2}\left(t\right)\right]\sim t^{A},$ (47)

where $A=A_{1}$ for $t\rightarrow 0$ and $A=A_{2}$ for $t\rightarrow\infty$.
The effect is called accelerating diffusion [68, 69, 70].

The expectation of TAMSD reads:

$\mathbb{E}\left[\delta_{\mathcal{A}}\left(\tau\right)\right]=\frac{1}{T-\tau}\left(p\frac{T^{A_{1}+1}-\left(T-\tau\right)^{A_{1}+1}-\tau^{A_{1}+1}}{A_{1}+1}+\left(1-p\right)\frac{T^{A_{2}+1}-\left(T-\tau\right)^{A_{2}+1}-\tau^{A_{2}+1}}{A_{2}+1}\right),$
(48)

Taking the first three components of Taylor expansion of $(T-\tau)^{A_{1}+1}$
and $(T-\tau)^{A_{2}+1}$ for $\tau/T<<1$ we obtain:

$\mathbb{E}\left[\delta_{\mathcal{A}}\left(\tau\right)\right]\sim p\left(\tau
T^{A_{1}-1}+\frac{\tau^{A_{1}+1}}{\left(A_{1}+1\right)T}-A_{1}\tau^{2}T^{A_{1}-2}\right)+\left(1-p\right)\left(\tau
T^{A_{2}-1}+\frac{\tau^{A_{2}+1}}{\left(A_{2}+1\right)T}-A_{2}\tau^{2}T^{A_{2}-2}\right).$
(49)

Thus, for the asymptotics at $\tau/T<<1$ (the inequality typical for
experiments [19]) we get:

$\mathbb{E}\left[\delta_{\mathcal{A}}\left(\tau\right)\right]\sim\frac{\tau}{T}\left(pT^{A_{1}}+(1-p)T^{A_{2}}\right),$
(50)

therefore for the trajectories long enough the $T$-dependence of the
expectation value is determined by larger exponent.

The EB parameter in this case is given by the following formula:

$EB_{\mathcal{A}}\left(\tau\right)=\frac{4\int_{0}^{\tau}\int_{0}^{T-\tau-
t_{2}}\left(p\left(\left(t_{1}+\tau\right)^{A_{1}}-\left(t_{1}+t_{2}\right)^{A_{1}}\right)+\left(1-p\right)\left(\left(t_{1}+\tau\right)^{A_{2}}-\left(t_{1}+t_{2}\right)^{A_{2}}\right)\right)^{2}dt_{1}dt_{2}}{\int_{0}^{T-\tau}\int_{0}^{T-\tau}\left(p\left(\left(t_{1}+\tau\right)^{A_{1}}-t_{1}^{A_{1}}\right)+\left(1-p\right)\left(\left(t_{1}+\tau\right)^{A_{2}}-t_{1}^{A_{2}}\right)\right)\left(p\left(\left(t_{2}+\tau\right)^{A_{1}}-t_{2}^{A_{1}}\right)+\left(1-p\right)\left(\left(t_{2}+\tau\right)^{A_{2}}-t_{2}^{A_{2}}\right)\right)dt_{1}dt_{2}}.$
(51)

The first hitting time in the barrier $b$ in the considered case is as
follows:

$f_{\tau_{b}}\left(t\right)=\frac{pA_{1}b}{\sqrt{2\pi}}\exp\left\\{-\frac{b^{2}}{2t^{A_{1}}}\right\\}t^{-1-\frac{A_{1}}{2}}+\frac{(1-p)A_{2}b}{\sqrt{2\pi}}\exp\left\\{-\frac{b^{2}}{2t^{A_{2}}}\right\\}t^{-1-\frac{A_{2}}{2}}.$
(52)

Notify that for $t\gg 1$ we have the following:
$\exp\left\\{-\frac{b^{2}}{2t^{A_{1}}}\right\\},\,\exp\left\\{-\frac{b^{2}}{2t^{A_{1}}}\right\\}\in\mathcal{O}\left(1\right)$
and $t^{-1-\frac{A_{2}}{2}}\in o\left(t^{-1-\frac{A_{1}}{2}}\right)$. Taking
into account these facts, we obtain the following asymptotic formula for (52)
in case $t\gg 1$:

$f_{\tau_{b}}\left(t\right)\sim t^{-1-A_{1}/2}.$ (53)

######  Example 5.2.

The second considered distribution is the beta distribution on the interval
$\left[A_{1},A_{2}\right]$, $0<A_{1}<A_{2}$ with pdf

$f_{\mathcal{A}}(a)=\frac{\left(a-A_{1}\right)^{\gamma-1}\left(A_{2}-a\right)^{\beta-1}}{\mathbb{B}\left(\gamma,\beta\right)\left(A_{2}-A_{1}\right)^{\gamma+\beta-1}}\mathbbm{1}_{\left[A_{1},A_{2}\right]},$
(54)

where $\gamma,\beta>0$ are beta distribution parameters and
$\mathbb{B}\left(\cdot,\cdot\right)$ is beta function. The moment generating
function of the beta distribution is given by:

$M_{\mathcal{A}}(s)=e^{A_{1}s}\prescript{}{1}{F_{1}\left(\alpha,\alpha+\beta,s\left(A_{2}-A_{1}\right)\right)},$
(55)

where $\prescript{}{1}{F_{1}\left(\cdot,\cdot,\cdot\right)}$ is a confluent
hypergeometric function (we refer to [71] for more details).

The pdf of SBMRE in this case is given by:

$p_{B_{\mathcal{A}}}\left(x,t\right)=\frac{1}{\left(A_{2}-A_{1}\right)^{\gamma+\beta-1}\mathbb{B}\left(\gamma,\beta\right)}\int_{A_{1}}^{A_{2}}\frac{\left(a-A_{1}\right)^{\gamma-1}\left(A_{2}-a\right)^{\beta-1}\exp\left\\{-\frac{x^{2}}{2t^{a}}\right\\}}{\sqrt{2\pi
t^{a}}}da,\quad x\in\mathbb{R}.$ (56)

In the considered case one can express the second moment of SBMRE as

$\mathbb{E}\left[B_{\mathcal{A}}^{2}\left(t\right)\right]=t^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t\right)\right)}.$
(57)

Using the asymptotics of the hypergeometric function [71] we obtain that (57)
behaves asymptotically as

$\mathbb{E}\left[B_{\mathcal{A}}^{2}\left(t\right)\right]\sim\frac{\Gamma\left(\gamma+\beta\right)t^{A_{1}}}{\Gamma\left(\gamma\right)\left(\left(A_{2}-A_{1}\right)\log\left(1/t\right)\right)^{\gamma}}$
(58)

for $t\ll 1$ and as

$\mathbb{E}\left[B_{\mathcal{A}}^{2}\left(t\right)\right]\sim\frac{\Gamma\left(\gamma+\beta\right)t^{A_{2}}}{\Gamma\left(\gamma\right)\left(\left(A_{2}-A_{1}\right)\log\left(t\right)\right)^{\beta}}$
(59)

for $t\gg 1$.

The expectation of TAMSD for SBMRE with beta distributed anomalous diffusion
exponent is as follows:

$\mathbb{E}\left[\delta_{\mathcal{A}}\left(\tau\right)\right]=\frac{1}{T-\tau}\int_{0}^{T-\tau}\left(e^{A_{1}\log\left(t+\tau\right)}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t+\tau\right)\right)}-e^{A_{1}\log
t}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log
t\right)}\right)dt.$ (60)

The EB parameter in this case takes the form:

$EB_{\mathcal{A}}\left(\tau\right)=\frac{\mathcal{N}_{\mathcal{A}}\left(\tau\right)}{\mathcal{D}_{\mathcal{A}}\left(\tau\right)},$
(61)

where

$\displaystyle\mathcal{N}_{\mathcal{A}}\left(\tau\right)=\frac{4}{\left(T-\tau\right)^{2}}\int_{0}^{\tau}\int_{0}^{T-\tau-
t_{2}}\left(\left(t_{1}+\tau\right)^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t_{1}+\tau\right)\right)}-\right.$
(62)
$\displaystyle\left.-\left(t_{1}+t_{2}\right)^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t_{1}+t_{2}\right)\right)}\right)^{2}dt_{1}dt_{2}$

and

$\displaystyle\mathcal{D}_{\mathcal{A}}\left(\tau\right)=\frac{1}{\left(T-\tau\right)^{2}}\int_{0}^{T-\tau}\int_{0}^{T-\tau}\left(\left(t_{1}+\tau\right)^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t_{1}+\tau\right)\right)}-t_{1}^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t_{1}\right)\right)}\right)$
(63)
$\displaystyle\left(\left(t_{2}+\tau\right)^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t_{2}+\tau\right)\right)}-t_{2}^{A_{1}}\prescript{}{1}{F_{1}\left(\gamma,\gamma+\beta,\left(A_{2}-A_{1}\right)\log\left(t_{2}\right)\right)}\right)dt_{1}dt_{2}.$

The expression for the first hitting time in barrier $b$ can be represented by
the following formula:

$f_{\tau_{b}}\left(t\right)=\frac{b}{\sqrt{2\pi}\left(A_{2}-A_{1}\right)^{\gamma+\beta-1}\mathbb{B}\left(\gamma,\beta\right)}\int_{A_{1}}^{A_{2}}a\left(a-A_{1}\right)^{\gamma-1}\left(A_{2}-a\right)^{\beta-1}\exp\left\\{-\frac{b^{2}}{2t^{a}}\right\\}t^{-1-\frac{a}{2}}da.$
(64)

To obtain the asymptotics of (64) for $t\gg 1$ let us note that
$\exp\left\\{-\frac{b^{2}}{2t^{a}}\right\\}\in\mathcal{O}\left(1\right)$ as a
function of $t$, $a$ is bounded by $A_{2}$, and
$t^{-1-\frac{a}{2}}\in\mathcal{O}\left(t^{-1-\frac{A_{1}}{2}}\right)$ for all
$a$. Then we calculate the integral
$\int_{A_{1}}^{A_{2}}\left(a-A_{1}\right)^{\gamma-1}\left(A_{2}-a\right)^{\beta-1}da$,
which is equal to a constant independent on $t$. Finally, for (64) we have the
following asymptotic formula:

$f_{\tau_{b}}\left(t\right)\sim t^{-1-A_{1}/2}.$ (65)

##  6 Numerical analysis

In this section, we present the results of numerical analysis demonstrating
main properties of the SBMRE (discussed theoretically in the previous section)
for two-point distribution and beta distribution of the anomalous diffusion
exponent. In addition, we compare the main characteristics of SBMRE with the
corresponding properties of SBM.

For all cases presented, we assume $A_{1}=0.5$ and $A_{2}=1.5$. These
parameters are selected to cover both subdiffusive and superdiffusive cases.
Additionally, for comparison purposes, we select the same values of $\alpha$
in the SBM case, namely $\alpha=0.5$ and $\alpha=1.5$. For the two-point
distribution of the anomalous diffusive exponent, we consider three $p$
values, namely $p=0.1$, $p=0.5$, and $p=0.9$. In such cases, we observe the
dominance of $A_{2}$ (for $p=0.1$), the balance between $A_{1}$ and $A_{2}$
(for $p=0.5$), and the dominance of $A_{1}$ (for $p=0.9$), respectively. In
order to enforce similar behavior for the beta distribution of the anomalous
diffusive exponent, we consider three scenarios:
$\left(\gamma,\beta\right)=\left(0.7,0.3\right)$,
$\left(\gamma,\beta\right)=\left(0.5,0.5\right)$, and
$\left(\gamma,\beta\right)=\left(0.3,0.7\right)$.

In Fig. 1 we present the exemplary trajectories for SBM for two considered
$\alpha$ values, while in Fig. 2 we demonstrate the trajectories of SBMRE for
two-point and beta distributed anomalous diffusive exponent. It is evident
that basing solely on these plots makes it difficult to distinguish between
the specific processes, thus prompting a deeper exploration of their
properties.

Figure 1:  Example trajectories of SBM with $\alpha=0.5$ and $\alpha=1.5$.
Figure 2:  Example trajectories of SBMRE with a two-point and beta
distribution of the anomalous diffusion exponent with different values of $p$,
$\gamma$ and $\beta$; $A_{1}=0.5$, $A_{2}=1.5$.

In Fig. 3, top panels, we present the pdfs depicting the position of SBMRE
with the two-point distribution of the anomalous diffusion exponent for
$t=0.5$, alongside with the pdfs representing the position of subdiffusive and
superdiffusive SBM. It is observed that for $p=0.1$ the pdf of SBMRE aligns
more closely with that of superdiffusive SBM, and as $p$ increases, it tends
towards the subdiffusive counterpart, consistent with intuition. In Fig. 3,
bottom panels, corresponding pdfs for $t=10$ are presented, where a transition
from subdiffusive to superdiffusive cases is apparent. In Fig. 7 we
demonstrate the similar plots for SBMRE with beta distributed anomalous
diffusion exponent.

Figure 3: Comparison of pdfs of the position of SBMRE with two-point
distribution of the anomalous diffusion exponent and SBM. Top panels: $t=0.5$.
Bottom panels: $t=10$.

In Fig. 4, top panels, we demonstrate short-time behavior of the second moment
of SBMRE with a two-point distribution of the anomalous diffusion exponent
exhibiting characteristics akin to both subdiffusive and superdiffusive SBM
(depending on $p$ parameter). Notably, we discern that at $p=0.1$, the SBMRE
aligns more closely with the superdiffusive scenario of SBM. As the parameter
$p$ increases, there is a noticeable shift toward proximity with the
subdiffusive scenario. In Fig. 4, bottom panels, a similar trend is observed
for long time periods, but notably, the transition is from subdiffusive to
superdiffusive SBM cases. The second moment for SBMRE with the beta-
distributed exponent and SBM is presented in Fig. 8. One can observe similar
behavior as in case of the two-point distribution.

In Fig. 9 (top panels) we depict the asymptotic formula (47) of the second
moment for SBMRE for two-point distribution case for short times together with
the exact formula presented in (46). Additionally, in Fig. 9 (bottom panels)
we depict the asymptotic formula (47) of the second moment for SBMRE for long
times together with the exact formula presented in (46). Similar comparison
for beta distribution of the anomalous diffusion exponent is presented in Fig.
10. As can be observed, in all cases the theoretical asymptotics perfectly
agree with the numerical analysis.

The expected value of TAMSD for two-point distribution of the anomalous
diffusion exponent is presented in Fig. 5. Moreover, in Fig. 11, we depict the
asymptotic formula (50) of the expected value of the TAMSD for SBMRE with a
two-point distribution of the anomalous diffusion exponent together with the
exact value presented in (48). Here we can see the perfect agreement between
theoretical formula and values received by numerical analysis. The expected
value of TAMSD for the beta distributed anomalous diffusion exponent is
visible in Fig. 12. The behavior of this statistic for beta distribution is
similar to the one observed for two-point distribution for appropriate cases
of parameter values.

Figure 4: Comparison of the second moments of SBMRE with two-point
distribution of the anomalous diffusion exponent and SBM. Top panels: short
times Bottom panels: long times.
Figure 5: Comparison of the expected values of TAMSD of SBMRE with two-point
distribution of the anomalous diffusion exponent and SBM. Here $T=10$.

Finally, in Fig. 6 we plot pdfs of the first hitting time to the barrier $b=1$
for SBMRE with the two-point distribution and compare the pdf shapes with
those for both subdiffusive and superdiffusive SBM. Domination of the lower
exponent at long times is clearly visible. For beta distributed anomalous
diffusion exponent the same comparison is presented in Fig. 13.

Figure 6: Comparison of pdf of the hitting time in barrier $b=1$ for SBMRE
with two-point distribution of the anomalous diffusion exponent and SBM.

##  7 Conclusions

In this paper, we investigate the probabilistic properties of SBMRE, scaled
Brownian motion with random anomalous diffusion exponent. This process
maintains the characteristics of SBM at the individual trajectory level,
albeit with randomly varying exponent across the trajectories. SBM is
recognized as one of the generic random processes with an anomalous diffusion
behavior, which finds various applications. The motivation for our study comes
from modern biophysical experiments, which clearly demonstrate that standard
anomalous diffusion processes with fixed parameters are insufficient to
describe the observed phenomena. Consequently, a natural modification of SBM
is to incorporate random parameters into classical anomalous diffusion
systems, a scenario previously discussed in the context of fractional Brownian
motion leading to FBM with a random Hurst exponent [55], and very recently in
doubly stochastic version of continuous time random walk [72] and Lévy-walk-
like Langevin dynamics with random parameters [73].

Here, we expand upon the methodology of anomalous diffusion processes with
random parameters by analyzing SBMRE. We established the mathematical
framework of this process by examining its pdf, and the q-th absolute moment
for the general distribution of the anomalous diffusion exponent.
Additionally, we presented the expected value of the TAMSD of SBMRE and
discussed its ergodicity breaking property. Furthermore, we explored the
properties of the first hitting time pdf and we proven that SBMRE is a
martingale, similarly as its stochastic exponential. As special cases of the
distributions of the anomalous diffusive exponent, we engage the two-point and
beta distributions. In such instances, we discuss the asymptotics of the
characteristics of SBMRE. Theoretical findings are validated through numerical
analysis, and compared with similar properties of SBM. Similarly to FBMRE, the
diffusion of SBMRE is determined by smaller exponent at short times while for
sufficiently large times the larger exponent is dominant, the effect called
accelerating diffusion. As for the expectation value of TAMSD routinely
measured in experiments, the dependence on the trajectory length is determined
by larger exponent if the trajectory is long enough. We also observed that the
long-time asymptotics of the first hitting time pdf of SBMRE is determined by
smaller exponent. These characteristics are essential for comprehensive
understanding the process considered and discerning its differences from SBM.
Moreover, they are significant for practical applications.

It is known that the characterization of the anomalous diffusion model from
the measurement of individual trajectory is a challenging task that requires
new approaches based, in particular, on machine learning or Bayesian inference
[74]. Recently, the latter one was implemented for SBM in [75] and combined
with the implementation of Bayesian inference for FBM presented in [76] to
demonstrate model selection and inference of the anomalous diffusion exponent
for these two models. Also, the statistical-based approach for identification
of such models was proposed in [44]. As a future study we plan to extend these
methods for anomalous diffusion models with random parameters, such as SBMRE
and FBMRE.

## Acknowledgments

The work of AW was supported by National Center of Science under Opus Grant
2020/37/B/HS4/00120 “Market risk model identification and validation using
novel statistical, probabilistic, and machine learning tools”.

## References

  * [1] J. P. Bouchaud and A. Georges.  Anomalous diffusion in disordered media: Statistical mechanisms, models and physical applications.  Physics Reports, 195(4–5):127–293, November 1990.
  * [2] R. M. Mazo.  Brownian Motion.  Oxford University Press, October 2008.
  * [3] R. L. Schilling and L. Partzsch.  Brownian Motion: An Introduction to Stochastic Processes.  DE GRUYTER, May 2012.
  * [4] A. N. Kolmogorov.  Wienersche Spiralen und einige andere interessante Kurven in Hilbertscen Raum.  C.R. (Doklady) Acad. Sci. URSS (NS), 26:115, 1940.
  * [5] B. B. Mandelbrot and J. W. Van Ness.  Fractional Brownian motions, fractional noises and applications.  SIAM review, 10(4):422–437, 1968.
  * [6] J. Beran.  Statistics for long-memory processes.  Chapman & Chal, 1994.
  * [7] F. J. Molz, H. H. Liu, and J. Szulga.  Fractional Brownian motion and fractional Gaussian noise in subsurface hydrology: A review, presentation of fundamental properties, and extensions.  Water Resources Research, 33(10):2273–2286, 1997.
  * [8] D. A. Benson, M. M. Meerschaert, and J. Revielle.  Fractional calculus in hydrologic modeling: A numerical perspective.  Advances in Water Resources, 51:479–497, 2013.
  * [9] A. Pashko.  Simulation of telecommunication traffic using statistical models of fractional Brownian motion.  In 2017 4th International Scientific-Practical Conference Problems of Infocommunications. Science and Technology (PIC S&T), pages 414–418, 2017.
  * [10] J. Frecon, G. Didier, N. Pustelnik, and P. Abry.  Non-Linear Wavelet Regression and Branch & Bound Optimization for the Full Identification of Bivariate Operator Fractional Brownian Motion.  IEEE Transactions on Signal Processing, 64(15):4040–4049, 2016\.
  * [11] Y. Chang and S. Chang.  A fast estimation algorithm on the Hurst parameter of discrete-time fractional Brownian motion.  IEEE Transactions on Signal Processing, 50(3):554–559, 2002.
  * [12] S. Liu and S. Chang.  Dimension estimation of discrete-time fractional Brownian motion with applications to image texture classification.  IEEE Transactions on Image Processing, 6(8):1176–1184, 1997.
  * [13] Z. Gelbaum and M. Titus.  Simulation of Fractional Brownian Surfaces via Spectral Synthesis on Manifolds.  IEEE Transactions on Image Processing, 23(10):4383–4388, 2014.
  * [14] S. Rostek and R. Schöbel.  A note on the use of fractional Brownian motion for financial modeling.  Economic Modelling, 30:30–35, 2013.
  * [15] K. Maraj, D. Szarek, G. Sikora, and A. Wyłomańska.  Time-averaged mean squared displacement ratio test for Gaussian processes with unknown diffusion coefficient.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 31(7):073120, 2021.
  * [16] W. Xiao, W. Zhang, X. Zhang, and Y. Wang.  Pricing currency options in a fractional Brownian motion with jumps.  Economic Modelling, 27(5):935–942, 2010.
  * [17] D. Ernst, M. Hellmann, J. Köhler, and M. Weiss.  Fractional Brownian motion in crowded fluids.  Soft Matter, 8:4886, 2012.
  * [18] F. Höfling and T. Franosch.  Anomalous transport in the crowded world of biological cells.  Rep. Prog. Phys., 76:046602, 2013.
  * [19] J. Metzler, R.and Jeon, A. G. Cherstvy, and E. Barkai.  Anomalous diffusion models and their properties: non-stationarity, non-ergodicity, and ageing at the centenary of single particle tracking.  Phys. Chem. Chem. Phys., 16(44):24128–24164, 2014.
  * [20] D. Szarek, K. Maraj-Zygmat, G. Sikora, D. Krapf, and A. Wylomanska.  Statistical test for anomalous diffusion based on empirical anomaly measure for Gaussian processes.  Computational Statistics & Data Analysis, 168:107401, 2022.
  * [21] D. Krapf, N. Lukat, E. Marinari, R. Metzler, G. Oshanin, C. Selhuber-Unkel, A. Squarcini, L. Stadler, M. Weiss, and X. Xu.  Spectral Content of a Single Non-Brownian Trajectory.  Phys. Rev. X, 9:011019, 2019.
  * [22] G. K. Batchelor.  Diffusion in a field of homogeneous turbulence: II. The relative motion of particles.  Mathematical Proceedings of the Cambridge Philosophical Society, 48(2):345–362, 1952.
  * [23] S. C. Lim and S. V. Muniandy.  Self-similar Gaussian processes for modeling anomalous diffusion.  Physical Review E, 66(2):021114, 2002.
  * [24] K. Debicki, Z. Michna, and T. Rolski.  On the supremum from Gaussian processes over infinite horizon.  Probability and Mathematical Statistics, 18:83, 1998.
  * [25] K. Debicki.  Asymptotics of the supremum of scaled Brownian motion.  Probability and Mathematical Statistics, 21, 2001.
  * [26] M. Magdziarz.  Lamperti transformation of scaled Brownian motion and related Langevin equations.  Communications in Nonlinear Science and Numerical Simulation, 83:105077, 2020.
  * [27] F. Thiel and I. M. Sokolov.  Scaled Brownian motion as a mean-field model for continuous-time random walks.  Physical Review E, 89(1):012115, 2014.
  * [28] J. Jeon, A. V. Chechkin, and R. Metzler.  Scaled Brownian motion: a paradoxical process with a time dependent diffusivity for the description of anomalous diffusion.  Phys. Chem. Chem. Phys., 16(30):15811–15817, 2014.
  * [29] H. Safdari, A. G. Cherstvy, A. V. Chechkin, F. Thiel, I. M. Sokolov, and R. Metzler.  Quantifying the non-ergodicity of scaled Brownian motion.  Journal of Physics A: Mathematical and Theoretical, 48(37):375002, 2015.
  * [30] V. Sposini, R. Metzler, and G. Oshanin.  Single-trajectory spectral analysis of scaled Brownian motion.  New Journal of Physics, 21(7):073043, 2019.
  * [31] A. Bodrova, A. V. Chechkin, A. G. Cherstvy, and R. Metzler.  Quantifying non-ergodic dynamics of force-free granular gases.  Physical Chemistry Chemical Physics, 17(34):21791–21798, 2015.
  * [32] A. S. Bodrova, A. V. Chechkin, A. G. Cherstvy, H. Safdari, I. M. Sokolov, and R. Metzler.  Underdamped scaled Brownian motion: (non-)existence of the overdamped limit in anomalous diffusion.  Scientific Reports, 6(1), 2016.
  * [33] A. S. Bodrova, A. V. Chechkin, and I. M. Sokolov.  Scaled Brownian motion with renewal resetting.  Physical Review E, 100(1):012120, 2019.
  * [34] A. S. Bodrova, A. V. Chechkin, and I. M. Sokolov.  Nonrenewal resetting of scaled Brownian motion.  Physical Review E, 100(1):012119, 2019.
  * [35] S. Lee, H.and Song, J. Kim, and J. Sung.  Survival Probability Dynamics of Scaled Brownian Motion: Effect of Nonstationary Property.  Bulletin of the Korean Chemical Society, 40(9):847–850, 2019.
  * [36] Y. Li, K. Suleiman, and Y. Xu.  Anomalous diffusion, non-Gaussianity, nonergodicity, and confinement in stochastic-scaled Brownian motion with diffusing diffusivity dynamics.  Physical Review E, 109(1):014139, 2024.
  * [37] A. G. Cherstvy, D. Vinod, E. Aghion, I. M. Sokolov, and R. Metzler.  Scaled geometric Brownian motion features sub- or superexponential ensemble-averaged, but linear time-averaged mean-squared displacements.  Physical Review E, 103(6):062127, June 2021.
  * [38] A. Valdés Gómez and F. J. Sevilla.  Fractional and scaled Brownian motion on the sphere: The effects of long-time correlations on navigation strategies.  Physical Review E, 108(5):054117, 2023.
  * [39] A. Molini, P. Talkner, G. G. Katul, and A. Porporato.  First passage time statistics of Brownian motion with purely time dependent drift and diffusion.  Physica A: Statistical Mechanics and its Applications, 390(11):1841–1852, 2011.
  * [40] K. E. Bassler, J. L. McCauley, and G. H. Gunaratne.  Nonstationary increments, scaling distributions, and variable diffusion processes in financial markets.  Proceedings of the National Academy of Sciences, 104(44):17287–17290, 2007.
  * [41] D. S. Novikov, J. H. Jensen, J. A. Helpern, and E. Fieremans.  Revealing mesoscopic structural universality with diffusion.  Proceedings of the National Academy of Sciences, 111(14):5088–5093, 2014.
  * [42] S. Thapa, S. Park, Y. Kim, J. Jeon, R. Metzler, and M. A. Lomholt.  Bayesian inference of scaled versus fractional Brownian motion.  Journal of Physics A: Mathematical and Theoretical, 55(19):194003, 2022.
  * [43] A. Grzesiek, J. Gajda, A. Wyłomańska, and S. Sundar.  Discriminating between scaled and fractional Brownian motion via p-variation statistics.  International Journal of Advances in Engineering Sciences and Applied Mathematics, 10(1):9–14, 2018.
  * [44] M. Balcerek, K. Burnecki, G. Sikora, and A. Wyłomańska.  Discriminating Gaussian processes via quadratic form statistics.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 31(6), 2021.
  * [45] A. A. Sadoon and Y. Wang.  Anomalous, non-Gaussian, viscoelastic, and age-dependent dynamics of histonelike nucleoid-structuring proteins in live escherichia coli.  Physical Review E, 98:042411, 2018.
  * [46] R. Benelli and M. Weiss.  From sub-to superdiffusion: Fractional Brownian motion of membraneless organelles in early C. elegans embryos.  New Journal of Physics, 23:063072, 2021.
  * [47] K. Speckner and M. Weiss.  Single-particle tracking reveals anti-persistent subdiffusion in cell extracts.  Entropy, 23(7):892, 2021.
  * [48] A. G. Cherstvy, S. Thapa, C. E. Wagner, and R. Metzler.  Non-Gaussian, non-ergodic, and non-Fickian diffusion of tracers in mucin hydrogels.  Soft Matter, 15(12):2526–2551, 2019.
  * [49] C. Beck and E. G. Cohen.  Superstatistics.  Physica A: Statistical mechanics and its applications, 322:267–275, 2003.
  * [50] C. Beck.  Superstatistical Brownian motion.  Progress of Theoretical Physics Supplement, 162:29–36, 2005.
  * [51] M. V. Chubynsky and G. W. Slater.  Diffusing Diffusivity: A Model for Anomalous, yet Brownian, Diffusion.  Physical Review Letters, 113(9):098302, August 2014.
  * [52] R. Metzler and A. V. Chechkin.  Non-Gaussian stochastic transport.  In B. J. West and S. Bianco, editors, Paolo Grigolini and 50 Years of Statistical Physics., pages 71–93. Cambridge Scholars Publishing, 2023\.
  * [53] J. Lévy-Véhel and R. F. Peltier.  Multifractional Brownian motion: definition and preliminary results.  Rapport de recherche de l’INRIA, 2645, 1995.
  * [54] A. Ayache and M. S. Taqqu.  Multifractional processes with random exponent.  Publicacions Matemàtiques, 49(2):459, 2005.
  * [55] M. Balcerek, K. Burnecki, S. Thapa, A. Wyłomańska, and A. V. Chechkin.  Fractional Brownian motion with random Hurst exponent: Accelerating diffusion and persistence transitions.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 32(9):093114, 2022.
  * [56] D. Han, N. Korabel, R. Chen, M. Johnston, A. Gavrilova, V. J. Allan, S. Fedotov, and T. A. Waigh.  Deciphering anomalous heterogeneous intracellular transport with neural networks.  ELife, 9:e52224, 2020.
  * [57] N. Korabel, D. Han, A. Taloni, G. Pagnini, S. Fedotov, V. Allan, and T. A. Waigh.  Local analysis of heterogeneous intracellular transport: Slow and fast moving endosomes.  Entropy, 23(8):958, 2021.
  * [58] M. Balcerek, A. Wyłomańska, K. Burnecki, R. Metzler, and D. Krapf.  Modelling intermittent anomalous diffusion with switching fractional Brownian motion.  New Journal of Physics, 25(10):103031, 2023.
  * [59] W. Wang, M. Balcerek, K. Burnecki, A. V. Chechkin, S. Janusonis, J. Slezak, T. Vojta, A. Wylomanska, and R. Metzler.  Memory-multi-fractional Brownian motion with continuous correlations.  Phys. Rev. Res., 5:L032025, 2023.
  * [60] M.A.F. dos Santos, L. Menon, and D. Cius.  Superstatistical approach of the anomalous exponent for scaled Brownian motion.  Chaos, Solitons & Fractals, 164:112740, 2022.
  * [61] S. M. Rytov, Y. A. Kravtsov, and V. I. Tatarskii.  Principles of statistical radiophysics 1.  Springer-Verlag New York Inc, United States, 1987.
  * [62] Y. He, S. Burov, R. Metzler, and E. Barkai.  Random Time-Scale Invariant Diffusion and Transport Coefficients.  Physical Review Letters, 101(5):058101, 2008.
  * [63] B. Øksendal.  Stochastic Differential Equations.  Springer Berlin Heidelberg, 2003.
  * [64] N.G. de Bruijn.  Asymptotic Methods in Analysis.  Bibliotheca mathematica. Dover Publications, 1981.
  * [65] J. Penot.  Analysis.  Springer International Publishing, 2016.
  * [66] M. Magdziarz.  Path properties of subdiffusion—a martingale approach.  Stochastic Models, 26(2):256–271, 2010.
  * [67] D. Revuz and M. Yor.  Continuous Martingales and Brownian Motion.  Springer Berlin Heidelberg, 1999.
  * [68] A. V. Chechkin, R. Gorenflo, and I. M. Sokolov.  Retarding subdiffusion and accelerating superdiffusion governed by distributed-order fractional diffusion equations.  Physical Review E, 66(4):046129, October 2002.
  * [69] A. V. Chechkin, V. Yu. Gonchar, R. Gorenflo, N. Korabel, and I. M. Sokolov.  Generalized fractional diffusion equations for accelerating subdiffusion and truncated lévy flights.  Physical Review E, 78(2):021111, August 2008.
  * [70] A. Chechkin, I. M. Sokolov, and J. Klafter.  Natural and Modified Forms of Distributed-Order Fractional Diffusion Equations, pages 107–127.  WORLD SCIENTIFIC, October 2011.
  * [71] R. Beals and R. Wong.  Special Functions: A Graduate Text.  Cambridge University Press, 2010.
  * [72] M. Arutkin and S. Reuveni.  Doubly stochastic continuous time random walk.  Phys. Rev. Res., 6:L012033, 2024.
  * [73] Y. Chen, X. Wang, and M. Ge.  Lévy-walk-like Langevin dynamics with random parameters.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 34(1):013109, 01 2024.
  * [74] G. Muñoz-Gil, G. Volpe, M.A. Garcia-March, and et al.  Objective comparison of methods to decode anomalous diffusion.  Nature Communications, 12:6253, 2021.
  * [75] S. Thapa, B. Pinchasik, and Y. Shokef.  Emergent clustering due to quorum sensing interactions in active matter.  New Journal of Physics, 26(2):023010.
  * [76] J. Krog, L. H. Jacobsen, F. W. Lund, D. Wüstner, and M. A. Lomholt.  Bayesian model selection with fractional brownian motion.  Journal of Statistical Mechanics: Theory and Experiment, 2018(9):093501, September 2018.

##  Appendix A Additional plots.

Figure 7: Comparison of pdfs of the position of SBMRE with beta distribution
of the anomalous diffusion exponent and SBM. Top panels: $t=0.5$. Bottom
panels: $t=10$.

Figure 8:  Comparison of second moments of SBMRE with beta distribution of the
anomalous diffusion exponent and SBM. Top panels: short times. Bottom panels:
long times.

Figure 9:  Comparison of second moments of SBMRE with two-point distribution
of the anomalous diffusion exponent and their asymptotics (see formula (47)).
Top panels: short times. Bottom panels: long times.

Figure 10:  Comparison of second moments of SBMRE with beta distribution of
the anomalous diffusion exponent and their asymptotics (see formulas (58),
(59)). Top panels: short times. Bottom panels: long times.  Figure 11:
Asymptotic (see formula (50)) of the expectation of the TAMSD of SBMRE with a
two-point distribution of the anomalous diffusion exponent. Here $T=10$.
Figure 12: Comparison of the expected values of TAMSD of SBMRE with beta
distribution of the anomalous diffusion exponent and SBM. Here $T=10$.  Figure
13: Comparison of pdf of the hitting time in barrier $b=1$ for SBMRE with beta
distribution of the anomalous diffusion exponent and SBM.


11affiliationtext: Imaging Biomarkers and Computer-Aided Diagnosis Laboratory,
Radiology and Imaging Sciences, National Institutes of Health Clinical Center,
Bethesda, MD 20892-1182, USA.

# DeepOrgan: Multi-level Deep Convolutional Networks for Automated Pancreas
Segmentation

Holger R. Roth<EMAIL_ADDRESS><EMAIL_ADDRESS>Le Lu Amal Farag Hoo-
Chang Shin Jiamin Liu Evrim Turkbey Ronald M. Summers

###### Abstract

Automatic organ segmentation is an important yet challenging problem for
medical image analysis. The pancreas is an abdominal organ with very high
anatomical variability. This inhibits previous segmentation methods from
achieving high accuracies, especially compared to other organs such as the
liver, heart or kidneys. In this paper, we present a probabilistic bottom-up
approach for pancreas segmentation in abdominal computed tomography (CT)
scans, using multi-level deep convolutional networks (ConvNets). We propose
and evaluate several variations of deep ConvNets in the context of
hierarchical, coarse-to-fine classification on image patches and regions, i.e.
superpixels. We first present a dense labeling of local image patches via
P-ConvNet and nearest neighbor fusion. Then we describe a regional ConvNet
($R_{1}{-}\mathrm{ConvNet}$) that samples a set of bounding boxes around each
image superpixel at different scales of contexts in a “zoom-out” fashion. Our
ConvNets learn to assign class probabilities for each superpixel region of
being pancreas. Last, we study a stacked $R_{2}{-}\mathrm{ConvNet}$ leveraging
the joint space of CT intensities and the $P{-}\mathrm{ConvNet}$ dense
probability maps. Both 3D Gaussian smoothing and 2D conditional random fields
are exploited as structured predictions for post-processing. We evaluate on CT
images of 82 patients in 4-fold cross-validation. We achieve a Dice Similarity
Coefficient of 83.6$\pm$6.3% in training and 71.8$\pm$10.7% in testing.

##  1 Introduction

Segmentation of the pancreas can be a prerequisite for computer aided
diagnosis (CADx) systems that provide quantitative organ volume analysis, e.g.
for diabetic patients. Accurate segmentation could also necessary for computer
aided detection (CADe) methods to detect pancreatic cancer. Automatic
segmentation of numerous organs in computed tomography (CT) scans is well
studied with good performance for organs such as liver, heart or kidneys,
where Dice Similarity Coefficients (DSC) of $>$90% are typically achieved Wang
et al. (2014), Chu et al. (2013), Wolz et al. (2013), Ling et al. (2008).
However, achieving high accuracies in automatic pancreas segmentation is still
a challenging task. The pancreas’ shape, size and location in the abdomen can
vary drastically between patients. Visceral fat around the pancreas can cause
large variations in contrast along its boundaries in CT (see Fig. 3). Previous
methods report only 46.6% to 69.1% DSCs Wang et al. (2014), Chu et al. (2013),
Wolz et al. (2013), Farag et al. (2014). Recently, the availability of large
annotated datasets and the accessibility of affordable parallel computing
resources via GPUs have made it feasible to train deep convolutional networks
(ConvNets) for image classification. Great advances in natural image
classification have been achieved Krizhevsky et al. (2012). However, deep
ConvNets for semantic image segmentation have not been well studied Mostajabi
et al. (2014). Studies that applied ConvNets to medical imaging applications
also show good promise on detection tasks Cireşan et al. (2013), Roth et al.
(2014). In this paper, we extend and exploit ConvNets for a challenging organ
segmentation problem.

##  2 Methods

We present a coarse-to-fine classification scheme with progressive pruning for
pancreas segmentation. Compared with previous top-down multi-atlas
registration and label fusion methods, our models approach the problem in a
bottom-up fashion: from dense labeling of image patches, to regions, and the
entire organ. Given an input abdomen CT, an initial set of superpixel regions
is generated by a coarse cascade process of random forests based pancreas
segmentation as proposed by Farag et al. (2014). These pre-segmented
superpixels serve as regional candidates with high sensitivity ($>$97%) but
low precision. The resulting initial DSC is $\sim$27% on average. Next, we
propose and evaluate several variations of ConvNets for segmentation
refinement (or pruning). A dense local image patch labeling using an axial-
coronal-sagittal viewed patch ($P{-}\mathrm{ConvNet}$) is employed in a
sliding window manner. This generates a per-location probability response map
$P$. A regional ConvNet ($R_{1}{-}\mathrm{ConvNet}$) samples a set of bounding
boxes covering each image superpixel at multiple spatial scales in a “zoom-
out” fashion Mostajabi et al. (2014), Girshick et al. (2014) and assigns
probabilities of being pancreatic tissue. This means that we not only look at
the close-up view of superpixels, but gradually add more contexts to each
candidate region. $R_{1}$-ConvNet operates directly on the CT intensity.
Finally, a stacked regional $R_{2}{-}\mathrm{ConvNet}$ is learned to leverage
the joint convolutional features of CT intensities and probability maps P.
Both 3D Gaussian smoothing and 2D conditional random fields for structured
prediction are exploited as post-processing. Our methods are evaluated on CT
scans of 82 patients in 4-fold cross-validation (rather than “leave-one-out”
evaluation Wang et al. (2014), Chu et al. (2013), Wolz et al. (2013)). We
propose several new ConvNet models and advance the current state-of-the-art
performance to a DSC of 71.8 in testing. To the best of our knowledge, this is
the highest DSC reported in the literature to date.

###  2.1 Candidate region generation

We describe a coarse-to-fine pancreas segmentation method employing multi-
level deep ConvNet models. Our hierarchical segmentation method decomposes any
input CT into a set of local image superpixels
$S=\left\\{S_{1},\ldots,S_{N}\right\\}$. After evaluation of several image
region generation methods Achanta et al. (2012), we chose entropy rate Liu et
al. (2011) to extract $N$ superpixels on axial slices. This process is based
on the criterion of DSCs given optimal superpixel labels, in part inspired by
the PASCAL semantic segmentation challenge Everingham et al. (2014). The
optimal superpixel labels achieve a DSC upper-bound and are used for
supervised learning below. Next, we use a two-level cascade of random forest
(RF) classifiers as in Farag et al. (2014). We only operate the RF labeling at
a low class-probability cut $>$0.5 which is sufficient to reject the vast
amount of non-pancreas superpixels. This retains a set of superpixels
$\\{S_{\mathrm{RF}}\\}$ with high recall ($>$97%) but low precision. After
initial candidate generation, over-segmentation is expected and observed with
low DSCs of $\sim$27%. The optimal superpixel labeling is limited by the
ability of superpixels to capture the true pancreas boundaries at the per-
pixel level with $DSC_{\mathrm{max}}=80.5\%$, but is still much above previous
state-of-the-art Wang et al. (2014), Chu et al. (2013), Wolz et al. (2013),
Farag et al. (2014). These superpixel labels are used for assessing ‘positive’
and ‘negative’ superpixel examples for training. Assigning image regions
drastically reduces the amount of ConvNet observations needed per CT volume
compared to a purely patch-based approach and leads to more balanced training
data sets. Our multi-level deep ConvNets will effectively prune the coarse
pancreas over-segmentation to increase the final DSC measurements.

###  2.2 Convolutional neural network (ConvNet) setup

We use ConvNets with an architecture for binary image classification. Five
layers of convolutional filters compute and aggregate image features. Other
layers of the ConvNets perform max-pooling operations or consist of fully-
connected neural networks. Our ConvNet ends with a final two-way layer with
softmax probability for ‘pancreas’ and ‘non-pancreas’ classification (see Fig.
1). The fully-connected layers are constrained using “DropOut” in order to
avoid over-fitting by acting as a regularizer in training Srivastava et al.
(2014). GPU acceleration allows efficient training (we use cuda-
convnet2111https://code.google.com/p/cuda-convnet2).

Figure 1: The proposed ConvNet architecture. The number of convolutional
filters and neural network connections for each layer are as shown. This
architecture is constant for all ConvNet variations presented in this paper
(apart from the number of input channels): $P{-}\mathrm{ConvNet}$,
$R_{1}{-}\mathrm{ConvNet}$, and $R_{2}{-}\mathrm{ConvNet}$.

Figure 2: The first layer of learned convolutional kernels using three
representations: a) 2.5D sliding-window patches ($P{-}\mathrm{ConvNet}$), b)
CT intensity superpixel regions ($R_{1}{-}\mathrm{ConvNet}$), and c) CT
intensity + $P_{0}$ map over superpixel regions ($R_{2}{-}\mathrm{ConvNet}$).

###  2.3 $P{-}\mathrm{ConvNet}$: Deep patch classification

We use a sliding window approach that extracts 2.5D image patches composed of
axial, coronal and sagittal planes within all voxels of the initial set of
superpixel regions $\\{S_{\mathrm{RF}}\\}$ (see Fig. 3). The resulting ConvNet
probabilities are denoted as $P_{0}$ hereafter. For efficiency reasons, we
extract patches every $n$ voxels and then apply nearest neighbor
interpolation. This seems sufficient due to the already high quality of
$P_{0}$ and the use of overlapping patches to estimate the values at skipped
voxels.

Figure 3: Axial CT slice of a manual (gold standard) segmentation of the
pancreas. From left to right, there are the ground-truth segmentation contours
(in red), RF based coarse segmentation $\\{S_{\mathrm{RF}}\\}$, and the deep
patch labeling result using $P{-}\mathrm{ConvNet}$.

###  2.4 $R{-}\mathrm{ConvNet}$: Deep region classification

We employ the region candidates as inputs. Each superpixel
$\in\\{S_{\mathrm{RF}}\\}$ will be observed at several scales $N_{s}$ with an
increasing amount of surrounding contexts (see Fig. 4). Multi-scale contexts
are important to disambiguate the complex anatomy in the abdomen. We explore
two approaches: $R_{1}{-}\mathrm{ConvNet}$ only looks at the CT intensity
images extracted from multi-scale superpixel regions, and a stacked
$R_{2}{-}\mathrm{ConvNet}$ integrates an additional channel of patch-level
response maps $P_{0}$ for each region as input. As a superpixel can have
irregular shapes, we warp each region into a regular square (similar to RCNN
Girshick et al. (2014)) as is required by most ConvNet implementations to
date. The ConvNets automatically train their convolutional filter kernels from
the available training data. Examples of trained first-layer convolutional
filters for $P{-}\mathrm{ConvNet}$, $R_{1}{-}\mathrm{ConvNet}$,
$R_{2}{-}\mathrm{ConvNet}$ are shown in Fig. 2. Deep ConvNets behave as
effective image feature extractors that summarize multi-scale image regions
for classification.

Figure 4: Region classification using $R{-}\mathrm{ConvNet}$ at different
scales: a) one-channel input based on the intensity image only, and b) two-
channel input with additional patch-based $P{-}\mathrm{ConvNet}$ response.

###  2.5 Data augmentation

Our ConvNet models ($R_{1}{-}\mathrm{ConvNet}$, $R_{2}{-}\mathrm{ConvNet}$)
sample the bounding boxes of each superpixel $\in\\{S_{\mathrm{RF}}\\}$ at
different scales $s$. During training, we randomly apply non-rigid
deformations $t$ to generate more data instances. The degree of deformation is
chosen so that the resulting warped images resemble plausible physical
variations of the medical images. This approach is commonly referred to as
data augmentation and can help avoid over-fitting Krizhevsky et al. (2012),
Cireşan et al. (2013). Each non-rigid training deformation $t$ is computed by
fitting a thin-plate-spline (TPS) to a regular grid of 2D control points
$\left\\{\omega_{i};i=1,2,\ldots,k\right\\}$. These control points are
randomly transformed within the sampling window and a deformed image is
generated using a radial basic function $\phi(r)$, where
$t(x)=\sum\nolimits^{k}_{i=1}c_{i}\phi\left(\left\|x-\omega_{i}\right\|\right)$
is the transformed location of $x$ and $\\{c_{i}\\}$ is a set of mapping
coefficients.

###  2.6 Cross-scale and 3D probability aggregation

At testing, we evaluate each superpixel at $N_{s}$ different scales. The
probability scores for each superpixel being pancreas are averaged across
scales: $p(x)=\frac{1}{N_{s}}\sum\nolimits_{i=1}^{N_{s}}p_{i}(x)$. Then the
resulting per-superpixel ConvNet classification values $\\{p_{1}(x)\\}$ and
$\\{p_{2}(x)\\}$ (according to $R_{1}{-}\mathrm{ConvNet}$ and
$R_{2}{-}\mathrm{ConvNet}$, respectively), are directly assigned to every
pixel or voxel residing within any superpixel $\in\\{S_{\mathrm{RF}}\\}$. This
process forms two per-voxel probability maps $P_{1}(x)$ and $P_{2}(x)$.
Subsequently, we perform 3D Gaussian filtering in order to average and smooth
the ConvNet probability scores across CT slices and within-slice neighboring
regions. 3D isotropic Gaussian filtering can be applied to any $P_{k}(x)$ with
$k=0,1,2$ to form smoothed $G(P_{k}(x))$. This is a simple way to propagate
the 2D slice-based probabilities to 3D by taking local 3D neighborhoods into
account. In this paper, we do not work on 3D supervoxels due to computational
efficiency222Supervoxel based regional ConvNets need at least one-order-of-
magnitude wider input layers and thus have significantly more parameters to
train. and generality issues. We also explore conditional random fields (CRF)
using an additional ConvNet trained between pairs of neighboring superpixels
in order to detect the pancreas edge (defined by pairs of superpixels having
the same or different object labels). This acts as the boundary term together
with the regional term given by $R_{2}{-}\mathrm{ConvNet}$ in order to perform
a min-cut/max-flow segmentation Boykov and Funka-Lea (2006). Here, the CRF is
implemented as a 2D graph with connections between directly neighboring
superpixels. The CRF weighting coefficient between the boundary and the unary
regional term is calibrated by grid-search.

##  3 Results & Discussion

####  3.0.1 Data:

Manual tracings of the pancreas for 82 contrast-enhanced abdominal CT volumes
were provided by an experienced radiologist. Our experiments are conducted
using 4-fold cross-validation in a random hard-split of 82 patients for
training and testing folds with 21, 21, 20, and 20 patients for each testing
fold. We report both training and testing segmentation accuracy results. Most
previous work Wang et al. (2014), Chu et al. (2013), Wolz et al. (2013) uses
leave-one-patient-out cross-validation protocols which are computationally
expensive (e.g., $\sim 15$ hours to process one case using a powerful
workstation Wang et al. (2014)) and may not scale up efficiently towards
larger patient populations. More patients (i.e. 20) per testing fold make the
results more representative for larger population groups.

####  3.0.2 Evaluation:

The ground truth superpixel labels are derived as described in Sec. 2.1. The
optimally achievable DSC for superpixel classification (if classified
perfectly) is 80.5%. Furthermore, the training data is artificially increased
by a factor $N_{s}\times~{}N_{t}$ using the data augmentation approach with
both scale and random TPS deformations at the $R{-}\mathrm{ConvNet}$ level
(Sec. 2.5). Here, we train on augmented data using $N_{s}=4$, $N_{t}=8$. In
testing we use $N_{s}=4$ (without deformation based data augmentation) and
$\sigma=3$ voxels (as 3D Gaussian filtering kernel width) to compute smoothed
probability maps $G(P(x))$. By tuning our implementation of Farag et al.
(2014) at a low operating point, the initial superpixel candidate labeling
achieves the average DSCs of only 26.1% in testing; but has a 97% sensitivity
covering all pancreas voxels. Fig. 5 shows the plots of average DSCs using the
proposed ConvNet approaches, as a function of $P_{k}(x)$ and $G(P_{k}(x))$ in
both training and testing for one fold of cross-validation. Simple Gaussian 3D
smoothing (Sec. 2.6) markedly improved the average DSCs in all cases. Maximum
average DSCs can be observed at $p_{0}=0.2$, $p_{1}=0.5$, and $p_{2}=0.6$ in
our training evaluation after 3D Gaussian smoothing for this fold. These
calibrated operation points are then fixed and used in testing cross-
validation to obtain the results in Table 3.0.2. Utilizing
$R_{2}{-}\mathrm{ConvNet}$ (stacked on $P{-}\mathrm{ConvNet}$) and Gaussian
smoothing ($G(P_{2}(x))$), we achieve a final average DSC of 71.8% in testing,
an improvement of 45.7% compared to the candidate region generation stage at
26.1%. $G(P_{0}(x))$ also performs well wiht 69.5% mean DSC and is more
efficient since only dense deep patch labeling is needed. Even though the
absolute difference in DSC between $G(P_{0}(x))$ and $G(P_{2}(x))$ is small,
the surface-to-surface distance improves significantly from 1.46$\pm$1.5mm to
0.94$\pm$0.6mm, (p$<$0.01). An example of pancreas segmentation at this
operation point is shown in Fig. 6. Training of a typical
$R{-}\mathrm{ConvNet}$ with $N\times~{}N_{s}\times N_{t}=\sim{850k}$
superpixel examples of size $64\times 64$ pixels (after warping) takes
$\sim$55 hours for 100 epochs on a modern GPU (Nvidia GTX Titan-Z). However,
execution run-time in testing is in the order of only 1 to 3 minutes per CT
volume, depending on the number of scales $N_{s}$. Candidate region generation
in Sec. 2.1 consumes another 5 minutes per case.

To the best of our knowledge, this work reports the highest average DSC with
71.8% in testing. Note that a direct comparison to previous methods is not
possible due to lack of publicly available benchmark datasets. We will share
our data and code implementation for future
comparisons333http://www.cc.nih.gov/about/SeniorStaff/ronald_summers.html444http://www.holgerroth.com/.
Previous state-of-the-art results are at $\sim$68% to $\sim$69% Wang et al.
(2014), Chu et al. (2013), Wolz et al. (2013), Farag et al. (2014). In
particular, DSC drops from 68% (150 patients) to 58% (50 patients) under the
leave-one-out protocol Wolz et al. (2013). Our results are based on a 4-fold
cross-validation. The performance degrades gracefully from training
(83.6$\pm$6.3%) to testing (71.8$\pm$10.7%) which demonstrates the good
generality of learned deep ConvNets on unseen data. This difference is
expected to diminish with more annotated datasets. Our methods also perform
with better stability (i.e., comparing 10.7% versus 18.6% Wang et al. (2014),
15.3% Chu et al. (2013) in the standard deviation of DSCs). Our maximum test
performance is 86.9% DSC with 10%, 30%, 50%, 70%, 80%, and 90% of cases being
above 81.4%, 77.6%, 74.2%, 69.4%, 65.2% and 58.9%, respectively. Only 2
outlier cases lie below 40% DSC (mainly caused by over-segmentation into other
organs). The remaining 80 testing cases are all above 50%. The minimal DSC
value of these outliers is 25.0% for $G(P_{2}(x))$. However Wang et al.
(2014), Chu et al. (2013), Wolz et al. (2013), Farag et al. (2014) all report
gross segmentation failure cases with DSC even below 10%. Lastly, the
variation $CRF(P_{2}(x))$ of enforcing $P_{2}(x)$ within a structured
prediction CRF model achieves only 68.2% $\pm$4.1%. This is probably due to
the already high quality of $G(P_{0})$ and $G(P_{2})$ in comparison.

Figure 5: Average DSCs as a function of un-smoothed $P_{k}(x),k=0,1,2$, and 3D
smoothed $G(P_{k}(x)),k=0,1,2$, ConvNet probability maps in training (left)
and testing (right) in one cross-validation fold

.

Table 1: 4-fold cross-validation: optimally achievable DSCs, our initial
candidate region labeling using $S_{RF}$, DSCs on $P(x)$ and using smoothed
$G(P(x))$, and a CRF model for structured prediction (best performance in
bold).

o l—c—c—c—c—c—c—c—c—c DSC (%) | Opt. | $S_{RF(x)}$ | $P_{0}(x)$ | $G(P_{0}(x))$ | $P_{1}(x)$ | $G(P_{1}(x))$ | $P_{2}(x)$ | $G(P_{2}(x))$ | $CRF(P_{2}(x))$
---|---|---|---|---|---|---|---|---|---
Mean | 80.5 | 26.1 | 60.9 | 69.5 | 56.8 | 62.9 | 64.9 | 71.8 | 68.2
Std | 3.6 | 7.1 | 10.4 | 9.3 | 11.4 | 16.1 | 8.1 | 10.7 | 4.1
Min | 70.9 | 14.2 | 22.9 | 35.3 | 1.3 | 0.0 | 33.1 | 25.0 | 59.6
Max | 85.9 | 45.8 | 80.1 | 84.4 | 77.4 | 87.3 | 77.9 | 86.9 | 74.2

Figure 6: Example of pancreas segmentation using the proposed
$R_{2}{-}\mathrm{ConvNet}$ approach in testing. a) The manual ground truth
annotation (in red outline); b) the $G(P_{2}(x))$ probability map; c) the
final segmentation (in green outline) at $p_{2}=0.6$ (DSC=82.7%).

##  4 Conclusion

We present a bottom-up, coarse-to-fine approach for pancreas segmentation in
abdominal CT scans. Multi-level deep ConvNets are employed on both image
patches and regions. We achieve the highest reported DSCs of 71.8$\pm$10.7% in
testing and 83.6$\pm$6.3% in training, at the computational cost of a few
minutes, not hours as in Wang et al. (2014), Chu et al. (2013), Wolz et al.
(2013). The proposed approach can be incorporated into multi-organ
segmentation frameworks by specifying more tissue types since ConvNet
naturally supports multi-class classifications Krizhevsky et al. (2012). Our
deep learning based organ segmentation approach could be generalizable to
other segmentation problems with large variations and pathologies, e.g.
tumors.

##### Acknowledgments:

This work was supported by the Intramural Research Program of the NIH Clinical
Center. The final publication will be available at Springer.

## References

  * Achanta et al. (2012) Achanta, R., A. Shaji, K. Smith, A. Lucchi, P. Fua, and S. Susstrunk (2012).  Slic superpixels compared to state-of-the-art superpixel methods.  Pattern Analysis and Machine Intelligence, IEEE Transactions on 34(11), 2274–2282.
  * Boykov and Funka-Lea (2006) Boykov, Y. and G. Funka-Lea (2006).  Graph cuts and efficient ND image segmentation.  IJCV 70(2), 109–131.
  * Chu et al. (2013) Chu, C., M. Oda, T. Kitasaka, K. Misawa, M. Fujiwara, Y. Hayashi, Y. Nimura, D. Rueckert, and K. Mori (2013).  Multi-organ segmentation based on spatially-divided probabilistic atlas from 3D abdominal CT images.  In MICCAI, pp. 165–172. Springer.
  * Cireşan et al. (2013) Cireşan, D. C., A. Giusti, L. M. Gambardella, and J. Schmidhuber (2013).  Mitosis detection in breast cancer histology images with deep neural networks.  In MICCAI, pp. 411–418. Springer.
  * Everingham et al. (2014) Everingham, M., S. A. Eslami, L. Van Gool, C. K. Williams, J. Winn, and A. Zisserman (2014).  The PASCAL visual object classes challenge: A retrospective.  IJCV 111(1), 98–136.
  * Farag et al. (2014) Farag, A., L. Lu, E. Turkbey, J. Liu, and R. M. Summers (2014).  A bottom-up approach for automatic pancreas segmentation in abdominal CT scans.  MICCAI Abdominal Imaging workshop, arXiv preprint arXiv:1407.8497.
  * Girshick et al. (2014) Girshick, R., J. Donahue, T. Darrell, and J. Malik (2014).  Rich feature hierarchies for accurate object detection and semantic segmentation.  In CVPR, pp. 580–587. IEEE.
  * Krizhevsky et al. (2012) Krizhevsky, A., I. Sutskever, and G. E. Hinton (2012).  Imagenet classification with deep convolutional neural networks.  In NIPS, pp. 1097–1105.
  * Ling et al. (2008) Ling, H., S. K. Zhou, Y. Zheng, B. Georgescu, M. Suehling, and D. Comaniciu (2008).  Hierarchical, learning-based automatic liver segmentation.  In CVPR, pp. 1–8. IEEE.
  * Liu et al. (2011) Liu, M.-Y., O. Tuzel, S. Ramalingam, and R. Chellappa (2011).  Entropy rate superpixel segmentation.  In CVPR, pp. 2097–2104. IEEE.
  * Mostajabi et al. (2014) Mostajabi, M., P. Yadollahpour, and G. Shakhnarovich (2014).  Feedforward semantic segmentation with zoom-out features.  arXiv preprint arXiv:1412.0774.
  * Roth et al. (2014) Roth, H. R., L. Lu, A. Seff, K. M. Cherry, J. Hoffman, S. Wang, J. Liu, E. Turkbey, and R. M. Summers (2014).  A new 2.5D representation for lymph node detection using random sets of deep convolutional neural network observations.  In MICCAI, pp. 520–527. Springer.
  * Srivastava et al. (2014) Srivastava, N., G. Hinton, A. Krizhevsky, I. Sutskever, and R. Salakhutdinov (2014).  Dropout: A simple way to prevent neural networks from overfitting.  JMLR 15(1), 1929–1958.
  * Wang et al. (2014) Wang, Z., K. K. Bhatia, B. Glocker, A. Marvao, T. Dawes, K. Misawa, K. Mori, and D. Rueckert (2014).  Geodesic patch-based segmentation.  In MICCAI, pp. 666–673. Springer.
  * Wolz et al. (2013) Wolz, R., C. Chu, K. Misawa, M. Fujiwara, K. Mori, and D. Rueckert (2013).  Automated abdominal multi-organ segmentation with subject-specific atlas generation.  TMI 32(9), 1723–1730.