Stochastic Modeling

class: center, middle
 
.title[Global observations of long-lived vortices in Earth’s oceans]
 
.author[Jonathan Lilly]
.institution[Theiss Research, La Jolla, California]
 
.date[May 12, 2021]
 
.note[Created with [{Liminal}](https://github.com/jonathanlilly/liminal) using [{Remark.js}](http://remarkjs.com/) + [{Markdown}](https://github.com/adam-p/markdown-here/wiki/Markdown-Cheatsheet) + [{KaTeX}](https://katex.org)]

---
class: center
##Some Long-Lived Planetary Vortices

[{Voyager 1 time lapse by Björn Jónsson and Ian Regan}](https://lightsinthedark.com/2010/11/13/sweet-sixteen-jupiter-in-motion)

---
class: center
##Long-Lived Vortices in Earth’s Oceans

Sea surface height slope magnitude from a global simulation by H. Simmons, University of Alaska Fairbanks.

---
class: center
##A Case Study from the Nordic Seas

A 15 km radius, apparently permanent vortex off the Norwegian coast is  maintained through vortex merger dynamics.

Nondimensionally the strongest such feature on Earth.

Many open questions remain regarding the dynamics and climalogical impacts of such vortex structures.

[{Trodahl, Isachsen, Nilsson, Lilly, and Kristensen (2020)}](https://doi.org/10.1175/JPO-D-20-0029.1)      
[{Bosse, Fer, Lilly, and Søiland (2019)}](https://doi.org/10.1038/s41598-019-49599-8)

---
class: center, middle
# Part I:  Eddy Detection in Lagrangian Data

---
class: center
## High-Tech Observations of the Ocean

---
class: center
## The Global Surface Drifter Dataset

Note: non-uniform sampling distribution, temporal variation,
        superposition of spatial scales

---
class: center
## Mean Surface Current Speed

Formed by binning in latitude and longitude, then averaging.
Easy to compute maps of low-order statistics: mean, variance, etc.

Clearly does not capture full richness of dataset.

What else can be done?

---
class: left
## Properties of the Surface Drifter Dataset
--

1.  Non-homogeneously distributed in space
--

2.  Non-uniformly sampled in time
--

3.  Irregular length or duration
--

4.  Non-stationary statistics
--

5.  Superposition of various physical processes
--

## Processes of Interest
--

1.  Low-frequency diffusive behavior
--

2.  Intermediate-frequency coherent eddy structures
--

3.  Higher-frequency wave motions
--

4.  Small-scale roughness
--

##Strategy
--

Construct particular analysis methods aimed at isolating various physical processes.  Here, focus on eddies (#2) and dispersion (#1).
--

Initially, work within simplified dynamics of numerical models.

---
class: center

###  A Numerical Simulation of an Unstable Current

A highly idealized version of the Gulf Stream.
Note formation of vortices or “coherent eddies”.

---
class: center

###  Identifying Vortices from Particle Trajectories

Vortices (loops) are identified using only particles (dots).

---
class: center

###  Extraction of Time-Varying Vortex Currents

Vortices can be identified, studied, and distinguished from waves.
From Lilly, Scott, and Olhede (2011).

---
class: left
##Elements of Lagrangian Vortex Extraction

This analyis method consists of four steps:

1.  *Extracting* oscillatory velocity features using a technique called *wavelet ridge analysis*. .cite[Lilly, Scott, and Olhede (2011)]
2.  *Thresholding* with a suitable error quantity to remove false positives. .cite[Lilly and Olhede (2012a)]
3.  *Linking* ellipse properties to spatially-integrated properties using an extended version of Stokes' theorem.  .cite[Lilly (2018)]
3.  *Testing* statistical confidence through comparison with a null hypothesis.  .cite[Lilly and Pérez-Brunius (2021b)]

---
class: left
## A Model for a Trajectory with Eddies

We conceptualize the motion of a particle trapped in an eddy as a time-varying ellipse.

<div>\[z_o(t)=x_o(t)+\mathrm{i}y_o(t)=\mathrm{e}^{\mathrm{i}\theta(t)}\left[a(t)\cos\phi(t) + \mathrm{i} b(t)\sin \phi(t)\right]\]</div>

We construct a type of  “filter” that can isolate such signals.

This lets us approximately split a time series $z(t)=z\_o(t)+z\_\epsilon(t)$ into an oscillatory $z\_o(t)$ and stochastic $z\_\epsilon(t)$ portion.

---
class: left
##Recovering Modulated Oscillations
The wavelet transform is a tool for recovering or estimating the properties of modulated oscillation immersed in background noise.

The idea is to project the time series onto an oscillatory test signal, or wavelet, to find a “best fit” frequencies at each moment. The fits are then chained together into a continuous curve called a ridge.

For a vector-valued signal `$\mathbf{x}_o(t)=\begin{bmatrix}x_1(t) & x_2(t) & \cdots & x_N(t)\end{bmatrix}^T$` in noise `$\mathbf{x}_\epsilon(t)$`, `$\mathbf{x}(t)=\mathbf{x}_o(t)+\mathbf{x}_\epsilon(t)$`,
define the wavelet transform
`\[\mathbf{w}_\mathbf{x}(t ,s) \equiv \int_{-\infty}^{\infty} \frac{1}{s}
\psi^*\left(\frac{\tau-t}{s}\right)\,\mathbf{x}(\tau)\,\mathrm{d} \tau\]`
(also a vector) and then find the wavelet ridges `$s(t)$` from
`\[\frac{\partial}{\partial s}\, \left\|\mathbf{w}_\mathbf{x}(t ,s)\right\| = 0,\quad\quad
\frac{\partial^2}{\partial s^ 2}\, \left\|\mathbf{w}_\mathbf{x}(t ,s)\right\| < 0.\]`
The oscillation is estimated simply by `$\widehat{\mathbf{x}_o}(t)\equiv\Re\left\{\mathbf{w}_\mathbf{x}(t,s(t))\right\}.$`

.cite[See Delprat et al. (1992), Lilly and Olhede (2009, 2010, 2012a). ]

---
class: center
##Example of Vortex Extraction
<img style="width:100%" src="../figures/bandwidth_example.png">

The ridge is the curve made by tracing out the maximum modulus of the wavelet transform.
---
class: center
##Example of Vortex Extraction
<img style="width:100%" src="../figures/bandwidth_looper.png">

The transform value along this curve is an estimate of the oscillatory signal component, visualized as ellipses.

Ridge analysis separates modulated elliptical signal from low-frequency meandering and higher-frequency variability.

Physically the ellipses quantify the properties of oceanic vortices.

---
class: left
##Analysis of Modulated Oscillations

It often happens that one wishes to analyze a signal of the form

`\[x(t)=a(t)\cos \phi(t) \]`

which is known as a *modulated oscillation* or *AM/FM* or *quasi-periodic* signal.

An extremely powerful method, the *analytic signal* lets an *instantaneous* amplitude, frequency, and bandwidth be uniquely associated with a given modulated oscillation.

These ‘instantaneous moments’ decompose the frequency-domain moments of the signal's Fourier spectrum across time.

For example, the power-weighted time average of the instantaneous frequency is the mean frequency of the spectrum.

This theory is primarily due to .cite[Gabor (1946)], with important contributions from .cite[Bedrosian (1963)], .cite[Vakman (1996)], .cite[Picinbono (1997)], and others.

---
class: left
##The Analytic Signal Method

The instantenous moments are defined by using the analytic signal
`\[x_+(t)\equiv x(t) +\mathrm{i}\mathcal{H} x(t) \equiv x(t)+\mathrm{i}\frac{1}{\pi}\,\int_{-\infty}^\infty\frac{x(u)}{t-u}\,\mathrm{d} u\]`
with `$\mathcal{H}$` being the Hilbert transform.

The amplitude is uniquely defined as `$a(t)=|x_+(t)|$` and frequency as `$\omega(t) = \frac{d}{dt}\arg x_+(t)$`.

**Uniqueness:** If you wish to assign the instantaneous moments using a linear filter, and if you wish your method to recover the correct results for a strictly periodic signal `$x(t)=a_o\cos(\omega_o t)$`—a property called *harmonic correspondence*— there is no alternative to using the analytic signal .cite[(Vakman, 1996)].

**Recoverability:**  With `$x(t)$` generated by `$x(t)=a(t)\cos \phi(t)$`, the generating amplitude `$a(t)$` and phase `$\phi(t)$` are exactly recoverable from `$x(t)$` provided `$a(t)$` is supported on strictly lower frequencies than is `$e^{i\phi(t)}$` .cite[(Bedrosian, 1963)].

---
class: left
##Generalization to 2D and 3D

Instantaneous moments generalize to multivariate oscillations.

A bivariate modulated oscillation `$\mathbf{x}=[x(t)\,\, y(t)]^T$` is the same as a time-varying ellipse .cite[(Lilly and Olhede, 2010a)].

A trivariate modulated oscillation `$\mathbf{x}=[x(t)\,\, y(t)\,\, z(t)]^T$` is equivalent to a modulated ellipse within a rotating, tilting plane .cite[(Lilly, 2011)].

The time-varying ellipse properties may be recovered from the vector-valued analytic signal `$\mathbf{x}_+(t)\equiv\mathbf{x}(t)+\mathrm{i}\mathcal{H}\mathbf{x}(t)$`.

---
class: left
##Recoverability Condition in 2D

Q.  Under what conditions are the time-varying properties of an ellipse exactly recoverable from the trajectory of one particle?

A.  When the ellipse geometry evolves more slowly than circulation of particles around the periphery.

Ellipse parameters for top half are exactly recoverable.

---
class: left
##Modulated Oscillations in Noise

With noise present, we write `$\mathbf{x}(t)=\mathbf{x}_o(t)+\mathbf{x}_\epsilon(t)$`, where `$\mathbf{x}_o(t)$` is an oscillation and `$\mathbf{x}_\epsilon(t)$` is noise.

In this case, applying the analytic signal method directly fails.

The wavelet transform can be used as basis for recovering or estimating the properties of modulated oscillation immersed in background noise.

Note that the analytic wavelet transform implicitly combines filtering, and creation of the analytic signal, into one step.

*Wavelet ridge analysis*, due to .cite[Delprat (1992)] and investigated further by .cite[Mallat (1999)],  allows the recovery of the oscillation.

This method lets us strip off the noise, resulting in an estimate of the analytic part of the oscillatory portion of the signal, `$\mathbf{x}_o(t)$`.

---
class: left
##A Bias Estimate

The wavelet ridge method has two sources of error, *variance* arising form the noise, and *bias* arising from the degree of modulation of the oscillatory signal.

An internal estimate of the bias, formed for an unknown modulated oscillation from the wavelet transform itself, was derived by .cite[Lilly and Olhede (2010b)].

The first-order departure from a pure sinusoid is called the *bandwidth*, which has terms like `$a'(t)/a(t)$` (amplitude modulation).

The *second-order* departure has terms like `$a''(t)/a(t)$`, which we call the *curvature*.  This turns out to control the bias of the estimate.

Wavelet ridge analysis was generalized to handle  modulated multivariate oscillations by .cite[Lilly and Olhede (2010a)].

A bias estimate for the multivariate case was derived by .cite[Lilly and Olhede (2012)], and is a multivariate analogue of the curvature.

Thresholding on the bias is critical for removing false positives!

---
class: center
##Lagrangian Analysis of a Nonlinear Eddy

An eddy at 24°N tracked by 512×256=131,072 particles, by J. Early.
Color is estimated enclosed vorticity in inferred ellipses.
---
class: center
##Estimated vs. Actual Enclosed Vorticity

Log10 histogram of estimated vs. true enclosed vorticity.

---
class: center
## Application to the Global Drifter Dataset

Apply multivariate ridge detection algorithm to 20K time series.

---
class: center
###  Ridge Detection in the Global Drifter Dataset

<img style="width:90%" src="../figures/alleddies.png"> 
Red = cyclonic rotation, blue = anticyclonic.

---
class: center
###  Removing False Positives

After comparison with a null hypothesis of red noise 
.cite[see Lilly and Pérez-Brunius (2021b)]

---
class: center, middle
# Part II: Stochastic Modelling

---
class: center
##Ocean Turbulence

An idealized model of oceanic turbulence by Jeffrey Early.
Shading is speed with particle trajectories shown in color.

---
class: center
##Modeling Lagrangian Trajectories

A simple stochastic model for trajectories in oceanic turbulence is the *Matérn process* or damped Fractional Brownian motion.

This three-parameter model lets us match (i) the kinetic energy, (ii) the dispersion, and (iii) the degree of small-scale roughness.

---
class: left
##The Matérn Process

A relatively little-known random process called the Matérn process .cite[(Matérn, 1960)] has a spectrum given by
`\[S_{zz}(\omega) = \frac{A^2}{(\omega^2+\lambda^2)^\alpha}\]`
where `$A$` sets the energy level, `$\alpha$` sets the high-frequency slope, and `$\lambda$` determines a low-frequency transition to a constant spectral value.

--
Recently, we investigated this process in detail, and showed that it is an excellent model for trajectories in ocean turbulence.  We speculate it will prove suitable for many other processes as well.

We showed the Matérn process is a generalization of *fractional Brownian motion* .cite[(Mandelbrot and Van Ness, 1968)] to include a damping term, `$\lambda$`, that is the essential ingredient for modeling dispersion.

.cite[Lilly, Sykulski, Early, and Olhede (2017). Fractional Brownian motion, the Matérn process, and stochastic modeling of turbulent dispersion. <a href="./papers/lilly17-npg.pdf">[3.2 Mb pdf] </a>]
---
class: left
## An Hierarchy of Stochastic Models

Let `$z(t)=u(t)+iv(t)$` represent a complex-valued velocity.

Then we have the following correspondences between stochastic differential equations (SDEs) and velocity spectra `$S_{zz}(\omega)$`.

Brownian motion:

`\[\frac{dz}{d t}  = A\,\frac{d w}{d t} \quad\longrightarrow\quad S_{zz}(\omega) = \frac{A^2}{\omega^2}\]`

Fractional Brownian motion (fBm):

`\[\frac{d^\alpha z}{d t^\alpha}  = A\,\frac{d w}{d t} \quad\longrightarrow\quad S_{zz}(\omega) = \frac{A^2}{\omega^{2\alpha}}\]`

Damped Fractional Brownian motion a.k.a. the Matérn process:

`\[\left[\frac{d}{d t} + \lambda \right]^\alpha z = A\,\frac{d w}{d t} \quad\longrightarrow\quad S_{zz}(\omega) = \frac{A^2}{(\omega^2+\lambda^2)^\alpha}\]`

(These should be written as the corresponding *integral* equations.)

---
class: left
## Another Hierarchy of Stochastic Models

Damped Brownian motion a.k.a. Ornstein-Uhlenbeck (OU) process:
`\[\frac{dz}{d t} + \lambda z = A\,\frac{d w}{d t} \label{maternsde}\quad\longrightarrow\quad S_{zz}(\omega) = \frac{A^2}{\omega^2+\lambda^2}\]`

Damped Brownian motion + spin a.k.a. complex OU (cOU) process:
`\[\frac{dz}{d t} - i\Omega z+ \lambda z = A\,\frac{d w}{d t} \quad\longrightarrow\quad S_{zz}(\omega) = \frac{A^2}{(\omega-\Omega)^2+\lambda^2}\]`

Damped fBm + spin a.k.a. Matérn + spin a.k.a. fractional cOU:
`\[\left[\frac{d}{d t} - i\Omega+ \lambda \right]^\alpha z = A\,\frac{d w}{d t} \quad\longrightarrow\quad S_{zz}(\omega) = \frac{A^2}{[(\omega-\Omega)^2+\lambda^2]^\alpha}\]`

Matérn + spin generalizes cOU to fractional orders, thus unifying OU+cOU+fBm into a single larger family.
---
class: left
##Slope and Roughness

The `$\alpha$` parameter sets the slope in fractional Brownian motion, or at high frequencies in the Matérn process, as `$S(\omega)\sim \omega^{-2\alpha}$`.

It also sets the degree of small-scale roughness in the trajectories.
In the above, we go from slopes of `$\omega^{-3}$` at left to `$\omega^{-1}$` at right.

Steeper spectral slope = more smooth.

Spectral slope can be formally linked to a mathematical measure of roughness called the *fractal dimension*, though in my opinion, fractal dimension seems to have little practical value.
---
class: center
##Slope and Self-Similarity

The `$\alpha$` parameter also sets *aspect ratio* of statistical self-similarity.

---
class: left
##Diffusivity is the Spectrum's Value at Zero

The velocity autocovariance and spectrum are defined as

`\[R_{zz} (\tau)\equiv\left\langle z(t+\tau)\,z^*(t)\right\rangle = \frac{1}{2\pi}\int_{-\infty}^\infty e^{i \omega \tau} S_{zz}(\omega) \, d\omega\]`

The corresponding diffusivity is given by
`\[\kappa \equiv \lim_{t\longrightarrow\infty}
 \frac{1}{4} \frac{d}{d t} \left\langle x^2(t) + y^2(t)\right\rangle
 =\frac{1}{4} \int_{-\infty}^\infty R_{zz}(\tau)\, d \tau =\frac{1}{4}S_{zz}(0)\]`

with the last equality following from `$\int_{-\infty}^\infty e^{-i \omega \tau} R_{zz}(\tau) \, d\tau=  S_{zz}(\omega)$`.

Matérn:  `$ \,\,\,\,S_{zz}(\omega) = \frac{A^2}{(\omega^2+\lambda^2)^\alpha} \quad\longrightarrow\quad \kappa =  \frac{1}{4}\frac{A^2}{\lambda^{2\alpha}}$`

fBm: `$ \quad\quad S_{zz}(\omega)= \frac{A^2}{\omega^{2\alpha}} \quad\quad\,\,\longrightarrow\quad \kappa =  \infty$`

Infinity means the diffusivity increases without bound.  Thus the Matérn, unlike fBm, can model diffusive processes.

---
class: center
##Comparison of Models

Top left: numerical simulations.  Top right: Matérn process.
Bottom left: white noise.  Bottom right: red noise (Brownian).

---
class: left
##Argument for Generality

We conjecture that the Matérn may be of general usefulness:

“Systems are often characterized by a pressure to grow—represented by a forcing—together with some drag or resistance on that growth, represented by a damping. After a sufficiently long time, the forcing and the damping equilibrate and one reaches a bounded state. This leads to the speculation that many time series that are well described as fBm over relatively short timescales may be better matched by the Matérn process over longer timescales. More generally, the Matérn process adds a third parameter (damping) to the two parameters (amplitude together with spectral slope or the Hurst parameter) of fBm, thus permitting a wider range of spectral forms to be accommodated. It is therefore reasonable to think that the Matérn process could be of broad interest in many areas in which fBm has already proven itself useful.”

.cite[Lilly, Sykulski, Early, and Olhede (2017)]

---
class: center,
## The Global Surface Drifter Dataset

Now we will apply the Matérn process to the global drifter dataset.

---
class: center,
## Spectra of Drifter Velocities

We see a low-frequency maximium and tidal maxima—all vertical lines—and a maximum along a curving line, the inertial frequency.

---
class: center,
##A Conceptual Model for the Spectrum

\begin{equation}
S(\omega) = \overset{\mathrm{background}}{\overbrace{\frac{A_b^2}{\left[\omega^2+\lambda_b^2\right]^\alpha}}}
+\overset{\mathrm{inertial\,\, oscillations}}{\overbrace{\frac{A_o^2}{\left(\omega-f_o\right)^2+\lambda_o^2}}}+
\overset{\mathrm{semidiurnal\,\,tide}}{\overbrace{\frac{A_s^2}{\left(\omega-f_s\right)^2+\lambda_s^2}}}
\end{equation}

This is just several versions of the Matérn spectrum added together.

---
class: left
## Details of the Stochastic Modeling

We have developed a method for stochastic modeling of trajectories in ocean turbulence, and inferring parameters from large datasets.

This is an example of *parametric spectral analysis*:  rather than trying to estimate the spectrum directly, one begins with a *conceptual model*, controlled by some parameters, and then finds the choice of parameters which provide the best match to the data.

The full method consists of several ingredients:

*  Create an appropriate stochastic model for particle trajectories.
.cite[Sykulski, Olhede, Lilly, and Danioux (2015)]
* Employ the  Matérn  process as a building block.  .cite[Lilly, Sykulski, Early, and Olhede (2017)]
* Parameter estimation is best done in the frequency domain.
.cite[Whittle (1953)]
* Parameter estimation must be adjusted to handle *complex-valued* time series. .cite[Sykulski, Olhede, Lilly, and Early (2017)]
* Parameter estimation must be corrected for subtle bias due to small sample size effects.  .cite[Sykulski, Olhede, Guillaumin, Lilly, and Early (2018), in review.]

---
class: center
## Spectra in Latitude Bands
<img style="width:75%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/drifterspectralinestransparency.png">

Averages in `$\pm 5 ^\circ $` bins.  Frequency is nondimensionalized as `$\omega/f_o$`.
Inertial peak broadens equatorward.

---
class: center
## Comparison with Spectral Fit

In general, spectral model compares very well with estimated spectra.  Inertial and semidiurnal peaks are well captured.

---
class: center
## Fit to Background vs. Cyclonic Side

“Background” fit to anticyclonic (inertial side), matches well the spectrum on the cyclonic (anti-inertial) side, apart from tides.

---
class: center
##Application to the Global Drifter Dataset

The nondimensionalized damping parameter from a parametric fit to the Matérn, only on the anti-inertial side.

Shorter ‘memory’ in more energetic regions.
We are not aware of any theory for this.
---
class: center
##Application to the Global Drifter Dataset

The first global map of the slope parameter `$\alpha$`, with `$|\omega|^{-2\alpha}$`.
Slopes vary from `$|\omega|^{-1}$` to `$|\omega|^{-3}$`, with `$|\omega|^{-2}$` over major currents.

Does not fit the conventional wisdom of only `$|\omega|^{-2n}$` slopes,
as in e.g. .cite[Berloff and McWilliams (2002b)]

---
class: center, middle
# Part III:  Alongtrack Eddy Detection

---
class: left
## Motivation

The ridge analysis is intended for signals consisting of modulated oscillations in noise.

Another type of signal that may occur is short-duration, isolated ‘bursts’ or  ‘impulses’—that is, events which may be represented as a wavelet or the temporal integral of a wavelet.

This suggests a model for a real, univariate signal of the form

`\begin{equation}\label{signalmodel}
x(t) = \sum_{n=1}^N \Re\left\{c_n \psi\left(\frac{t-t_n}{\rho_n}\right)\right\} +x_\epsilon(t)
\end{equation}`

consisting of superpositions of amplitude scaled, stretched, and phase-shifted versions of some basis signal `$\psi(t)$`, called the *element*.

It is assumed that the different realizations of the element are sufficiently separated in time and frequency such that they do not interfere, in a way that will be made precise later.

---
class: left
## Ingredients

Element analysis consists of several steps:

1.  Choice of element
1.  Event detection
2.  Rejection of statistically insignificant events
3.  Rejection of non-isolated events

Developed in .cite[Lilly, J. M. (2017). Element analysis: a wavelet-based method for analyzing time-localized events in noisy time series. Proceedings of the Royal Society of London, Series A, 473 (2200): 20160776, 1–28. [{10.7 Mb Pdf}](http://rspa.royalsocietypublishing.org/content/473/2200/20160776)]
---
class: center
## Localized Wavelets

The method is suitable for events that themselves resemble wavelets, or the temporal integral of a wavelet (the `$\beta=0$` case).

---
class: left
## Signal Model

Our signal model, using the generalized Morse wavelets, is
`\begin{equation}\label{morseelementmodel}
x(t) = \sum_{n=1}^N \Re \left\{ c_n \psi_{\mu,\gamma}\left(\frac{t-t_n}{\rho_n}\right)\right\} +x_\epsilon(t)
\end{equation}`
where `$\mu$` is the most suitable value of the order or `$\beta$` parameter.

The `$n$`th event is then characterized by coefficient `$c_n$`, time `$t_n$`, and scale `$\rho_n$`.  We wish to *estimate* these unknown quantities.

We then take the wavelet transform with an order `$\beta$` wavelet in the same `$\gamma$` family, leading to

`\begin{equation}\label{transformofelementmodel}
w_{\beta,\gamma}(\tau,s)=\frac{1}{2}\sum_{n=1}^N c_n\int_{-\infty}^{\infty} \frac{1}{s}
\psi_{\beta,\gamma}^*\left(\frac{t-\tau}{s}\right)\psi_{\mu,\gamma}\left(\frac{t-t_n}{\rho_n}\right)\,d t+\varepsilon_{\beta,\gamma}(\tau,s).
\end{equation}`

Owing to the properties of the generalized Morse wavelets, the integral in the above equation has a closed-form expression in terms of a rescaled `$\psi_{\beta+\mu,\gamma}(t)$` wavelet.

---
class: left
## Event Detection

We then find *isolated maxima* of the transform, that is, time-scale `$(\tau,s)$` points at which
`\begin{multline}
\quad\quad\quad\quad\frac{\partial}{\partial \tau}\left|w_{\beta,\gamma}(\tau,s)\right| =\frac{\partial}{\partial s}\left|w_{\beta,\gamma}(\tau,s)\right|=0, \\\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\frac{\partial^2}{\partial \tau^2}\left|w_{\beta,\gamma}(\tau,s)\right|<0,\quad\frac{\partial^2}{\partial s^2}\left|w_{\beta,\gamma}(\tau,s)\right|<0.\label{maxconditions}
\end{multline}`

Because we have a simple expression for `$w_{\beta,\gamma}(\tau,s)$`, we can relate the event properties  `$c_n$` `$t_n$`, and `$\rho_n$` to the value of `$w_{\beta,\gamma}(\tau,s)$` at a maxima.

This lets us work backwards from the maxima points to estimates of the event properties.

---
class: center
##An Example of Events in White Noise

The events are based on the `$\psi_{2,2}(t)$` wavelet in this case.
Grey dots are all maxima, black are significant and isolated.

---
class: left
##Assessing Statistical Significance

This is the key step for large-scale applications.

The obvious approach is to simulate many noise realizations suitable for the time series of interest, and then apply the method.

This works great for one time series but not for 20,000.

Fortunately, there is a workaround.

It turns out that all we need to know is the statistics of a 5-vector, containing this point, the next point, the previous point, the point above, and the point below.

Actually, given the assumption of power-law Gaussian noise, having a spectrum of the form `$A^2/\omega^{2\alpha}$`, the covariance statistics of this 5-vector are known analytically.

This means that all one needs to do is simulation realizations of this 5-vector in order to determine the probability that a point at a given scale will be a maximum of a certain magnitude, thus establishing statistical significance.

---
class: center
##Event Distributions for the Example

Color = 5-vector realizations, gray dots = all detected events
Gray squares = values of true events
Black circles = estimated values of significant events.
Lines are significance curves for event detection rates.

---
class: left
##Assessing Isolation

This is done  through the identification of *regions of influence* around isolated maxima, related to but distinct from concentration regions.

Analytic expressions for these regions of influence may be found.

An event is rejected if it is within the region of influence (subject to a threshold) of another, larger-amplitude event.

This turns out to be mostly important for rejecting artifacts due to closely-spaced events.

---
class: center
## A Real-World Application

Measurement tracks for a sea surface height satellite.
One of the main windows into the ocean circulation.

---
class: center
## A Real-World Application

From an short altimeter track segment crossing the Labrador Sea.
---
class: left
# Conclusions

We discussed three techniques for studying coherent vortices and their impacts on Earth:

1.  Wavelet ridge analysis for studying oscillatory features
2.  Stochastic modeling
3. Element analysis for the detection of localized events

These methods offer promising avenues for research into large-scale vortex and turbulence phenomenon on Earth, and possibly other planets as well!

---
class: center, middle
# Thank you!