Home | About | Journals | Submit | Contact Us | Français |

**|**HHS Author Manuscripts**|**PMC4465221

Formats

Article sections

- Abstract
- 1. Introduction
- 2. Tempered fractional diffusion
- 3. Tempered fractional Calculus
- 4. Tempered fractional Brownian motion
- 5. Tempered fractional difference operator
- 6. ARTFIMA model
- 7. Turbulence
- 8. Summary
- References

Authors

Related links

J Comput Phys. Author manuscript; available in PMC 2016 July 15.

Published in final edited form as:

J Comput Phys. 2015 July 15; 293: 14–28.

doi: 10.1016/j.jcp.2014.04.024PMCID: PMC4465221

NIHMSID: NIHMS590724

MARK M. MEERSCHAERT, Department of Statistics and Probability, Michigan State University, East Lansing MI 48823 ; Email: ude.usm.tts@debucm URL: http://www.stt.msu.edu/users/mcubed/

Fractional derivatives and integrals are convolutions with a power law. Multiplying by an exponential factor leads to tempered fractional derivatives and integrals. Tempered fractional diffusion equations, where the usual second derivative in space is replaced by a tempered fractional derivative, govern the limits of random walk models with an exponentially tempered power law jump distribution. The limiting tempered stable probability densities exhibit semi-heavy tails, which are commonly observed in finance. Tempered power law waiting times lead to tempered fractional time derivatives, which have proven useful in geophysics. The tempered fractional derivative or integral of a Brownian motion, called a tempered fractional Brownian motion, can exhibit semi-long range dependence. The increments of this process, called tempered fractional Gaussian noise, provide a useful new stochastic model for wind speed data. A tempered difference forms the basis for numerical methods to solve tempered fractional diffusion equations, and it also provides a useful new correlation model in time series.

Fractional derivatives were invented by Leibnitz soon after the more familiar integer order derivatives [39, 54], but have only recently become popular in applications. They are now used to model a wide variety of problems in physics [25, 36, 39, 52, 53, 57, 68], finance [23, 28, 44, 38, 63, 64], biology [4, 2, 22, 27, 37], and hydrology [1, 7, 8, 15, 17, 65]. Fractional derivatives can be most easily understood in terms of their connection to probability [47, 48]. Einstein [20] explained the connection between random walks, Brownian motion, and the diffusion equation ${\partial}_{t}p(x,t)={\partial}_{x}^{2}p(x,t)$. Sokolov and Klafter [66] review the modern theory of anomalous diffusion, where the integer order derivatives in the diffusion equation are replaced by their fractional analogues: ${\partial}_{t}^{\beta}p(x,t)={\partial}_{x}^{\alpha}p(x,t)$. A fractional space derivative of order *α <* 2 corresponds to heavy tailed power law particle jumps $\mathbb{P}[J>x]\approx {x}^{-\alpha}$ (the famous Lévy flight), while a fractional time derivative of order *β <* 1 models heavy tailed power law waiting times $\mathbb{P}[W>t]\approx {t}^{-\beta}$ between jumps. Hence fractional space derivatives model anomalous super-diffusion, where a plume of particles spreads faster than the traditional diffusion equation predicts, and fractional time derivatives model anomalous sub-diffusion.

The goal of this paper is to describe a new variation on the fractional calculus, where power laws are tempered by an exponential factor. This exponential tempering has both mathematical and practical advantages. Mantegna and Stanley [41] proposed a truncated Lévy flight to capture the natural cutoff in real physical systems. Koponen [30] introduced the tempered Lévy flight as a smoother alternative, without a sharp cutoff. Cartea and del-Castillo-Negrete [12] developed the tempered fractional diffusion equation that governs the probability densities of the tempered Lévy flight. Unlike the truncated model, tempered Lévy flights offer a complete set of statistical physics and numerical analysis tools. Random walks with exponentially tempered power law jumps converge to a tempered stable motion [13]. Probability densities of the tempered stable motion solve a tempered fractional diffusion equation that describes the particle plume shape [3], just like the original Einstein model for traditional diffusion. Tempered fractional derivatives are approximated by tempered fractional difference quotients, and this facilitates finite difference schemes for solving tempered fractional diffusion equations [3]. The tempered diffusion model has already proven useful in applications to geophysics [46, 73, 74] and finance [10, 11]. In finance, the tempered stable process models price fluctuations with *semi-heavy tails*, resembling a pure power law at moderate time scales, but converging to a Gaussian at long time scales [5]. Since the anomalous diffusion eventually relaxes into a traditional diffusion profile at late time, this model is also called *transient anomalous diffusion* [74].

Kolmogorov [29] invented a new stochastic model for turbulence in the inertial range. Mandelbrot and Van Ness [40] pointed out that this stochastic process is the fractional derivative of a Brownian motion, and coined the name *fractional Brownian motion*. Fractional Brownian motion can exhibit a very useful property called *long range dependence*, where correlations fall off like a power law with time lag. A new variation on this model, called *tempered fractional Brownian motion*, is the tempered fractional derivative of a Brownian motion [49]. Tempered fractional Brownian motion can exhibit *semi-long range dependence*, with correlations that fall off like a power law at moderate time scales, but then eventually become short-range dependent at long time scales. This extends the Kolmogorov model for turbulence to also include low frequencies, and in fact tempered fractional Brownian motion provides a time-domain stochastic process model for the famous Davenport spectrum of wind speed [6, 16, 26, 55], which is used to design electric power generation facilities.

We begin by recalling the connection between random walks, Brownian motion, and the diffusion equation (see [48] for complete details). Given a random walk of mean zero particle jumps *S*(*n*) = *X*_{1} + . . . + *X _{n}*, the Central Limit Theorem implies that

Tempered fractional diffusion applies an exponential tempering factor to the particle jump density. Consider a random walk *S ^{ε}*(

$${f}_{\epsilon}\left(x\right)={C}_{\epsilon}^{-1}{x}^{-\alpha -1}{e}^{-\lambda x}{1}_{[\epsilon ,\infty )}\left(x\right)\phantom{\rule{1em}{0ex}}\text{where}\phantom{\rule{1em}{0ex}}{C}_{\epsilon}={\int}_{\epsilon}^{\infty}{x}^{-\alpha -1}{e}^{-\lambda x}dx$$

(1)

using the incomplete gamma function, and define the Poisson jump rate

$${\lambda}_{\epsilon}=D\frac{\alpha}{\Gamma (1-\alpha )}{C}_{\epsilon}$$

(2)

for any *ε* > 0. To ease notation, we begin with the case of positive jumps.

**Theorem 2.1. ***Suppose* 0 < α < 1. *Given a random walk*
${S}^{\epsilon}\left(n\right)={X}_{1}^{\epsilon}+\cdots +{X}_{n}^{\epsilon}$ with independent jumps, each having probability density function (1)*, and an independent Poisson process*
${N}_{t}^{\epsilon}$
*with rate*
(2)*, as ε* → 0 *we have the convergence*

$${S}^{\epsilon}\left({N}_{t}^{\epsilon}\right)\Rightarrow A\left(t\right)$$

(3)

where the limit is a tempered stable process whose probability density function p(*x, t*) *has Fourier transform*

$$\widehat{p}(k,t)={e}^{-tD[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}]}$$

(4)

*for any t* ≥ 0.

*Proof.* The Poisson random variable satisfies $\mathbb{P}({N}_{t}^{\epsilon}=n)={e}^{-{\lambda}_{e}t}{\left({\lambda}_{\epsilon}t\right)}^{n}\u2215n!$ for *n* ≥0, and then

$${P}_{\epsilon}(x,t)=\mathbb{P}({S}^{\epsilon}\left({N}_{t}^{\epsilon}\right)\le x)=\sum _{n=0}^{\infty}\mathbb{P}({S}^{\epsilon}\left(n\right)\le x\mid {N}_{t}^{\epsilon}=n)\mathbb{P}({N}_{t}^{\epsilon}=n)$$

by the law of total probability. Apply the Fourier-Stieltjes transform $\widehat{f}\left(k\right)=\int {e}^{-ikx}F\left(dx\right)$ where *f*(*x*) = *F*′(*x*), nothing that ${\widehat{f}}_{\epsilon}{\left(k\right)}^{n}$ is the Fourier transform of probability distribution of *S ^{ε}*(

$${\widehat{p}}_{\epsilon}(k,t)=\sum _{n=0}^{\infty}{\widehat{f}}_{\epsilon}{\left(k\right)}^{n}{e}^{-{\lambda}_{\epsilon}t}\frac{\left({\lambda}_{\epsilon}t\right)}{n!}={e}^{-{\lambda}_{\epsilon}t(1-{\widehat{f}}_{\epsilon}\left(k\right))}$$

(5)

using the Taylor series for the complex exponential function. Rewrite this in the form

$$\begin{array}{cc}\hfill {\widehat{p}}_{\epsilon}(k,t)=& \mathrm{exp}\left[t{\lambda}_{\epsilon}{\int}_{\epsilon}^{\infty}({e}^{-ikx}-1){f}_{\epsilon}\left(x\right)dx\right]\hfill \\ \hfill =& \mathrm{exp}\left[tD\frac{\alpha}{\Gamma (1-\alpha )}{\int}_{\epsilon}^{\infty}({e}^{-ikx}-1){x}^{-\alpha -1}{e}^{-\lambda x}dx\right].\hfill \end{array}$$

A straightforward computation using integration by parts and the formula for the gamma probability density [48, Eq. (3.14)] shows that

$${\left(ik\right)}^{\alpha}=\frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}(1-{e}^{-ikx}){x}^{-\alpha -1}dx\phantom{\rule{1em}{0ex}}\text{for}\phantom{\rule{thickmathspace}{0ex}}0<\alpha <1.$$

(6)

Use this formula twice to compute

$$\begin{array}{cc}\hfill \frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}({e}^{-ikx}-1){e}^{-\lambda x}{x}^{-\alpha -1}dx=& \frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}({e}^{-(ik+\lambda )x}-1){x}^{-\alpha -1}dx-\frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}({e}^{-\lambda x}-1){x}^{-\alpha -1}dx\hfill \\ \hfill =& -[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}]\hfill \end{array}$$

(7)

when 0 < *α* < 1. Then as *ε* → 0, we get

$${\widehat{p}}_{\epsilon}(k,t)\to \widehat{p}(k,t)={e}^{-tD[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}]}$$

for any *t* ≥ 0.

It follows from Theorem 2.1 that (*k,t*) = *e*^{–tD[(λ+ik)α–λα]} is the Fourier transform of the *tempered Lévy flight A*(*t*), a stochastic process with tempered power law jumps. We define the *tempered fractional derivative*
${\partial}_{x}^{\alpha ,\lambda}f\left(x\right)$ as the function with Fourier transform $[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}]\widehat{f}\left(k\right)$ when 0 < *α* < 1. Note that

$${\partial}_{t}\widehat{p}(k,t)=-D[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}]\widehat{p}(k,t)$$

and invert the Fourier transform to see that the probability density *p*(*x, t*) of the *tempered Lévy flight A*(*t*) is the point source solution to the *tempered fractional diffusion equation*

$${\partial}_{t}p(x,t)=-D{\partial}_{x}^{\alpha ,\lambda}p(x,t).$$

We will also define the negative tempered fractional derivative ${\partial}_{-x}^{\alpha ,\lambda}f\left(x\right)$ as the function with Fourier transform [(λ − *ik*)* ^{α}* λ

$${f}_{\epsilon}\left(x\right)={C}_{\epsilon}^{-1}[p{x}^{-\alpha -1}{e}^{-\lambda x}{1}_{[\epsilon ,\infty )}\left(x\right)+q{\mid x\mid}^{-\alpha -1}{e}^{-\lambda \mid x\mid}{1}_{(-\infty ,-\epsilon ]}\left(x\right)]$$

(8)

for *p, q ≥* 0 with *p* + *q* = 1 leads to a tempered stable limit process *A*(*t*) having both positive and negative jumps. A simple extension of Theorem 2.1 shows that (3) holds, and now the limit density *p*(*x, t*) solves the tempered fractional diffusion equation

$${\partial}_{t}p(x,t)=-pD{\partial}_{x}^{\alpha ,\lambda}p(x,t)-qD{\partial}_{-x}^{\alpha ,\lambda}p(x,t).$$

To write the tempered fractional derivative in real space, use the formula (7) to see that

$$[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}]\widehat{f}\left(k\right)=\frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}\left(\widehat{f}\left(k\right)-{e}^{-iky}\widehat{f}\left(k\right)\right){e}^{-\lambda y}{y}^{-\alpha -1}dy$$

and then use the shift property of the Fourier transform ∫ *e ^{−ikx}f*(

$${\partial}_{x}^{\alpha ,\lambda}f\left(x\right)=\frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}(f\left(x\right)-f(x-y)){e}^{-\lambda y}{y}^{-\alpha -1}dy$$

(9)

for 0 < *α* < 1. A similar argument shows that

$${\partial}_{-x}^{\alpha ,\lambda}f\left(x\right)=\frac{\alpha}{\Gamma (1-\alpha )}{\int}_{0}^{\infty}(f\left(x\right)-f(x-y)){e}^{-\lambda y}{y}^{-\alpha -1}dy.$$

(10)

**Theorem 2.2. ***Suppose* 1 < *α* < 2. *Given the mean-centered random walk*
${S}^{\epsilon}\left(n\right)={X}_{1}^{\epsilon}+\cdots +{X}_{n}^{\epsilon}$ with independent jumps, each having probability density function (1)
*and mean*
${\mu}_{\epsilon}=\alpha {C}_{\epsilon}^{1\u2215\alpha}\u2215(\alpha -1)$, *and an independent Poisson process*
${N}_{t}^{\epsilon}$
*with rate*
(2), *as ε →* 0 *we have the convergence*

$${S}^{\epsilon}\left({N}_{t}^{\epsilon}\right)-t{\lambda}_{\epsilon}{\mu}_{\epsilon}\Rightarrow A\left(t\right)$$

(11)

where the limit is a mean zero tempered stable Lévy motion whose probability density function p(*x, t*) *has Fourier transform*

$$\widehat{p}(k,t)={e}^{-tD[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}-ik\alpha {\lambda}^{\alpha -1}]}$$

(12)

*for any t* ≥ 0.

*Proof.* An easy computation shows that the jump variable *X ^{ε}* with probability density (1) has mean ${\mu}_{\epsilon}=\int x{f}_{\epsilon}\left(x\right)dx=\alpha {C}_{\epsilon}^{1\u2215\alpha}\u2215(\alpha -1)$. Then it follows from (5) that the density of ${S}^{\epsilon}\left({N}_{t}^{\epsilon}\right)-t{\lambda}_{\epsilon}{\mu}_{\epsilon}$ has Fourier transform

$${\widehat{p}}_{\epsilon}(k,t)={e}^{-{\lambda}_{\epsilon}t(1-{\widehat{f}}_{\epsilon}\left(k\right))}{e}^{-ikt{\lambda}_{\epsilon}{\mu}_{\epsilon}}$$

for any *t* ≥ 0. Rewrite this in the form

$$\begin{array}{cc}\hfill {\widehat{p}}_{\epsilon}(k,t)=& \mathrm{exp}\left[-t{\lambda}_{\epsilon}{\int}_{\epsilon}^{\infty}(1-{e}^{-ikx}+ikx){f}_{\epsilon}\left(x\right)dx\right]\hfill \\ \hfill =& \mathrm{exp}\left[tD\frac{\alpha (\alpha -1)}{\Gamma (2-\alpha )}{\int}_{\epsilon}^{\infty}({e}^{-ikx}-1-ikx){x}^{-\alpha -1}{e}^{-\lambda x}dx\right]\hfill \end{array}$$

using the fact that Γ(2 − *α*) = (1 *− α*)Γ(1 *α*). Use (6) and another integration by parts [48, p. 58] to see that

$${\left(ik\right)}^{\alpha}=\frac{\alpha (\alpha -1)}{\Gamma (2-\alpha )}{\int}_{0}^{\infty}({e}^{-ikx}-1+ikx){x}^{-\alpha -1}dx\phantom{\rule{1em}{0ex}}\text{for}\phantom{\rule{thickmathspace}{0ex}}1<\alpha <2.$$

(13)

Let *C* = *α*(*α −* 1)/Γ(2 − *α*) and compute

$$\begin{array}{cc}\hfill C{\int}_{0}^{\infty}({e}^{-iky}-1+iky){e}^{-\lambda y}{y}^{-\alpha -1}dy=& C{\int}_{0}^{\infty}({e}^{-(\lambda +ik)y}-1+(\lambda +ik)y){y}^{-\alpha -1}dy\hfill \\ \hfill -& C{\int}_{0}^{\infty}({e}^{-\lambda y}-1+\lambda y){y}^{-\alpha -1}dy-C{\int}_{0}^{\infty}({e}^{-\lambda y}-1)iky\phantom{\rule{thinmathspace}{0ex}}{y}^{-\alpha -1}dy\hfill \\ \hfill =& {(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}-ik\alpha {\lambda}^{\alpha -1}\hfill \end{array}$$

(14)

using (13). Then as *ε* → 0, we get

$${\widehat{p}}_{\epsilon}(k,t)\to \widehat{p}(k,t)={e}^{-tD[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}-ik\alpha {\lambda}^{\alpha -1}]}$$

for any *t* ≥ 0.

It follows from Theorem 2.2 that (*k, t*) = *e*^{−tD[(λ+ik)α−λα−ikαλα−1]} is the Fourier transform of the tempered Lévy flight *A*(*t*). We define the *tempered fractional derivative*
${\partial}_{x}^{\alpha ,\lambda}f\left(x\right)$ as the function with Fourier transform [(*λ* + *ik*)* ^{α}* − λ

$${\partial}_{t}p(x,t)=D{\partial}_{x}^{\alpha ,\lambda}p(x,t).$$

(15)

Figure 1 shows a typical particle path, illustrating the effect of tempering to cool the long particle jumps.

To write this tempered fractional derivative in real space, use the formula (14) to see that

$$[{(\lambda +ik)}^{\alpha}-{\lambda}^{\alpha}-ik\alpha {\lambda}^{\alpha -1}]\widehat{f}\left(k\right)=\frac{\alpha (\alpha -1)}{\Gamma (2-\alpha )}\widehat{f}\left(k\right){\int}_{0}^{\infty}({e}^{-iky}-1+iky){e}^{-\lambda y}{y}^{-\alpha -1}dy$$

and invert the Fourier transform to get

$${\partial}_{x}^{\alpha ,\lambda}f\left(x\right)=\frac{\alpha (\alpha -1)}{\Gamma (2-\alpha )}{\int}_{0}^{\infty}(f(x-y)-f\left(x\right)+y{f}^{\prime}\left(x\right)){e}^{-\lambda y}{y}^{-\alpha -1}dy$$

(16)

for 1 *<* α < 2.

The negative tempered fractional derivative ${\partial}_{-x}^{\alpha ,\lambda}f\left(x\right)$ is the function with Fourier transform [(λ − *ik*)* ^{α} −λ^{α}* +

$${\partial}_{-x}^{\alpha ,\lambda}f\left(x\right)=\frac{\alpha (\alpha -1)}{\Gamma (2-\alpha )}{\int}_{0}^{\infty}(f(x+y)-f\left(x\right)-y{f}^{\prime}\left(x\right)){e}^{-\lambda y}{y}^{-\alpha -1}dy.$$

(17)

The particle jump density (8) leads to a tempered stable limit process *A*(*t*) having both positive and negative jumps. A simple extension of Theorem 2.2 shows that (11) holds, and now the limit density *p*(*x, t*) solves the tempered fractional diffusion equation

$${\partial}_{t}p(x,t)=pD{\partial}_{x}^{\alpha ,\lambda}p(x,t)+qD{\partial}_{-x}^{\alpha ,\lambda}p(x,t).$$

If *p* = *q* then this produces a symmetric profile, with *semi-heavy tails*. That is, the tails of the probability density *p*(*x, t*) resemble a power law *p*(*x, t*) |*x*|^{−α−1} for large |*x*| and moderate *t* > 0, but as *t* → ∞ the tails relax to a Gaussian profile. These semi- heavy tails are typical in applications to finance [10, 11]. Figure 2 fits a symmetric tempered stable density function to macroeconomic data on inflation rates. The data shows a sharper peak and a heavier tail than the best fitting Gaussian curve.

Symmetric tempered stable model fit (solid line) with *α* = 1.1, λ = 12, *D* = 0.1 to one-step bilinear ARMA forecast errors for annual inflation rates. The best fitting normal density (dotted line) fails to capture the sharp data peak.

The tempered space-fractional diffusion equation models transient super-diffusion, where a particle plume initially spreads faster than the traditional diffusion equation predicts, but later relaxes to a typical diffusion profile. The corresponding model for transient sub-diffusion assumes tempered power law waiting times between particle jumps. Suppose that the *n* th particle jump *S*(*n*) = *X*_{1} + · · · + *X _{n}* occurs at a random time

A general result from CTRW limit theory [45, Theorem 2.1] shows that the CTRW converges at late time to *A*(*E _{t}*), a tempered Lévy flight with the time variable replaced by an independent

$$\stackrel{~}{h}(u,s)={\int}_{0}^{\infty}{e}^{-st}h(u,t)dt={s}^{-1}{\psi}_{D}\left(s\right){e}^{-u{\psi}_{D}\left(s\right)},$$

see [45, Theorem 3.1] for details. Since *x* = *A*(*u*) has density function *p*(*x, u*) and *u* = *E _{t}* has density function

$$q(x,t)={\int}_{0}^{\infty}p(x,u)h(u,t)du.$$

Recall that the tempered Lévy density has Fourier transform (*k, u*) = *e ^{−uψA}*

$$\stackrel{\u2012}{q}(k,s)={\int}_{0}^{\infty}\widehat{p}(k,u)\stackrel{~}{h}(u,s)du=\frac{{s}^{-1}{\psi}_{D}\left(s\right)}{{\psi}_{D}\left(s\right)+{\psi}_{A}\left(k\right)},$$

then rearrange and invert to see that the limit density solves the space-time tempered fractional diffusion equation

$${\partial}_{t}^{\beta ,\eta}p(x,t)=pD{\partial}_{x}^{\alpha}p(x,t)+qD{\partial}_{-x}^{\alpha}p(x,t)+\delta \left(x\right)f\left(t\right)$$

(18)

where the forcing term *f*(*t*) = * _{D}*(

Given any λ > 0, we define the positive tempered fractional integral of a suitable function *f*(*x*) by

$${\mathbb{I}}_{+}^{\alpha ,\lambda}f\left(x\right)=\frac{1}{\Gamma \left(\alpha \right)}{\int}_{-\infty}^{x}f\left(u\right){(x-u)}^{\alpha -1}{e}^{-\lambda (x-u)}du,$$

(19)

and the negative tempered fractional integral by

$${\mathbb{I}}_{-}^{\alpha ,\lambda}f\left(x\right)=\frac{1}{\Gamma \left(\alpha \right)}{\int}_{x}^{\infty}f\left(u\right){(u-x)}^{\alpha -1}{e}^{-\lambda (u-x)}du.$$

(20)

If λ = 0, these formulae reduce to the well-known Riemann-Liouville fractional integrals [48, 61]. Hence we will call these operators the *Riemann-Liouville tempered fractional integrals*. Since these operators are convolutions, their Fourier transforms are products. Since *g*(*x*) = *x*^{α−1}*e ^{−λx}*

For functions in the fractional Sobolev space

$${W}^{\alpha ,2}\left(\mathbb{R}\right)\u2254\{f\in {L}^{2}\left(\mathbb{R}\right):{\int}_{\mathbb{R}}{({\lambda}^{2}+{k}^{2})}^{\alpha}{\mid \widehat{f}\left(k\right)\mid}^{2}dk<\infty \},$$

we define the *Riemann-Liouville tempered fractional derivatives*
${\mathbb{D}}_{\pm}^{\alpha ,\lambda}f\left(x\right)$ to be the functions with Fourier transform (λ ± *ik*)* ^{α}*(

$${\mathbb{D}}_{\pm}^{\alpha ,\lambda}f\left(x\right)={\partial}_{\pm x}^{\alpha ,\lambda}f\left(x\right)+{\lambda}^{\alpha}f\left(x\right)\phantom{\rule{1em}{0ex}}\text{for}\phantom{\rule{thickmathspace}{0ex}}0<\alpha <1.$$

(21)

Since ${\partial}_{\pm x}^{\alpha ,\lambda}f\left(x\right)$ has Fourier transform [(*λ ± ik*)* ^{α}* − λ

$${\mathbb{D}}_{\pm}^{\alpha ,\lambda}f\left(x\right)={\partial}_{\pm x}^{\alpha ,\lambda}f\left(x\right)+{\lambda}^{\alpha}f\left(x\right)\pm \alpha {\lambda}^{\alpha -1}{f}^{\prime}\left(x\right)\phantom{\rule{1em}{0ex}}\text{for}\phantom{\rule{thickmathspace}{0ex}}1<\alpha <2.$$

(22)

Using the shift property of the Fourier transform, we can also relate the Riemann-Liouville tempered fractional derivative to the (untempered) fractional derivative ${\partial}_{x}^{\alpha}f\left(x\right)$ from Section 2. Since ∫ *e ^{−ikx}e^{λx}f*(

$${\mathbb{D}}_{+}^{\alpha ,\lambda}f\left(x\right)={e}^{-\lambda x}{\partial}_{x}^{\alpha}\left[{e}^{\lambda x}f\left(x\right)\right]\phantom{\rule{1em}{0ex}}\text{for any}\phantom{\rule{thickmathspace}{0ex}}\alpha >0.$$

(23)

Using the inner product $\langle f,g\rangle =\int f\left(x\right)\overline{g\left(x\right)}dx$ on ${L}^{2}\left(\mathbb{R}\right)$, along with the Plancherel Theorem *f, g* = *, ĝ,* it follows that

$$\langle {\mathbb{D}}_{+}^{\alpha ,\lambda}f,g\rangle =\int {(\lambda +ik)}^{\alpha}\widehat{f}\left(k\right)\overline{\widehat{g}\left(k\right)}dk=\int \widehat{f}\left(k\right)\overline{{(\lambda -ik)}^{\alpha}\widehat{g}\left(k\right)}dk=\langle f,{\mathbb{D}}_{-}^{\alpha ,\lambda}g\rangle $$

so that the positive and negative Riemann-Liouville tempered fractional derivatives are adjoint operators in the Hilbert space ${L}^{2}\left(\mathbb{R}\right)$. Essentially the same argument shows that the positive and negative Riemann-Liouville tempered fractional integrals are also ${L}^{2}\left(\mathbb{R}\right)$ adjoints. The semigroup property

$${\mathbb{D}}_{\pm}^{\alpha ,\lambda}{\mathbb{D}}_{\pm}^{\beta ,\lambda}f={\mathbb{D}}_{\pm}^{\alpha +\beta ,\lambda}f\phantom{\rule{1em}{0ex}}\phantom{\rule{1em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{1em}{0ex}}{\mathbb{I}}_{\pm}^{\alpha ,\lambda}{\mathbb{I}}_{\pm}^{\beta ,\lambda}f={\mathbb{I}}_{\pm}^{\alpha +\beta ,\lambda}f$$

follows by another easy Fourier transform argument.

The space of rapidly decreasing functions $\mathcal{S}\left(\mathbb{R}\right)$ consists of the infinitely differentiable functions $g:\mathbb{R}\to \mathbb{R}$ such that

$$\underset{x\in \mathbb{R}}{\mathrm{sup}}\mid {x}^{n}{g}^{\left(m\right)}\left(x\right)\mid <\infty ,$$

where *n, m* are non-negative integers, and *g*^{(m)} is the derivative of order *m*. The space of tempered distributions ${\mathcal{S}}^{\prime}\left(\mathbb{R}\right)$ consists of continuous linear functionals on $\mathcal{S}\left(\mathbb{R}\right)$. The Fourier transform maps ${\mathcal{S}}^{\prime}\left(\mathbb{R}\right)$ into itself. If $f:\mathbb{R}\to \mathbb{R}$ is polynomial growth, so that ∫ |*f*(*x*)|(1 + |*x*|)* ^{−p}dx* < ∞ for some

$${\mathbb{D}}_{+}^{\alpha ,\lambda}f\left(x\right)={\mathbb{D}}_{+}^{n,\lambda}{\mathbb{I}}_{+}^{n=\alpha ,\lambda}f\left(x\right)$$

where ${\mathbb{D}}_{+}^{n,\lambda}f\left(x\right)={e}^{-\lambda x}{\left[{e}^{\lambda x}f\left(x\right)\right]}^{\left(n\right)}$ is a tempered derivative of integer order, defined by (23). There are many technical issues behind the fractional calculus of tempered distributions [61, Chapter 8], and it would be interesting to extend that theory for tempered fractional calculus.

As noted in Section 2, a random walk *S*(*n*) = *X*_{1} + · · · + *X _{n}* with mean zero, finite variance particle jumps converges to a Brownian motion

$${B}_{H,\lambda}\left(t\right)\u2254\int \left[{e}^{-\lambda {(t-x)}_{+}}{(t-x)}_{+}^{H-1\u22152}-{e}^{-\lambda {(0-x)}_{+}}{(0-x)}_{+}^{H-1\u22152}\right]B\left(dx\right)$$

(24)

where (*x*)_{+} = **1**_{[0,∞)}(*x*), and 0^{0} = 0. The Gaussian stochastic process (24) has mean zero, and its increments are stationary, but not independent. Using the fact that *B*(*c dx*) *c*^{1/2}*B*(*dx*) (same distribution), it is not hard to check that tempered fractional Brownian motion has a nice scaling property: *B _{H,λ}*(

Stochastic integrals with respect to Brownian motion can also be constructed using white noise theory. Heuristically, one writes ∫ *f*(*x*)*B*(*dx*) = *f*(*x*)*W* (*x*)*dx* where the white noise *W* (*x*) = *B′* (*x*) is the (weak) derivative of the Brownian motion. This integral does not exist path-wise, since the paths of a Brownian motion are not differentiable. However, we can define *W* (*f*) := ∫ *f*(*x*)*W* (*x*)*dx* = *f, W* as a random distribution (generalized function) on ${L}^{2}\left(\mathbb{R}\right)$ by specifying that this object has the same (Gaussian) distribution and covariance function as ∫ *f*(*x*)*B*(*dx*) [31]. Taking *f*(*x*) = **1**_{[0,t]}(*x*), it follows that $B\left(t\right)={\int}_{0}^{t}W\left(x\right)dx$, and in this sense, the white noise *W* (*x*) is the (distributional) derivative of *B*(*x*). Suppose that *H* > 1/2, and compute

$${\mathbb{I}}_{-}^{H-1\u22152,\lambda}{1}_{[0,t]}\left(x\right)-\lambda {\mathbb{I}}_{-}^{H+1\u22152,\lambda}{1}_{[0,t]}\left(x\right)={e}^{-\lambda {(t-x)}_{+}}{(t-x)}_{+}^{H-1\u22152}-{e}^{-\lambda {(0,-x)}_{+}}{(0-x)}_{+}^{H-1\u22152}$$

for any *t* > 0. Using the fact that positive and negative Riemann-Liouville tempered fractional integrals are ${L}^{2}\left(\mathbb{R}\right)$ adjoints, we can therefore use white noise theory to write

$${B}_{H,\lambda}\left(t\right)=\langle {1}_{[0,t]}\left(x\right),{\mathbb{I}}_{+}^{H-1\u22152,\lambda}W\left(x\right)-\lambda {\mathbb{I}}_{+}^{H+1\u22152,\lambda}W\left(x\right)\rangle $$

which shows that tempered fractional Brownian motion with *H* > 1/2 can be obtained by integrating a linear combination of tempered fractional integrals of a white noise. The same is true for 0 < *H* < 1/2, if we interpret the tempered fractional integral of negative order as a tempered fractional derivative, see [51] for more details. Now if we take λ = 0 and 0 < *H* < 1, we recover the fact that fractional Brownian motion is the fractional integral of a white noise [40, 56].

Fractional Brownian motion with 1/2 < *H* < 1 is a popular model for *long range dependence*. The increments *X _{n}* =

Tempered fractional derivatives can also be defined as the limit of a tempered fractional difference quotient.

**Theorem 5.1. ***If f and its derivatives up to order n* > 1 + α exist and are absolutely integrable, then the Riemann-Liouville tempered fractional derivative

$${\mathbb{D}}_{+}^{\alpha ,\lambda}f\left(x\right)=\underset{h\to 0}{\mathrm{lim}}{h}^{-\alpha}{\Delta}_{\lambda}^{\alpha}f\left(x\right)$$

(25)

*where*

$${\Delta}_{\lambda}^{\alpha}f\left(x\right)=\sum _{j=0}^{\infty}{w}_{j}{e}^{-\lambda jh}f(x-jh)\phantom{\rule{1em}{0ex}}with\phantom{\rule{1em}{0ex}}{w}_{j}\u2254{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right)=\frac{{(-1)}^{j}\Gamma (1+\alpha )}{j!\Gamma (1+\alpha -j)}.$$

*Proof.* The proof is similar to [48, Proposition 2.1]. The fractional binomial formula

$$\sum _{k=0}^{\infty}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill k\hfill \end{array}\right){z}^{k}={(1+z)}^{\alpha}$$

holds for any *α* > 0 and any complex |*z*| ≤ 1. Since *e ^{−ikjh}*(

$$\begin{array}{cc}\hfill & {h}^{-\alpha}\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){e}^{-\lambda jh}{e}^{-ikjh}\widehat{f}\left(k\right)\hfill \\ \hfill & ={h}^{-\alpha}\sum _{j=0}^{\infty}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){(-1)}^{j}{e}^{-(\lambda +ik)jh}\widehat{f}\left(k\right)\hfill \\ \hfill & ={h}^{-\alpha}{(1-{e}^{-(\lambda +ik)h})}^{\alpha}f\left(k\right)\to {(\lambda +ik)}^{\alpha}\widehat{f}\left(k\right)\hfill \end{array}$$

as *h* → 0, using a Taylor expansion of the complex exponential function. Since *f* and its derivatives up to order *n* > 1 + *α* exist and are absolutely integrable, we have (*k*) ≤ *C*(1 + |*k*|* ^{n}*) for some

The tempered fractional difference quotient can be applied to develop simple, effective numerical schemes to solve the tempered fractional diffusion equation. For physical reasons [7, 8, 17, 65], the most important case in applications is 1 < *α* < 2. Consider the tempered fractional diffusion equation with drift

$${\partial}_{t}p(x,t)=-v\left(x\right){\partial}_{x}p(x,t)+D\left(x\right){\partial}_{x}^{\alpha ,\lambda}p(x,t)+f(x,t)$$

(26)

where 1 < *α* ≤ 2, *v*(*x*) ≥ 0, and *D*(*x*) ≥ 0 on a finite domain *a* ≤ *x* ≤ *b*, with initial condition *p*(*x,* 0) = *p*_{0}(*x*), boundary conditions *p*(*a, t*) = 0 and *p*(*b, t*) = *B*(*t*), and forcing function *f*(*x, t*) for all *t* ≥ 0. Using Theorem 5.1 along with (22) and the fact that *h ^{−α}(1 − e^{−λh}*)

$${\partial}_{x}^{\alpha ,\lambda}f\left(x\right)+\alpha {\lambda}^{\alpha -1}{f}^{\prime}\left(x\right)=\underset{h\to 0}{\mathrm{lim}}{h}^{-\alpha}\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){e}^{-\lambda jh}f(x-jh)-{h}^{-\alpha}{(1-{e}^{-\lambda h})}^{\alpha}f\left(x\right).$$

To write the details of our numerical code, set *t _{i}* =

$${p}_{i+1,j}-{p}_{ij}=\left[-{v}_{j}\frac{{p}_{i+1,j}-{p}_{i+1,j-1}}{h}+{D}_{j}{h}^{-\alpha}\sum _{m=0}^{j}{g}_{m}{p}_{i+1,j-m}+{f}_{ij}\right]\Delta t$$

using the weights

$${g}_{j}=\{\begin{array}{cc}{w}_{j}{e}^{-jh\lambda}\hfill & \text{for}\phantom{\rule{thickmathspace}{0ex}}j\ge 1,\hfill \\ {w}_{j}{e}^{-jh\lambda}-{(1-{e}^{-h\lambda})}^{\alpha}\hfill & \text{for}\phantom{\rule{thickmathspace}{0ex}}j=0.\hfill \end{array}\phantom{\}}$$

However, this implicit scheme is unstable. The proof is very similar to [43, Proposition 2.3]. The explicit Euler scheme (replace *i* + 1 by *i* on the right-hand side) is also unstable. The proof is essentially identical to [43, Proposition 2.1].

In order to obtain a stable convergent numerical method, it is necessary to use a shifted finite difference approximation

$${\partial}_{x}^{\alpha ,\lambda}f\left(x\right)+\alpha {\lambda}^{\alpha -1}{f}^{\prime}\left(x\right)=\underset{h\to 0}{\mathrm{lim}}{h}^{-\alpha}\sum _{j=0}^{\infty}{g}_{j}f(x-(j-1)h),$$

(27)

where the exponentially tempered Grünwald weights are defined by

$${g}_{j}=\{\begin{array}{cc}{w}_{j}{e}^{-(j-1)h\lambda}\hfill & \text{for}\phantom{\rule{thickmathspace}{0ex}}j\ne 1,\hfill \\ {w}_{1}-{e}^{h\lambda}{(1-{e}^{-h\lambda})}^{\alpha}\hfill & \text{for}\phantom{\rule{thickmathspace}{0ex}}j=1.\hfill \end{array}\phantom{\}}$$

(28)

The proof that (56) holds is similar to Theorem 5.1, see [3, Proposition 3] for details. Now consider the implicit Euler scheme

$${p}_{i+1,j}-{p}_{ij}=\left[-{v}_{j}\frac{{p}_{i+1,j}-{p}_{i+1,j-1}}{h}+{D}_{j}{h}^{-\alpha}\sum _{m=0}^{j+1}{g}_{m}{p}_{i+1,j-m+1}+{f}_{ij}\right]\Delta t.$$

(29)

**Theorem 5.2. ***The implicit Euler scheme*
(29)
*for the tempered fractional diffusion equation*
(26)
*using the weights*
(28)
*is consistent and unconditionally stable.*

*Proof.* The proof is similar to [43, Theorem 2.7]. The main difference is that the Grünwald weights *w _{j}* are replaced by the tempered weights

$$\sum _{j=0}^{\infty}{g}_{j}={e}^{h\lambda}\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){e}^{-jh\lambda}-{e}^{h\lambda}{(1-{e}^{-h\lambda})}^{\alpha}=0,$$

and hence this numerical method mass-preserving. Since 1 < *α* < 2 we have *w*_{1} = −*α* and *w _{j}* > 0 for all

$$\sum _{j=0}^{K}{g}_{j}<0$$

for any finite number *K*. The implicit Euler scheme can be written as a linear system of equations *AP _{i}*

Note that the shifted finite difference approximation (27) is essentially a convolution with the exponentially tempered Grünwald weights. An application of Stirling's formula shows that *w _{m}* ~ −

**Example 5.3.** Consider the tempered fractional diffusion equation (26) with *v*(*x*) = 0, *x* [0,1], *α* = 1.5, λ = 0.5, initial condition *u*(*x*, 0) = *x ^{βe−λx}*/Γ(

$$q(x,t)={e}^{-\lambda x-1}\frac{\Gamma (1+\beta -\alpha )}{\Gamma (\beta +1)}\left(\frac{(1-\alpha ){\lambda}^{\alpha}{x}^{\alpha +\beta}}{\Gamma (\beta +1)}+\frac{\alpha {\lambda}^{\alpha -1}{x}^{\alpha +\beta -1}}{\Gamma \left(\beta \right)}-\frac{2{x}^{\beta}}{\Gamma (1+\beta -\alpha )}\right)$$

and boundary conditions *u*(0*, t*) = 0 and *u*(1*, t*) = *B*(*t*) := *e ^{−λ−t}/*Γ(

$${\stackrel{\u2012}{A}}_{ij}=\{\begin{array}{cc}1+{v}_{i}{h}^{-1}\Delta t-{D}_{i}{h}^{-\alpha}\Delta t{g}_{1}\hfill & \text{when}\phantom{\rule{thickmathspace}{0ex}}j=i,\hfill \\ -{v}_{i}{h}^{-1}\Delta t-{D}_{i}{h}^{-\alpha}\Delta t{g}_{2}\hfill & \text{when}\phantom{\rule{thickmathspace}{0ex}}j=i-1,\hfill \\ -{D}_{i}{h}^{-\alpha}\Delta t{g}_{i-j+1}\hfill & \text{when}\phantom{\rule{thickmathspace}{0ex}}j<i-1\phantom{\rule{thickmathspace}{0ex}}\text{or}\phantom{\rule{thickmathspace}{0ex}}j=i+1,\hfill \\ 0\hfill & \text{otherwise}.\hfill \end{array}\phantom{\}}$$

Then the vector of interior points * _{i}* = [

The implicit Euler method in Theorem 5.2 has error *O*(Δ*x* + Δ*t*), consistent with Table 1. A Crank-Nicolson method for solving the tempered fractional diffusion equation (26) with 1 < *α* ≤ 2 is developed in [3]. The method is stable and consistent, and second order convergent in time. However, because the finite difference approximation (27) is only first order accurate, the method is *O*(Δ*t*^{2} + Δ*x*). An application of Richardson extrapolation yields *O*(Δ*t*^{2} + Δ*x*^{2}), see [3] for more details. It is also possible to construct particle tracking solutions based on the Monte Carlo simulation scheme in [3, Section 4], see [46, 74]. A wide variety of effective numerical methods have been developed to solve fractional diffusion equations, see for example [18, 33, 34, 35, 21, 59, 67, 71]. It would certainly be interesting to extend these methods to tempered fractional diffusion equations. Even in the untempered case, many interesting open problems remain, including a rigorous proof that fractional boundary value problems are well-posed, see [19] for a nice discussion.

Another useful application of the fractional difference operator arises in time series analysis. Given a sequence of independent Gaussian random variables *Z _{n}* with mean zero and variance

$${X}_{n}={\Delta}^{-\alpha}{Z}_{n}=\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill -\alpha \hfill \\ \hfill j\hfill \end{array}\right){Z}_{n-j}$$

where typically −1/2 < *α* < 1/2. Since Δ* ^{α}*Δ

$${n}^{-H}({X}_{1}+\cdots +{X}_{\left[nt\right]})\Rightarrow {B}_{H}\left(t\right)\phantom{\rule{1em}{0ex}}\phantom{\rule{1em}{0ex}}\text{as}\phantom{\rule{thickmathspace}{0ex}}n\to \infty $$

for all *t* > 0, see [69, Theorem 4.6.1]. The fractional noise can also exhibit long range dependence. The covariance between *X _{j}* and

Fractional differencing can be combined with the standard ARMA(*p, q*) model for time series, to yield the ARFIMA(*p, α, q*) time series model. The ARMA(*p, q*) model is

$${X}_{t}-\sum _{j=1}^{p}{\varphi}_{j}{X}_{t=j}={Z}_{t}+\sum _{i=1}^{q}{\theta}_{i}{Z}_{t-i}$$

and we say that *X _{t}* follows an ARFIMA(

$${f}_{X}\left(k\right)=\frac{1}{2\pi}\sum _{j=-\infty}^{+\infty}{e}^{ikh}\gamma \left(h\right).$$

The spectral density of the noise sequence *Z _{n}* is a constant

$${f}_{X}\left(k\right)=\frac{{\mid \u03f4\left({e}^{-ik}\right)\mid}^{2}}{{\mid \Phi \left({e}^{-ik}\right)\mid}^{2}}{f}_{Z}\left(k\right)$$

using the complex absolute value. Indeed, for any time invariant linear filter *X _{t}* = Σ

$${\Delta}^{\alpha}=\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){B}_{j}={(1-B)}^{\alpha}$$

in operator notation, a similar argument shows that the spectral density of an ARFIMA(*p, α, q*) time series is given by

$${f}_{X}\left(k\right)=\frac{{\sigma}^{2}}{2\pi}\frac{{\mid \u03f4\left({e}^{-ik}\right)\mid}^{2}}{{\mid \Phi \left({e}^{-ik}\right)\mid}^{2}}{\mid 1-{e}^{-ik}\mid}^{-2\alpha}.$$

As *k* → 0, the fractional term |1 − e^{−ik}|^{−2α} ~ |*k*|^{−2}* ^{α}* blows up at a power law rate. This spectral density is the hallmark of long range dependence.

The tempered fractional difference

$${X}_{n}={\Delta}_{\lambda}^{\alpha}{Z}_{n}=\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){e}^{-\lambda j}{Z}_{n-j}$$

can provide a more flexible model for (semi-)long range dependence. We say that *X _{t}* follows an ARTFIMA(

$${\Delta}_{\lambda}^{\alpha}=\sum _{j=0}^{\infty}{(-1)}^{j}\left(\begin{array}{c}\hfill \alpha \hfill \\ \hfill j\hfill \end{array}\right){e}^{-\lambda j}{B}^{j}={(1-{e}^{-\lambda}B)}^{\alpha}$$

in operator notation, it follows that the spectral density of an ARTFIMA(*p, α, λ, q*) time series is given by

$${f}_{X}\left(k\right)=\frac{{\sigma}^{2}}{2\pi}\frac{{\mid \u03f4\left({e}^{-ik}\right)\mid}^{2}}{{\mid \Phi \left({e}^{-ik}\right)\mid}^{2}}{\mid 1-{e}^{-(\lambda +ik)}\mid}^{-2\alpha}.$$

(30)

For small values of the tempering parameter λ ≈ 0, this spectral density grows like a power law as *k* decreases, but then remains bounded as *k* 0 because of tempering. More details on ARTFIMA(*p, α, λ, q*) time series will be included in a forthcoming paper by the authors.

Kolmogorov [29] invented fractional Brownian motion as a model for turbulence in the inertial range. The main idea is that the spectral density grows like a power law as the frequency decreases, and a scaling argument leads to *H* = 1/3. This anti-persistent model has been verified in many real world experiments, for a range of frequencies. Davenport [16] modified this model in his study of wind speed at low altitudes, recognizing that the observed spectral density cannot actually diverge at low frequencies, outside the inertial range. The Davenport model for wind turbulence is *s _{t}* =

$$4800D{V}_{10}\frac{{x}^{2}}{{(1+{x}^{2})}^{4\u22153}}$$

(31)

where *V*_{10} is the mean velocity (m/sec) at an altitude of 10 meters, *D* is the corresponding drag coefficient, and *x* = 1200*k*/*V*_{10} where *k* is the wave number, see also Li and Kareem [32].

The spectral density of a stochastic process is the Fourier transform of its covariance function. The increment process *X _{t}* =

$$h\left(k\right)=\frac{C{\mid {e}^{-ik}-1\mid}^{2}}{{({\lambda}^{2}+{k}^{2})}^{H+1\u22152}}$$

(32)

for all real *k*, where *C* > 0 is some constant [49, Remark 4.3]. Since *e ^{−ik}* − ~ −

$$h\left(k\right)\approx C\frac{{k}^{2}}{{({\lambda}^{2}+{k}^{2})}^{H+1\u22152}},$$

which shows that the stochastic process *X _{t}* has the Davenport spectrum with

Figure 4 plots the spectral density of this tempered fractional Gaussian noise, along with the corresponding fractional Gaussian noise (the case λ = 0). Here we take λ = 0.01 to represent an average wind speed of 12 m/sec. The low frequency range corresponds to large eddy production. In the high frequency range, energy is dissipated via shear stress. In the intermediate (inertial) range, energy is transferred from smaller to larger eddies, and this is where the Kolmogorov model pertains. In that model, the spectral density is proportional to |*k*|^{1−2H} for moderate frequencies. For more details, see for example Beaupuits et al. [6].

In the inertial range, it can be seen from Figure 4 that the spectrum of tempered fractional Gaussian noise is almost indistinguishable from that of (untempered) fractional Gaussian noise. However, there is a significant difference at low frequencies. The spectral density of fractional Gaussian noise diverges to infinity, which is not physical. However, the spectral density of tempered fractional Gaussian noise remains bounded, providing a more reasonable model for turbulence, consistent with the Davenport spectrum.

The Davenport model for wind turbulence is widely used to design windmill farms, offshore oil platforms, and large antennae [6, 16, 26, 55], but to date the model provides only a spectrum for wind speed data. Tempered fractional Gaussian noise, the increments of a tempered fractional Brownian motion, provides a stochastic process model in the time domain, which can also be useful in applications. For example, it can be used to simulate realistic wind speed data. It may also be useful to apply tempered fractional Brownian motion with *H* ≈ 1/3 as a model for turbulence in a more general setting. In the inertial range, it agrees well with the original Kolmogorov model (fractional Brownian motion), but it would seem to provide a more reasonable physical model at low frequencies.

The ARTFIMA process from Section 6 is the discrete time analogue of a tempered fractional Gaussian noise. Hence it could also provide a useful model for turbulence, when applied for example to a time series of turbulent velocities. From (30) we see that the ARTFIMA(0*, α, λ,* 0) spectral density is proportional to |1 − *e ^{−}^{(λ+ik)}|^{−2α}* and so, for small values of the tempering parameter λ ≈ 0, the spectral density is approximately |

Fractional calculus is a very useful analytical toolbox, with many diverse applications to all areas of science and engineering. This paper is promoting a new variation called tempered fractional calculus, as a more flexible alternative with considerable promise for practical applications. A fractional derivative (or integral) is a (distributional) convolution with a power law. A tempered fractional derivative multiplies that power law kernel by an exponential factor. Fractional diffusion equations govern random walks with a power law probability density for the particle jumps. Multiplying that power law by an exponential tempering factor leads to a tempered fractional diffusion equation. Point source solutions to those equations are tempered stable probability densities, with semi-heavy tails that transition from power law to Gaussian. The emerging theory of tempered fractional calculus provides a sound mathematical basis for applications. Tempered fractional Brownian motion, the tempered fractional derivative or integral of a Brownian motion, extends the popular fractional Brownian motion, and provides a flexible model for semi-long range dependence that transitions from long to short range correlations. The tempered fractional difference operator can provide the basis for effective numerical schemes to solve tempered fractional differential equations. It is also useful as a model for time series analysis, where ARTFIMA(*p, α, λ, q*) models incorporate semi-long range dependence in a natural way. Tempered fractional Brownian motion provides a new stochastic process model for turbulence, consistent with the Davenport spectrum that extends the Kolmogorov theory of turbulence beyond the inertial range, and the ARTFIMA(*p, α, λ, q*) model provides a useful discrete time analogue.

This work was partially supported by NSF grant DMS-1025486, a Fujian Governmental Scholarship, and the National Natural Science Foundation of China (Grant No. 11101344).

**Publisher's Disclaimer: **This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

MARK M. MEERSCHAERT, Department of Statistics and Probability, Michigan State University, East Lansing MI 48823 ; Email: ude.usm.tts@debucm URL: http://www.stt.msu.edu/users/mcubed/

FARZAD SABZIKAR, Department of Statistics and Probability, Michigan State University, East Lansing MI 48823 ; Email: ude.usm.tts@2akizbas.

JINGHUA CHEN, School of Sciences, Jimei University, Xiamen, China, 361500 ; Email: moc.361@zdzdhjc..

1. Baeumer B, Meerschaert MM, Benson DA, Wheatcraft SW. Subordinated advection-dispersion equation for contaminant transport. Water Resour. Res. 2001;37:1543–1550.

2. Baeumer B, Kovács M, Meerschaert MM. Fractional reproduction-dispersal equations and heavy tail dispersal kernels. Bull. Math. Biol. 2007;69:2281–2297. [PubMed]

3. Baeumer B, Meerschaert MM. Tempered stable Lévy motion and transient super-diffusion. J. Comput. Appl. Math. 2010;233:243–2448.

4. Barkai E, Garini Y, Metzler R. Strange kinetics of single molecules in living cells. Phys. Today. 2012;65:29.

5. Barndorff-Nielsen OE. Processes of normal inverse Gaussian type. Finance Stoch. 1998;2:41–68.

6. Pérez Beaupuits JP, Otárola A, Rantakyrö FT, Rivera RC, Radford SJE, Nyman L-Å. Analysis of wind data gathered at Chajnantor. ALMA Memo. 2004;497

7. Benson D, Wheatcraft S, Meerschaert M. Application of a fractional advection-dispersion equation. Water Resour. Res. 2000;36:1403–1412.

8. Benson D, Schumer R, Meerschaert M, Wheatcraft S. Fractional dispersion, Lévy motions, and the MADE tracer tests. Transport in Porous Media. 2001;42:211–240.

9. Brockwell PJ, Davis RA. Time Series: Theory and Methods. 2nd Ed. Springer-Verlag; New York: 1991.

10. Carr P, Geman H, Madan DB, Yor M. The fine structure of asset returns: An empirical investigation. J. Business. 2002;75:303–325.

11. Carr P, Geman H, Madan DB, Yor M. Stochastic volatility for Lévy processes. Math. Finance. 2003;13:345–382.

12. Cartea Á, del-Castillo-Negrete D. Fluid limit of the continuous-time random walk with general Lévy jump distribution functions. Phys. Rev. E. 2007;76:041105. [PubMed]

13. Chakrabarty A, Meerschaert MM. Tempered stable laws as random walk limits. Statist. Probab. Lett. 2011;81:989–997.

14. Cohen S, Rosiński J. Gaussian approximation of multivariate Lévy processes with applications to simulation of tempered stable processes. Bernoulli. 2007;13:195–210.

15. Cushman JH, Ginn TR. Fractional advection-dispersion equation: A classical mass balance with convolution-Fickian flux. Water Resour. Res. 2000;36:3763–3766.

16. Davenport AG. The spectrum of horizontal gustiness near the ground in high winds. Q. J. Royal Meteor. Soc. 1961;87:194–211.

17. Deng Z-Q, Bengtsson L, Singh VP. Parameter estimation for fractional dispersion model for rivers. Environ. Fluid Mech. 2006;6:451–475.

18. Diethelm K, Ford NJ, Freed AD, Luchko Yu. Algorithms for the fractional calculus: A selection of numerical methods. Comput. Methods Appl. Mech. Engrg. 2005;194:743–773.

19. Du Q, Gunzburger M, Lehoucq RB, Zhou K. Analysis and approximation of nonlocal diffusion problems with volume constraints. SIAM Review. 2012;54:667–696.

20. Einstein A. On the movement of small particles suspended in a stationary liquid demanded by the molecular kinetic theory of heat. Ann. Phys. 1905;17:549–560.

21. Ervin VJ, Roop JP. Variational solution of fractional advection dispersion equations on bounded domains in Rd. Numer. Meth. P.D.E. 2007;23:256–281.

22. Fedotov S, Iomin A. Migration and proliferation dichotomy in tumor-cell invasion. Phys Rev Lett. 2007;98:8101. [PubMed]

23. Gorenflo R, Mainardi F, Scalas E, Raberto M. Trends Math. Birkhäuser; Basel: 2001. Fractional calculus and continuous-time finance. III. The diffusion limit. Mathematical finance (Konstanz, 2000) pp. 171–180.

24. Herrmann R. Fractional Calculus an Introduction for Physicists. World Scientific; Singapore: 2011.

25. Hilfer R. Applications of Fractional Calculus in Physics. World Scientific; 2000.

26. Jang J-J, Jyh-Shinn G. Analysis of maximum wind force for offshore structure design. J. Marine Sci. Tech. 1999;7:43–51.

27. Jeon J-H, Martinez-Seara Monne H, Javanainen M, Metzler R. Anomalous diffusion of phospholipids and cholesterols in a lipid bilayer and its origins. Phys. Rev. Lett. 2012;109:188103. [PubMed]

28. Jurlewicz A, Wy lomańska A, Żebrowski P. Coupled continuous-time random walk approach to the Rachev-Rüschendorf model for financial data. Physica A. 2009;388:407–418.

29. Kolmogorov AN. Wiener spiral and some other interesting curves in Hilbert space. Dokl. Akad. Nauk SSSR. 1940;26:115–118.

30. Koponen I. Analytic approach to the problem of convergence of truncated Lévy flights towards the Gaussian stochastic process. Phys. Rev. E. 1995;52:1197–1199. [PubMed]

31. Kuo HH. White noise distribution theory. CRC Press; Boca Raton, Florida: 1996.

32. Li Y, Kareem A. ARMA systems in wind engineering. Probabilistic Engin. Mech. 1990;5:49–59.

33. Liu F, Ahn V, Turner I, Zhuang P. Numerical simulation for solute transport in fractal porous media. ANZIAM J. 2004;45:C461–C473.

34. Liu F, Ahn V, Turner I. Numerical solution of the space fractional Fokker-Planck equation. J. Comput. Appl. Math. 2004;166:209–219.

35. Lynch VE, Carreras BA, del-Castillo-Negrete D, Ferreira-Mejias KM, Hicks HR. Numerical methods for the solution of partial differential equations of fractional order. J. Comput. Phys. 2003;192:406–421.

36. Magdziarz M, Weron A, Weron K. Fractional Fokker-Planck dynamics: Stochastic representation and computer simulation. Phys. Rev. E. 2007;75:016708. [PubMed]

37. Magin RL. Fractional Calculus in Bioengineering. Begell House; 2006.

38. Mainardi F, Raberto M, Gorenflo R, Scalas E. Fractional Calculus and continuous-time finance II: the waiting-time distribution. Physica A. 2000;287:468–481.

39. Mainardi F. Fractional Calculus and Waves in Linear Viscoelasticity: An Introduction to Mathematical Models. World Scientific; 2010.

40. Mandelbrot B, Van Ness J. Fractional Brownian motion, fractional noises and applications. SIAM Review. 1968;10:422–437.

41. Mantegna RN, Stanley HE. Stochastic process with ultraslow convergence to a Gaussian: The truncated Lévy flight. Phys. Rev. Lett. 1994;73:2946–2949. [PubMed]

42. Meerschaert MM, Scheffler H-P. Limit theorems for continuous time random walks with infinite mean waiting times. J. Applied Probab. 2004;41:623–638.

43. Meerschaert MM, Tadjeran C. Finite difference approximations for fractional advection-dispersion flow equations. J. Comput. Appl. Math. 2004;172:65–77.

44. Meerschaert MM, Scalas E. Coupled continuous time random walks in finance. Physica A. 2006;370:114–118.

45. Meerschaert MM, Scheffler H-P. Triangular array limits for continuous time random walks. Stoch. Proc. Appl. 2008;118:1606–1633.

46. Meerschaert MM, Zhang Y, Baeumer B. Tempered anomalous diffusions in heterogeneous systems. Geophys. Res. Lett. 2008;35:L17403–L17407.

47. Meerschaert MM. In: Fractional calculus, anomalous diffusion, and probability, Fractional Dynamics. Metzler R, Klafter J, editors. World Scientific; Singapore: 2012. pp. 265–284.

48. Meerschaert MM, Sikorskii A. Stochastic Models for Fractional Calculus. De Gruyter; Berlin: 2012.

49. Meerschaert MM, Sabzikar F. Tempered fractional Brownian motion. Statist. Prob. Lett. 2013;83:2269–2275.

50. Meerschaert MM, Dogan M, Van Dam RL, Hyndman DW, Benson DA. Hydraulic conductivity fields: Gaussian or not? Water Resour. Res. 2013;49:4730–4737. [PMC free article] [PubMed]

51. Meerschaert MM, Sabzikar F. Stochastic integration for tempered fractional Brownian motion. Stoch. Proc. Appl. 2014;124:2363–2387. [PMC free article] [PubMed]

52. Metzler R, Klafter J. The random walk's guide to anomalous diffusion: A fractional dynamics approach. Phys. Rep. 2000;339:1–77.

53. Metzler R, Klafter J. The restaurant at the end of the random walk: recent developments in the description of anomalous transport by fractional dynamics. J. Physics A. 2004;37:R161–R208.

54. Miller K, Ross B. An Introduction to the Fractional Calculus and Fractional Differential Equations. Wiley and Sons; New York: 1993.

55. Norton DJ. Mobile offshore platform wind loads. Proc. 13th Offshore Techn. Conf., OTC 4123. 1981;4:77–88.

56. Pipiras V, Taqqu M. Integration questions related to fractional Brownian motion. Prob. Th. Rel. Fields. 2000;118:251–291.

57. Piryatinska A, Saichev AI, Woyczynski WA. Models of anomalous diffusion: Subdiffusive and superdiffusive cases. Physica A. 2005;349:375–424.

58. Ramos-Fernandez G, Mateos JL, Miramontes O, Cocho G, Larralde H, Ayala-Orozco B. Lévy walk patterns in the foraging movements of spider monkeys (Ateles geoffroyi). Behav. Ecol. Sociobiol. 2003;55:223–230.

59. Roop JP. Computational aspects of FEM approximation of fractional advection dispersion equations on bounded domains in R^{2}. J. Comput. Appl. Math. 2006;193:243–268.

60. Rosiński J. Tempering stable processes. Stoch. Proc. Appl. 2007;117:677–707.

61. Samko SG, Kilbas AA, Marichev OI. Fractional Integrals and Derivatives. Gordon and Breach. 1993

62. Samorodnitsky G, Taqqu M. Stable non-Gaussian Random Processes. Chapman and Hall; New York: 1994.

63. Scalas E, Gorenflo R, Mainardi F. Fractional calculus and continuous-time finance. Physica A. 2000;284:376–384.

64. Scalas E. The Complex Networks of Economic Interactions. Springer; Berlin: 2006. Five years of continuous-time random walks in econophysics. pp. 3–16.

65. Schumer R, Benson DA, Meerschaert MM, Wheatcraft SW. Eulerian derivation of the fractional advection-dispersion equation. J. Contaminant Hydrol. 2001;48:69–88. [PubMed]

66. Sokolov IM, Klafter J. Anomalous diffusion spreads its wings. Physics World. 2005;18:29–32.

67. Tadjeran C, Meerschaert MM, Scheffler H-P. A second order accurate numerical approximation for the fractional diffusion equation. J. Comput. Phys. 2006;213:205–213.

68. Tarasov VE. Fractional vector calculus and fractional Maxwell's equations. Ann. Phys. 2008;323:2756–2778.

69. Whitt W. Stochastic-Process Limits. Springer; New York: 2002.

70. Yosida K. Functional Analysis. Sixth edition Springer; 1980.

71. Yuste SB, Acedo L. An explicit finite difference method and a new von Neumann type stability analysis for fractional diffusion equations. SIAM J. Numer. Anal. 2005;42:1862–1874.

72. Zaslavsky G. Fractional kinetic equation for Hamiltonian chaos. Physica D. 1994;76:110–122.

73. Zhang Y, Meerschaert MM. Gaussian setting time for solute transport in fluvial systems. Water Resour. Res. 2011;47:W08601.

74. Zhang Y, Meerschaert MM, Packman AI. Linking fluvial bed sediment transport across scales. Geophys. Res. Lett. 2012;39:L20404.

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |