# Front Propagation from Radiative Sources

Theodore Steele<EMAIL_ADDRESS>Department of Applied Mathematics and
Theoretical Physics, University of Cambridge  Kinwah Wu<EMAIL_ADDRESS>Mullard Space Science Laboratory, University College London

###### Abstract

Fronts are regions of transition from one state to another in a medium. They
are present in many areas of science and applied mathematics, and modelling
them and their evolution is often an effective way of treating the underlying
phenomena responsible for them. In this paper, we propose a new approach to
modelling front propagation, which characterises the evolution of structures
surrounding radiative sources. This approach is generic and has a wide range
of applications, particularly when dealing with the propagation of phase or
state transitions in media surrounding radiation emitting objects. As an
illustration, we show an application in modelling the propagation of
ionisation fronts around early stars during the cosmological Epoch of
Reionisation (EoR) and show that the results are consistent with those of
existing equations but provide much richer sources of information.

##  I Introduction

A front is a transitional region between two states. They can be infinitesimal
in thickness or may take the form of a wide gradient. They may also be either
static or dynamic, depending upon the nature of the system and states being
modelled. The study of the dynamics of fronts is important to many areas of
sciences, in particular physics and applied mathematics, including the study
of the effects of ionising radiation, of chemical catalysis by radiation, of
the spread of fires, and of pressure shock waves.

A number of algorithms have been developed to model the propagation of fronts
in various scenarios, many of which rely upon finding solutions to the Eikonal
equation of motion. Examples of Eikonal based models are the level-set method
and its variants [1, 2, 3, 4], which model the propagation of fronts from a
set of initial conditions. A simple yet useful class of algorithms in the
restricted level-set methods are known as the fast-marching algorithms [5, 6],
in which a monotonically evolving front is defined only by its local
conditions, again through treatment of the Eikonal equation. This
simplification greatly reduces the computational cost. A variant of the fast
marching algorithms, which attempts to retain their computational advantages
over other level set methods but allow for non-monotonically evolving fronts
[7] was also developed to give the flexibility to account for non-uniform
interaction around the fronts.

However, these models would not be idea for describing the propagation of
fronts driven by linear radiation emitted from a spatially finite source.
Examples of scenarios such models could describe is the propagation of
ionisation fronts around radiative astrophysical objects such as stars [see
e.g. 8] or photochemical reactions in a medium [e.g. 9] irradiated by a light
source.

In this paper, we present a new model of front propagation which was
specifically developed to describe hydrogen ionisation around radiative
astrophysical objects, but which can be applied to any scenario in which a
source drives a front radially outwards. This model divides the space around
the object into a polar grid of cells, each of which contains a number of
particles and, over a series of time steps, models the progression of the
front as it transmits through these cells. The precision of the model
ultimately depends upon the sizes of the cells and time steps which, when the
model is applied numerically, will be limited only by the computational
resources available. A similar method was developed in [10] for the study of
astrophysical ionisation in media of homogeneous densities, but their study
differed in its aims and applications from ours.

In Sec II, we introduce the model and describe its workings in a general
scenario. In Sec. III, we apply the model to the case of a Population-III star
during the cosmological Epoch of Reionisation (EoR) and compare the results we
obtain to those from a common analytic approximation, showing that the results
are not only consistent but provide us with richer information than those pre-
existing methods. Finally, in Sec. IV, we highlight the capabilities of this
model and its potential for modelling front propagation in complex media that
researchers encounter in a variety of disciplines both within and beyond
physics.

##  II The Radiative Front Model

Our model begins with a very generic assumption: radiation is being emitted
from a point source and travelling radially away from that source. As such,
the model would be applicable to any system which features a source
transmitting radiation into a surrounding medium which could be expected to
illicit an alteration in that medium with minimal scattering.

The space around the source is divided into cells arranged into a polar grid,
with the source at the origin. Taking the example of such a grid in 3D space,
the volume of each cell is given by

$V_{i,\theta,\phi}=\frac{4}{3}\frac{\pi}{N_{\theta}N_{\phi}}\left[\;\\!(r_{i+1,\theta,\phi})^{3}-(r_{i,\theta,\phi})^{3}\;\\!\right]\
,$ (1)

where $N_{\theta}$ is the number of cells defined in the $\theta$ dimension,
$N_{\phi}$ is the number of cells defined in the $\phi$ dimension, and
$r_{i,\theta,\phi}$ is the inner radius of the cell in question. The medium
under consideration is divided into particles distributed between the cells in
accordance with the nature of the system being modelled. The number of
particles in a given cell with coordinates $\left(t,i,\theta,\phi\right)$ that
are in the initial state will be labelled $n_{\mathrm{I};t,i,\theta,\phi}$,
and the number in the affected state will be labelled
$n_{\mathrm{II};t,i,\theta,\phi}$.

Once the initial grid has been established, the radiating object is treated as
being activated. Assuming the speed of sound of its radiation is a known
function, which we label as $v$, the maximum distance that radiation could
have travelled radially from the source is computed at each time step as
$v(t-1)\Delta t$, where $\Delta t$ is the length of each time step, such that
$(t-1)\Delta t$ is the time that has elapsed by the beginning of that time
step. Every cell that is within this distance is then treated as being in
contact with the radiation.

The medium in question has an absorbativity parameter associated with it. In
the case of a physical or chemical reaction this will be a function, $f$, of
the particle cross-section of the substance and the particle number in a given
cell. We can therefore calculate the probability of each unit of radiation
affecting a given unit of medium in a given cell as

$\Pi_{t,i,\theta,\phi}=1-e^{f_{t,i,\theta,\phi}}\ ,$ (2)

such that the number of affected medium particles per time step will be given
by

$\Gamma_{t,i,\theta,\phi}=J_{t,i,\theta,\phi}\;\\!\Pi_{t,i,\theta,\phi}\ ,$
(3)

where $J_{t,i,\theta,\phi}$ is the flux of radiation entering that cell.

In order to calculate $J_{t,i,\theta,\phi}$, we will need to not only account
for the angular distribution of radiation but also for the absorption in cells
radially preceding the cell in question. We define the radial absorption
function as

$s_{t,i,\theta,\phi}=\sum_{j=1}^{i-1}f_{t,j,\theta,\phi}$ (4)

such that $s_{t,i,\theta,\phi}$ gives the fraction of radiation that reaches
the $i$th radial cell in a given angular bin. Thus, the flux of radiation
entering the cell will be given by

$J_{t,i,\theta,\phi}=J_{t,1,\theta,\phi}\ e^{s_{t,i,\theta,\phi}}\ .$ (5)

Having calculated the influx of radiation into a given cell at a given time
step, we recalculate the number of affected and unaffected particles as

$\displaystyle n_{\mathrm{I};t,i,\theta,\phi}$
$\displaystyle=n_{\mathrm{I};,t-1,i,\theta,\phi}-\Gamma_{t,i,\theta,\phi}\ ,$
(6) $\displaystyle n_{\mathrm{II};t,i,\theta,\phi}$
$\displaystyle=n_{\mathrm{II};,t-1,i,\theta,\phi}+\Gamma_{t,i,\theta,\phi}\ ,$
(7)

in a manner that the overall particle count in a given cell is preserved. Any
other effects which might affect the particle counts, such as the stochastic
processes that can reverse the radiation’s effects, can be accounted for by
adding additional terms into these equations. Having recalculated the particle
numbers in each cell at a given time step, the absorption fraction can either
be calculated again as a function of that cell’s parameters at that time step
or more simply recalculated as

$f_{t+1,i,\theta,\phi}=f_{t,i,\theta,\phi}\
\left(\;\\!\frac{n_{\mathrm{I};t,i,\theta,\phi}}{n_{\mathrm{I};0,i,\theta,\phi}}\;\\!\right)\
,$ (8)

accounting for the alteration in the number of reactive particles in the
medium at a given time step by taking their ratio with the initial number, and
the radial absorption function would then be recalculated in kind.

In effect, at every time step and every cell, we are calculating the amount of
radiation reaching that cell from the source, taking into account absorption
by preceding cells in the same angular bin, calculating the reaction rate in
the cell due to that radiation, and preparing our map for a new set of
calculations to take place at the next time step, taking into account the
ensuing changes in the medium. The mechanisms of this model is illustrated in
the flow chart in Fig. 1.

Figure 1: A flow chart showing the basic underlying mechanism of the radiative
front model. $r$ is the radial distance to the inner edge of the cell in
question, $t$ is the time step coordinate, $\Delta t$ is the real time
duration of each time step, and $v$ is the speed of the radiation causing the
front, such that the cell in question can be considered to be in contact with
the radiation in any timestep for which $v(t-1)\Delta t>r$. In the case of a
background which allows for non-trivial transport of the medium through which
the front is propagating, transport terms will be included in the cell
parameter calculations to represent this. In the event that the radiation is
travelling at less than the speed of sound in that medium, an additional
series of statements would need to be added within the loop over the map to
find which cells are within range of medium transport and recalculate those
cell’s parameters even when they are out of range of the radiation.

##  III Example Application - Ionisation in Astrophysical Media

###  III.1 Background

The earliest stars were composed almost entirely of hydrogen and helium, and
they are known as Population-III stars. [11, 12]. They were large and bright
and emitted predominantly ultraviolet radiation. They were also short lived,
with lifespans generally not more than a few million years. They are believed
to have first come into existence during the transition from what are called
the cosmological dark ages, during which the baryonic component of the
Universe was primarily composed of neutral hydrogen (HI), to the EoR, during
which that HI was predominantly converted into the ionised hydrogen (HII) that
composes most of the baryonic matter of the Universe today [13, 14], beginning
at a redshift111In cosmology, redshift is a measure of the expansion of the
Universe; at a redshift of $z=15$, objects in the Universe will have been, on
average, $1/15$th as far apart as they are now. of approximately $z=15$ [15].

An understanding of the ionisation dynamics around individual objects as
Population-III stars is essential for high precision studies of the EoR, in
particular, those being carried out by the 21cm surveys [16, 17, 18] though
there is no precise model able to describe in detail the evolution, spatial
extent and morphological structures of the fronts these sources create which
determine the tomographic features to be observed in the 21cm spectral lines
[19].

When ionising light travels through a medium containing neutral atoms, each
photon will have a probability of ionising an atom of

$\Pi=1-e^{-\tau}\ ,$ (9)

where

$\tau=\sigma d$ (10)

is the optical depth of the path traversed for interaction cross-section
$\sigma$ and the column density of the path $d=\rho_{\mathrm{HI}}l$ for HI
number density $\rho_{\mathrm{HI}}$ and distance $l$. In the case of neutral
hydrogen (HI) being converted to ionised hydrogen (HII), $\sigma=6.3\times
10^{-18}\leavevmode\nobreak\ \mathrm{cm^{-2}}$.

$\Gamma_{\mathrm{B}}=n_{\mathrm{e}}n_{\mathrm{p}}\alpha\leavevmode\nobreak\ ,$
(11)

for electron number density $n_{\mathrm{e}}$, proton number density
$n_{\mathrm{p}}$, and B-recombination coefficient [20, 21, 22]

$\alpha\approx\frac{2.6\times 10^{-13}}{(T/10^{4}{\rm K})^{4/5}}$ (12)

with units of volume per unit time, where $T$ is the average temperature of
the medium.

Due to B-recombinations, there will always be a theoretical limit to how large
a photoionised bubble can be for a given radiation source. That limit, given
under the assumption that the system remains radiating and unchanging except
for the interaction of the radiation and the medium for all of time, is given
by the Strömgren radius [23]:

$r_{\mathrm{S}}=\left(\frac{3\phi}{4\pi
n_{\mathrm{e}}n_{\mathrm{p}}\alpha}\right)^{\frac{1}{3}}\leavevmode\nobreak\
,$ (13)

where $\phi$ is the total flux of photoionising particles being emitted. This
can be trivially rearranged to see that it is defined as the radius of a
sphere around a radiative object at which the flux of radiation is equal to
the number of B-recombinations; assuming that each photon can ionise only one
particle, this will of course be the limiting radius of an ionisation bubble.
In the event that the ionising particles are able to ionise more than one atom
each, this could be corrected by simply multiplying $\phi$ by the number of
atoms each particle could ionise.

Of course, Eq. (13) is not an entirely accurate model even for the idealised
assumptions we are making; since the flux of particles will decrease with
distance from the source, the number density of electrons and protons will
also decrease, such that $r_{\mathrm{S}}$ should more realistically be phrased
as an integral over spacetime with variable number densities. Since Eq. (13)
would generally be calculated under the assumption that within the ionised
region the number densities of ionised products equals the initial number
density of neutral hydrogen and, if we were to proceed with such a more
detailed analogue of the Strömgren radius, we would expect that they would
decrease with distance, we would find that such a calculation would yield a
slightly higher result. Therefore, we would expect the radius calculated with
Eq. (13) to actually be slightly smaller than the maximum radius reached by
the ionisation front of a star. Effectively, Eq. (13) assumes a hard cutoff at
which complete ionisation gives way to complete neutrality, while failing to
take into account the existence of an intermediating front at all.

###  III.2 Implementation

We wish to model the physics described in Sec. III.1 with the model we
described in Sec. II. Specifically, we wish to model the evolution of the
ionisation front around a single Population-III star in a neutral background
medium during the EoR by treating the star as a source of radiation at the
origin of a polar coordinate system and studying the effects of the radiation
on the matter in the cell’s that compose the spherical polar map. The speed of
the radiation is the speed of light such that cells are considered to be in
contact with the source when $r_{i,\theta,\phi}\leq c(t-1)\Delta t$ where $c$
is the speed of light and we are assuming a negligible difference between the
speed of light in a vacuum and the speed of sound of the light in the medium.

The absorption function of a given cell is given by its optical depth

$\tau_{t,i,\theta,\phi}=\sigma\;\\!\rho_{\mathrm{HI};t,i,\theta,\phi}l\leavevmode\nobreak\
,$ (14)

for cross-section $\sigma$ and cell depth $l=r_{i+1}-r_{i}$. From the
definition of $\tau$, we have that $e^{-\tau}$ is the fraction of radiation
that transmits through the medium without scattering.

The radial absorption function is then given by a sum over the optical depths
of the preceding cells in the same angular bin,

$s_{t,i,\theta,\phi}=\sum_{j=1}^{i-1}\;\\!\tau_{t,j,\theta,\phi}\ .$ (15)

In order to calculate the flux of radiation reaching each cell at each time
step, we must not only account for the division of the star’s overall flux
throughout the angular components of the map, but also for the absorption that
takes place before they reach a given cell. Thus, the radiation intensity at a
given cell will be a function of the overall stellar luminosity222Luminosity
is defined as being the overall radiation flux emitted by a light emitting
source. $L$, the number of angular components $N_{\theta}$ and $N_{\phi}$, and
the radial optical depth $S_{t,i,\theta,\phi}$.

Thus, at a given cell which we label cell $i$, we have that the ionising
radiation intensity is given by

$J_{t,i,\theta,\phi}=\frac{Le^{-s_{t,i,\theta,\phi}}}{N_{\theta}N_{\phi}}\Delta
t\leavevmode\nobreak\ .$ (16)

The ionisation rate within that cell is then given by

$\Gamma_{t,i,\theta,\phi}=\frac{JN}{h\nu}(1-e^{-\tau_{t,i,\theta,\phi}})\ ,$
(17)

where $N$ is the number of atoms each photon is capable of ionising, given by
$N=\lfloor h\nu/E_{\mathrm{i}}\rfloor$ where $E_{\mathrm{i}}$ is the energy
needed to ionise a ground state atom. In the event that multiple types of
radiation are emitted, we can generalise Eq. (17) to become a sum over the
different wavelengths emitted or, in the event of a continuous distribution of
radiation wavelengths, this could be further generalised to

$\Gamma_{t,i,\theta,\phi}=\int{\rm d}j\
\frac{J_{j}N_{j}}{h\nu_{j}}\;\\!\left(1-e^{-\tau_{j,t,i,\theta,\phi}}\right)\
,$ (18)

where $N_{j}=\lfloor h\nu_{j}/E_{\mathrm{i}}\rfloor$. However, for this
simplified study we content ourselves to treat the star’s radiation as being
purely composed of $13.6\mathrm{eV}$ ultraviolet light; given the heavily
ultraviolet spectrum of a Population-III star, this is not a particularly
inaccurate treatment and will allow accurate results to a level of precision
that is satisfactory for the present study.

Given that HI has only one electron, we treat
$n_{\mathrm{e}}=n_{\mathrm{p}}=n_{\mathrm{HII}}$, such that Eq. (11) can be
implemented as

$\Gamma_{\mathrm{B};t,i,\theta,\phi}={(n_{\mathrm{HII};t,i,\theta,\phi})}^{2}\;\\!\alpha_{t,i,\theta,\phi}\;\\!V_{i,\theta,\phi}\;\\!\Delta
t\ .$ (19)

As a first approximation, we treat the temperature of all cells as being
constant, at $10^{4}\mathrm{K}$. This is an appropriate temperature of such
ionised media333The temperature will be heavily dependent on the chemical
composition of the medium; since this paper concerns itself with demonstrating
the validity of the proposed method rather than a precise study of
astrophysical ionisation, we content ourselves to use a constant temperature
of $T=10^{4}\leavevmode\nobreak\ \mathrm{K}$ in both our analytic calculations
with Eq. (13) and our implementation of the radiative front model..
Furthermore, we treat the astrophysical medium as being entirely composed of
hydrogen, such that at all times we will have $n_{\mathrm{e}}=n_{\mathrm{p}}$
and $N=\lfloor\;\\!h\nu/13.6\;\\!{\rm eV}\;\\!\rfloor$. In a detailed study
relevant to 21cm surveys, we would want to use a more detailed analysis with
variable temperatures and a medium containing helium, but for the purposes of
demonstrating the radiative front model and its application to the evolution
of photoionised fronts, these approximations are valid, particularly when one
bears in mind that the Strömgren radius given in Eq. (13) assumes a purely
hydrogen background and is itself a function of temperature, which we can also
set to a constant $10^{4}\leavevmode\nobreak\ \mathrm{K}$.

Thus, at each time step we are altering each cell’s HI and HII count by

$\displaystyle n_{\mathrm{HI};t,i,\theta,\phi}$
$\displaystyle=n_{\mathrm{HI};t-1,i,\theta,\phi}-\Gamma_{t,i,\theta,\phi}+\Gamma_{\mathrm{B};t,i,\theta,\phi}\leavevmode\nobreak\
,$ (20) $\displaystyle n_{\mathrm{HII};t,i,\theta,\phi}$
$\displaystyle=n_{\mathrm{HII};t-1,i,\theta,\phi}+\Gamma_{t,i,\theta,\phi}-\Gamma_{\mathrm{B};t,i,\theta,\phi}\leavevmode\nobreak\
,$ (21)

in a manner that naturally preserves the overall number of particles per cell,
before recalculating the optical depths as

$\displaystyle\tau_{t+1,i,\theta,\phi}$
$\displaystyle=\sigma\rho_{\mathrm{HI};t,i,\theta,\phi}l\leavevmode\nobreak\
,$ (22) $\displaystyle s_{t+1,i,\theta,\phi}$
$\displaystyle=\sum_{j=1}^{i-1}\tau_{t+1,j,\theta,\phi}\leavevmode\nobreak\ ,$
(23)

in preparation for the next time step.

###  III.3 Application

We demonstrate the model by taking the example of a $100\leavevmode\nobreak\
\mathrm{M_{\odot}}$ Population-III star in a homogeneous background composed
entirely of HI with a number density of $\rho_{\mathrm{HI}}=5\times
10^{-6}\;\\!(1+z)^{3}\mathrm{\leavevmode\nobreak\ cm^{-3}}$ [24], where the
bracketed term accounts for the reduction in density with the expansion of the
Universe. Such a star would have a lifespan of approximately
$2.7\leavevmode\nobreak\ \mathrm{Myr}$ and would emit $2.4\times
10^{51}\mathrm{\leavevmode\nobreak\ eV\leavevmode\nobreak\ s}^{-1}$ of
predominantly ultraviolet radiation [12], which for the sake of simplicity we
treat as being purely $13.6\leavevmode\nobreak\ \mathrm{eV}$ light. With these
parameters, we generated homogeneous density maps with densities equal to the
average densities at $z=15$, $z=10$, and $z=5$.

Before using the model to obtain new results, we would like to confirm that it
is consistent with existing ones. In order to do this, we allow the stars to
exist for unrealistically long lifespans and compare the maximum radii of
their ionised bubbles to the Strömgren radii calculated from Eq. (13). In Fig.
2, we plot the maximum extent of the front against the Strömgren radii
calculated from Eq. (13) at the same three redshifts, where the maximum extent
of the front is defined in the algorithm as being the inner radius of the
cells at which
$1-n_{\mathrm{HI};t,i,\theta,\phi}/n_{\mathrm{HI};0,i,\theta,\phi}>0$ and
$1-n_{\mathrm{HI};t,i+1,\theta,\phi}/n_{\mathrm{HI};0,i+1,\theta,\phi}=0$. In
these cases, we adjusted the run time of the model to account for the
respective times it took for the front to reach an equilibrium at the
different redshifts, changing the length of each time step accordingly, as
well as altering the radius of each cell. In the $z=15$ case, we used cell
radii of $2.5\leavevmode\nobreak\ \mathrm{pc}$ and time steps of
$450\leavevmode\nobreak\ \mathrm{years}$, in the $z=10$ case we used a cell
radius of $5.2\leavevmode\nobreak\ \mathrm{pc}$ and time steps of
$1000\leavevmode\nobreak\ \mathrm{years}$, and in the $z=5$ case we used cell
radii of $18\leavevmode\nobreak\ \mathrm{pc}$ and time steps of
$5000\leavevmode\nobreak\ \mathrm{years}$. Due to the isotropy of the map, we
used $N_{\theta}=N_{\phi}=1$ in order to maximise computational efficiency.

We can see that the curves end up slightly exceeding the Strömgren radius, as
expected. We also see that in higher densities, ionisation fronts tend to
reach this equilibrium position sooner; in the case of $z=5$, for example, it
took approximately $500\leavevmode\nobreak\ \mathrm{Myr}$ for the front to
stop expanding, while in the case of $z=15$ it too approximately
$50\leavevmode\nobreak\ \mathrm{Myr}$. This is in keeping with the notion
that, in higher density, regions, the higher rate of B-recombinations will
bring about a total radiation absorption in a smaller time frame as well as a
smaller volume.

Figure 2: The maximum extent of the ionisation front, defined as being the
radius of the cell at which
$1-n_{\mathrm{HI},t,i,\theta,\phi}/n_{\mathrm{HI,0},t,i,\theta,\phi}>0$ and
$1-n_{\mathrm{HI},t,i+1,\theta,\phi}/n_{\mathrm{HI,0},t,i+1,\theta,\phi}=0$,
over periods of time far greater than the realistic lifespan of a Population-
III star assuming a constant density as given by the average densities at
$z=15$ (top left panel), $z=10$ (top right panel), and $z=5$ (bottom panel),
showing that, as expected, the ionisation fronts slightly exceed the Strömgren
radii calculated from Eq. (13), plotted alongside the model’s results as grey
lines.

Having confirmed that the model is consistent with expected results, we may
use it to generate new results by studying the composition of the photoionised
region and the front that borders it in more detail. We therefore set the
stars to exist for their realistic lifespan of $2.7\leavevmode\nobreak\
\mathrm{Myr}$ and ran the model with time steps of $2000\leavevmode\nobreak\
\mathrm{years}$ and cell radii of, in the case of $z=15$,
$17.14\leavevmode\nobreak\ \mathrm{pc}$, in the case of $z=10$,
$8.57\leavevmode\nobreak\ \mathrm{pc}$, and in the case of $z=5$,
$5\leavevmode\nobreak\ \mathrm{pc}$, once again setting
$N_{\theta}=N_{\phi}=1$ for efficiency since we are assuming a homogeneous
cosmological background.

Running the model thus, we obtain the results shown in Fig. 3. These figures
show the evolution of the angularly symmetric front with time, defined through
the ratio
$1-n_{\mathrm{HI}}/n_{\mathrm{HI,0}}=n_{\mathrm{HII}}/n_{\mathrm{HII,0}}$, the
fraction of the medium which has been ionised. As can clearly be seen,
photoionised bubbles in higher density media as are presented by the higher
redshift samples evolve more slowly and with a thinner front. Indeed, not only
do the higher density samples present smaller bubbles and fronts, but the
ratio of the width of the front to the width of the photoionised region is
much lower, indicating that the radiation is able to penetrate partially
ionised cells much less effectively due to their higher count of HI particles.
We also find that the front during the earlier phases of expansion is much
thinner than this and only becomes thicker with a more apparent gradient as
the ionised bubble becomes larger, which is what one might intuitively expect
as the radiation becomes more dispersed and so begins to ionise the media it
is in contact with at a lower rate.

Figure 3: The evolution in time of the photoionisation front created by a
Population-III star emitting $2.4\times 10^{51}\mathrm{\leavevmode\nobreak\
eV\leavevmode\nobreak\ s}^{-1}$ of purely $13.6\leavevmode\nobreak\
\mathrm{eV}$ ultraviolet light consistently from its creation to its
destruction at an age of $2.7\mathrm{Myr}$ at $z=15$ (top left (panel), $z=15$
(top right panel), and $z=5$ (bottom panel), under the approximation of a
constant background density at each redshift. The hydrogen to the outside of
the front is almost entirely neutral, that inside it is almost entirely
ionised with a small amount of constant recombination preventing full
ionisation from ever occurring, and the region within the front is defined as
being the region that has a reduced ionised to neutral hydrogen ratio
$1-n_{\mathrm{p}}/n_{\mathrm{H}}$ of between $0$ and $1$. As one would expect,
the photoionised bubble evolves more quickly and to a larger size in eras with
a lower matter density and the front at a given time is wider, as the lower
optical depth of the medium allows a greater depth of space to be being
ionised simultaneously. Cells in which $n_{\mathrm{HI}}=n_{\mathrm{HI,0}}$ at
a given time step are not plotted and therefore make up the region beyond the
blue curve of the outer layer of the front.

Taking another example, that of an inhomogeneous background density, once
again with a supersonic ionisation front and the same stellar parameters, we
define a density map in which the background is primarily defined as being the
average density at $z=5$ with a region of overdensity at $\theta=3\pi/2$ which
decreases exponentially around that centre until it vanishes at $\theta=0$ and
$\theta=\pi$. Allowing the star to live for the unrealistic timeframe examined
for the results shown in Fig. 2 and studying again the maximum extent of the
front, we find that a slice in time and the homogeneous angular dimension
gives a final ionisation front at $2.7\leavevmode\nobreak\ \mathrm{Myr}$ which
varies with density, as shown in Fig. 4. This serves to highlight the models
ability to study supersonic radiative transport in inhomogeneous media, such
as in a region of space containing nebulae and other sources of over- and
underdensity.

Figure 4: A slice in $\phi$ of the maximum extent of the ionisation front,
defined as being the radius of the cell at which
$1-n_{\mathrm{HI};t,i,\theta,\phi}/n_{\mathrm{HI,0};t,i,\theta,\phi}>0$ and
$1-n_{\mathrm{HI};t,i+1,\theta,\phi}/n_{\mathrm{HI,0};t,i+1,\theta,\phi}=0$,
after a period of time of $500\leavevmode\nobreak\ \mathrm{Myr}$ with a
density contrast defined at $r\geq 9\leavevmode\nobreak\ \mathrm{kPc}$ around
$\theta=3\pi/2$ and exponentially decreasing around that angle before
vanishing at $\theta=0$ and $\theta=\pi$.

##  IV Discussion

We have proposed a new method for modelling the evolution of fronts around
radiative objects and have examined the example application of the ionisation
bubble around a star in a neutral medium. We have shown that this method
agrees with existing methods for estimating basic parameters of the front
while providing significantly more information than those methods, allowing us
to precisely model its evolution in time as a function of all relevant
physical parameters. We propose that this method could find many uses in
varied fields, from modelling the effects of ionising radiation as examined
here, to modelling chemical photocatalysis, to studying more varied phenomenon
such as the propagation of ant expeditions as they explore the territory
around their nests. We intend to continue the development of this model and
use it for astrophysical applications, specifically using it to study ionising
objects during the cosmological epoch of reionisation and how the information
this model provides can be applied to upcoming studies of 21cm spectroscopy
and we encourage members of other academic groups to apply this method with
appropriate modifications and parameters to their own fields.

The current formulation of the Radiative Front Model relies upon a static
medium through which the front is propagating; this works for a supersonic
front, but would run into issues when describing subsonic fronts as they are
much more likely to encounter alterations due to fluid transport caused by the
change of state they herald. Furthermore, the model assumes that no radiation
is scattered to non-radiative angles.

The subsonic front issue could be corrected for by introducing a set of medium
transport equations into the model which allow the particle number to change
between cells in a manner that preserves the overall particle number over the
set of cells.

The scattering of radiation could be accounted for in one of a number of ways;
the simplest would be to incorporate a scalar field into the grid which
affects each cell as a function of the conditions of its surrounding cells,
effectively modelling the non-radiative scattering by taking an average of its
expected effects and applying them directly to each cell. Alternatively, the
radiation could be modelled coming out of each cell in which non-radiative
scattering is occurring directly before being mapped onto the initial polar
coordinate system of cells; however, this method would be significantly more
computationally expensive.

Another possible modification which could improve both the precision and
computational efficiency of the model would involve dividing each cell into
smaller subcells which could be selectively explored depending upon the
circumstances; for example, regions well within or well outside the front
could be treated as single units, while transitional front cells could be
divided up to get a smoother description of the front as it evolves. This
could be performed by assigning regions that have been transitioned enough to
be considered no longer a part of the front but a part of the altered region
with a logical statement and then grouping all cells in contact with each
other with that logical statement together as one larger cell, while taking
all cells that are in contact with the radiation but are still considered to
be a part of the front into smaller cells. However, for the purposes of this
paper which was intended to describe and demonstrate our new method, we
restricted ourselves to the simpler example of using fixed cell sizes.

We leave it to a future paper to explore these potential extensions of the
model in more detail and limit ourselves to the study of entirely radiative,
supersonic fronts in this paper.

## References

  * Sethian [1996] J. A. Sethian,  _Level Set Methods and Fast MArching Methods_ (Cambridge University Press, Cambridge, 1996).
  * Gibou  _et al._ [2018] F. Gibou, R. Fedkiw, and S. Osher, A review of level-set methods and some recent applications, Journal of Computational Physics 353, 82 (2018).
  * Osher [1993] S. Osher, A level set formulation for the solution of the dirichlet problem for hamilton–jacobi equations, SIAM Journal on Mathematical Analysis 24, 1145 (1993), https://doi.org/10.1137/0524066 .
  * Osher and Fedkiw [2003] S. Osher and R. Fedkiw,  _The Level Set Methods and Dynamic Implicit Surfaces_ (Springer-Verlag, New York, 2003).
  * Sethian [1996] J. A. Sethian, A fast marching level set method for monotonically advancing fronts., Proceedings of the National Academy of Sciences 93, 1591 (1996), https://www.pnas.org/doi/pdf/10.1073/pnas.93.4.1591 .
  * Sethian [1999] J. A. Sethian, Fast marching methods, SIAM Review 41, 199 (1999).
  * Tcheng and Nave [2016] A. Tcheng and J.-C. Nave, A fast-marching algorithm for nonmonotonically evolving fronts, SIAM Journal on Scientific Computing 38, A2307 (2016).
  * Goldsworthy [1961] F. A. Goldsworthy, Ionization Fronts in Interstellar Gas and the Expansion of HII Regions, Philosophical Transactions of the Royal Society of London Series A 253, 277 (1961).
  * Ducharme  _et al._ [1990] S. Ducharme, J. Hautala, and P. C. Taylor, Photodarkening profiles and kinetics in chalcogenide glasses, Phys. Rev. B 41, 12250 (1990).
  * Altay and Wise [2015] G. Altay and J. Wise, Rabacus: A Python Package for Analytic Cosmological Radiative Transfer Calculations, Astron. Comput. 10, 73 (2015), arXiv:1502.02798 [astro-ph.CO] .
  * Glover [2005] S. Glover, The Formation Of The First Stars In The Universe, Space Science Review 117, 445 (2005), arXiv:astro-ph/0409737 .
  * Wise [2011] J. Wise, First Light, in  _New Horizons in Astronomy_ (2011) p. 13.
  * Loeb and Barkana [2001] A. Loeb and R. Barkana, The Reionization of the Universe by the first stars and quasars, Ann. Rev. Astron. Astrophys. 39, 19 (2001), arXiv:astro-ph/0010467 .
  * Liu and Bromm [2020] B. Liu and V. Bromm, When did Population III star formation end?, MNRAS 497, 2839 (2020), arXiv:2006.15260 [astro-ph.GA] .
  * Koopmans  _et al._ [2015] L. V. E. Koopmans _et al._ , The Cosmic Dawn and Epoch of Reionization with the Square Kilometre Array, PoS AASKA14, 001 (2015), arXiv:1505.07568 [astro-ph.CO] .
  * Di Matteo  _et al._ [2004] T. Di Matteo, B. Ciardi, and F. Miniati, The 21 Centimeter emission from the reionization epoch: Extended and point source foregrounds, MNRAS 355, 1053 (2004), arXiv:astro-ph/0402322 .
  * Pritchard and Loeb [2012] J. R. Pritchard and A. Loeb, 21-cm cosmology in the 21st century, Rep. Prog. Phys. 75, 086901 (2012), arXiv:1109.6012 .
  * Smoot, George F. and Debono, Ivan [2017] Smoot, George F. and Debono, Ivan, 21 cm intensity mapping with the five hundred metre aperture spherical telescope, A&A 597, A136 (2017).
  * Mirocha  _et al._ [2018] J. Mirocha, R. H. Mebane, S. R. Furlanetto, K. Singal, and D. Trinh, Unique signatures of Population III stars in the global 21-cm signal, MNRAS 478, 5591 (2018), arXiv:1710.02530 [astro-ph.GA] .
  * Osterbrock and Ferland [2006] D. E. Osterbrock and G. J. Ferland,  _Astrophysics of gaseous nebulae and active galactic nuclei_ (University Science Book, Sausalito, CA, 2006).
  * Seaton [1959] M. J. Seaton, The solution of capture-cascade equations for hydrogen, MNRAS 119 (1959).
  * Pengelly [1964] R. M. Pengelly, Recombination spectra, I, MNRAS 127 (1964).
  * Strömgren [1939] B. Strömgren, The Physical State of Interstellar Hydrogen., Astrophys. J. 89, 526 (1939).
  * Crighton  _et al._ [2015] N. H. M. Crighton _et al._ , The neutral hydrogen cosmological mass density at z = 5, MNRAS 452, 217 (2015), arXiv:1506.02037 [astro-ph.CO] .


# Evaluating Deep Learning Models and Adversarial Attacks on Accelerometer-
Based Gesture Authentication

Elliu Huang111Department of Computer Science, San Jose State University Fabio
Di Troia11footnotemark: 1 Mark Stamp11footnotemark: 1<EMAIL_ADDRESS>
###### Abstract

Gesture-based authentication has emerged as a non-intrusive, effective means
of authenticating users on mobile devices. Typically, such authentication
techniques have relied on classical machine learning techniques, but recently,
deep learning techniques have been applied this problem. Although prior
research has shown that deep learning models are vulnerable to adversarial
attacks, relatively little research has been done in the adversarial domain
for behavioral biometrics. In this research, we collect tri-axial
accelerometer gesture data (TAGD) from 46 users and perform classification
experiments with both classical machine learning and deep learning models.
Specifically, we train and test support vector machines (SVM) and
convolutional neural networks (CNN). We then consider a realistic adversarial
attack, where we assume the attacker has access to real users’ TAGD data, but
not the authentication model. We use a deep convolutional generative
adversarial network (DC-GAN) to create adversarial samples, and we show that
our deep learning model is surprisingly robust to such an attack scenario.

##  1 Introduction

With the ubiquity of technology, authentication has become an essential part
of everyday life. Passwords and PINs are the most common forms of
authentication, but biometrics are also popular. Biometric authentication
includes physiological (e.g., facial recognition and fingerprint) and
behavioral (e.g., gait and keystroke dynamics) approaches [5].

While physiological biometric authentication has proven to be highly
effective, sensors and equipment required for such approaches are usually
costly. Additionally, attackers can sometimes bypass such a system if they
have access to a copy of the required features [31]. On the other hand,
behavioral biometrics not only have the potential to be cost effective, they
may also be more secure, at least in cases where attackers have difficulty
imitating the relevant features. Furthermore, the non-intrusive nature of
behavioral biometrics may be considered desirable, in comparison to
physiological biometric authentication.

Gesture-based authentication is a relatively recent behavioral biometric that
has achieved promising results. There are various techniques for analyzing
gestures, including acceleration, angular motion, 3D motion, and a mix of the
three. Several machine learning techniques, including those we discuss in
Section 2, have been applied to the gesture-based authentication problem.

In this research, we explore the effectiveness of deep learning techniques on
gesture-based authentication. Our research is based on a new dataset that we
have collected. Given that our tri-axial accelerometer gesture data (TAGD) are
time series, we consider two time series classification (TSC) techniques:
support vector machines (SVM) and one-dimensional convolutional neural
networks (1D-CNN). When combined with feature extraction techniques, our SVM
model provides for rudimentary analysis of our TAGD data, as well as a basis
for comparison to our 1D-CNN model. We also generate adversarial samples using
generative adversarial networks (GAN) and use these samples to explore the
robustness of our 1D-CNN model against a realistic adversarial attack.

The remainder of this paper is structured as follows. Section 2 discusses
relevant work in the field of gesture-based authentication and adversarial
attacks on biometric authentication. Section 3 provides an overview of the
dataset, including specific steps of the data collection process and data
preprocessing techniques. Machine learning techniques, including generating
adversarial samples, and adversarial strategies are introduced in Section 4.
We present our experimental results for our classification models and an
adversarial attack in Section 5. Our conclusion and a brief discussion of
future research directions are provided in Section 6.

##  2 Related Work

Relative to the vast research literature on behavioral biometrics, there is
comparatively little work on gesture-based authentication. In this section, we
provide an overview of research on gesture-based authentication and
adversarial attacks on such security systems.

Most gesture-based authentication techniques can be categorized into two main
methods, namely, touchscreen and motion gestures [11]. There are, however,
studies that combine both into a single authentication system [8].

Touchscreen-based gesture authentication methods typically analyze touch
dynamics, including various inputs recorded from a touchscreen interface such
as finger size and pressure. One approach consisted of collecting finger
behavior and position data and authenticated users via SVMs [3]. Another study
employed particle swarm optimization to find patterns in touchscreen dynamics
[24].

Motion gestures generally rely on accelerometer and gyroscope data to analyze
the acceleration and angular motion of the mobile device. Prior research in
this domain has applied dynamic time warping (DTW) [22], SVMs [23], and hidden
Markov models (HMM) [15] to authenticate users. One gesture-based approach
employed a more sophisticated method that involved the “leap motion”
controller that collects 3D motion data and applied similarity thresholds to
authenticate users [19]. Another approach analyzed full-body and hand-gestures
in 3D space using two-stream CNNs [32].

Adversarial attacks on gesture-based authentication is not a well-researched
area, but there are a handful of relevant studies in the general field of
behavioral biometrics. One study analyzed behavioral mouse dynamics and found
that deep learning authentication models were susceptible to adversarial
attacks [28]. Adversarial samples have been generated using a fast gradient
sign method (FGSM) to create perturbations in the data, with a gated recurrent
unit (GRU) then used to generate adversarial samples. Another study [1]
analyzed the resilience of continuous touch-based authentication systems
(TCAS) to adversarial attacks. These researchers found that their TCAS trained
with the help of generative adversarial networks (GAN) had a lower false
acceptance rate than that of vanilla TCAS. The paper [2] reports on
experiments with randomization attacks on gesture-based security systems that
use SVMs, and finds that their models are highly vulnerable to adversarial
attacks.

Similarly, adversarial learning on time series classification has not seen
much research. Most adversarial attacks involve small perturbations of the
original data using state-of-the-art FGSM or basic iterative method (BIM) in
order to “trick” a classification or regression model. A study of adversarial
attacks on multivariate time series regression found that three of the most
popular deep learning models—CNNs, GRUs, and long short-term memory
(LSTM)—were highly susceptible to such attacks [25]. Another study also used
FGSM and BIM to create small perturbations in time series data, which
significantly lowered the classification accuracy of deep learning models for
vehicle sensors and electricity consumption data [13].

##  3 Dataset

In this section, we give an overview of the data collected and specify the
steps in the data collection process. We also include a discussion of data
preprocessing and the feature engineering techniques that we have employed.

###  3.1 Data Collection

In this research, we collect users’ tri-axial accelerometer gesture data
(TAGD) while the user holds a smartphone and writes their “signature” in the
air, similar to the process in [18]. The accelerometer sensor of the Physics
Toolbox Sensor Suite [10] was used to collect data. A screenshot of the
application is shown in Figure 1. The tri-axial acceleration is represented as
separate curves—the red curve represents acceleration along the x-axis, the
green curve represents acceleration along the y-axis, the blue curve
represents acceleration along the z-axis, and the white curve represents the
total magnitude of acceleration. Data is collected at a frequency of 100 Hz,
i.e., 100 data points per second. We performed data collection solely on the
Apple iOS platforms in order to reduce the possibility of smartphone type
being a confounding variable. We note that iPhone models varied from iPhone 8
to iPhone X.

|  |
---|---|---
(a) Accelerometer sensor |  | (b) App settings
Figure 1: Screenshots of Physics Toolbox Sensor Suite app

After installing the app, he user performs the following steps to collect
signature data.

  1. 1.

Tap the red button to start recording accelerometer data

  2. 2.

Move the smartphone in the air to draw a signature

  3. 3.

Tap the red button again to stop recording data

  4. 4.

Upload data in the form of a CSV file into Google Drive

These steps were repeated 50 times to collect 50 signatures for each user.

In total, we collected 50 signatures from each of 46 different users. Users
typically chose their initials as their signature, but they were free to
create their own unique signature. Generally, the time to write each signature
varied between 3 and 7 seconds, and the entire data collection process for one
individual required about 20 minutes. Our dataset is freely available for use
by other researchers at [17]. A sample of our TAGD data is shown in Figure 2.

Figure 2: Sample tri-axial acceleration time sequence

###  3.2 Data Preprocesssing

The raw tri-axial data is of variable length. As discussed below, we resize
all signatures to the same size, as most traditional classification techniques
require input data of a fixed size. We have applied several feature
engineering techniques to the resulting time series.

####  3.2.1 Feature Engineering

As noted above, each signature is a temporal sequence of tri-axial
accelerometer data. First, we extract statistical features based on the the
acceleration for each axis, ignoring the sequential nature of the data. The
resulting distributions vary—we compute the following statistical measures of
shape, center, and spread for each of the three axes.

Length ($L$)

— The number of data points in the signature.

Mean ($\mu$)

— The center of the distribution is the mean

${\mu}=\frac{1}{n}\sum_{i=1}^{n}x_{i}=\frac{1}{n}(x_{1}+x_{2}+\cdots+x_{i})$

Median ($m$)

— The median is another measure of the center of the distribution.

Standard deviation ($\sigma$)

— The standard deviation

$\sigma=\mathchoice{{\hbox{$\displaystyle\sqrt{\frac{\displaystyle\sum_{i=1}^{n}(x_{i}-\mu)}{n-1}\,}$}\lower
0.4pt\hbox{\vrule
height=20.72234pt,depth=-16.57794pt}}}{{\hbox{$\textstyle\sqrt{\frac{\displaystyle\sum_{i=1}^{n}(x_{i}-\mu)}{n-1}\,}$}\lower
0.4pt\hbox{\vrule
height=19.75568pt,depth=-15.80461pt}}}{{\hbox{$\scriptstyle\sqrt{\frac{\displaystyle\sum_{i=1}^{n}(x_{i}-\mu)}{n-1}\,}$}\lower
0.4pt\hbox{\vrule
height=19.11122pt,depth=-15.28905pt}}}{{\hbox{$\scriptscriptstyle\sqrt{\frac{\displaystyle\sum_{i=1}^{n}(x_{i}-\mu)}{n-1}\,}$}\lower
0.4pt\hbox{\vrule height=19.11122pt,depth=-15.28905pt}}}$

measures the variability in the signature data.

Kurtosis ($k$)

— The kurtosis is computed as

$k=\frac{\displaystyle\sum_{i=1}^{n}(x_{i}-\mu)^{4}}{n\sigma^{4}}$

and it measures the weight of the tails relative to the center of the
distribution and provides additional information related to the signature
motion.

Skewness ($s$)

— The symmetry of the distribution, which is computed as

$s=\frac{\displaystyle\sum_{i=1}^{n}(x_{i}-\mu)^{3}}{n\sigma^{3}}$

can help us understand the “smoothness” of motion in a signature.

All of these features provide some information about underlying patterns in
users’ signatures.

For each signature, we calculate a feature vector of the form

$(L,\mu_{x},\mu_{y},\mu_{z},m_{x},m_{y},m_{z},\sigma_{x},\sigma_{y},\sigma_{z},k_{x},k_{y},k_{z},s_{x},s_{y},s_{z})$

consisting of the measures of shape, center, and spread of the distribution,
as discussed above. This feature vector of 16 elements is utilized only in our
SVM experiments, below.

####  3.2.2 Time Series Resampling

Since 1D-CNNs and GANs require feature vectors of fixed length, we use tslearn
[30], we resize all the TAGDs to length 400. The resizing function in tslearn
interpolates for arrays less than the target size. We chose to resample all
the time series to length 400 since the median length of the sequences is 380,
while the mean length is 400.

##  4 Implementation

In this section, we discuss the machine learning techniques that form the
basis of our experiments. These techniques, namely, support vector machines,
convolutional neural networks, and generative adversarial networks, will serve
as the basis for our classification experiments in Section 5. We also
introduce our strategies for adversarial attacks.

###  4.1 Support Vector Machines

Support vector machines (SVMs) are one of the most popular supervised learning
techniques for classification and regression. SVMs attempt to find the optimal
separating hyperplane between two labeled sets of training data [27]. However,
a dataset need not be linearly separable, in which case we can employ the
“kernel trick.” As depicted in Figure 3, the kernel trick maps the input data
into a higher-dimensional space where it is more likely to be linearly
separable. The kernel trick, together with “soft margin” calculations that
allow for classification errors, makes an SVM an extremely powerful and
flexible tool in the field of machine learning.

| ${{\displaystyle\phi}\atop{\displaystyle\implies}}$ |
---|---|---
Figure 3: The kernel trick

A classification problem with a small training sample size and high
dimensionality is prone to overfitting [9]. Feature selection techniques can
help to prevent this problem by discarding features, with minimal loss—or even
improvements—in performance. In our SVM experiments, we used support vector
machine recursive feature elimination (SVM-RFE) for feature selection. SVM-RFE
consists of eliminating the least significant feature (based on linear SVM
weights), then training a model on the reduced feature set. This process is
repeated until the desired number of features is reached, or the performance
degrades beyond acceptable limits.

###  4.2 1D Convolutional Neural Networks

Typically, convolutional neural networks (CNN) are associated with feature
extraction and classification for images, which generally involves two-
dimensional convolutional neural networks (2D-CNN). In such a model, 2D data
(e.g., images) are fed into a CNN and classified via a final fully-connected
layer.

In this paper, we do not consider images; instead, we have temporal sequences
of fixed length. Such data is suitable for one-dimensional convolutional
neural networks (1D-CNN). While not as common as 2D-CNNs, 1D-CNNs are used for
signal processing and sequence classification, with numerous applications in
biomedical and civil engineering [21].

The architecture of a 1D-CNN is analogous to that of a 2D-CNN, with the key
differences being the dimensionality of the input data and the convolution
operation. 2D-CNNs typically use a rectangular kernel that slides from left to
right, top to bottom. In contrast, 1D-CNNs employ a kernel that spans some
number of variables and slides along a vector.

###  4.3 Adversarial Strategy

We also consider how our deep learning models perform under adversarial
attacks. While several studies analyze adversarial attacks involving small
perturbations of the original data [20], we explore a scenario where we assume
the intruder has access to a real users’ gesture data. We test both poisoning
and evasion attacks using learning models to generative adversarial samples;
specifically, we use a type of generative adversarial network (GAN) to produce
adversarial samples.

####  4.3.1 Deep Convolutional Generative Adversarial Networks

Two competing neural networks are trained in a GAN—a generative network and a
discriminative network—with the generative network creating fake data that is
designed to defeat the discriminative network. The two networks are trained
simultaneously following a game-theoretic approach. In this way, both networks
improve, with the ultimate objective being a model (discriminative,
generative, or both) that is stronger than it would have been if it was
trained only on the real training data; see [12] for additional details.

An overview of GAN structure is depicted in Figure 4. Among many other uses,
GANs have been used to generate realistic adversarial samples.

Figure 4: Overview of GAN structure

We use deep convolutional GANs (DC-GAN) to replicate our time series data in a
form that will serve as adversarial samples [4]. The generator and
discriminator models are both based on 1D-CNN models. The generator
essentially performs the functions of a convolutional layer in reverse—the
input is an arbitrary sequence of values and it uses transposed convolution
layers to shape the data into a desired form. The discriminator can be based
on a traditional convolutional neural network. The fully-connected layer
outputs a value between $\hbox{}-1$ and $\hbox{}+1$, with a negative output
indicating a fake sample and positive output indicating a real sample. The
generator and discriminator are connected by a loss function, which provides
feedback to both models. Over several epochs of training, the generator should
becomes better at generating adversarial samples while the discriminator
should become better at distinguishing between real and fake samples.

In our GAN model, the generator has an input consisting of a sequence of 100
values sampled from a normal distribution. After the data passes through three
transposed convolution layers, the sequence of 100 values are transformed into
a sequence of the same size as the TAGD, that is, $400\times 3$. The
discriminator model closely follows the architecture of the 1D-CNN
classification model outlined above, with the major difference being the
binary output of the fully-connected layer. The generator and discriminator
are connected by a binary cross entropy loss function. In our experiments, we
vary the number of training epochs to see how effective the adversarial
samples are in breaking down our model.

In Figure 5, we have examples of TAGD generated using DC-GANs with different
training epochs. The number of training epochs is directly related to how well
the the DC-GAN model can replicate data. As we can see in Figure 5 (a) and
(b), the data is quite random doesn’t resemble the real TAGD in 2, whereas
Figures 5 (c) and (d) appear closer to the real TAGD. As we train the DC-GAN
more epochs, the adversarial samples resemble the real data more.

|
---|---
(a) Acceleration sequence after 10 epochs | (b) Acceleration sequence after 25 epochs
|
(c) Acceleration sequence after 50 epochs | (d) Acceleration sequence after 100 epochs
Figure 5: DC-GAN generated acceleration sequences

####  4.3.2 Adversarial Attack

Similar to [26], we “poison” our training dataset with adversarial samples
generated from DC-GANs, meaning that we mix in adversarial samples with our
real training dataset. Then, we train our 1D-CNN on the poisoned dataset and
try to classify real data. The accuracy of classifying with the poisoned
training dataset suggests how well our 1D-CNN can survive a poisoning attack
while simultaneously indicating how well our DC-GAN can generate adversarial
samples.

##  5 Experiments and Results

In this section, we present and analyze the results of the experiments
outlined in the previous section. For our first experiment, we provide the
results of a multiclass classification problem using SVMs, based on the
statistical features discussed in Section 3.2.1. Then we apply a deep learning
technique, 1D-CNNs, in the multiclass classification problem. Both of these
techniques, SVMs and 1D-CNNs, produce strong results. Then we move on to
adversarial learning where we use GANs to generate adversarial samples. With
the GAN-generated adversarial samples, we show that our deep learning model is
robust under a poisoning type of adversarial attack.

We use several metrics in our multiclass classification problem. The most
basic metric we consider is the accuracy, which is calculated from a 46 by 46
confusion matrix. Of course, higher accuracy indicates a more successful
model. We also measure the false acceptance rate (FAR), which is the rate at
which a different user is classified as the actual user, and the false reject
rate (FRR), which is the rate at which real users are mis-classified as other
users. Since we are dealing with multiclass classification, we calculate FAR
and FRR for each individual user and report the average for the 46 different
users [14]. Intuitively, lower FAR and FRR indicates greater success in the
classification model.

In our SVM experiments, we use SVM and RFE libraries from scikit-learn. In our
1D-CNN and GAN experiments, we implemented Keras libraries [6] to develop our
models. For all the experiments, we use an 80-20 train-test split, i.e. 80% of
the data is used to train the model, while the remaining 20% is used for
testing.

###  5.1 SVM Results

We first considered rudimentary experiments with various kernels and values of
parameters and found that a linear kernel with regularization parameter
$C=1000$ worked the best for multiclass classification. As a result, for all
experiments in this section, we use an SVM with linear kernel and $C=1000$.

Here, we train our SVM model on the feature vector described in Section 3.2.1
and use SVM-RFE to select the strongest features. We analyze the relationship
between the number of features and the FAR, FPR, and accuracy.

Using all 16 features discussed in Section 3.2.1, the SVM achieved a 95%
classification accuracy and 0.0014 and 0.057 FAR and FRR, respectively. These
results (and more) are summarized in Figure 6.

$L,\mu,m,\sigma,k,s$$\mu,m,\sigma,k,s$$\mu,m,\sigma,s$$\mu,m,\sigma(y,z),s(y)$$\mu(x,z),m,s(y)$$0.00$$0.20$$0.40$$0.60$$0.80$$1.00$$0.95$$0.94$$0.92$$0.89$$0.76$FeaturesAccuracy$0.00$$0.05$$0.10$$0.15$$0.20$$0.25$$0.0014$$0.0013$$0.0018$$0.0024$$0.0053$FAR
and FRR$0.057$$0.059$$0.083$$0.110$$0.220$AccuracyFARFRR Figure 6: Linear SVM
results for selected combinations of features

As we eliminate features based on the rankings determined by SVM-RFE, the
classification accuracy in Figure 6 (c) generally decreases, dropping to about
89% with 9 of the 16 features. This is still quite strong, considering the
number of classes. As we can see from Figure 6, the accuracy generally
decreases slightly as we eliminate more features, although the accuracy does
not drop below 90% until we have eliminated 7 features. Similarly, the FAR is
exceptionally low even as features are eliminated, staying below 1% for every
combination of features. However, the FRR starts at around 5% and increases to
more than 20% once we have eliminated 10 features.

###  5.2 1D-CNN Results

We also experiment with a 1D-CNN as our classification model, where the model
follows the architecture in [7]. The 1D-CNN model is trained on temporal
sequences of fixed length, as described in Section 4. The model performs one-
dimensional convolutions along the time axis with RELU activation functions
after each convolution layer. The first convolution layer produces 128
filters, while the second convolution layer produces 256 filters. A dropout
layer follows the convolution layers to prevent overfitting [16]. Then we have
a 1D max-pooling layer to downsample the data and highlight the key features.
The output of the max-pooling layer is passed through a flatten layer, which
is then passed through two fully-connected layers. The last fully-connected
layer produces the value that corresponds to the classification. This model is
illustrated in Figure 7.

Figure 7: 1D-CNN architecture

After fine-tuning hyper-parameters (dropout rate and number of filters), we
moved on to experiment with kernel size and stride length. We found the
accuracy hovered around 89% to 91%. According to [29], the performance of our
1D-CNN should be more receptive to changes in hyper-parameters involving the
convolution layers.

We tested different combinations of kernels and stride lengths in the
convolution layers. These results are summarized in Table 1. We see that most
of the results are fairly strong, ranging from a low accuracy of 90% to a high
of 94%. Generally, as the kernel size increases, the accuracy increases.
Similarly, as the stride length increases, the accuracy tends to increase.
This is most likely due to the fact that larger kernels and stride lengths
produce more refined features for the fully-connected layers.

Table 1: Stride length versus kernel size stride | kernel
---|---
3 | 5 | 10 | 25
1 | 0.90457 | 0.89391 | 0.91276 | 0.90870
3 | 0.92696 | 0.93348 | 0.93630 | 0.93587
6 | 0.93087 | 0.94000 | 0.94109 | 0.94022

###  5.3 Adversarial Results

We trained a GAN to generate adversarial samples, and using these samples to
determine whether our deep learning model is robust under adversarial attack.
First, we determine how close our GAN-generated data is to real data through a
simulated poisoning attack. Then we test how well the adversarial samples can
evade an authenticator model.

As outlined in Section 4.3.2, we poison our training dataset with an
increasing percentage of fake data. All real samples are always included, so
any changes in classification accuracy should be caused by our adversarial
(fake) samples. We trained our DC-GAN for different numbers of epochs and
different numbers of adversarial samples. Generally, a higher number of epochs
results in better imitations of the original data. Increasing the number of
adversarial samples in the training dataset gauges how well our 1D-CNN can
resist large-scale poisoning.

The results of our poisoning attacks are given in Table 2. The accuracy
remains relatively high at more than 90%, even for high training epochs and up
to a 1:1 ratio of real to fake data. Classification accuracy does decrease
slightly when there are more adversarial samples in the training dataset, but
the loss in accuracy is not large. These results show that our 1D-CNN is
highly resistant to poisoning attacks.

Table 2: Adversarial attack results (1840 real samples) $\mbox{adversarial}_{\vphantom{N}}\atop\mbox{samples}^{\vphantom{N}}$ | epochs
---|---
10 | 25 | 50 | 100
0100 | 0.93609 | 0.93130 | 0.94761 | 0.94565
0250 | 0.94587 | 0.94848 | 0.93783 | 0.94196
0600 | 0.94087 | 0.93978 | 0.94002 | 0.94152
1840 | 0.94000 | 0.94152 | 0.92239 | 0.93435

We conjecture that the reason for the limited success of our adversarial
attack is that it is exceedingly difficult to generate realistic signatures of
the type considered in this paper. While the real signatures in our dataset
vary wildly, those that we generated using DC-GAN appear to exhibit more
homogeneity. We plan to investigate this issue further in future work.

##  6 Conclusion and Future Work

Previous research has shown that SVMs are a viable techniques for
accelerometer-based gesture authentication [18]. In this paper, we expanded on
and improved upon previous work. First, we refined the feature selection
process with SVM-RFE to select the best features, while maintaining a high
classification accuracy. Then, we used deep learning models, specifically
1D-CNNs, for classification. We obtained strong results, with greater than 90%
classification accuracy, slightly surpassing the accuracy of our SVM model.
Lastly, we experimented with adversarial attacks on our 1D-CNN model, namely,
poisoning and evasion attacks. These simulate realistic attacks, assuming an
intruder has access to the real data, but not the model itself. Our results
indicate that our 1D-CNN is robust under such attacks, achieving greater than
90% accuracy for poisoning attacks and near perfect accuracy for evasion
attacks.

For future work, additional machine learning techniques could be considered.
For example, long-short term memory (LSTM) models could be used instead of
1D-CNNs, since LSTMs generally perform well on sequential data. Additionally,
as mentioned in Section 5.3, we would like to explore alternative methods of
generating adversarial samples to determine whether we can improve on the
limited adversarial attack results obtained with DC-GANs.

## References

  * [1] Mohit Agrawal, Pragyan Mehrotra, Rajesh Kumar, and Rajiv Ratn Shah.  Defending touch-based continuous authentication systems from active adversaries using generative adversarial networks.  arXiv preprint arXiv:2106.07867, 2021.
  * [2] Mohammad Al-Rubaie and J Morris Chang.  Reconstruction attacks against mobile-based continuous authentication systems in the cloud.  IEEE Transactions on Information Forensics and Security, 11(12):2648–2663, 2016.
  * [3] Ala Abdulhakim Alariki and Azizah Abdul Manaf.  Touch gesture authentication framework for touch screen mobile devices.  Journal of Theoretical & Applied Information Technology, 62(2), 2014.
  * [4] Sukarna Barua, Sarah Monazam Erfani, and James Bailey.  Fcc-gan: A fully connected and convolutional net architecture for gans.  arXiv preprint arXiv:1905.02417, 2019.
  * [5] Debnath Bhattacharyya, Rahul Ranjan, Farkhod Alisherov, Minkyu Choi, et al.  Biometric authentication: A review.  International Journal of u-and e-Service, Science and Technology, 2(3):13–28, 2009.
  * [6] Jason Brownlee.  Deep learning with Python: develop deep learning models on Theano and TensorFlow using Keras.  Machine Learning Mastery, 2016.
  * [7] Jason Brownlee.  1d convolutional neural network models for human activity recognition, Aug 2020.
  * [8] Attaullah Buriro, Bruno Crispo, Filippo Delfrari, and Konrad Wrona.  Hold and sign: A novel behavioral biometrics for smartphone user authentication.  In 2016 IEEE security and privacy workshops (SPW), pages 276–285. IEEE, 2016.
  * [9] Xue-wen Chen and Jong Cheol Jeong.  Enhanced recursive feature elimination.  In Sixth International Conference on Machine Learning and Applications (ICMLA 2007), pages 429–435. IEEE, 2007.
  * [10] Chrystian Vieyra.  Physics toolbox sensor suite, 2016.
  * [11] Gradeigh D Clark and Janne Lindqvist.  Engineering gesture-based authentication systems.  IEEE Pervasive Computing, 14(1):18–25, 2015.
  * [12] Antonia Creswell, Tom White, Vincent Dumoulin, Kai Arulkumaran, Biswa Sengupta, and Anil A Bharath.  Generative adversarial networks: An overview.  IEEE Signal Processing Magazine, 35(1):53–65, 2018.
  * [13] Hassan Ismail Fawaz, Germain Forestier, Jonathan Weber, Lhassane Idoumghar, and Pierre-Alain Muller.  Adversarial attacks on deep neural networks for time series classification.  In 2019 International Joint Conference on Neural Networks (IJCNN), pages 1–8. IEEE, 2019.
  * [14] Margherita Grandini, Enrico Bagli, and Giorgio Visani.  Metrics for multi-class classification: an overview.  arXiv preprint arXiv:2008.05756, 2020.
  * [15] Dennis Guse.  Gesture-based user authentication on mobile devices using accelerometer and gyroscope.  Master’s thesis, Technische Universität Berlin, 2017.
  * [16] Geoffrey E Hinton, Nitish Srivastava, Alex Krizhevsky, Ilya Sutskever, and Ruslan R Salakhutdinov.  Improving neural networks by preventing co-adaptation of feature detectors.  arXiv preprint arXiv:1207.0580, 2012.
  * [17] Elliu Huang.  Gesture dataset, 2021.  Available from the authors upon request.
  * [18] Elliu Huang, Fabio Di Troia, Mark Stamp, and Preethi Sundaravaradhan.  A new dataset for smartphone gesture-based authentication.  In ICISSP, pages 771–780, 2021.
  * [19] Satoru Imura and Hiroshi Hosobe.  A hand gesture-based method for biometric authentication.  In International Conference on Human-Computer Interaction, pages 554–566. Springer, 2018.
  * [20] Hassan Ismail Fawaz, Germain Forestier, Jonathan Weber, Lhassane Idoumghar, and Pierre-Alain Muller.  Adversarial attacks on deep neural networks for time series classification.  arXiv e-prints, pages arXiv–1903, 2019.
  * [21] Serkan Kiranyaz, Onur Avci, Osama Abdeljaber, Turker Ince, Moncef Gabbouj, and Daniel J Inman.  1d convolutional neural networks and applications: A survey.  Mechanical systems and signal processing, 151:107398, 2021.
  * [22] Jiayang Liu, Lin Zhong, Jehan Wickramasuriya, and Venu Vasudevan.  uwave: Accelerometer-based personalized gesture recognition and its applications.  Pervasive and Mobile Computing, 5(6):657–675, 2009.
  * [23] Duo Lu, Kai Xu, and Dijiang Huang.  A data driven in-air-handwriting biometric authentication system.  In 2017 IEEE International Joint Conference on Biometrics (IJCB), pages 531–537. IEEE, 2017.
  * [24] Yuxin Meng, Duncan S Wong, Roman Schlegel, et al.  Touch gestures based biometric authentication scheme for touchscreen mobile phones.  In International conference on information security and cryptology, pages 331–350. Springer, 2012.
  * [25] Gautam Raj Mode and Khaza Anuarul Hoque.  Adversarial examples in deep learning for multivariate time series regression.  In 2020 IEEE Applied Imagery Pattern Recognition Workshop (AIPR), pages 1–10. IEEE, 2020.
  * [26] Luis Muñoz-González, Bjarne Pfitzner, Matteo Russo, Javier Carnerero-Cano, and Emil C Lupu.  Poisoning attacks with generative adversarial nets.  arXiv preprint arXiv:1906.07773, 2019.
  * [27] Mark Stamp.  Introduction to machine learning with applications in information security.  CRC Press, Taylor & Francis Group, Boca Raton, FL, 2018.
  * [28] Yi Xiang Marcus Tan, Alfonso Iacovazzi, Ivan Homoliak, Yuval Elovici, and Alexander Binder.  Adversarial attacks on remote user authentication using behavioural mouse dynamics.  In 2019 International Joint Conference on Neural Networks (IJCNN), pages 1–10. IEEE, 2019.
  * [29] Wensi Tang, Guodong Long, Lu Liu, Tianyi Zhou, Jing Jiang, and Michael Blumenstein.  Rethinking 1d-cnn for time series classification: A stronger baseline.  arXiv preprint arXiv:2002.10061, 2020.
  * [30] Romain Tavenard.  tslearn documentation, 2021.
  * [31] Cong Wu, Kun He, Jing Chen, Ziming Zhao, and Ruiying Du.  Liveness is not enough: Enhancing fingerprint authentication with behavioral biometrics to defeat puppet attacks.  In 29th $\\{$USENIX$\\}$ Security Symposium ($\\{$USENIX$\\}$ Security 20), pages 2219–2236, 2020.
  * [32] Jonathan Wu, Prakash Ishwar, and Janusz Konrad.  Two-stream cnns for gesture-based verification and identification: Learning user style.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops, pages 42–50, 2016.