††thanks: These two authors contributed equally††thanks: These two authors
contributed equally

# Laser pulse-droplet interaction enables the deformation and fragmentation of
droplet array

D. Chaitanya Kumar Rao  Department of Mechanical Engineering, Indian Institute
of Science, Bangalore 560012, India  Department of Physics, University of
Gothenburg, SE-412 96 Gothenburg, Sweden  Awanish Pratap Singh  Saptarshi Basu
<EMAIL_ADDRESS>Department of Mechanical Engineering, Indian Institute of
Science, Bangalore 560012, India

###### Abstract

Droplet-droplet interactions is ubiquitous in various applications ranging
from medical diagnostics to enhancing and optimizing liquid jet propulsion. We
employ an experimental technique where the laser pulse interacts with a
micron-sized droplet and causes optical breakdown. The synergy of a nanosecond
laser pulse and an isolated spherical droplet is accurately controlled and
manipulated to influence the deformation and fragmentation of an array of
droplets. We elucidate how the fluid dynamic response (such as drop-drop and
shock-drop interactions) of an arrangement of droplets can be regulated and
optimally shaped by laser pulse energy and its interplay with the optical
density of liquid target. A new butterfly type breakup is revealed, which is
found to result in controlled and efficient fragmentation of the outer
droplets in an array. The spatio-temporal characteristics of a laser-induced
breakdown dictate how shock wave and central droplet fragments can influence
outer droplets. The incident laser energy and pulse width employed in this
work are representative of diverse industrial applications such as surface
cleaning, nano-lithography, microelectronics, and medical procedures such as
intraocular microsurgery.

Droplet array, drop-drop interactions, laser-induced breakup, shock-drop
interaction

##  I Introduction

Effective liquid atomization is crucial in a broad spectrum of industrial
applications such as combustion engines, spray drying, and medical nebulizers
[1, 2, 3]. For instance, the fuel atomization characteristics (such as drop-
drop interactions) in automotive and gas turbine engines are essential in
ascertaining combustion stability, efficiency, and pollutant emissions. To
evaluate these vital parameters, understanding the spatio-temporal
distribution of droplets that originate from the breakup process is essential
[4, 5, 1]. Considering that droplets constitute the sub-grid level of a spray,
it is important to scrutinize the deformation mechanisms and atomization
pertaining to a single isolated droplet [6, 7, 8].

Fragmentation of liquid can be achieved by a focused laser pulse. The focusing
of a high-intensity laser pulse within a fluid beyond a specific threshold
($\sim$ GW/cm2) leads to optical breakdown [9, 10, 11]. The breakdown is
observed as a luminous spark with intense sound, accompanied by the creation
of charged particles while the absorption of laser pulse by the liquid
increases. The breakdown leads to the formation of high-temperature plasma of
the order of 10${}^{6}\>{}^{o}$K [11]. The rapid expansion of plasma is
followed by the generation of a spherical shock wave that grows with time.

Laser-induced breakdown (LIB) inside medium such as liquid droplets provide a
prospect of spatially isolated energy deposition that cannot be accomplished
with any other optical method. The development of well-established ultra-short
pulse lasers has made it feasible to achieve spatial and temporal accuracy in
wide ranging applications such as nano-lithography [12, 13], surface cleaning
[14, 15, 16, 17], combustion [18, 19], and microsurgery [20, 21].

Laser shock cleaning process based on LIB has also demonstrated the ability to
remove nanoscale particles ($\sim$ 10 nm) from silicon wafers, and it has
applications in semiconductor industry [17]. These micro/nano-particle
contaminations alter precision manufacturing quality in the semiconductor
industry and high-threshold optical components. However, the conventional
cleaning techniques do not function efficiently at micro-and nano-length
scales [17]. Laser-induced breakup of a micron-sized droplet can produce a
high-speed jet with velocities of the order of 1000 m/s. It was demonstrated
that the impingement of the atomized droplets on the contaminated substrates
could efficiently remove the nanoparticles under hydrodynamic drag forces
[22]. The distinctive attribute of LIB allows its applicability in biological
and biomedical applications [23, 24, 20, 21]. The laser-induced optical
breakdown has been employed in sophisticated non-intrusive intraocular
microsurgery with nanosecond to femtosecond laser pulse [25, 20, 21].

The focusing of a laser pulse within a droplet can result in cavitation,
deformation, or fragmentation of a freely falling [26, 13] or acoustically
levitated [27, 28] droplets. Gelderblom et al. [26] investigated the
deformation while Klien et al. [13] studied the fragmentation of freely
falling droplets with a nanosecond pulse laser. The significance of shape of
laser pulse and droplet properties on the droplet deformation and breakup was
examined. Avila and Ohl [27] investigated the deformation, breakup, and shock
dynamics corresponding to acoustically levitated droplets of different sizes.
Based on the droplet size and laser energy, distinct breakup modes (from
course fragmentation to rapid atomization) have been reported. Although a
number of studies have been conducted on the interaction of a laser-pulse with
a single droplet, the literature lacks any study on laser-pulse interaction
with an array leading to drop-drop or shock-drop interactions. Moreover, the
direct interaction of the laser pulse with a droplet is not always beneficial,
especially when the target surfaces are photosensitive (e.g. eye).

In this paper, we analyze the breakup of an array of acoustically levitated
droplets with a nanosecond pulse laser. We report the breakup dynamics of
water and Diesel droplets upon their interaction with a laser pulse and their
subsequent collision (droplet fragments and shock wave) with the neighboring
droplets. Furthermore, we have attempted to elucidate questions such as how
the fragmentation modes (based on laser energy and droplet properties) of
central droplets affect the degree of deformation and breakup of neighboring
droplets? Or what role mechanisms such as ligament-mediated breakup play in
the narrow or broader size distribution of the secondary droplets?

##  II Materials and Methods

The droplet arrays are vertically arranged with three different droplet sizes
ranging from 200 $\mu$m to 600 $\mu$m. The mono-disperse droplets with
distinct physico-chemical properties (water and Diesel) are levitated on the
pressure nodes of a single-axis levitator (tec5) with 100 kHz frequency. Water
is used in the present work since it has been established that the optical
breakdown in distilled water is comparable to the breakdown thresholds in
ocular medium [29, 30]. Furthermore, the mechanical properties of water are
similar to those of an aqueous tumor and the vitreous substance [31].
Similarly, Diesel is used as a representative practical liquid that is widely
used not only in internal combustion and gas turbine engines but also in the
manufacturing, construction, and oil and gas industries [5, 4].

Figure 1: Schematic representing (A) Side view (displaying acoustically
levitated droplets) and (B) top view (displaying optical arrangement and
visualization method) of the experimental setup. (C) Representation of shadow
side and illumination side with respect to laser direction and droplet
location. (D-J) Breakdown zones and corresponding breakup modes with respect
to laser direction and droplet location.

A Q-switched Nd: YAG laser with a wavelength of 532 nm (InnoLas Spitlight 400)
was used to create the laser-induced breakdown to fragment the central
droplet. The laser beam with a Gaussian spatial profile was expanded using a
Galilean beam expander (approximately 25.5 mm) and focused via a laser-grade
achromatic doublet lens ($f=$ 150 mm). The laser delivered a pulse duration of
6 ns with incident energies ranging from 5 mJ to 80 mJ and was used in a
single pulse mode. The beam waist was calculated to be $\sim$ 3 $\mu$m
following Singh et al. [18]. Laser pulses were aimed precisely at the central
droplet that was placed near the focal point ( $f\pm$ 1600 $\mu$m). The
droplet location was chosen such that different modes of breakdown are
observed. The laser pulse energies in front of and behind the breakdown
location were measured using a pyroelectric sensor.

Figure 1A & 1B display the schematic diagram of the experimental apparatus. A
high-speed camera (Photron SA5) and a high-speed laser for illumination
(CAVILUX® Smart UHS, 640 nm) were utilized to investigate the droplet breakup
and the expansion of shock waves generated by the breakdown process. The high-
speed images were recorded at 40 kHz and the single shot images (with laser
illumination) were acquired at 50 ns interval to capture the shock waves and
the subsequent evolution of the droplet fragments. The spatial resolution of
recorded images was 9.98 $\mu$m/pixel. A digital delay generator (BNC 745T)
was used to synchronize the high-speed camera, laser for optical breakdown,
and laser for illumination. To achieve high repeatability, the experiments
were performed at least 10 times.

##  III Results and Discussion

The Opto-hydrodynamic phenomena induced during laser pulse-droplet interaction
depend on energy absorption characteristics (strength and location of the LIB)
inside the droplet. The features of the droplet dynamics during the
interaction mainly rely on the energy distribution of the incident laser
(transmission, absorption, and reflection) and the droplet properties (optical
density and breakdown threshold). Therefore, hydrodynamic phenomena
(deformation, fragmentation, and propulsion of fragments) can be regulated by
altering either laser energy or droplet properties.

When the laser pulse is focused inside the droplet, some part of the energy is
transmitted, and the droplet absorbs the rest. A part of absorbed energy
creates plasma and causes the droplet to expand and fragment concurrently, and
the rest of the energy is taken away by the shock wave. We controlled the
location and mode of optical breakdown, which dictates the direction of sheet
or propulsion of fragments. It is found that the mode of optical breakdown
changes with fluid properties. Irrespective of the mode of breakdown, the
breakup of the target droplet is effective when most of the absorbed energy is
utilized to break the central droplet. The breakup of the surrounding droplets
depends on the combined effect of shock wave and central droplet fragments,
and it is efficient when there is a trade-off between these two parameters.
Various time scales relevant to the laser-induced breakup of droplet array are
1. Pulse width (6 ns), 2. shock wave timescale ($\sim$ 10 $\mu$s), 3. breakup
timescale of center droplet ($\sim$ 100 $\mu$s), and 4. breakup timescale of
outer droplets ($\sim$ 1 ms). The morphology of both the inner and outer
droplets remains nearly unaffected through the laser pulse with time in the 6
ns scale.

In the upcoming sections, we discuss how the laser energy, mode of breakdown,
and droplet properties influence the breakup of the central and surrounding
droplets.

###  III.1 Modes of atomization

Fig. 2A and 2B show the fragmentation of the center droplet for 5 mJ laser
energy corresponding to 200 $\mu$m and 600 $\mu$m droplets, respectively. A
sequence of breakup events for a 400 $\mu$m droplet for the same incident
energy is shown in Fig. 2C. Irrespective of the droplet size, the fragments
propel in only one direction (opposite to laser beam direction). A clear
distinction in the intensity of breakup is apparent for different-sized
droplets (see Fig. 2A-C). The difference in the breakup strength can be
attributed to the deposited energy per unit volume available to the droplets
for the same incident energy. The influence of acoustic pressure on the
trajectories of fragments produced by the atomization process is assumed to be
negligible due to the short timescales associated with the dynamics of
fragments ($<$ 100 $\mu$s).

With the increment in incident energy (From 5 mJ to 80 mJ), the breakup
pattern of the center droplet progressively changes, and the orientation of
fragments begins to shift from one direction to all directions (see Fig. 2C at
$\Delta$t = 25 $\mu$s). It is witnessed that for 5 mJ and 10 mJ laser
energies, a single breakdown occurs towards the shadow side of the droplet
(see Fig. 1C). The breakdown prompts the droplet to open from shadow side and
results in propulsion of fragments in the opposite direction (illumination
side). In contrast, for 40 mJ and 80 mJ energies, multiple breakdowns occur
(equally probable in shadow side, droplet center, and illumination side),
causing the droplet to open from multiple locations, which in turn disperses
the fragments in all directions. In the case of 20 mJ, there is a likelihood
of both single and multiple breakdowns, and therefore, it appears to be the
transition energy from single-mode to multi-mode breakdown.

Figure 2: Atomization of an array of water droplets corresponding to droplets
of (A) 200 $\mu$m (5 mJ), (B) 600 $\mu$m (5 mJ), and (C) 400 $\mu$m ( 5 mJ to
80 mJ). (D) High magnification images of the breakup of upper droplet
corresponding to 5 mJ and 10 mJ pulse energy. (E) Normalized secondary droplet
size distribution of water corresponding to 5 mJ pulse energy ( Ds and D0 are
the secondary and parent droplet diameters respectively) (F) Comparison of the
secondary droplet size distribution for different pulse energies (5 mJ, 10 mJ,
and 20 mJ) for the first stage of droplet array breakup. Images shown in (A)
and (B) are captured using high-speed imaging while images in (C) correspond
to single-shot and high-speed imaging. The scale bar represents 1 mm.

For moderately low incident energy (5 mJ), the fragments from the center
droplet interact with the surrounding droplets at an angle, $\alpha\approx
35.5^{o}$. Here, $\alpha$ denotes the divergence angle of the fragments (see
Fig. 2B). The uncertainty in the measurement of $\alpha$ is below 2.5%. The
impact of fragments on the surrounding droplets induces minor breakup ensued
by shape oscillations. These shape oscillations do not cause a further breakup
of the surrounding droplets (upper and lower droplets). It is apparent from
Fig. 2C and high magnification images displayed in Fig. 2D that only a small
portion of the central droplet fragments interact with surrounding droplets
while the larger mass of the fragments propels towards the illumination side.
As soon as the fragments strike the surrounding droplets, a thin sheet forms
that quickly collapses and, in turn, traps air (Fig. 2D). Meanwhile, the thin
sheet develops a rim and consequently leads to the creation and growth of
ligaments or threads. These ligaments become unstable and eventually
experience breakup due to Rayleigh-Plateau instability [32, 33, 34]. The
experimental ligament breakup time ($\sim$ 75 $\mu$s) is consistent with the
capillary timescale ($\tau_{c}=\sqrt{\frac{\rho\times d_{l}{{}^{3}}}{\sigma}}$
= 69 $\mu$s), where $\rho$ and $\sigma$ are the density and surface tension of
liquid, and dl is the ligament diameter. When the laser energy is increased to
10 mJ, $\alpha$ is observed to increase ($\sim 41.5^{o}$). Two stages of
breakup are observed; the first stage represents thin sheet formation followed
by ligament development and breakup, while the second stage signifies sheet
stretching and resulting ligament growth and breakup. Similar to 5 mJ, air
entrapment is also witnessed in the case of 10 mJ due to the initial impact of
fragments onto the surrounding droplets. However, multiple air bubbles were
seen to be trapped inside the surrounding droplets ($\Delta$t = 2250 $\mu$s)
that consequently undergo coalescence ($\Delta$t = 5025 $\mu$s). With an
additional increase in laser energy (20 mJ), $\alpha$ is observed to decrease
($\sim 23^{o}$ ) due to breakdown mode transition (Fig. 2C and Movie S1).
Therefore, the fragments now collide with the surrounding droplets more evenly
than lower laser energies (5 mJ and 10 mJ). The breakup associated with 20 mJ
energy is also accompanied by two stages, although no air entrapment is
noticed in the residual surrounding droplets. Since the sheet formed due to
the impact of fragments keeps expanding, the likelihood of air entrapment
inside the surrounding droplet is marginal. For the laser energies greater
than 20 mJ, no air entrapment is witnessed.

At higher laser energies (40 mJ and 80 mJ), the fragments collide more
symmetrically and evenly with the surrounding droplets since $\alpha\approx
0^{o}$ due to multi-mode breakdown. However, an additional stage of breakup
(third stage) is witnessed, where the expanding sheet becomes substantially
thin and subsequently collapses and undergoes further breakup due to
significant shape oscillations. This additional breakup stage is associated
with the airflow induced by shock wave (will be discussed in the upcoming
section).

Figure 3: (A) High-speed images of breakup in Diesel array droplets
corresponding to 400 $\mu$m droplet ( 5 mJ to 80 mJ). (B) High magnification
images of breakup of upper droplet corresponding to 5 mJ, 20 mJ, and 80 mJ
pulse energy. (C) Butterfly type breakup (strobe and laser-illuminated) in
Diesel droplets via laser-induced breakdown of center droplet. (D) Normalized
secondary droplet size distribution of Diesel corresponding to 5 mJ pulse
energy (where Ds and D0 are the secondary and parent droplet diameters
respectively). (E) Secondary droplet size distribution of Diesel compared to
water droplets for different pulse energies (5 mJ, 10 mJ, and 20 mJ) for the
first stage of droplet array breakup. (F) Secondary droplet size distribution
for multiple stages of droplet breakup for pulse energies of 40 mJ and 80 mJ.
The scale bar represents 1 mm.

The secondary droplets originating from the surrounding droplets were observed
to follow gamma distribution, which is consistent with literature pertinent to
secondary droplet formation via ligament-mediated breakup [35, 1]. Moreover,
Villermaux [1] has shown that compared to distributions such as Poisson or
log-normal distribution, gamma distributions are more effective for fitting
drop-size distribution data. The size distribution of secondary droplets for 5
mJ laser energy is shown in Fig. 2E, and in the present study, a similar
distribution is observed for all the laser energies. For a single (5 mJ) and
two-stage breakup events (10 mJ and 20 mJ), it is noted that the size of the
secondary droplets is dictated by whether the breakup is dominated by the
impact of fragments or ligament-mediated pinch-off. The secondary droplets
formed primarily due to the impact of fragments are usually smaller than the
droplets generated via ligament-mediated breakup (Fig. 2F). Compared to all
other laser energies, it is found that the breakup of the center droplet is
more efficient for 20 mJ, and a further increase in energy results in over-
killing (i.e. pulse-energy exceeds beyond what is required for efficient
droplet fragmentation).

The breakup modes discussed above (for water) are compared with a practical
oil such that its breakdown threshold ($\sim 10^{10}$ W/cm2) is sufficiently
lower than water ($\sim 10^{11}$ W/cm${}^{2})$. In this work, we have used
Diesel as a representative transportation fuel, which has wide applications in
several combustion engines and where drop-drop interactions are dominant.
Unlike water, due to Diesel’s relatively lower breakdown threshold, multiple
breakdowns occur even at low laser energy (5 mJ), and interestingly, a new
type of breakup is observed, which resembles a butterfly (Fig. 3A, $\Delta$
t=50 $\mu$s). In the butterfly type breakup, the ring represents the body of
the butterfly (see Fig. 3A, 40 mJ, 50 $\mu$s) while the wings resemble a thin
sheet (Fig. 3C and Movie S2). A laser-illuminated high-speed image of the
breakup event indicates that the sheet is fragmented and is composed of a
network of holes. For 5 mJ to 20 mJ, the ring radially expands and strikes the
surrounding droplets head-on. However, for 40 mJ and 80 mJ, the ring collides
with the edge of surrounding droplets due to potential over-killing of
breakdown zone. This over-killing causes an imbalance in the energy deposition
at multiple locations, which dictates the direction of the ring. Although not
as striking as in Diesel, this butterfly-like breakup is also seen for water
droplets at the laser energies of 40 mJ and 80 mJ (Fig. 2C).

Fig. 3B shows the high magnification images of upper droplets corresponding to
5 mJ, 20 mJ, and 80 mJ laser energies. In contrast to water, initial thin
sheet formation does not occur in Diesel, probably due to its relatively lower
surface tension and higher viscosity (Fig 3B, 5mJ). Soon after the impact of
fragments, the droplet deforms and stretches, and subsequently undergoes
Rayleigh-Plateau breakup. The experimental ligament breakup time ($\sim 25\mu
s$) matches well with the capillary time ($\sim 20\mu s$). It is evident that,
unlike water, the breakup in Diesel is essentially due to ligament-mediated
breakup (for 5 mJ to 20 mJ), resulting in a comparatively larger size of
secondary droplets. For 20 mJ, the surrounding droplet significantly deforms,
flattens, and ultimately fragments into several secondary droplets (via
ligament-mediated breakup) following shape oscillations.

Beyond 20 mJ laser energy, the influence of shock wave induced airflow becomes
predominant. Accordingly, instead of a ligament-mediated breakup, the airflow
ruptures the stretching sheet resulting in the multiple sized secondary
droplets (Fig. 3B). Similar to water, the secondary droplets follow gamma size
distribution (Fig. 3D) since the breakup is principally ligament-mediated,
especially for 5 mJ to 20 mJ (Fig. 3E). A clear distinction in the trend of
droplet size distribution is noticed between water and Diesel for laser
energies ranging from 5 mJ to 20 mJ (Fig. 3E). In the case of water, as the
laser energy is raised, the secondary droplet size also increases due to the
dominance of ligament-mediated breakup. In contrast, for Diesel, an increase
in laser energy leads to a reduction in the secondary droplet size as the
ligament-mediated breakup becomes less imperative. At high laser energy (for
40 mJ and 80 mJ), the trend for secondary droplet sizes becomes identical for
both water and Diesel due to the dominance of shock-induced airflow (Fig. 3F).

It is important to note that the breakup dynamics (such as ligament-mediated
breakup and sheet breakup) are significantly dominated by the relative
influence of energy deposited in the center droplet, the energy for the
breakup of center droplet, and shock-induced airflow. In the following
section, we will address causes of ligament-mediated breakup, sheet breakup,
over-killing of breakdown zone, and the influence of shock wave on the breakup
of the surrounding droplets.

###  III.2 Dynamics of Shock wave

Figure 4: Single-shot images corresponding to evolution of early dynamics of
the shock wave for water and Diesel center droplet via laser-induced breakdown
for (A) 5mJ, (B) 20mJ, and (C) 80 mJ pulse energy. (D) Interaction of shock
wave with water and Diesel upper droplet for 80 mJ pulse energy. (E)
Comparison of experimental and theoretical evolution of shock wave radius with
time for water (single-mode breakdown) and Diesel (multi-mode breakdown)
droplets. Distribution of transmitted energy, shock wave energy, and energy
for breakup for different laser pulse energies corresponding to (F) Water and
(G) Diesel. The scale bar represents 1 mm.

The focusing of laser pulse inside the droplet results in optical breakdown.
The laser-induced breakdown triggers the radial outward movement of the plasma
and shock wave (multiple shock waves in case of multi-mode breakdown). The
plasma and shock wave move together up to a few hundred nanoseconds
(approximately between 200-400 ns, depending on the pulse energy) (Fig. 4A).
Subsequently, the shock wave detaches from the plasma and takes away most of
the absorbed energy. Therefore, the energy left with plasma causes the breakup
of the central droplet, and the generated fragments move radially outward
behind the shock wave.

The shock wave appears at approximately 350 ns following the laser-droplet
interaction (Fig. 4A). The shock wave is visible in the images along with the
vapor plume and the mist cloud (two-phase mixture of vapor and fragments of
drop, noticeable as a gray to black haze). For a short period (up to 1
$\mu$s), the fragments propel just behind the shock wave. The white spot in
the center is the high-temperature plasma visible for at least 5 $\mu$s. For a
water droplet at 5 mJ energy, a single shock wave evolves with time due to a
single-mode breakdown (Fig. 4A). For the same energy, multi-mode breakdown
occurs in a Diesel droplet, resulting in multiple shock waves and initiation
of breakup from both the shadow and illumination side (1 $\mu$s in Fig. 4A).
This multi-mode breakdown results in the earlier mentioned butterfly type
breakup. It is also evident that the shock wave generated for water is
stronger compared to Diesel (Fig. 4A). This signifies that the energy for
breakup for water is significantly lower than that of Diesel for low laser
energy (5 mJ and 10 mJ). For 20 mJ laser energy, compared to water, LIB in
Diesel results in a higher shock strength (and therefore higher shock
velocity) and lower energy for breakup (lower fragments velocity) (Figs. 4B
and S1, and Movies S3 and S4). A similar trend follows for higher laser
energies (40 mJ and 80 mJ), where the shock strength is more significant for
Diesel compared to water (Fig. 4C).

As discussed earlier, at high laser energy, the shock wave interaction with
surrounding droplets is more dominant. When an incident shock wave interacts
with a convex surface (surrounding droplet), Mach reflection occurs when the
angle of incidence is below a critical value (known as detachment angle) [36].
A triple point comprising an incident shock wave, reflected shock wave, and a
Mach stem can be seen in Fig. 4D ($\Delta$t = 1.1 $\mu$s). The evolution of
the reflected shock can be seen between $\Delta$t = 1.1 and 1.5 $\mu$s in Fig.
4D). Mach-Mach collision takes place when the incident shock passes the
surrounding droplet. Due to the Mach-Mach collision (1.5 $\mu$s), a locally
supersonic flow is created, resulting in standing shock waves (2 $\mu$s). The
shock-droplet interaction in this study is similar to shock interaction with a
solid sphere considering no stripping takes place. Shock-induced stripping is
observed for Diesel surrounding droplets, which can be attributed to a high
shock strength and lower surface tension (7.5 $\mu$s).

For the breakup of central droplet, only the energy remaining with plasma is
valuable. However, for the surrounding droplets, both shock wave and central
droplet fragments play a significant role. Therefore, it is essential to
determine the shock wave energy for the breakup of both central and
surrounding droplets. In this regard, we have quantified the shock wave energy
using a model presented by Jones [37]. A similar model has been used in
various other studies to calculate the shock wave energy after LIB [38, 39,
40]. Fig. 4E shows the comparison of experimental data (shock radius) with the
Jones model. It can be seen from the figure that the Jones model closely
matches the experimental data for single (in water) and multi-mode (in Diesel)
type breakdowns. The shock wave decelerates very quickly (from Ma$\approx$10
to Ma$\approx$2) within a period of 2 $\mu$s. The experimental uncertainty due
to the laser’s pulse-to-pulse variation is measured to be less than 3%. Fig.
4F and 4G show the energy balance (transmitted energy, shock wave energy, and
energy for breakup) of the optical breakdown in water and Diesel,
respectively. In this study, we have neglected the energy losses due to
radiation, reflection, and scattering as it has been reported that their
combined contribution is less than one percent of the incident energy [31,
41].

In energy balance, the transmitted energy depends on the absorption
characteristics of the fluid. It is observed that Diesel absorption is higher
compared to water at lower pulse energy due to its relatively lower breakdown
threshold. Therefore, the center droplet breakup is more efficient at a lower
energy (5 mJ) in Diesel (73 %) compared to water (49 %) due to a higher
probability of breakdown. For water, central droplet breakup is more efficient
(64 %) at 20mJ laser energy. As the energy increases (higher than 5 mJ in
Diesel and 20 mJ in water), it over-kills the breakdown zone. It is noticed
that both absorption characteristics and pulse energy determine the mode of
the breakup. If the absorption is higher, like Diesel, the laser pulse-droplet
interaction will result in multi-mode breakdown even at lower energies (see
Fig. 1G-J), resulting in a butterfly-type breakup. However, if the breakdown
threshold is relatively higher (like water), it may result in a single
breakdown at three possible locations, i.e., shadow, center, and illumination
side of the droplet (see Fig. 1D-F). The location of the LIB determines the
direction of movement of the fragments. If the breakdown occurs in the
illumination side of the droplet, the fragments will move towards the shadow
side (see Fig. 1C). When LIB takes place in the center of the droplet, the
fragments move equally in all directions. And when it occurs in the shadow
side of droplet, the fragments move towards the illumination side. It is worth
noting that we can control and manipulate the breakdown within the droplet,
hence the movement of the droplet fragments. However, once the breakdown zone
is over-killed, it causes multiple breakdowns (hence multiple shocks) within
the droplet and leads to the butterfly type breakup, regardless of the fluid.

Due to this over-kill, the multiple shocks merge within a few hundred
nanoseconds and travel as a single intense shock (see Fig. 4B). The potent
shock wave takes away most of the absorbed energy, and only limited energy is
available to break the center droplet. Since over-killing starts in Diesel
from significantly less energy ($>$5 mJ), the shock velocity is faster
compared to water; whereas, the fragment velocity is much less than the water
(see Figs. 4D and S1). Furthermore, it is observed that irrespective of the
droplet properties, after the over-killing of breakdown is reached, shock wave
energy increases linearly while energy for breakup decreases linearly (Fig. 4F
and 4G). In contrast, the variation in energy for shock wave and breakup is a
nonlinear process when there is no over-killing.

The breakup of the surrounding droplets is dependent on the combined effect of
shock-induced flow and fragments of the center droplet. It is observed that
the breakup is more efficient at higher energies irrespective of the over-
killing of the breakdown zone. Therefore, to understand the role of shock-
induced flow on the breakup of surrounding droplets, we created a breakdown in
significantly small-sized droplets ($\sim$ 20 $\mu$m) such that the influence
of fragments on the surrounding droplets is negligible. The creation of
breakdown in a small-sized droplet resembles laser-induced breakdown in air
[18, 19, 42].

When a significantly smaller central droplet ($\sim$20 $\mu$m) is suspended
while keeping the surrounding droplets large enough (600 $\mu$m), the airflow
behind the shock wave interacts and deforms the surrounding droplets,
indicating that the shock wave has a significant influence on their breakup
(Fig. S2). It can be deduced that the shock wave itself can contribute to the
substantial deformation and breakup of the surrounding droplets at high laser
energy ($>$ 40 mJ) due to significant shape oscillations. This breakup process
by shock wave does not expose the droplets to laser irradiation, which has
great significance when objects are photosensitive. See supplementary material
(Figs. S2 and S3) for further details.

##  IV Conclusion

In summary, by employing a non-intrusive technique, we show that the
deformation, breakup, and shock dynamics corresponding to a droplet array can
be accurately controlled and manipulated. The change in morphology of droplet
array originating from an interplay of incident laser pulse energy (ranging
from 5 mJ to 80 mJ) and liquid properties (water and Diesel) is studied in
detail. For water, the optimum breakup of the array (both central and
surrounding droplets) is achieved at 20 mJ laser pulse energy, beyond which
over-killing of the breakdown zone occurs. In the case of Diesel, over-killing
begins to take place after 5 mJ energy, and the breakup is effective at
different energies for the central (5 mJ) and surrounding droplets (80 mJ). It
is observed that irrespective of the droplet properties, after the over-
killing of breakdown is reached, shock wave energy increases linearly while
energy for breakup decreases linearly with incident laser energy. In contrast,
the variation in energy for shock wave and breakup is a nonlinear process when
there is no over-killing. In addition, we report a new butterfly type breakup,
which occurs due to multi-mode breakdown. This mode of breakup is found to be
an efficient way to deform and fragment the surrounding droplets. Since the
present method offers a direct visualization and control of breakup events
corresponding to droplet array, it can be further extended to study different
photosensitive materials/surfaces (e.g., eye) where a direct interaction of
laser with the target is not desirable. The findings of this work on the
interaction of laser with the droplet array have strong bearings on various
applications such as droplet-droplet interaction in combustion engines,
intraocular microsurgery, surface cleaning, to name a few.

## References

  * Villermaux [2007] E. Villermaux, Fragmentation, Annual Review of Fluid Mechanics 39, 419 (2007).
  * Yarin [2006] A. Yarin, DROP IMPACT DYNAMICS: Splashing, spreading, receding, bouncing…, Annual Review of Fluid Mechanics 38, 159 (2006).
  * Villermaux [2020] E. Villermaux, Fragmentation versus cohesion, Journal of Fluid Mechanics 898, 10.1017/jfm.2020.366 (2020).
  * Lefebvre [1988] A. Lefebvre,  _Atomization and Sprays_ (CRC Press, 1988).
  * Lefebvre and Ballal [2010] A. H. Lefebvre and D. R. Ballal,  _Gas Turbine Combustion_ (CRC Press, 2010).
  * Lyu  _et al._ [2021] S. Lyu, H. Tan, Y. Wakata, X. Yang, C. K. Law, D. Lohse, and C. Sun, On explosive boiling of a multicomponent leidenfrost drop, Proceedings of the National Academy of Sciences 118, e2016107118 (2021).
  * Rao  _et al._ [2017] D. C. K. Rao, S. Karmakar, and S. Basu, Atomization characteristics and instabilities in the combustion of multi-component fuel droplets with high volatility differential, Scientific Reports 7, 10.1038/s41598-017-09663-7 (2017).
  * Rao  _et al._ [2018] D. C. K. Rao, S. Karmakar, and S. Basu, Bubble dynamics and atomization mechanisms in burning multi-component droplets, Physics of Fluids 30, 067101 (2018).
  * Meyerand and Haught [1963] R. G. Meyerand and A. F. Haught, Gas breakdown at optical frequencies, Phys. Rev. Lett. 11, 401 (1963).
  * RAMSDEN and SAVIC [1964] S. A. RAMSDEN and P. SAVIC, A radiative detonation model for the development of a laser-induced spark in air, Nature 203, 1217 (1964).
  * Morgan [1975] C. G. Morgan, Laser-induced breakdown of gases, Reports on Progress in Physics 38, 621 (1975).
  * Torretti  _et al._ [2020] F. Torretti, J. Sheil, R. Schupp, M. M. Basko, M. Bayraktar, R. A. Meijer, S. Witte, W. Ubachs, R. Hoekstra, O. O. Versolato, A. J. Neukirch, and J. Colgan, Prominent radiative contributions from multiply-excited states in laser-produced tin plasma for nanolithography, Nature Communications 11, 10.1038/s41467-020-15678-y (2020).
  * Klein  _et al._ [2020] A. L. Klein, D. Kurilovich, H. Lhuissier, O. O. Versolato, D. Lohse, E. Villermaux, and H. Gelderblom, Drop fragmentation by laser-pulse impact, Journal of Fluid Mechanics 893, 10.1017/jfm.2020.197 (2020).
  * Ohl  _et al._ [2006] C.-D. Ohl, M. Arora, R. Dijkink, V. Janve, and D. Lohse, Surface cleaning from laser-induced cavitation bubbles, Applied Physics Letters 89, 074102 (2006).
  * Jang  _et al._ [2009] D. Jang, J. H. Oh, J.-M. Lee, and D. Kim, Enhanced efficiency of laser shock cleaning process by geometrical confinement of laser-induced plasma, Journal of Applied Physics 106, 014913 (2009).
  * Jang  _et al._ [2011] D. Jang, J.-G. Park, and D. Kim, Enhancement of airborne shock wave by laser-induced breakdown of liquid column in laser shock cleaning, Journal of Applied Physics 109, 073101 (2011).
  * Lee and Watkins [2001] J. M. Lee and K. G. Watkins, Removal of small particles on silicon wafer by laser-induced airborne plasma shock waves, Journal of Applied Physics 89, 6496 (2001).
  * Singh  _et al._ [2019] A. P. Singh, U. P. Padhi, R. Joarder, and A. Roy, Spatio-temporal effect of the breakdown zone in the laser-initiated ignition of atomized ethyl alcohol-air mixture, Applied Energy 247, 140 (2019).
  * Bradley  _et al._ [2004] D. Bradley, C. Sheppard, I. Suardjaja, and R. Woolley, Fundamentals of high-energy spark ignition with lasers, Combustion and Flame 138, 55 (2004).
  * Vogel and Venugopalan [2003] A. Vogel and V. Venugopalan, Mechanisms of pulsed laser ablation of biological tissues, ChemInform 34, 10.1002/chin.200318293 (2003).
  * Vogel  _et al._ [2005] A. Vogel, J. Noack, G. Hüttman, and G. Paltauf, Mechanisms of femtosecond laser nanosurgery of cells and tissues, Applied Physics B 81, 1015 (2005).
  * Ahn  _et al._ [2013] D. Ahn, J. Ha, and D. Kim, Development of an opto-hydrodynamic process to remove nanoparticles from solid surfaces, Applied Surface Science 265, 630 (2013).
  * Vogel [1997] A. Vogel, Nonlinear absorption: intraocular microsurgery and laser lithotripsy, Physics in Medicine and Biology 42, 895 (1997).
  * Palczewska  _et al._ [2020] G. Palczewska, J. Boguslawski, P. Stremplewski, L. Kornaszewski, J. Zhang, Z. Dong, X.-X. Liang, E. Gratton, A. Vogel, M. Wojtkowski, and K. Palczewski, Noninvasive two-photon optical biopsy of retinal fluorophores, Proceedings of the National Academy of Sciences 117, 22532 (2020).
  * Vogel  _et al._ [1995] A. Vogel, U. Parlitz, and S. Busch, Experimental investigations and numerical modeling of shock wave emission and cavitation bubble generation by picosecond and nanosecond optical breakdown, in  _Laser-Tissue Interaction VI_, edited by S. L. Jacques (SPIE, 1995).
  * Gelderblom  _et al._ [2016] H. Gelderblom, H. Lhuissier, A. L. Klein, W. Bouwhuis, D. Lohse, E. Villermaux, and J. H. Snoeijer, Drop deformation by laser-pulse impact, Journal of Fluid Mechanics 794, 676 (2016).
  * Avila and Ohl [2016] S. R. G. Avila and C.-D. Ohl, Fragmentation of acoustically levitating droplets by laser-induced cavitation bubbles, Journal of Fluid Mechanics 805, 551 (2016).
  * Zeng  _et al._ [2018] Q. Zeng, S. R. Gonzalez-Avila, S. T. Voorde, and C.-D. Ohl, Jetting of viscous droplets from cavitation-induced rayleigh–taylor instability, Journal of Fluid Mechanics 846, 916 (2018).
  * Docchio  _et al._ [1986] F. Docchio, C. Sacchi, and J. Marshall, Experimental investigation of optical breakdown thresholds in ocular media under single pulse irradiation with different pulse durations, Lasers in Ophthalmology 1, 83 (1986).
  * Vogel  _et al._ [2008] A. Vogel, N. Linz, S. Freidank, and G. Paltauf, Femtosecond-laser-induced nanocavitation in water: Implications for optical breakdown threshold and cell surgery, Physical Review Letters 100, 10.1103/physrevlett.100.038102 (2008).
  * Vogel  _et al._ [1999] A. Vogel, J. Noack, K. Nahen, D. Theisen, S. Busch, U. Parlitz, D. Hammer, G. Noojin, B. Rockwell, and R. Birngruber, Energy balance of optical breakdown in water at nanosecond to femtosecond time scales, Applied Physics B 68, 271 (1999).
  * Rao and Basu [2020a] D. C. K. Rao and S. Basu, Atomization modes for levitating emulsified droplets undergoing phase change, Experiments in Fluids 61, 10.1007/s00348-019-2873-4 (2020a).
  * Rao and Basu [2020b] D. C. K. Rao and S. Basu, Phenomenology of disruptive breakup mechanism of a levitated evaporating emulsion droplet, Experimental Thermal and Fluid Science 115, 110086 (2020b).
  * Rao and Karmakar [2018] D. C. K. Rao and S. Karmakar, Crown formation and atomization in burning multi-component fuel droplets, Experimental Thermal and Fluid Science 98, 303 (2018).
  * Kooij  _et al._ [2018] S. Kooij, R. Sijs, M. M. Denn, E. Villermaux, and D. Bonn, What determines the drop size in sprays?, Physical Review X 8, 10.1103/physrevx.8.031019 (2018).
  * Sembian  _et al._ [2016] S. Sembian, M. Liverts, N. Tillmark, and N. Apazidis, Plane shock wave interaction with a cylindrical water column, Physics of Fluids 28, 056102 (2016).
  * Jones [1968] D. L. Jones, Intermediate strength blast wave, Physics of Fluids 11, 1664 (1968).
  * Gebel  _et al._ [2015] G. C. Gebel, T. Mosbach, W. Meier, and M. Aigner, Laser-induced blast waves in air and their effect on monodisperse droplet chains of ethanol and kerosene, Shock Waves 25, 415 (2015).
  * An  _et al._ [2017] B. An, Z. Wang, L. Yang, G. Wu, J. Zhu, and X. Li, Experimental investigation of the shock loss and temporal evolution of hot plume resulting from dual-pulse laser-induced breakdown in quiescent air, Journal of Applied Physics 122, 193301 (2017).
  * Padhi  _et al._ [2020] U. P. Padhi, A. P. Singh, and R. Joarder, Experimental and numerical investigations of double pulse laser energy deposition in air, International Journal of Heat and Fluid Flow 82, 108563 (2020).
  * Joarder  _et al._ [2013] R. Joarder, G. Gebel, and T. Mosbach, Two-dimensional numerical simulation of a decaying laser spark in air with radiation loss, International Journal of Heat and Mass Transfer 63, 284 (2013).
  * Bak  _et al._ [2015] M. S. Bak, L. Wermer, and S. kyun Im, Schlieren imaging investigation of successive laser-induced breakdowns in atmospheric-pressure air, Journal of Physics D: Applied Physics 48, 485203 (2015).
  * Lamb [1993] H. Lamb,  _Hydrodynamics_ (Cambridge university press, 1993).

##  Appendix A Relative variation of shock and fragments velocity

Fig. S1 shows the relative dominance of shock velocity over fragments velocity
for water and diesel droplets at different laser energies. In case of water,
for 5 mJ laser energy, shock velocity is much more dominant compared to the
velocity of fragments, which indicates that the shock takes away most of the
absorbed energy while significantly low energy is left for the breakup of
droplets. A marginal increment in fragment velocity is observed for 10 mJ
energy. At 20 mJ energy, the breakup of central droplet is most efficient
since larger energy is available for the breakup of the droplet. Similarly,
the breakup of the surrounding droplets is also efficient at 20 mJ because the
fragments and the shock wave propel at similar velocity. As the laser energy
is further increased to 40 mJ and 80 mJ, over-killing of breakdown zone occurs
and the shock wave takes away most of the energy.

In case of diesel, since over-killing occurs at energy greater than 5 mJ, a
further increase in energy leads to a progressive decrease in fragment
velocity because of lower energy left for the breakup of centre droplet (Fig.
S1).

Figure S1: The variation of ratio of fragment velocity to shock velocity with
the incident laser energy.

##  Appendix B Deformation and breakup of droplet due to shape oscillations

When an external force is applied to a spherical droplet that is not disturbed
by its weight, it oscillates in different resonance modes. The resonance
frequency of the corresponding nth resonance mode depends on the physical
properties of the liquid and the droplet radius. For a small amplitude
oscillation of a spherical droplet immersed in another fluid, a general
theoretical expression of the resonance frequency is given by Lamb [43],

$\omega_{n}^{2}=\frac{n(n+1)(n-1)(n+2)}{(n+1)\rho+n\rho^{{}^{\prime}}}\frac{\sigma}{R^{3}}$
(1)

where $\rho$ and $\rho^{{}^{\prime}}$ are the density of the droplet and the
surrounding medium, respectively. $\sigma$ is the surface tension and $R$ is
the radius of the droplet.

The oscillation period obtained from equation B1 for Diesel and water are 968
$\mu$s and 629 $\mu$s, respectively, and it closely matches with experimental
results. See supplementary material (Figs. S2 and S3) for further details.

###  B.1 Multi-mode resonance

Figure S2: Deformation and breakup of surrounding droplets via shape
oscillations due to shock-induced airflow for Diesel and water droplets.

When the shock wave interacts with the surrounding droplet, initially, the
Diesel droplet deforms into a disk and begins to resonate ($\Delta$t = 625
$\mu$s in Fig. S2). Subsequently, necking starts to occur on the deformed
droplet, resulting in three secondary droplets. Similar deformation and
breakup also occurs for water droplets (Fig. S2).

Droplet deformation begins at 10 $\mu$s (in Diesel) and 15 $\mu$s (in water)
after the laser pulse is triggered. Fig. S2 (at 1000 $\mu$s in Diesel and 675
$\mu$s in water) shows the point at which the droplet crosses the relaxation
point after completing oscillation, indicating an oscillation period of 990
$\mu$s in Diesel and 660 $\mu$s in water. The oscillation period obtained from
equation B1 (n=3) for Diesel and water are 968 $\mu$s and 629 $\mu$s,
respectively, and it closely matches with experimental results.

###  B.2 Single mode resonance

Single mode breakup is observed in Diesel; however, the probability of its
occurrence is considerably low (Fig. S3). Due to the influence of shock
induced airflow the droplet first deforms into a disk. However, the deformed
droplet does not retract; instead, bag formation occurs in the opposite
direction to that of airflow (1350 $\mu$s). If the bag grows and experiences
rupture, a ligament is evolved which undergoes Rayleigh-Plateau breakup (3500
$\mu$s).

Figure S3: Bag formation and rupture in a diesel surrounding droplet due to
single-mode shape oscillations

##  Appendix C Supplementary movies

Movie S1. Atomization of an array of water droplets corresponding to 400
$\mu$m droplets for 20 mJ laser energy.
Movie S2. Atomization of an array of Diesel droplets corresponding to 400
$\mu$m droplets for 20 mJ laser energy.
Movie S3. Evolution of shock wave and fragments of water droplet via laser-
induced breakdown for 20mJ energy.
Movie S4. Evolution of shock wave and fragments of Diesel droplet via laser-
induced breakdown for 20mJ energy.


# Autoregressive based Drift Detection Method

Mansour Zoubeirou A Mayaki
Université Côte d’Azur
CNRS, Inria, I3S
France
email@email
&Michel Riveill
Université Côte d’Azur
CNRS, Inria, I3S
France

###### Abstract

In the classic machine learning framework, models are trained on historical
data and used to predict future values. It is assumed that the data
distribution does not change over time (stationarity). However, in real-world
scenarios, the data generation process changes over time and the model has to
adapt to the new incoming data. This phenomenon is known as concept drift and
leads to a decrease in the predictive model’s performance. In this study, we
propose a new concept drift detection method based on autoregressive models
called ADDM. This method can be integrated into any machine learning algorithm
from deep neural networks to simple linear regression model. Our results show
that this new concept drift detection method outperforms the state-of-the-art
drift detection methods, both on synthetic data sets and real-world data sets.
Our approach is theoretically guaranteed as well as empirical and effective
for the detection of various concept drifts. In addition to the drift
detector, we proposed a new method of concept drift adaptation based on the
severity of the drift.

_K_ eywords Concept drift detection $\cdot$ Data streams $\cdot$ Auto-
regressive model $\cdot$ Machine learning $\cdot$ Deep neural networks

##  1 Introduction

Thanks to progress in the field of big data and data analysis, machine
learning models and more particularly those based on deep neural networks
(Deep Learning) are nowadays experiencing phenomenal success. Since the
beginning of the 2010s, neural networks have been developing at high speed and
the fields of application are multiplying in all business sectors. In the
machine learning framework, models are trained on historical data and used to
predict future values. In this framework, we assume that future incoming data
streams are stationary, i.e., the data generating process does not change over
time. However, this assumption does not hold in most real-world applications
[1]. For example the statistical properties of a streaming data can change
over time due to seasonality or random events. This phenomenon is known in the
machine learning community as concept drift. In the presence of concept drift,
the model’s predictions become less accurate over time.

Machine learning models should therefore take in account concept drift and
update their weights at the right time. Detecting a concept drift is one of
the main challenges when learning with streaming data because of the high
speed and their large size sets which are not able to fit in the main memory
[2]. To deal with concept drift, many algorithms and methods (ADWIN, DDM,
KSWIN, PageHinkley) have been proposed in the literature. Most of these
algorithms detect concept drifts by tracking the changes in the model’s error
rate or using a distance function to measure the dissimilarity of the input
data distribution between some timestamps. These methods are very sensitive to
changes leading to large numbers of detected drifts and false alarms [1].
Moreover, most of these algorithms require full and immediate access to
ground-truth labels which is an unrealistic assumption in most real-world
applications.

To accurately detect concept drifts in stream data, we propose to integrate an
autogressive time series model inside the machine learning loop by considering
the model’s error as a time series. We used a self-exciting threshold auto-
regressive (SETAR) model [3] as the base autoregressive model. SETAR is a
nonlinear time series model and a special case of regime switching models in
which different models apply to different intervals of values of some key
features. Our method has two components: a machine learning model for the
learning task and a SETAR model that detects the changes in the learning
model’s error rate distribution. We call the new concept drift detection
method ADDM. This approach can be used with any type of predictive model
(Logistic regression, random forest etc.). Our results show that the new
concept drift detection method outperforms all state-of-the-art methods on six
(6) synthetic data sets and five (5) real-world data sets. ADDM is more
accurate and has a very low false alarms rate. A low false alarms rate is very
important in real-world application because retraining a machine model is
time-consuming and resource-intensive. This method also has some theoretically
guarantees as the parameters of the change detection component (SETAR model)
are estimated using ordinary least squared (OLS) [4]. Another advantage of our
method is that we can construct confidence intervals for the detected drift
points using statistical inference and subsampling [5]. In addition to the
drift detector, we proposed a new method of concept drift adaptation based on
the severity of the drift. The main idea is to aggregate the old and the new
models using an estimate of the dissimilarity between the old concept and the
new one as weights. The higher is the severity, the less relevant is the old
model.

The rest of the paper is outlined as follows. The Related Works section
discusses the notions related to concept drift and the other studies (or
articles) related to concept drifts detection methods/algorithms. The third
section is dedicated to the theoretical definition of the self-exciting
threshold auto-regressive model and the description of ADDM concept drift
detection method. In section Experimental Data sets, we describe the data sets
used for our experiments. The hyperparameters optimization section describes
in detail the models architecture, the performance metrics and the drift
detection algorithms hyper-parameters optimization. The results are presented
and discussed in section Results and discussions.

##  2 Related Works

In statistical (or machine) learning domain, concept drift occurs when the
statistical properties of the targeted variable $y$ varies arbitrary over time
due to a change in the input data X distribution. Concept drifts can be
categorised in three groups according to their sources [6]. The first type of
concept drift called virtual drift occurs when the data distribution changes
but does not affect the decision boundaries: $P_{t}(X)\neq P_{t+1}(X)$ while
$P_{t}(y|X)=P_{t+1}(y|X)$. Virtual concept drift is not well studied in
machine learning community because it does not affect the model’s outputs. The
second source of concept drift called actual drift happens when the drift
changes the target variable. Thus the a posterior probability of the data
changes in time while it distribution remains unchanged: $P_{t}(y|X)\neq
P_{t+1}(y|X)$, $P_{t}(X)=P_{t+1}(X)$. The last type of concept drift result
from the mixture of the two first sources: $P_{t}(y|X)\neq P_{t+1}(y|X)$ and
$P_{t}(X)\neq P_{t+1}(X)$. In practice it is very difficult to separate these
sources of concept drifts when learning with stream data. The drift detection
algorithms just try to detect the changes that affect the model’s output
without focusing on their sources.

###  2.1 Concept Drift Understanding

Concept drift understanding answers three main questions: When did the drift
occur, How severe is the change and Where are the drifts regions [6]. The when
refers to the fact that any concept drift detection algorithms should be able
to detect the timestamps where the data distribution changes significantly.
Recalling the definition of concept drift, when a drift occurs at time $t$, an
alarm signal is triggered and it also indicates that the learning system
should adapt to a new concept. Another important question is how much did the
data distribution change at the drift points (severity of the drifts). The
severity of concept drift quantifies the dissimilarity between the new concept
and the previous concept . The severity is defined as
$\Delta=\delta(P_{t}(X,y),P_{t+1}(X,y))$ where $\delta$ is a function that
measures the discrepancy of two data distributions when there is a drift at
timestamp $t$ [6]. The greater the value of $\Delta$, the larger the severity
of the concept drift. The severity gives an idea of how the learners should
adapt to the new concept. If the $\Delta$ is low we may just need to update
the learners without changing many parameters. In contrast if the drift is
severe, we may need to retrain a whole new model. The last question is to
identify where the drift regions (new concepts) are located. The drift regions
of concept drift are the sub-regions where the new concept and the previous
concept are located. These sub-regions are identified by finding parts of the
features space where $P_{t}(X,y)$ and $P_{t+1}(X,y)$ are statistically
different. In ensemble learning scenarios, detecting concept drift regions can
help predicting instances in stable regions. Moreover, when learning with an
artificial neural network we can use the knowledge of the concept drift
regions to put weights on the features.

Like all the state-of-the-art concept drift detection algorithms, ADDM
responds to the first question. It detects with high precision the drift
points by monitoring the model’s error rate. Contrary to the state-of-the-art
algorithms, ADDM can compute confidence intervals of each detected drift using
statistical hypothesis testing and subsampling methods.

###  2.2 Concept Drift Detection Methods in the Literature

Many algorithms and methods have been proposed in the literature to detect
concept drifts in stream data. These methods can be classified into three
categories in terms of the test statistics they apply [6]. The first category
called error rate-based drift detection algorithms refers to all the methods
that track changes in the online error rate of base models. The algorithms
trigger an alarm when there is a statistically significant increase or
decrease of the error rate at some timestamps. Our ADDM method belongs to this
category. The second category is data distribution-based drift detection
algorithms which use a distance function or a metric to quantify the
dissimilarity between the distribution of data before and after the suspected
drift timestamp [6]. These algorithm detect drifts directly from the input
data and try to detect the time and the location of the drifts. The last
category called multiple hypothesis test drift detection methods use multiple
hypothesis tests to detect concept drift. The most popular state-of-the-art
drift detection algorithms (ADWIN, DDM, KSWIN, PageHinkley) are error rate-
based. The ADaptive WINdowing (ADWIN) [7] is an adaptive sliding window
algorithm for concept drift detection in stream data. ADWIN require the user
to specify a sensitivity hyperparameter $\alpha\in(0,1)$ which allows the
algorithm to adjust to the input data. A drift is detected when two sub-
windows of a recent window of observations exhibit an absolute difference in
means larger than $\alpha$. The Drift Detection Method (DDM) is a concept
drift detection method based on the PAC learning model premise [8]. If the
algorithm detects an increase in the error rate higher than a calculated
threshold, an alarm is triggered, either change is detected or the algorithm
will warn the user that change may occur in the near future. The Page-Hinkley
(PH) concept drift detector detects changes by computing the observed values
and their mean up to the current moment [9]. The algorithm detects a concept
drift if the observed mean at some instant is greater than a threshold value
$\lambda$. The Kolmogorov-Smirnov Windowing (KSWIN) concept drift detection
method is based on the Kolmogorov-Smirnov (KS) statistical test [10]. Other
versions of these algorithms have been proposed by other authors: Learning
with Local Drift Detection (LLDD) [11], Early Drift Detection Method (EDDM)
[12], Heoffding’s inequality based Drift Detection Method (HDDM) [13], Dynamic
Extreme Learning Machine (DELM) [14].

Baier et al. [1] used neural network uncertainty instead of the model’s error
rate to detect concept drift. The authors proposed to use the Monte Carlo
Dropout technique to capture the model’s uncertainty. A drift is detected when
the model’s uncertainty increases or decreases significantly at some
timestamp(s). Their algorithm called Uncertainty Drift Detection (UDD) is
based on ADWIN algorithm. The main advantage of UDD method is that in contrast
of the majority of error based drift detection algorithms, it does not require
full and immediate access to ground-truth labels which is an unrealistic
assumption in most real-world use cases [1].

Yan et al. [15] proposed an algorithm based on Hoeffding’s inequality to
monitor the error rate and detect concept drift. The main idea of their method
is to use Hoeffding’s concentration inequality to examine the consistency of
the predictive error. This algorithm relies on the theorem that if the data
distribution is stationary if the difference between the predictive error at
time $t$ denoted $p_{t}$ as and its lower bound $p_{bayes}$ goes to 0 when the
number of training instances increases [16]. The authors used Hoeffding’s
Inequality to estimate the desired upper bound above which the error rate is
considered highly unstable. If the predictive error difference
$\Delta_{t}=p_{t}-p_{bayes}$ is very high after learning a large enough number
of instances, it means that their is a concept drift and the data distribution
has changed. According to there results, their method gives better results
than the state-of-the-art methods. Wang et al. [17] proposed a new concept
drift detection method for class imbalanced problem called DDM-OCI. The new
method is inspired from the DDM [8], instead of monitoring the model’s error
rate, the authors used the recall of the minority class to detect changes in
the data distribution. Their results show that, DDM-OCI responds to new
concepts faster than the model applying DDM. Greco and Tania [18] proposed a
real-time unsupervised per-label drift detection methodology based on
embedding distribution distances in deep learning models . Their method
exploits the inner representations assigned by a deep learning model to new
unseen data to detect drifts.

##  3 Our approach for learning Under Concept Drift

###  3.1 Self-exciting Threshold Autoregressive (SETAR) Models

The self-exciting threshold autoregressive model is a nonlinear time series
model proposed by Tong in 1978 [4]. This model is special case of regime
switching models in which different models apply to different intervals of
values of some key variable. This model has certain properties such as limit
cycles, amplitude dependent frequencies, and jump phenomena that can not be
captured by classic linear time series models [3]. Formally, a SETAR model
with k regimes can be written mathematically as follows [4]:

$Y_{t}=\phi_{0}^{(i)}+\phi_{1}^{(i)}Y_{t-1}+\dots+\phi_{p_{i}}^{(i)}Y_{t-p}+\epsilon_{t}^{i}$
(1) $r_{i-1}\leq Y_{t-d}<r_{i}$ (2)

Where $i=1\dots k$, $1<d\leq\max(p_{i})$ a positive integer, $Y_{t}$ is a time
series and $Y_{t-d}$ the threshold variable. The thresholds values are
$-\infty<r_{0}<r_{1}<\dots<r_{k}<+\infty$; for each regime $i$, the error term
${\epsilon_{t}^{i}}$ is a sequence of martingale differences satisfying:

$E(\epsilon_{t}^{i}|F_{t-1})=0,\hskip
8.5359ptsup_{t}=E(|\epsilon_{t}^{i}|^{\delta}|F_{t-1}|)<\infty$ $\hskip
8.5359pta.s.\hskip 8.5359ptfor\hskip 8.5359ptsome\hskip 8.5359pt\delta>2$

With $F_{t-1}$ a $\sigma$ field generated by
$\\{\epsilon_{t-j}^{i}|j=1,\dots;i=1,\dots,k\\}$. Such a process partitions
the one-dimensional Euclidean space into k regimes and follows a linear auto-
regressive model in each regime [4]. A two regime SETAR model can be written
as follows:

$Y_{t}=(\phi_{0}+\phi_{1}Y_{t-1}+\dots+\phi_{p}Y_{t-p})\mathbbm{1}_{Y_{t-d}\leq
r}+(\beta_{0}+\beta_{1}Y_{t-1}+\dots+\beta_{p}Y_{t-p})\mathbbm{1}_{Y_{t-d}>r}+\epsilon_{t}$
(3)

Where $p$ denotes the autoregressive level , $Y_{t-d}$ the threshold variable
and $r$ the threshold parameter. In principle, we would like $\epsilon_{t}$ to
be conditionally heteroskedastic, but for formal theory, we assume that
$\epsilon_{t}$ is $iid$ $(0,\sigma^{2})$. This model is called self-exciting
because the threshold variable is a function of the past values of the
endogenous variable $Y_{t}$. Since the SETAR model is a locally linear model,
ordinary least squares (OLS) techniques can be used to estimate its parameters
[19]. Under the assumption that the error $\epsilon_{t}$ is $iid$
$N(0,\sigma^{2})$, OLS is equivalent to maximum likelihood estimation. The
model in equation 1 can be rewritten as follows:

$Y_{t}=\Phi^{{}^{\prime}}X_{t}(r)+\epsilon_{t}$ (4)

$X_{t}=(1,Y_{t-1},\dots,Y_{t-p})^{{}^{\prime}}$

$X_{t}(r)=(X_{t}\mathbbm{1}_{Y_{t-d}\leq r},X_{t}\mathbbm{1}_{Y_{t-d}>r})$

Where $\Phi=(\phi,\beta)$, $\phi=(\phi_{0},\phi_{1},\dots,\phi_{p})$ and
$(\beta=\beta_{0},\beta_{1},\dots,\beta_{p})$.

The parameters of interest are $\Phi$ and $r$. For a given threshold value
$r$, the OLS estimate of $\Phi$ is [19]

$\hat{\Phi(r)}=\left(\sum_{t=1}^{T}X_{t}(r)^{{}^{\prime}}X_{t}(r)\right)^{{}^{\prime}}\left(\sum_{t=1}^{T}X_{t}(r)^{{}^{\prime}}Y_{t}\right)$
(5)

With estimated residuals
$\hat{\epsilon_{t}}(r)=Y_{t}-X_{t}(r)^{{}^{\prime}}\hat{\Phi(r)}$ and their
estimated variance

$\hat{\sigma_{T}^{2}}(r)=\frac{1}{T}\sum_{t=1}^{T}\hat{\epsilon_{t}}^{2}$ (6)

The estimation task is now reduced to finding the threshold values $r$ that
minimizes estimated residuals variance $\hat{\sigma_{T}^{2}}(r)$ which depends
exclusively on $r$. The other parameters are then computed by using equation
$\eqref{eq3}$. The threshold parameter is computed by minimizing equation
$\eqref{threshold}$ as follows:

$\hat{r}=\arg\min_{r\in R}\hat{\sigma_{T}^{2}}(r)$ (7)

Where $R=\\{y_{t-d},\text{for all}\hskip 8.5359ptd+1\leq t\leq T\\}$. $R$ has
a finite number of elements $[T-(d+1)]$.

Since the model’s parameters are estimated using ordinary least squared (OLS)
method, we have some theoretical guarantees of its convergence [4]. We can
also construct confidence intervals for the detected drift points by using
statistical inference and sub-sampling[5].

###  3.2 ADDM for Concept Drift Detection

Our ADDM method belongs to the category of concept drift detection methods
based on error rate monitoring. The model triggers an alarm when there is a
statistically significant increase or decrease of the error rate at some
timestamps. As proved by Gama et al.[16], if the data distribution is
stationary, the error rate should converge to its minimum value. Therefore if
the error rate becomes very high after predicting a large enough sample, it
may indicate that the data distribution has changed and the model is no longer
fit to the data. ADDM has two components: a machine learning model for the
prediction task and a SETAR model that detects the changes in the learning
model’s error rate $Y_{t}$. When new data instances arrive, they are predicted
at the time of arrival with the deep learning model. The prediction errors are
then computed and used as target variable $Y_{t}$ for the SETAR model (see
Fig.1). Referring to the definition of the SETAR model in equation (1), our
model forecasts the error rate using a linear combination of its past values
assuming that the behavior of the error rate changes once it enters a
different regime or concept. The threshold values estimated by equation (7)
correspond to the concept drift points in the data. We can also construct
confidence intervals for each detected drift using statistical inference and
sub-sampling [5].

###  3.3 Model Updating or Concept drift adaptation

After detecting concept drifts, we need to update the model so that it adapts
to the new data distribution. In the literature there are mainly three groups
of drift adaptation methods. The first method consist of retraining a whole
new model when a concept drift is detected. The second strategy is ensemble
method witch consist of aggregating a new model trained on the samples from
the new distribution with the old models. This strategy can save significant
effort to retrain a new model for recurring concepts. Ensemble methods
comprise a set of base classifiers that may have different types or different
parameters [6]. The last drift adaptation strategy consists of partially
updating the model when the underlying data distribution changes. This
strategy is more efficient than retraining an entire model when the drift only
occurs in well located regions (decision tree algorithm are suited in this
case because trees have the ability to examine and adapt to each sub-region
separately). This approach can be difficult to use in case of deep learning
models because these models are considered as black boxes and we don’t know
witch parameters to update.

We propose a new concept drift adaptation method that consist of aggregating
the old model with a new model trained on the most recent samples using the
severity of the concept drift. The main idea is aggregate the old and new
model using an estimate of the dissimilarity (severity of the drift) between
the old concept and the new one as weight. The higher is the severity, the
less relevant is the old model. The severity denoted $w_{t}$ is also used as
the weight of the new model in the final model. We first compute the third
quantile $Q_{3}$ in each regime and compute $w_{t}$ as follows:

$w_{t}=\frac{\max(Q_{3}^{0},Q_{3}^{t})}{Q_{3}^{0}+Q_{3}^{t}}$ (8)

Where $Q_{3}^{0}$ is the third quantile of the error rate in the old concept
and $Q_{3}^{t}$ in the new concept. We used the third quantile $Q_{3}$ just to
make sure that we have a good estimate of the error rate in each regime. One
can use any other quantile or aggregating metrics (mean, variance etc.). The
quantiles are more suited because they are less sensitive to extreme values.
The term $\max(Q_{3}^{0},Q_{3}^{t})$ is used to ensure that the model learned
under the new concept always gets the highest weight during aggregation. The
final model is defined as follows:

$M_{\hat{\phi}}^{t}\leftarrow M_{\hat{\phi}}^{0}\cdot(1-w_{t})+w_{t}\cdot
M_{new}$ (9)

Where $M_{\hat{\phi}}^{0}$ is the old model and $M_{new}$ the new model
trained with the most recent data. The new model $M_{new}$ can be learned on a
subset (or a window) containing the most recent data or on the whole data set.
The main advantage of our method is that it takes into account the severity of
the drift when updating the model. If the drift is very severe, the new model
has a much more important role than the old one and the influence of the old
model may fade. Another advantage is that, it’s very flexible and can be used
with any kind of model. For example when learning with artificial neural
networks, we can average the old and the new model parameters using $w_{t}$ or
just average their outputs.

ADDM algorithm is defined as follows:

Algorithm 1 ADDM algorithm

1:Input: Training data $D_{tr}$; Validation data $D_{val}$; Data stream
$D_{s}$

2:Train deep learning model: $M_{\hat{\phi}}^{0}\leftarrow M_{\phi}\cdot
fit(D_{tr},D_{val})$

3:Compute validation error:
$\hat{\epsilon}_{val}\leftarrow(y_{val}-\hat{y}_{val})^{2}$

4:Fixe a time window $w$

5:repeat

6: Receive incoming data instances $x_{t-w}$

7: Predict values: $\hat{y}_{t-w}\leftarrow M_{\hat{\phi}}^{0}\cdot
predict(x_{t-w})$

8: Compute error: $\hat{\epsilon}_{t-w}\leftarrow(y_{t-w}-\hat{y}_{t-w})^{2}$

9: Learn Setar model with $\hat{\epsilon}_{t-w}\cup\hat{\epsilon}_{val}$

10: if change is detected then

11: Compute drift severity:
$w_{t}\leftarrow\frac{\max(Q_{3}^{0},Q_{3}^{t})}{Q_{3}^{0}+Q_{3}^{t}}$

12: Get most recent labeled data $D_{recent}$

13: Train new model $M_{new}\leftarrow M_{\phi}\cdot fit(D_{recent})$

14: $M_{\hat{\phi}}^{t}\leftarrow M_{\hat{\phi}}^{0}\cdot(1-w_{t})+w_{t}\cdot
M_{new}$

15: end if

16:until $D_{s}$ ends

##  4 Experimental Data sets

In order to evaluate our method’s capabilities, we compared its performance to
those of seven (7) state-of-the-art methods on six (6) synthetic data sets
with artificial concept drifts and five (5) real-world data sets. The
synthetic data sets were simulated using the python scikit-multiflow package
[20]. The Friedman multi-variate regression data set [21] consist of teen
features each generated from a uniform distribution from the interval $[0,1]$.
The Friedman data set is commonly used to test concept drift detection methods
. In our experiments, we simulated three different versions of the Friedman
data sets with different types of drifts. The Brieman regression data set is
inspired by Baeir et al. [1]. The data set contains teen features, simulated
from uniform distribution. The Mixed data set was inspired by Gama et al. [8]
and has 6 attributes. Four attributes are relevant for classification: two
boolean attributes and two numeric attributes uniformly distributed from 0 to
1 [20]. The Agrawal stream generator was first introduced by Agrawal et al.
[22]. The generator generates a stream data set of nine features, six numeric
and three categorical [20] for binary classification task. We generated two
data sets from this generator with different types of concept drift. At the
end we have six synthetic data sets among which three regression data sets and
three classification data sets. In each data set, artificial concept drifts
were introduced by modifying the distribution of some features.

In addition to these synthetic data sets, we tested our method on five real-
world data sets. Note that all the following data sets are publicly available
on UCI Machine Learning Repository website [23]. The Panama electricity data
set contains historical records of Panama’s electricity demand and weather
measures from January 2015 until June 2020 [24]. The data set contains
historical electricity load, calendar information related to holidays (and
school period) and Weather variables, such as temperature, relative humidity,
precipitation, and wind speed, from three main cities in Panama [24]. The goal
is to predict the electricity demand using all available features. In this
data set, concept drift is present due to seasonal weather changes which
affects the electricity demand. The Italian air quality data set contains the
responses of a gas multi-sensor devices deployed on the field in the Italian
main cities. Hourly responses averages are recorded along with gas
concentrations references from a certified analyzer [25]. The data were
recorded from March 2004 to February 2005. They recorded some air quality
measures such as the hourly averaged concentrations for CO, Non Metanic
Hydrocarbons, Benzene, Total Nitrogen Oxides (NOx) and Nitrogen Dioxide (NO2).
The goal is to predict the benzene concentration (C6H6(GT)), which is a proxy
for air pollution. Like in the electricity data set, concept drift is present
due to seasonal weather changes. According to the World Air Quality Index
project, the air is very polluted in the Italian main cities between November
and February [26]. The NSW data set contains data from the Australian New
South Wales electricity market [8]. In this market, prices are flexible and
are affected by demand and supply of the market. The data set contains 45.312
instances and nine features dated from 7 May 1996 to 5 December 1998. The goal
is to predict if the electricity price goes up or down each 30 minutes.
Concept drift is present due to seasonal weather changes which affects the
electricity demand and its price. The gas sensor array drift data set contains
measurements from 16 chemical sensors utilized in simulations for drift
compensation in a discrimination task of 6 gases at various levels of
concentrations. The data set was gathered over a period of 36 months in a gas
delivery platform facility situated [27]. The goal is to achieve good
classification performance over time. Concept drifts are present in the data
due to sensor aging and external alterations. The Beijing Multi-Site air
quality data set contains 6 main air pollutants and 6 relevant meteorological
variables at multiple sites in Beijing. Each variable is measured hourly from
March 1st, 2013 to February 28th, 2017. The goal is to predict the PM2.5
variable which is a proxy of air quality measure. Concept drifts are present
in this data set due to seasonal weather changes.

##  5 Experimental Design

###  5.1 Performance Metrics

In the case of synthetic data sets we know exactly where the drifts occurred
so we can compare the detector’s outputs to them. For the synthetic data sets,
we use the following metrics: detection accuracy, True positive, false
positive (or false alarms) and the mean time to detection (MTD) . Contrary to
synthetic data sets, real-world data sets don’t have specified concept drift
points. It is therefore very difficult to evaluate concept drift detection
algorithms on them. In this case, we can’t use metrics like accuracy or true
positive rate to compare the detectors. To compare ADDM to the state-of-the-
art algorithms on real-world data sets, we use the mean squared error (MSE)
loss of the learning model in case of regression task and the cross entropy
loss for classification tasks. For each data set, if a drift is detected, a
new model is learned from scratch and evaluated on a subset of the most recent
data. The final performance of the detector is computed by averaging its
losses on all the detected regions. For each detection algorithm, we also take
in account the number of detections because in real-world applications a
detector that gives a large number of alarms is not optimal.

###  5.2 Hyper-parameters Optimization

In this article, we have used deep learning based models as the backbone of
the prediction step of the ADDM method, but any type of machine learning model
can be used (logistic regression, random forest, etc.). For the synthetic data
sets, our used simple multi-layer perceptron (MLP) neural networks with two
hidden layers. In case of real-world data sets, we used long short-term memory
(LSTM) neural networks architecture (see Table 3). Each model is trained and
validated on a subset of the data set where there is no concept drift. The
learned model is then used to predict values for new incoming data where
concept drift is suspected to be present. We then compute the error rate of
the model on new data and try to find if there are significant changes at some
timestamps. When learning with artificial neural networks, instead of
monitoring the error rate which requires total access to the true labels, we
can use the model’s uncertainty to detect concept drift as done by Baeir et
al. [1]. Using Monte Carlo Dropout, we can compute and monitor the model’s
uncertainty and detect concept drifts.

The state-of-the-art drift detection methods and algorithms ( DDM, ADWIN,
PageHinkley, HDDM,KSWIN) used in this study have some hyper-parameters that
should be well chosen carefully so the algorithm can adjust to the data set.
For each method/algorithm, we determined the optimal hyper-parameters by using
a subset of the data set that we called experimental set. For each synthetic
data set, we took a subset containing one concept drift as the experimental
set. Each algorithm is executed on the experimental sets to find the best
hyper-parameters. These hyper-parameters are then used in the final
experiments. We compared these state-of-the-art methods to our ADDM drift
detection method. As described in section 3, the SETAR model requires the user
to set some hyper-parameters. The main parameters are the time delay for the
threshold variable $d$, the auto regressive level $p$. In our study the
parameters values are: $d=2$ and $p=5$. We aim to evaluate and give a fair
comparison among the detectors concerning the performances of real concept
drift detection. The optimal hyper-parameters are listed in table 4.

##  6 Results and Discussions

In this section, we present the results of the state-of-the-art drift
detection methods and ADDM method on the experimental data sets. For each data
set, we compared ADDM to seven state-of-the-art concept drift detection
methods.

###  6.1 Results on Synthetic Data sets

Table 2 shows the experimental results of the drift detection algorithms on
the synthetic data sets. Recall that in case of synthetic data sets we know
the exact drift points so we can compare them to the detector’s outputs. These
results show that ADDM outperforms all other methods in terms of true
positives (TP) and false alarms (FA). It has a very low false alarm rate.
Despite the parameter optimization , note that the state-of-the-art methods
detect a very large number of drifts. This illustrates these algorithms’
problem of high reactivity leading to a large number of false positive drift
detection [1]. Contrary to state-of-the-art methods, ADDM is less sensitive to
small variation and only detects statistically significant drifts. This is
very important in real-world application because retraining a machine model is
time-consuming and resource-intensive. Fig. 1 and Fig. 2 show the results of
all the algorithms applied respectively to the Mixed and the  Brieman 2d
planes data sets. On these data sets, ADDM accurately detects all the drifts
unlike the other algorithms. The Mixed data set contains three concept drift
points and ADDM is the only method capable of accurately detecting all the
drift points (horizontal line with red plus (+) markers in Fig. 1). The
Brieman data set contains six (6) drift points, ADDM accurately detected five
of them (horizontal line with red plus (+) markers in Fig. 2) when none of the
other methods detects more than one drift.

Figure 1: Detected drifts on the Mixed data set. The vertical lines in green
are the actual drifts points. The horizontal lines show the detected drifts by
each detector. The horizontal line with red plus (+) markers shows our
method’s detections.

###  6.2 Results on Real-world Data sets

Table 1 shows the results of the drift detection algorithms on real-world data
sets. Recall that when comparing ADDM to state-of-the-art algorithms on real-
world data sets, we use the mean squared error loss of the learning model in
case of regression task and the cross entropy loss for classification tasks.
In order to have a better view of the detector’s performances, we also listed
the number of times the model was retrained ($nb\\_train$). The number of
retraining is equal to the number of detected drifts and gives us an idea of
how sensitive the detector is. As expected, the state-of-the-art detectors
lead to a large number of retraining. By combining the loss and the number of
retraining, in almost all cases ADDM outperforms the state-of-the-art
detectors. In the rare cases where other methods have outperformed ADDM, the
improvement is very small and the retraining is at least two times that of
ADDM. For example, on the Gas sensor drift data set, the KSWIN algorithm
retrained the model 42 times (loss=1.69) while ADDM retrained the model only 7
times (loss=1.89).

Table 1: Performance of detectors on real-world data sets. The model’s loss
(the lower the best) and the number of retrainings $nb\\_train$ (the lower the
less computationally expensive).

Method |  | Metric | Panama electricity | Italian air quality | Beijing air quality | Gas sensor drift | NSW data set
---|---|---|---|---|---|---|---
ADDM | loss | 0.020 | 0.007 | 0.006 | 1.89 | 0.290
nb_retrain | 7 | 5 | 5 | 7 | 5
ADWIN | loss | 0.032 | 0.02 | 0.016 | 1.727 | 0.234
nb_retrain | 7 | 0 | 195 | 21 | 12
DDM | loss | 0.030 | 0.014 | 0.014 | 1.714 | 0.253
nb_retrain | 18 | 4 | 145 | 41 | 17
KSWIN | loss | 0.034 | 0.011 | 0.016 | 1.69 | 0.255
nb_retrain | 10 | 6 | 132 | 42 | 10
PH | loss | 0.044 | 0.016 | 0.016 | 1.724 | 0.239
nb_retrain | 12 | 7 | 119 | 37 | 11
EDDM | loss | 0.026 | 0.008 | 0.004 | 1.748 | 0.276
nb_retrain | 10 | 3 | 78 | 33 | 8
HDDM_A | loss | 0.058 | 0.007 | 0.003 | 1.694 | 0.260
nb_retrain | 5 | 4 | 69 | 23 | 9
HDDM_W | loss | 0.059 | 0.009 | 0.004 | 1.696 | 0.288
nb_retrain | 5 | 3 | 58 | 27 | 8

##  Conclusion

Detecting concept drift is important in real-world applications as it leads to
a decrease in machine learning models performance. The traditional concept
drift detection methods are very sensitive to changes and leads to a large
number of false alarms. These methods also often require full access to the
true labels. In this paper, we propose a method that combines a machine
learning algorithm with autoregressive time series models to detect concept
drift in stream data. The main idea is to consider the error rate of a machine
learning model as a time series and model them with an autoregressive time
series model. We compared ADDM to seven (7) state-of-the-art concept drift
detection algorithms on six (6) synthetic data sets and five (5) real-world
data sets. The results show that it outperforms all the state-of-the-art
algorithms in terms of accuracy and has a very low false alarm rate. In
addition to the drift detection method, we proposed a new method of concept
drift adaptation based on the severity of the drift. The main idea is to
aggregate the old and new model using an estimate of the dissimilarity between
the old concept and the new one as weights. The higher is the severity, the
less relevant is the old model. In future works we aim to use auto-regressive
models to detect concept drifts using directly the input data instead of the
error rate or the model’s uncertainty. This can be very useful in cases where
the true labels are not available.

## References

  * [1] Lucas Baier, Tim Schlör, Jakob Schöffer, and Niklas Kühl.  Detecting concept drift with neural network model uncertainty.  arXiv preprint arXiv:2107.01873, 2021.
  * [2] Shweta Kadam.  A survey on classification of concept drift with stream data.  2019\.
  * [3] Howell Tong and Keng S Lim.  Threshold autoregression, limit cycles and cyclical data.  In Exploration Of A Nonlinear World: An Appreciation of Howell Tong’s Contributions to Statistics, pages 9–56. World Scientific, 2009.
  * [4] Ruey S Tsay.  Testing and modeling threshold autoregressive processes.  Journal of the American statistical association, 84(405):231–240, 1989.
  * [5] Jesus Gonzalo and Michael Wolf.  Subsampling inference in threshold autoregressive models.  Journal of Econometrics, 127(2):201–224, 2005.
  * [6] Jie Lu, Anjin Liu, Fan Dong, Feng Gu, Joao Gama, and Guangquan Zhang.  Learning under concept drift: A review.  IEEE Transactions on Knowledge and Data Engineering, 31(12):2346–2363, 2018.
  * [7] Albert Bifet and Ricard Gavalda.  Learning from time-changing data with adaptive windowing.  In Proceedings of the 2007 SIAM international conference on data mining, pages 443–448. SIAM, 2007.
  * [8] Joao Gama, Pedro Medas, Gladys Castillo, and Pedro Rodrigues.  Learning with drift detection.  In Brazilian symposium on artificial intelligence, pages 286–295. Springer, 2004.
  * [9] Ewan S Page.  Continuous inspection schemes.  Biometrika, 41(1/2):100–115, 1954.
  * [10] Christoph Raab, Moritz Heusinger, and Frank-Michael Schleif.  Reactive soft prototype computing for concept drift streams.  Neurocomputing, 416:340–351, 2020.
  * [11] Joao Gama and Gladys Castillo.  Learning with local drift detection.  In International conference on advanced data mining and applications, pages 42–55. Springer, 2006.
  * [12] Manuel Baena-Garcıa, José del Campo-Ávila, Raúl Fidalgo, Albert Bifet, R Gavalda, and Rafael Morales-Bueno.  Early drift detection method.  In Fourth international workshop on knowledge discovery from data streams, volume 6, pages 77–86, 2006.
  * [13] Isvani Frias-Blanco, José del Campo-Ávila, Gonzalo Ramos-Jimenez, Rafael Morales-Bueno, Agustin Ortiz-Diaz, and Yaile Caballero-Mota.  Online and non-parametric drift detection methods based on hoeffding’s bounds.  IEEE Transactions on Knowledge and Data Engineering, 27(3):810–823, 2014.
  * [14] Shuliang Xu and Junhong Wang.  Dynamic extreme learning machine for data stream classification.  Neurocomputing, 238:433–449, 2017.
  * [15] Myuu Myuu Wai Yan.  Accurate detecting concept drift in evolving data streams.  ICT Express, 6(4):332–338, 2020.
  * [16] Joao Gama, Raquel Sebastiao, and Pedro Pereira Rodrigues.  On evaluating stream learning algorithms.  Machine learning, 90(3):317–346, 2013.
  * [17] Shuo Wang, Leandro L. Minku, Davide Ghezzi, Daniele Caltabiano, Peter Tino, and Xin Yao.  Concept drift detection for online class imbalance learning.  In The 2013 International Joint Conference on Neural Networks (IJCNN), pages 1–10, 2013.
  * [18] Salvatore Greco and Tania Cerquitelli.  Drift lens: Real-time unsupervised concept drift detection by evaluating per-label embedding distributions.  In 2021 International Conference on Data Mining Workshops (ICDMW), pages 341–349, 2021.
  * [19] Bruce E Hansen.  Inference in tar models.  Unpublished working paper. Chestnut Hill, MA: Boston College Department of Economics (December), 1996.
  * [20] Jacob Montiel, Jesse Read, Albert Bifet, and Talel Abdessalem.  Scikit-multiflow: A multi-output streaming framework.  Journal of Machine Learning Research, 19(72):1–5, 2018.
  * [21] Jerome H Friedman.  Multivariate adaptive regression splines.  The annals of statistics, pages 1–67, 1991.
  * [22] Rakesh Agrawal, Tomasz Imielinski, and Arun Swami.  Database mining: A performance perspective.  IEEE transactions on knowledge and data engineering, 5(6):914–925, 1993.
  * [23] Dheeru Dua and Casey Graff.  UCI machine learning repository, 2017.
  * [24] Ernesto Aguilar Madrid and Nuno Antonio.  Short-term electricity load forecasting with machine learning.  Information, 12(2):50, 2021.
  * [25] Saverio De Vito, Ettore Massera, Marco Piga, Luca Martinotto, and Girolamo Di Francia.  On field calibration of an electronic nose for benzene estimation in an urban pollution monitoring scenario.  Sensors and Actuators B: Chemical, 129(2):750–757, 2008.
  * [26] About the world air quality index project.
  * [27] Alexander Vergara, Shankar Vembu, Tuba Ayhan, Margaret A Ryan, Margie L Homer, and Ramón Huerta.  Chemical gas sensor drift compensation using classifier ensembles.  Sensors and Actuators B: Chemical, 166:320–329, 2012.

##  7 Supplementary Material

Table 2: Detector performance on synthetic data sets. The columns show the
false alarms (FA), true positives (TP) and Mean Time to Detection (MTD) of the
detector. The MTD is expressed in seconds

Method |  | Metric | Agrawal 32 | Agrawal 3213 | Mixed | Friedman | friedman no return | Brieman 2d planes
---|---|---|---|---|---|---|---|---
ADDM | TP | 1 | 2 | 3 | 3 | 5 | 5
FA | 0 | 2 | 0 | 0 | 1 | 1
MTD | 1.5 | 3.657 | 0.368 | 21.823 | 11.912 | 20.711
ADWIN | TP | 0 | 4 | 2 | 1 | 0 | 0
FA | 7 | 15 | 11 | 7 | 10 | 2
MTD | 0.342 | 0.657 | 0.128 | 1.182 | 0.935 | 1.289
DDM | TP | 0 | 2 | 1 | 0 | 1 | 1
FA | 17 | 21 | 15 | 6 | 16 | 5
MTD | 0.372 | 0.612 | 0.126 | 1.435 | 0.916 | 1.348
KSWIN | TP | 0 | 0 | 0 | 0 | 1 | 1
FA | 18 | 21 | 19 | 4 | 9 | 6
MTD | 0.398 | 0.690 | 0.178 | 1.383 | 1.103 | 1.287
PH | TP | 0 | 1 | 1 | 0 | 0 | 1
FA | 15 | 19 | 17 | 3 | 11 | 7
MTD | 0.411 | 0.613 | 0.138 | 1.229 | 0.980 | 1.360
EDDM | TP | 0 | 5 | 3 | 0 | 1 | 3
FA | 14 | 22 | 17 | 4 | 17 | 3
MTD | 0.355 | 0.611 | 0.134 | 1.336 | 1.044 | 1.247
HDDM_A | TP | 0 | 2 | 3 | 0 | 2 | 3
FA | 17 | 23 | 15 | 5 | 17 | 3
MTD | 0.384 | 0.819 | 0.140 | 1.373 | 1.012 | 1.388
HDDM_W | TP | 0 | 1 | 3 | 0 | 0 | 3
FA | 17 | 19 | 18 | 6 | 16 | 3
MTD | 0.403 | 0.630 | 0.133 | 1.343 | 1.236 | 1.324

Figure 2: Detected drifts on the Brieman 2d planes data set. The vertical
lines in green are the actual drifts points. The horizontal lines show the
detected drifts by each detector. The horizontal line with red plus (+)
markers shows our method’s detections. Table 3: Neural network architecture
for each data set.

Data set | architecture | Hyper-parameters | Epochs
---|---|---|---
Friedman | MLP | (30,15,1) | 50
Friedman no return | MLP | (30,15,1) | 50
Brieman 2d planes | MLP | (30,15,1) | 50
Agrawal 32 | MLP | (20,10,1) | 50
Agrawal 3213 | MLP | (20,10,1) | 50
Mixed | MLP | (20,10,1) | 50
Panama electricity | LSTM | (120,60,30) | 50
Italian air quality | LSTM | (64,32,15) | 50
NSW electricity | LSTM | (120,60,30) | 50
Gas sensor drift | MLP | (120,120,6) | 50
Beijing air quality | LSTM | (120,60,30) | 50

Table 4: Description of each data set. Drifts refers to the number of drifts
(only available for synthetic data sets).The remaining columns show, the
optimal hyper-parameter according to the detector.

Data set | Samples | Features | Target | Drifts | ADWIN | Page-Hinkley | KSWIN
---|---|---|---|---|---|---|---
Friedman | 20000 | 3 | continuous | 3 | $1\times 10^{-4}$ | $1\times 10^{-6}$ | 0.001
Friedman no return | 20000 | 3 | continuous | 6 | $1\times 10^{-6}$ | $1\times 10^{-6}$ | 0.0034
Brieman 2d planes | 20000 | 11 | continuous | 6 | $1\times 10^{-3}$ | $1\times 10^{-6}$ | 0.0032
Agrawal 32 | 20000 | 3 | 2 classes | 1 | $5\times 10^{-6}$ | 50 | 0.0034
Agrawal 3213 | 20000 | 3 | 2 classes | 4 | $4.2\times 10^{-5}$ | 50 | 0.0029
Mixed | 20000 | 5 | continuous | 3 | 0.0441 | $1\times 10^{-6}$ | 0.0059
Panama electricity | 48048 | 17 | continuous | - | $1\times 10^{-5}$ | $1\times 10^{-10}$ | $1\times 10^{-10}$
Italian air quality | 8991 | 14 | continuous | - | $3.41\times 10^{-4}$ | $1\times 10^{-10}$ | $1\times 10^{-10}$
NSW electricity | 45312 | 9 | 2 classes | - | $1\times 10^{-5}$ | $1\times 10^{-10}$ | $1\times 10^{-10}$
Gas sensor drift | 13910 | 129 | 6 classes | - | $4\times 10^{-15}$ | 0.05 | $1\times 10^{-10}$
Beijing air quality | 420768 | 18 | continuous | - | $6\times 10^{-5}$ | $1\times 10^{-10}$ | $1\times 10^{-10}$