# Unified Framework of Forced Magnetic Reconnection and Alfvén Resonance

D. Urbanski  A. Tenerani  F. L. Waelbroeck  Department of Physics, The
University of Texas at Austin  Institute for Fusion Studies

###### Abstract

A unified linear theory that includes forced reconnection as a particular case
of Alfvén resonance is presented. We consider a generalized Taylor problem in
which a sheared magnetic field is subject to a time-dependent boundary
perturbation oscillating at frequency $\omega_{0}$. By analyzing the
asymptotic time response of the system, the theory demonstrates that the
Alfvén resonance is due to the residues at the resonant poles, in the complex
frequency plane, introduced by the boundary perturbation. Alfvén resonance
transitions towards forced reconnection, described by the constant-psi regime
for (normalized) times $t\gg S^{1/3}$, when the forcing frequency of the
boundary perturbation is $\omega_{0}\ll S^{-1/3}$, allowing the coupling of
the Alfvén resonances across the neutral line with the reconnecting mode, as
originally suggested in [1]. Additionally, it is shown that even if forced
reconnection develops for finite, albeit small, frequencies, the reconnection
rate and reconnected flux are strongly reduced for frequencies $\omega_{0}\gg
S^{-3/5}$.

††preprint: APS/123-QED

##  I Introduction

The dynamical formation of boundary layers has long been a problem of
fundamental interest to understand plasma heating in weakly-collisional
magnetized plasmas such as those found in space and in most laboratory
devices. Alfvén resonance and forced magnetic reconnection are two energy
conversion processes mediated by the formation of boundary layers and, as
such, have received much attention in the past [2, 3, 4, 5, 6, 7, 8, 9, 10,
11].

The problem of forced magnetic reconnection, known as the “Taylor problem”,
was first investigated analytically by Hahm and Kulsrud [2] in slab geometry,
for an initial sheared magnetic field subject to a boundary perturbation. They
demonstrated that the system evolves from an initial configuration that
preserves the topology of the magnetic field, towards a final equilibrium with
magnetic islands. This process can be described by four linear stages, with
nonlinear island dynamics thereafter. The first two stages essentially
correspond to the ideal evolution leading to the growth of a surface current
at the neutral line that increases linearly with time. If resistivity is
retained, the magnetic flux starts to reconnect at the neutral line,
increasing algebraically with time. The system transitions to the resistive
regime proper in the third stage, at around $t\sim S^{1/3}$, where $S$ is the
Lundquist number based on the magnetic field shear length and $t$ is the time
normalized to the corresponding Alfvén time. The long-time behavior represents
the fourth stage and is described by the constant-psi regime. The constant-psi
regime occurs over a time scale $t\sim S^{3/5}$, during which the reconnected
flux leads to the formation of magnetic islands due to the resistive
dissipation at the developed boundary layer, by ultimately reaching a
stationary equilibrium.

Alfvén resonance has also been studied for a long time for its role in plasma
heating via resonant absorption, with applications to the solar corona and
tokamak plasmas [12, 13, 14]. Alfvén resonance occurs when the wave is
polarized in the plane containing the mean sheared magnetic field and its
gradient. In the presence of an inhomogeneous magnetic field, normal-mode
analysis predicts a continuum spectrum of singular modes exhibiting a
logarithmic singularity at the location where $\omega=k_{\parallel}v_{a}(x)$,
$v_{a}(x)$ being the Alfvén speed and $k_{\parallel}$ the field-aligned wave
vector [7, 15]. When resistivity is included, such a singularity also implies
that if the plasma is driven externally, a steady state is achieved where the
wave energy is accumulated and absorbed at the resonant layer [8].

In 1999, Uberoi and Zweibel [1] pointed out similarities between forced
reconnection and Alfvén resonance, drawing from earlier studies of Alfvén
resonance at the magnetic field neutral line [16]. Indeed, the current density
grows initially linearly with time for Alfvén resonance and forced
reconnection, and resistive effects become important over a time scale $t\sim
S^{1/3}$ in both processes [8, 2]. Furthermore, they show that both processes
are described by the same governing equation in the limit of zero frequency.
Since the Alfvén resonance develops an inner layer $\delta$ that scales as
$\delta\sim S^{-1/3}$ [9], it was reasonably argued that the transition to
forced reconnection should occur for values of the frequency of the injected
wave $\omega\sim S^{-1/3}$. However, an explicit time-dependent solution to
the boundary value problem, demonstrating that Alfvén resonance and forced
magnetic reconnection are in fact intrinsically related, has not yet been
determined.

The purpose of this paper is to provide a unified theory that includes forced
reconnection as a particular case of Alfvén resonance. We solve the time
dependent, boundary value problem for both the early and long time stages of
the evolution of the system, by assuming a perturbation of arbitrary frequency
$\omega_{0}$ away from the neutral line along the lines of the Taylor problem.

The paper is organized as follows: in section 2 we present the model
equations; in section 3 and 4 we derive the most general solution to the time-
dependent boundary value problem (in Laplace and Fourier space) for arbitrary
boundary conditions. The solution is derived by matching the so-called outer
solution (section 3) to the inner layer solution (section 4) by applying the
asymptotic matching technique; in section 5 and 6 we consider the special case
of a standing wave of arbitrary frequency as boundary condition. We derive,
analytically or by inverting our solutions to time and configuration space
with Mathematica [17], the explicit space-time solution in the early-time
stage (section 5), and in the time-asymptotic stage (section 6). In section 6,
we discuss how the known solutions for forced reconnection and Alfvén
resonance are both recovered by analyzing in the complex plane the dominant
contributions to the stream function, and we compare our solutions with linear
2D simulations described in the appendix. Finally, we discuss our results and
provide a summary in section 7.

##  II Model equations and boundary conditions

We consider the two-dimensional resistive magnetohydrodynamic (MHD) equations
for an inviscid and incompressible plasma in slab geometry. In this model, the
unperturbed state is represented by a stationary plasma with uniform density
$\rho$ and anti-symmetric sheared magnetic field given by

$\displaystyle\bm{B_{0}}(x)={B}_{0}\frac{x}{a}{\bf\hat{y}},$ (1)

as is customary for the Taylor problem [2, 3, 6]. An inhomogeneous pressure or
out-of-plane magnetic field are understood to maintain the equilibrium. The
system is then perturbed at $t=0$ by a time-dependent displacement of the
boundary at $x=\pm a$, represented by a forcing function $\Xi(y,t)$. The
resulting perturbed magnetic and velocity fields are defined in terms of the
flux ($\psi$) and stream ($\phi$) functions given, respectively, by

$\bm{B}=\nabla\psi\times{\bf\hat{z}},$ (2)
$\bm{v}=\nabla\phi\times{\bf\hat{z}}.$ (3)

To solve the time dependent, initial value problem for $\psi$ and $\phi$, we
perform both a Fourier transform with respect to the $y$ variable, and a
Laplace transform in time. The stream and flux functions, $\phi$ and $\psi$,
are zero at $t=0$ and the driver is provided by the boundary condition. We
obtain the following set of linearized equations,

$g(\hat{\phi}^{\prime\prime}-k^{2}\hat{\phi})=ikx(\hat{\psi}^{\prime\prime}-k^{2}\hat{\psi}),$
(4)
$g\hat{\psi}=\frac{1}{S}(\hat{\psi}^{\prime\prime}-k^{2}\hat{\psi})+ikx\hat{\phi},$
(5)

where prime denotes the derivative with respect to $x$, and $g$ and $k$ are
the Laplace and Fourier variables defined with the following conventions,

$\hat{f}(g)=\mathcal{L}(f)=\int_{0}^{\infty}f(t)e^{-gt}dt,$ (6)

$\bar{f}(k)=\mathcal{F}(f)=\frac{1}{2\pi}\int_{-\infty}^{\infty}f(y)e^{-iky}dy,$
(7)

for an arbitrary function $f$. For ease of notation, hereafter we will denote
with hats Laplace transforms (therefore functions of $g$), which may or may
not depend on the Fourier variable $k$. In case of ambiguity, we will indicate
explicitly the variables upon which a given function depends. Here, length
scales have been normalized to the slab width $a$, and time to the Alfvén time
defined by $\tau_{a}=a/v_{a}$, where $v_{a}=B_{0}/\sqrt{4\pi\rho}$ is the
Alfvén speed. The Lundquist number $S=av_{a}/\eta$ is the ratio of the
resistive times scale $\tau_{r}=a^{2}/\eta$ to $\tau_{a}$, $\eta$ being the
magnetic diffusivity.

We look for solutions to Eqs. (4)-(5) in the domain $x=[-1,1]$, subject to a
displacement of the boundaries $\hat{\Xi}_{\pm}(k,g)$. The time-dependent
solution for a given Fourier mode will then be obtained by taking the inverse
Laplace transform. The form of the boundary displacement (the forcing) will be
left general for much of the derivation, but it is assumed periodic in the
azimuthal $y$ direction and symmetric about the midplane $x=0$. The general
boundary conditions that complete our problem are therefore represented by

$\hat{\psi}(\pm 1,k,g)=\hat{\Xi}(k,g),$ (8) $\hat{\phi}(\pm
1,k,g)=\pm\frac{g}{ik}\hat{\Xi}(k,g),$ (9)

with Eq. (9) following from the ideal limit of Eq. (5).

##  III The Time Dependent Boundary Layer Problem

Essential to the theories of both forced magnetic reconnection and Alfvén
resonance is the formation of boundary layers about their resonant surfaces.
We assume that boundary layers develop for $|x|\ll 1$. This means that the
regime we are exploring for the Alfvén resonance is that of a sub-Alfvénic
resonant frequency $\omega_{0}\tau_{a}\ll 1$. Therefore, for $S\gg 1$, the
plasma can be assumed to be governed by ideal MHD (or $S=\infty$) except in
the immediate vicinity of the boundary layer, which contains large gradients.
In the boundary layer, or inner layer, resistive MHD must be considered,
however, because of the large gradients, one can assume that
$\partial^{2}/\partial x^{2}\gg k^{2}$. In this way, equations (4) and (5) can
be solved by breaking the domain into two regions: an ideal “outer region”
away from the boundary layer, and an “inner layer” which contains all the
physics of the problem but is treated as one dimensional. These regions are
then asymptotically matched for a fully general solution valid in the entire
domain.

###  III.1 The Outer Region

The region away from the boundary layer where gradients are small composes
much of the plasma and is treated as both ideal and in steady state. Equations
(4) and (5) become, under these conditions,

$\displaystyle\hat{\psi}^{\prime\prime}=k^{2}\hat{\psi},$ (10)
$\displaystyle\hat{\phi}=\frac{g}{ikx}\hat{\psi}.$ (11)

The general solution to Eq. (10) that satisfies the boundary condition given
by Eq. (8) and is even in $x$ is given by [2]

$\hat{\psi}^{out}=\hat{\psi}_{0}(k,g)\left[\cosh(kx)-\frac{\sinh(k|x|)}{\tanh(k)}\right]+\hat{\Xi}(k,g)\frac{\sinh(k|x|)}{\sinh(k)},$
(12)

where $\hat{\psi}_{0}(k,g)=\hat{\psi}^{out}(x=0)$. The outer solution for
$\hat{\phi}$ follows directly from Eq. (11). The solution (12) consists of two
contributions: a contribution satisfying homogeneous boundary conditions (no
forcing) that takes a time dependent value of amplitude $\hat{\psi}_{0}(k,g)$
at the rational surface $x=0$ and an ideal contribution describing a
“screened” plasma response that vanishes at $x=0$ that satisfies the boundary
condition. As can be seen from Eq. (12), $\hat{\psi}^{out}$ exhibits a jump in
its derivative which we define here as

$\displaystyle\Delta\Psi^{\prime}\equiv-\frac{2k}{\tanh(k)}\hat{\psi}_{0}+\frac{2k}{\sinh(k)}\hat{\Xi}.$
(13)

The first term on the right-hand-side of Eq. (13) contains the
$\Delta^{\prime}$ parameter that determines the stability to tearing mode of
an unforced current sheet equilibrium [18], which is negative in this case
(thus, the equilibrium defined by Eq. (1) is stable to tearing mode). The
second term represents the effect of the forcing through the function
$\hat{\Xi}$.

Close to the origin, the outer solutions for the flux and stream functions,
under series expansion, tend to

$\hat{\psi}^{out}(x\rightarrow
0)=\hat{\psi}_{0}+\frac{1}{2}\Delta\Psi^{\prime}|x|,$ (14)
$\displaystyle\hat{\phi}^{out}(x\rightarrow
0)=\frac{g}{ikx}\hat{\psi}_{0}\left(1+\frac{1}{2}\frac{\Delta\Psi^{\prime}}{\hat{\psi}_{0}}|x|\right).$
(15)

The limit of the outer solution for $x\rightarrow 0$ expressed by Eq. (14) and
(15), provides the boundary condition for the inner layer solution as it tends
away from the resonance. We anticipate that the inner layer equations can be
solved for $\hat{\phi}$ by applying a Fourier transform with respect to the
variable $x$ [19, 20]. With that in mind, it follows from Eq. (15) that the
asymptotic behavior of the Fourier transform of ${\phi}^{out}(x)$, denoted as
$\tilde{\phi}^{out}(\theta)$, at large $\theta$ is given by

$\displaystyle\tilde{\phi}^{out}(\theta\rightarrow\infty)=-\frac{g}{k}\frac{\hat{\psi_{0}}}{2}\frac{1}{\theta}\left(|\theta|+\frac{1}{\pi}\Delta_{g}\right),$
(16)

where from now on we use the tilde $\tilde{[\ ]}$ to denote functions in
$\theta$ space, and

$\Delta_{g}\equiv\frac{\Delta\Psi^{\prime}}{\hat{\psi}_{0}}$ (17)

is the layer response. The asymptotic form in Eq. (16) must be matched with
the inner layer solution.

###  III.2 The Inner Layer

In the inner layer, where gradients are large and resistive effects must be
included, equations (4) and (5) are approximated to

$\displaystyle g\hat{\phi}^{\prime\prime}=ikx\hat{\psi}^{\prime\prime},$ (18)
$\displaystyle\hat{\psi}=\frac{1}{gS}\hat{\psi}^{\prime\prime}+\frac{ikx}{g}\hat{\phi},$
(19)

respectively. Just like equations (4)-(5), the system of equations (18)-(19)
admit solutions for $\hat{\psi}$ and $\hat{\phi}$ with a defined parity. In
the following we seek solutions with $\hat{\psi}$ even in $x$ and $\hat{\phi}$
odd, and we will show only the $x\geq 0$ part of the solutions. The solution
in the entire domain can then be obtained by extension, according to the
appropriate parity.

To make the above set of equations analytically tractable, we Fourier
transform in the $x$ direction. This approach allows one to reduce the
original fourth order system of differential equations with respect to the
variable $x$, to a second order system in $\theta$ that is analytically
solvable [21]. In this way, the differential equations in terms of $\theta$
are

$\displaystyle\theta^{2}\tilde{\phi}=-\frac{k}{g}\frac{\partial}{\partial\theta}(\theta^{2}\tilde{\psi}),$
(20)
$\displaystyle\left(1+\frac{\theta^{2}}{gS}\right)\tilde{\psi}=-\frac{k}{g}\frac{\partial}{\partial\theta}\tilde{\phi}.$
(21)

From the definition $J_{z}=-\nabla^{2}\psi\approx\theta^{2}\tilde{\psi}$ and
equation (21), we derive the following expression for $\tilde{J}_{z}$ in terms
of $\tilde{\phi}$,

$\tilde{J}_{z}=-\frac{(Sk)^{1/3}}{Q}\frac{\theta^{2}}{(1+\frac{\theta^{2}}{(Sk)^{2/3}Q})}\frac{d}{d\theta}\tilde{\phi},$
(22)

where we have introduced the normalized Laplace variable

$\displaystyle Q=\frac{gS^{1/3}}{k^{2/3}}.$ (23)

Next, by dividing equation (20) by $\theta^{2}$, differentiating it with
respect to $\theta$, and making use of Eq. (22), we find an expression
entirely in terms of the current density

$\displaystyle\theta^{2}\frac{d}{d\theta}\left(\frac{(Sk)^{2/3}}{\theta^{2}}\frac{d}{d\theta}\tilde{J_{z}}\right)=Q^{2}\left(1+\frac{\theta^{2}}{(Sk)^{2/3}Q}\right)\tilde{J_{z}}.$
(24)

Finally, with the change of variable

$z=\frac{\sqrt{Q}\theta^{2}}{(Sk)^{2/3}},$ (25)

and with the ansatz $\tilde{J}_{z}(z)=e^{-\frac{z}{2}}\tilde{f}(z)$ (based on
the expected behavior of the solution by taking the limit $z\gg 1$ of Eq.
(24)), Eq. (24) takes the form of Kummer’s confluent hypergeometric equation
for the function $\tilde{f}(z)$:

$z\tilde{f}^{\prime\prime}+\left(-\frac{1}{2}-z\right)\tilde{f}^{\prime}-\frac{(Q^{3/2}-1)}{4}\tilde{f}=0.$
(26)

The solution to Kummer’s equation that is well behaved for
$z\rightarrow\infty$ [21] is the Tricomi confluent hypergeometric function
$U(a,b,z)$ [22]. In terms of this function, the solution to Eq. (24) for the
current density takes the form

$\tilde{J}_{z}=Ae^{-\frac{z}{2}}U\left(\frac{(Q^{3/2}-1)}{4},-\frac{1}{2},z\right),$
(27)

where $A$ is an amplitude to be determined via asymptotic matching. By using
the recurrence properties of the exponential and of the function $U(a,b,z)$,
together with equation (20) upon substitution of
$\tilde{J}_{z}=\theta^{2}\tilde{\psi}$, we determine the inner layer solution
for $\tilde{\phi}$ as:

$\tilde{\phi}^{in}=\frac{2A}{(Sk)^{2/3}Q^{1/4}}\frac{e^{-\frac{z}{2}}}{\sqrt{z}}\left[U\left(\frac{(Q^{3/2}-1)}{4},\frac{1}{2},z\right)-\frac{1}{2}U\left(\frac{(Q^{3/2}-1)}{4},-\frac{1}{2},z\right)\right].$
(28)

In the limit of small $\theta$ (or large $x$), the solution to be
asymptotically matched with Eq. (16) is

$\tilde{\phi}^{in}_{\theta\rightarrow
0}=\frac{2\sqrt{\pi}A}{(Sk)^{2/3}Q^{1/4}}\frac{1}{\theta}\left(\frac{Q^{5/4}}{(Q^{3/2}+1)}\frac{(Sk)^{1/3}}{\Gamma\left(\frac{Q^{3/2}+1}{4}\right)}-\frac{2}{\Gamma\left(\frac{Q^{3/2}-1}{4}\right)}\theta\right)$
(29)

with $\Gamma$ the Euler Gamma function.

##  IV The Matched Inner Layer Solution

The most general matched inner layer solution is determined by Eq.
$\eqref{eq:Phi Inner}$ with the following amplitude $A$ determined by
asymptotic matching,

$A=\frac{\hat{\psi}_{0}(k,g)}{8\sqrt{\pi}}Q^{5/4}(Sk)^{1/3}\Gamma\left(\frac{(Q^{3/2}-1)}{4}\right),$
(30)

with

$\hat{\psi}_{0}(k,g)=\hat{\Xi}(k,g)\frac{\frac{2k}{\sinh[k]}}{\Delta_{g}+\frac{2k}{\tanh[k]}}\equiv\hat{\Xi}(k,g)\frac{\Delta_{f}}{\Delta_{g}-\Delta_{0}},$
(31)

and

$\Delta_{g}=-2\pi(Sk)^{1/3}\frac{Q^{5/4}}{(Q^{3}-1)}\frac{\Gamma\left(\frac{(Q^{3/2}-1)}{4}+1\right)}{\Gamma\left(\frac{(Q^{3/2}-1)}{4}+\frac{1}{2}\right)}.$
(32)

Here $\Delta_{f}=\frac{2k}{\sinh{k}}$ and $\Delta_{0}=-\frac{2k}{\tanh{k}}$,
with $\Delta_{0}$ being the conventional jump of the derivative of
$\hat{\psi}^{out}$ at the origin (the $\Delta^{\prime}$ parameter of the
tearing mode).

Equation (28), together with Eqs. (30)-(32), completely describes the stream
function for both forced magnetic reconnection and Alfvén resonance, and it
includes the full time evolution from the early stages, $t\ll S^{1/3}$, when
the evolution is essentially ideal and the boundary layer is being formed, all
the way to the long-time behavior determined by resistivity.

We note that Eq. (28) includes as special cases the eigenfunctions that
describe the ideal kink, the resistive kink (also referred to as non constant-
psi regime for reconnection in slab geometry) and the ordinary tearing mode
(constant-psi regimes) under appropriate asymptotic expansions [23, 21, 24,
20]111We note a typo on equation (31) of [20].. The important difference is
that here we are solving for the time-dependent problem rather than the
eigenvalue problem and, therefore, $g$ is a variable; additionally, we include
the effect of a general forcing function at the boundary through the function
$\hat{\Xi}$ by allowing not only (forced) reconnection, but also Alfvén
resonances.

The solution given in Eq. (28) will now be considered for two distinct stages
of the dynamical evolution of the system. The early-time ideal regime
corresponds to the limit $Q\gg 1$ ($t\ll S^{1/3}$). The time-asymptotic
resistive regime, when resistivity has allowed for the formation of a finite
width boundary layer, corresponds to the limit $Q\ll 1$ ($t\gg S^{1/3}$).
These two stages, which forced reconnection and Alfvén resonance have in
common, are naturally recovered through two asymptotic expansions of the
general solution (28) with respect to the parameter $Q$. The transition
between the ideal and the resistive regimes occurs for $Q\sim 1$, or for times
$t\sim S^{1/3}$. Thus, with Eq. (28) we recover the expected time scales.

To make progress, and to provide a specific solution to discuss how forced
reconnection and Alfvén resonances are recovered from Eq. (28) in the two
aforementioned stages, we will consider a specific boundary forcing given by a
standing wave of the following form,

$\Xi(y,t)=\Xi_{0}\cos(\omega_{0}t)\cos(k_{0}y),$ (33)

of small amplitude $\Xi_{0}$. The Fourier and Laplace transformed boundary
condition is, therefore,

$\hat{\Xi}(k,g)=\frac{\Xi_{0}}{2}\frac{g}{(g^{2}+\omega_{0}^{2})}\left[\delta(k-k_{0})+\delta(k+k_{0})\right].$
(34)

Even though we specialize our results to the forcing in (33), the validity of
our discussion holds in general.

##  V The Ideal Regime

The ideal solution can be recovered from the general solution given by Eq.
(28) by taking the limit $Q\gg 1$ (or $g\gg S^{-1/3}$) and, simultaneously,
$z/Q^{3/2}\ll 1$. The leading order contributions, after imposing the
appropriate parity for the stream function, yield

$\displaystyle\tilde{\phi}(\theta,k,g)=\frac{\hat{\psi}_{0}}{2}\frac{e^{-\frac{g}{|k|}\theta}}{\theta}.$
(35)

Under this limit, $\Delta_{0}\ll\Delta_{g}$ and is thus neglected. Then,
$\hat{\psi}_{0}$ takes the form

$\hat{\psi}_{0}=\hat{\Xi}\frac{\Delta_{f}}{\Delta_{g}}$

with

$\Delta_{g}=-\pi\frac{k}{g}.$

The inverse Fourier transform with respect to $\theta$ of Eq. (35) is:

$\hat{\phi}(x,k,g)=i\hat{\psi}_{0}\arctan\left(\frac{|k|}{g}x\right).$ (36)

The solution in Eq. (36) just recovered can be derived also by solving Eq.
(20)-(21) by imposing $S=\infty$ (e.g., [6]). The branch cuts and poles of
this solution, assuming the forcing given by the standing wave in Eq. (34), is
reported in Figure 1. The poles introduced by the boundary condition
$\hat{\Xi}(k,g)$, in green, provide the long-time finite frequency response
which drives the Alfvén resonance. The orange branch cut arises from the
complex arctangent contribution and is responsible for the forced reconnection
solution when $\omega_{0}\rightarrow 0$.

For the standing wave boundary condition, the ideal stream function can be
inverted analytically by following the Bromwich contour represented in Figure
1. The solution in terms of the variables $x$, $y$, and $t$ is given by

$\phi=-\Phi\left\\{\cos(\omega_{0}t)\left[CI(z_{+})-CI(z_{-})+\ln\left(\frac{z_{-}}{z_{+}}\right)\right]-\frac{2\sin(k_{0}xt)}{\omega_{0}t}+\sin(\omega_{0}t)\left[SI(z_{+})+SI(z_{-})\right]\right\\},$
(37)

where

$z_{\pm}=(k_{0}x\pm\omega_{0})t,\quad\Phi=\frac{\Xi_{0}}{\pi}\frac{\omega_{0}\sin[k_{0}y]}{\sinh[k_{0}]},$
(38)

and $SI(z_{\pm})$, $CI(z_{\pm})$ are the sine and cosine integral functions,
respectively. We see that for $\omega_{0}\neq 0$ we recover the expected
logarithmic singularity at the resonant site, where the phase velocity of the
driven boundary perturbation, $v_{ph}=\omega_{0}/k_{0}$, is equal to the local
Alfvén speed. The logarithmic contribution to the stream function originates
from the residues at the poles $g=\pm\omega_{0}$, while the sine and cosine
integrals arise from evaluating the jump of the Laplace integral across the
branch cut. For long times, the contribution of the logarithmic response to
the resonant poles dominates over that of the sine and cosine integrals, which
initially provide much of the oscillatory behavior of the solution. In the
limit $\omega_{0}\rightarrow 0$, Eq. (37) describes forced reconnection at the
neutral line and recovers the results of [2].

Figure 1: The branches and poles of the ideal solution (36) for the standing
wave boundary forcing (34) are shown here in the complex plane. In black is
the Bromwich contour for the inverse Laplace transform. The poles introduced
by a finite frequency forcing determine the resonant response and are
represented in green. The branch cut of the complex $\arctan$ is represented
in orange and it contributes to the ideal stage of the forced reconnection
solution.

The ideal current density can be found either from $\eqref{eq:Current Density
From Phi}$ in the limit of infinite $S$ or from the asymptotic expansion of
(27). In either case, the ideal current density for a generic boundary
condition is

$\tilde{J}_{z}=\frac{\psi_{0}}{2}\frac{k}{g}\left(\frac{g}{|k|}\theta+1\right)e^{-\frac{g}{|k|}\theta},$
(39)

and in the case of the standing wave boundary condition the ideal current
density can be inverted to obtain

$J_{z}=-\Xi_{0}\frac{2}{\pi}\frac{k_{0}^{4}\cos(k_{0}y)}{\sinh(k_{0})}\frac{x^{2}t}{(k_{0}^{2}x^{2}-\omega_{0}^{2})}\left[\cos(k_{0}xt)+\frac{(k_{0}^{2}x^{2}-3\omega_{0}^{2})}{(k_{0}^{2}x^{2}-\omega_{0}^{2})}\frac{\sin(k_{0}xt)}{k_{0}xt}+\frac{2\omega_{0}^{3}}{k_{0}^{2}x^{2}t}\frac{\sin(\omega_{0}t)}{(k_{0}^{2}x^{2}-\omega_{0}^{2})}\right].$
(40)

For $t\ll 1$ Eq. (40) is given by

$J_{z}=-\Xi_{0}\frac{4}{\pi}\frac{k_{0}^{2}\cos(k_{0}y)}{\sinh(k_{0})}t\left[1+\frac{1}{6}(k_{0}^{2}x^{2}+\omega_{0}^{2})t^{2}\right]+\mathcal{O}(t^{4}).$
(41)

As can be seen from Eq. (41), the current density’s amplitude increases
linearly with time regardless of the forcing frequency to lowest order. As
time goes by, the more general Eq. (40) recovers the result of [2] for
$\omega_{0}=0$, where the amplitude of $J_{z}$ increases linearly with time at
$x=0$, and it is localized in a region whose width decreases as $1/t$. Outside
the origin, the current density’s amplitude still increases linearly with time
but its spatial behavior is highly oscillatory. For $\omega_{0}\neq 0$, the
current density at the origin is a purely oscillating function. For long
times, however, it grows faster than linear near the resonant point where it
builds up by eventually leading to a singularity.

##  VI The Resistive Regime

We separate now the discussion of the time-asymptotic behavior, corresponding
to $Q\ll 1$, into two cases: the case of small or zero frequency,
$\omega_{0}\ll S^{-1/3}$, for which we expect to describe forced magnetic
reconnection, and the case of high frequency $\omega_{0}\gg S^{-1/3}$.

###  VI.1 Forced Reconnection

In the limit $Q\ll 1$ (or $g\ll S^{-1/3}k^{2/3}$) Eq. (28) with Eq. (30)
approximates to

$\tilde{\phi}=-\frac{\Gamma[\frac{3}{4}]}{2\pi}\psi_{0}\frac{g}{k}\left(\frac{1}{|k|}\sqrt{\frac{g}{S}}\theta^{2}\right)^{1/4}K_{\frac{1}{4}}\left(\frac{1}{|k|}\sqrt{\frac{g}{S}}\frac{\theta^{2}}{2}\right),$
(42)

where $K_{\nu}$ is the modified Bessel function of the second kind. This is
colloquially known as the constant-psi regime, or phase D of [2]. In this
limit, both $\Delta_{g}$ and $\Delta_{0}$ are retained in $\hat{\psi}_{0}$
with $\Delta_{g}$ approximated to

$\Delta_{g}=\sqrt{2}(\frac{g^{5}S^{3}}{k^{2}})^{1/4}\Gamma(3/4)^{2}.$ (43)

From the definition of $\Delta_{g}$, we observe the time scale $t\sim S^{3/5}$
after which reconnection reaches a steady state [2, 3]. The flux function at
the origin $\hat{\psi}_{0}$, given below, is consistent with that found by [2,
3], albeit for different boundary conditions,

$\displaystyle\hat{\psi}_{0}$ $\displaystyle=$
$\displaystyle\frac{\frac{2k}{\sinh(k)}\hat{\Xi}}{\sqrt{2}(\frac{g^{5}S^{3}}{k^{2}})^{1/4}\Gamma(3/4)^{2}+2\frac{k}{\tanh(k)}}.$
(44)

Figure 2: The branches and poles of the constant-psi solution (45) are shown
here in the complex plane. In green are the resonant poles of the boundary
condition, in blue the singular contributions of the $g^{5/4}$ root, and in
red is the branch of the Modified Bessel and Struve functions which extends
from $g=-\infty$ to $g=0$.

The stream function in Eq. (42) can be inverted from Fourier space to
configuration space yielding, for the standing wave boundary condition, the
following result,

$\hat{\phi}=\Xi_{0}\sqrt{\frac{\pi}{2}}\Gamma(\frac{3}{4})\frac{\sin(k_{0}y)}{\sinh(k_{0})}\frac{g^{2}}{(g^{2}+\omega_{0}^{2})}\left(k_{0}\sqrt{\frac{S}{g}}\right)^{3/4}\sqrt{x}\left(\frac{I_{1/4}\left(k_{0}\sqrt{\frac{S}{g}}\frac{x^{2}}{2}\right)-L_{1/4}\left(k_{0}\sqrt{\frac{S}{g}}\frac{x^{2}}{2}\right)}{\sqrt{2}(\frac{g^{5}S^{3}}{k_{0}^{2}})^{1/4}\Gamma(3/4)^{2}+2\frac{k_{0}}{\tanh(k_{0})}}\right),$
(45)

where $I_{\nu}$ is the modified Bessel function of the first kind and
$L_{\nu}$ is the modified Struve function. The branch and poles of this
solution with the Bromwich contour are sketched in Figure 2. The green poles
represent once more the resonant contributions of the forcing that merge with
the branch cut, represented in red, in the limit of $\omega_{0}\rightarrow 0$.
The blue poles are the roots of the denominator in the last term in Eq. (45).
It is noted that we recover the same branches and poles as [2], but in their
case the resonant poles of the boundary overlap with the branch cut at $g=0$.

To find the time dependent stream function in this regime, we use Mathematica
to numerically take the inverse Laplace transform of Eq. (45) about the
Bromwich contour just discussed. In the limit $\omega_{0}\rightarrow 0$, the
contributions to the integral are the residues of Eq. (45) evaluated at the
$g^{5/4}$ poles (the blue poles in Fig. 2) and the jump across the branch
extending from $g=-\infty$ to $g=0$. For small but non zero frequencies
$\omega_{0}$, we must include also the contribution of the residues at $g=\pm
i\omega_{0}$, which add an oscillatory contribution to the solution that is
non-negligible for $t\rightarrow\infty$.

The resulting stream function $\phi(x,y,t)$ is compared with results from a
numerical simulation in Figure 3, where dashed lines correspond to the
theoretical solutions and solid lines correspond to results from ab-initio
numerical simulations. Details on the numerical code can be found in the
Appendix. The constant-psi solution is represented by the purple, cyan and
green colors for $\omega_{0}=0$, $\omega_{0}=0.0002$, and $\omega_{0}=0.002$,
respectively, where we have fixed $y=\pi/2$ and $t=3000$, which is
sufficiently long compared to the resistive time scale $t\sim S^{1/3}$. The
chosen frequencies all satisfy the condition $\omega_{0}\ll S^{-1/3}$ and, as
can be seen, the solution corresponding to the constant-psi regime matches
very well with the simulation results for each frequency.

Figure 3: The time-asymptotic solution for the stream function $\phi(x,y,t)$
(dashed lines) compared with numerical simulations (solid lines). The dashed
green, cyan, and purple lines are the inverse Laplace transformed constant-psi
solutions, Eq. (45), for $\omega_{0}=0.002$, $\omega_{0}=0.0002$, and
$\omega_{0}=0$, respectively. These solutions are best seen in the inset plot,
showing an enlarged picture of the region around the neutral line. The dashed
blue and red lines are the inverse Laplace transformed Alfvén resonance
solutions, Eq. (51), for $\omega_{0}=0.1$ and $\omega_{0}=0.05$, respectively.
The constant-psi solution for $\omega_{0}=0.1$ has been also plotted (dotted
line) to demonstrate the breakdown of the constant-psi solution as the system
transitions towards Alfvén resonance for increasing frequencies.

The current density is found by taking the inverse Fourier transform of (22)
in the limit $Q\ll\theta^{2}$

$\hat{J}_{z}=-(Sk)\int^{\infty}_{-\infty}e^{ix\theta}\frac{d}{d\theta}\tilde{\phi},$
(46)

which is identically

$\hat{J}_{z}=-i(Sk)x\hat{\phi}(x,k,g),$ (47)

from the definition of the Fourier transform. By using the appropriate
$\hat{\phi}$ for this limit for constant-psi (from equation (45) in Fourier
$k$ space) $\hat{J_{z}}$ is given by

$\hat{J}_{z}=\Xi_{0}\sqrt{\frac{\pi}{2}}\Gamma(\frac{3}{4})\frac{k_{0}\cos(k_{0}y)}{\sinh(k_{0})}\frac{g^{2}S}{(g^{2}+\omega_{0}^{2})}\left(k_{0}\sqrt{\frac{S}{g}}x^{2}\right)^{3/4}\left(\frac{I_{1/4}\left(k_{0}\sqrt{\frac{S}{g}}\frac{x^{2}}{2}\right)-L_{1/4}\left(k_{0}\sqrt{\frac{S}{g}}\frac{x^{2}}{2}\right)}{\sqrt{2}(\frac{g^{5}S^{3}}{k_{0}^{2}})^{1/4}\Gamma(3/4)^{2}+2\frac{k_{0}}{\tanh(k_{0})}}\right).$
(48)

We note that the current density has the same branches and poles as the stream
function in Eq. $\eqref{eq:Phi Const g}$, outlined in Figure 2. The current
$J_{z}$ can thus be found as a function of time by computing the inverse
Laplace transform over the same contour as we did for $\phi$.

In principle, Eq. (48) describes also the current density at the origin under
an appropriate limit by expanding the Modified Bessel and Struve functions,
$I_{1/4}$ and $L_{1/4}$ respectively, for small arguments. For simplicity, to
find the current density at $x=0$ required to calculate the reconnection rate,
we follow the method adopted in other works by employing the relation

$\displaystyle J_{z}(0,y,t)=S\frac{\partial\psi_{0}}{\partial
t}=S\mathcal{L}^{-1}(g\hat{\psi}_{0}),$ (49)

obtainable from (46) by imposing $x=0$ in the exponential, making use of the
relationship of $\tilde{\phi}$ and $\tilde{\psi}$, Eq. (19), and, finally, by
applying the limit $Q\ll\theta^{2}$. We have also assumed here that
reconnection has been developed sufficiently such that $\psi$ is approximately
constant across the layer, as is the case for the constant-psi approximation.
For the boundary condition employed in this paper, Eq. (49) yields,

$\hat{J}_{z}(0,y,g)=-\Xi_{0}\frac{2k_{0}\cos(k_{0}y)}{\sinh(k_{0})}\frac{g^{2}S}{(g^{2}+\omega_{0}^{2})}\left(\frac{1}{\sqrt{2}(\frac{g^{5}S^{3}}{k_{0}^{2}})^{1/4}\Gamma(3/4)^{2}+2\frac{k_{0}}{\tanh(k_{0})}}\right).$
(50)

We now evaluate numerically its inverse Laplace transform via the Bromwich
contour in Figure 2 (using Mathematica) to calculate the reconnection rate
$R(t)=J_{z}(0,t)/S$ at a given position $y$.

The time evolution of the reconnection rate $R(t)$ in the cases of
$\omega_{0}=0$ and $\omega_{0}=0.002$ is plotted in Figure 4, left and middle
panels, where we compare $R$ from the 2D simulation (orange dashed line) with
the ideal (blue) and constant-psi (green) solutions. The reconnection rate
behaves as expected for the zero frequency case, (left panel): it initially
grows linearly with time as described by the ideal solution and for $t\gg
S^{1/3}$ it follows the constant-psi solution. After times much longer than
the reconnection time scale $t\sim S^{3/5}$, the reconnection rate approaches
zero as the flux function at the origin reaches a steady state. The
reconnection rate has a first peak in the intermediate regime, between the
ideal and constant-psi regimes, when the resistive layer forms at $t\sim
S^{1/3}$. In the case of a small but non zero frequency (middle panel), the
reconnection rate follows an analogous evolution but, over a longer time
scale, instead of reaching a steady state the reconnection rate oscillates in
response to the boundary driver. The fact that the small but non zero
frequency case also gives rise to reconnection isn’t surprising, as discussed
in the Introduction, and was analyzed previously by [25] for the eigenvalue
problem and by [26] in the study of the interaction of resonant magnetic
perturbations with rotating plasmas. In a later section we discuss in more
details the effects of finite frequency on reconnection.

Figure 4:  The left and middle panels show the reconnection rate as a function
of time from numerical simulations compared with the two asymptotic solutions:
the ideal solution, Eq. (40) divided by $S$ in blue, and the inverse Laplace
transform of the constant-psi solution, Eq. (49) in green. We show results for
$\omega_{0}=0$ (left panel) and $\omega_{0}=0.002$ (middle panel). In the
right panel, the simulated current density for $\omega_{0}=0.1$ has been
plotted at the resonant site. It is compared with the same ideal solution
(blue color) and with the resistive solution solution given by the inverse
Fourier transform of Eq. (54) (green color). The parameters used here are
$S=10^{6}$, $k_{0}=1$, $y=\pi/2$.

###  VI.2 Transition to Alfvén Resonance

In the previous subsection, we analyzed the constant-psi regime of Eq. (28)
for small frequencies ($\omega_{0}\ll S^{-1/3}$). As discussed, when
$\omega_{0}\ll S^{-1/3}$ the solution corresponding to the constant-psi regime
matches very well with the simulation results. However, as the frequency is
increased beyond $\omega_{0}\sim S^{-1/3}$, the agreement diverges and the
constant-psi solution breaks down. This is shown for example by the dotted
blue line in Figure 3, that represents the constant-psi solution but for a
high frequency case, $\omega_{0}=0.1$; as can be seen, the dotted line does
not reproduce simulation results, which are represented by the solid blue
line. Therefore, for large frequencies we must consider a different limit to
describe Alfvén resonance.

For large frequencies, $\omega_{0}\gg S^{-1/3}$, the residues at the resonant
poles of the boundary driver, represented in green in Figure 2, contribute to
the time-asymptotic response of the system. The contributions from both the
jump at the branch cut and the residues at the poles located in the region
$\Re({g})<0$ have vanishing contributions for $t\rightarrow\infty$.
Additionally, if the frequency is large enough such that $\omega_{0}\gg
S^{-1/3}k^{2/3}$, one cannot use the asymptotic expansion of Eq. (28) and (30)
for $Q\ll 1$ to evaluate the residues at $g=\pm i\omega_{0}$, as it was done
for small frequency. Instead, one should use the asymptotic expansion for
$Q\gg 1$ to calculate the residues at those poles, except this time by
retaining resistive corrections. As gradients become strong around the
resonant point, the effects of resistivity will be to resolve and smooth out
the logarithmic singularity of the ideal solution. Performing the same
expansion for $Q\gg 1$ as with the ideal case but by retaining higher-order
terms in $z/Q^{3/2}\ll 1$, we find that the stream function can be
approximated by

$\displaystyle\tilde{\phi}(\theta,y,g)=\frac{\hat{\psi}_{0}}{2}\frac{1}{\theta}e^{-\frac{g}{|k|}\theta\left(1+\frac{1}{6}\frac{\theta^{2}}{gS}\right)+\mathcal{O}\left(\frac{\theta^{5}}{gS^{2}}\right)}\left[1+\mathcal{O}\left(\frac{\theta^{2}}{gS}\right)\right],$
(51)

which recovers the ideal solution in the limit of infinite $S$, as well as the
same exponential dependence found in [10] in the context of normal mode
analysis. As just discussed, in the time asymptotic limit, the dominant
contribution is given by the residues evaluated at the resonant poles of the
boundary condition, which in our particular case give

$\tilde{\phi}(\theta,y,t)=i\frac{\Xi_{0}}{\pi}\frac{\omega_{0}\sin(k_{0}y)}{\sinh(k_{0})}\sin\left(\omega_{0}(t-\frac{\theta}{k_{0}})\right)\frac{e^{-\frac{1}{6}\frac{\theta^{3}}{Sk_{0}}}}{\theta}.$
(52)

To find the spatio-temporal stream function we finally used Mathematica to
evaluate the inverse Fourier transform from $\theta$ to $x$ of Eq. (52). The
resulting solution for $\omega_{0}=0.1$ and $\omega_{0}=0.05$ is reported in
Figure 3, dashed blue and red lines, and compared with simulation results
represented by the solid blue and red lines. It is apparent that Eq. (51),
when Fourier-inverted, captures the long-time behavior of Alfvén resonances.
We see that the contribution of the resistivity is both to remove the
logarithmic singularity and to provide a finite width to the resonant boundary
layer. The boundary layer, as can be seen by inspection of Eq. (52), scales as
$\delta\sim S^{-1/3}$.

The current density is evaluated from the expansion of (27) while maintaining
higher orders in the exponential,

$\tilde{J}_{z}(\theta,k,g)=\frac{\psi_{0}}{2}\frac{k}{g}\left(\frac{g}{|k|}\theta+1\right)e^{-\frac{g}{|k|}\theta\left(1+\frac{1}{6}\frac{\theta^{2}}{gS}\right)}.$
(53)

Again we evaluate the residues at the resonant poles of the boundary condition
to find the long time behavior, yielding the following solution to be
numerically inverted from $\theta$ to configuration space,

$\tilde{J}_{z}(\theta,y,t)=-\frac{\Xi_{0}}{\pi}\frac{k_{0}\cos(k_{0}y)}{\sinh(k_{0})}e^{-\frac{1}{6}\frac{\theta^{3}}{Sk_{0}}}\left(\cos(\zeta)-\sin(\zeta)\frac{\omega_{0}}{k_{0}}\theta\right),$
(54)

where $\zeta=\omega_{0}(t-\frac{\theta}{k_{0}})$. The current density at the
resonant site, found through numerical inversion of (54), has been plotted in
the right frame of Figure 4. Just like with the reconnection rate for small
frequencies, we see very good agreement between the ideal solution and the
numerical simulation up to nearly $t\sim S^{1/3}$. Over longer time scales,
the system reaches a steady-state due to resistive dissipation, and the
solution in Eq. (54) provides the correct time-asymptotic behavior where an
oscillating steady-state is reached.

Figures 3 and 4 together demonstrate the transition from Alfvén resonance to
forced magnetic reconnection as $\omega_{0}$ varies from $\omega_{0}\gg
S^{1/3}$ to $\omega_{0}\ll S^{1/3}$. As the frequency approaches zero, the
resonant sites shifts towards the neutral line, since Alfvén resonance
necessarily occurs where the phase velocity of the driver is equal to the
local Alfvén speed, $\omega_{0}/k=v_{a}(x)$. When the resonant layers, which
are located on either side of the neutral line, overlap across the neutral
line as a consequence of their finite widths, the resonance couples with the
reconnection mode. The subsequent section discusses the effect of finite
frequency on such a coupling.

###  VI.3 Discussion: scalings laws and effect of finite frequency on
reconnection

Figure 5: Plots of the reconnection rate and reconnected flux as functions of
$\omega_{0}$ evaluated at the time of maximum reconnection rate in the
constant-psi regime, for $S=10^{8}$. They are calculated from the inverse
Laplace transforms of Eq. (49), divided by $S$, and Eq. (44), which are
applicable for $\omega_{0}\ll S^{-1/3}$. The red (left most) line in each plot
indicates $\omega_{0}=S^{-3/5}$. The orange (middle) line indicates the
location of half maximum half width of $\psi_{0}$. The green (right most) line
indicates $\omega_{0}=S^{-1/3}.$ In the frequency range beyond the half-width
of $\psi_{0}$ but below $\omega_{0}=S^{-1/3}$, we observe a power law scaling
$R\sim\omega_{0}^{-0.25}$ and $\psi_{0}\sim\omega_{0}^{-1.23}$.

To further understand the effect of finite frequency on forced reconnection,
in this section we analyze two sets of scalings predicted by the theory. The
first set is the scaling of the width $\delta$ of the boundary layer for
forced reconnection and Alfvén resonance as functions of the Lundquist number
$S$. The second set is the scaling of the constant-psi reconnection rate and
reconnected flux at the neutral line as functions of the driving frequency
$\omega_{0}$.

Through inspection of the current density equation in the constant-psi regime
which describes forced reconnection, Eq. (48), we find the theoretical scaling
law for $\delta$. For this regime we assume the scaling $g\sim S^{-3/5}$,
found from the poles in the denominator. Upon substitution of $g$ for $S$ in
Eq. (48) we find the argument of the Modified Bessel and Struve functions to
be of order unity when $x\sim S^{-2/5}$. This indicates a layer width scaling
$\delta\sim S^{-2/5}$, which is the known scaling for this process [2].

For Alfvén resonance we perform the same analysis of its equivalent current
density equation, applicable for large driving frequencies $\omega_{0}\gg
S^{-1/3}$, Eq. (54). For this scaling we need only the argument of the
exponential to determine the inner layer width dependence on $S$. We see that
the argument is of order unity when $\theta\sim S^{1/3}$ which is indicative
of a layer width scaling $\delta\sim S^{-1/3}$, as expected [9, 8, 10].

Although it is not possible to provide a theoretical solution for the
intermediate regime when $Q\sim 1$ (which characterizes the region at which
our approximations mismatch with simulation in Figure 4), we can nevertheless
use the general solution given in Eq. (28) to infer the expected $\delta$.
From the definition of the parameter $Q$, Eq. (23), the scaling $g\sim
S^{-1/3}$ naturally emerges. Replacing $g\sim S^{-1/3}$ in Eq. (28) we find
the argument of the Tricomi confluent hypergeometric functions of Eq. (28),
$z=\sqrt{Q}\theta/(Sk)^{2/3}$, is of order 1 for $\theta\sim S^{1/3}$ or
$\delta\sim S^{-1/3}$. The regime $Q\sim 1$, which corresponds to $t\sim
S^{1/3}$, applies to both Alfvén resonance and forced reconnection. In
particular, it shows that forced reconnection starts off in the non-constant-
psi regime, as first discussed in [2], before continuing on through the
constant-psi. Thus, we see that when resistivity becomes important at $t\sim
S^{1/3}$ both forced reconnection and Alfvén resonance exhibit the same layer
width, allowing for their coupling, when the forcing frequency is sufficiently
small.

We make now analytical predictions of how the reconnection rate and the
reconnected flux scale with the forcing frequency $\omega_{0}$. We conduct
this analysis under the constant-psi regime, i.e., $R$ is calculated by using
Eq. (49) divided by $S$, and (44) is used for $\psi_{0}$, that in this
approximation represents the reconnected flux at $x=0$.

For frequencies $\omega_{0}\ll S^{-3/5}$, the system has time to evolve from
the ideal to the constant-psi regime almost unaffected by the oscillating
driver, and we expect that for such low frequencies both $R$ and $\psi_{0}$
are unaffected by $\omega_{0}$. When the frequency approaches the inverse of
the constant-psi time scale, which is around the time when the peak in
reconnection occurs (see Fig. 4), we expect the reconnection process to be
affected by the finite frequency driver. Now, in the frequency range
$S^{-3/5}\ll\omega_{0}\ll S^{-1/3}$ one can neglect the term
$\Delta_{0}=2k_{0}/\tanh(k_{0})$ in the denominator of Eqs. (44) and (49).
Additionally, as demonstrated previously, as $\omega_{0}$ increases the
contribution to the long time solution is provided by residues at $g=\pm
i\omega_{0}$. By analyzing under this light the residues of the inverse
Laplace transforms of Eqs. (49) and (44) at $g=\pm i\omega_{0}$, it can be
shown that the expected scalings are $R\sim\omega_{0}^{-1/4}$ and
$\psi_{0}\sim\omega_{0}^{-5/4}$. This would suggest that although
$\omega_{0}\sim S^{-1/3}$ determines the transition from forced reconnection
to Alfvén resonance, the reconnection rate is reduced in the range of
frequencies $S^{-3/5}\ll\omega_{0}\ll S^{-1/3}$. For frequencies larger than
$\omega_{0}\gg S^{-1/3}$ reconnection is completely suppressed by the Alfvén
resonance.

To measure the effect of finite frequency on reconnection, in Figure 5 we show
the plot of the normalized, maximum reconnection rate $R$ and corresponding
reconnected flux function $\psi_{0}$ as functions of $\omega_{0}$ in the left
and right panel, respectively, for $S=10^{8}$. As explained above, the
constant-psi regime was assumed in making these plots. Therefore, we used Eqs.
(49), divided by $S$, for $R$ and (44) for $\psi_{0}$, and calculated their
inverse Laplace transforms at the time for which the reconnection rate is at
its maximum. With the vertical red, orange and green (left, central, and
right) lines, we mark $\omega_{0}=S^{-3/5}$, the half-width of $\psi_{0}$, and
$\omega_{0}=S^{-1/3}$. As can be seen, the predicted scalings of
$R\sim\omega_{0}^{-1/4}$ and $\psi_{0}\sim\omega_{0}^{-5/4}$ are very closely
matched in Figure 5, in which we found power law scalings of
$R\sim\omega_{0}^{-0.25}$ and $\psi_{0}\sim\omega_{0}^{-1.23}$ for the
reconnection rate and reconnected flux, respectively, shown by the fits with
the dashed line. That is we find, as expected, a strong suppressing effect on
reconnection for frequencies less than the transition frequency
$\omega_{0}\sim S^{-1/3}$ and higher than $\omega_{0}\sim S^{-3/5}$.

##  VII Summary and conclusions

In this paper we have derived and analyzed the general solution to the
resistive, linearized MHD equations (Eqs. 4-5) describing the evolution in
time of a plasma with a sheared magnetic field (Eq. 1) subjected to an
oscillatory boundary perturbation turned on at time $t=0$ (Eqs. 8-9). Through
asymptotic matching, we determined the spatial dependence (Eq. (28)) and the
amplitude (Eq. (30)) of the general solution for the perturbed stream function
in the inner layer. Such a solution shows explicitly that the key parameter to
describe different regimes in the evolution of the system, as well as the
transition from forced reconnection to Alfvén resonance, is
$Q=gS^{1/3}k^{-2/3}$. We have then discussed two asymptotic limits of Eq.
(28), $Q\gg 1$ and $Q\ll 1$, and we demonstrated that our theory describes
Alfvén resonance when the frequency of the forcing is $\omega_{0}\gg
S^{-1/3}$, and forced reconnection for small to zero frequency, $\omega_{0}\ll
S^{-1/3}$. Spatio-temporal solutions have been derived for a standing wave
boundary condition, given in Eq. (33), which we compared with linear
simulations.

In the ideal limit, defined by $Q\gg 1$ and $z/Q^{3/2}\ll 1$ (up to first
order), our general solution recovers that of forced reconnection for
$\omega_{0}=0$; for $\omega_{0}\neq 0$, our solution describes Alfvén
resonance. In particular, the known logarithmic singularity at the resonance
point where $\omega_{0}/k=v_{a}(x)$ is due to the residues at the poles of the
boundary driver, and contribute to the time-asymptotic response of the system,
as can be seen from the Bromwich contour in Figure 1.

The non-ideal time-asymptotic response of the system corresponds to $Q\ll 1$.
We have provided the corresponding asymptotic expansion of the inner layer
stream function for $Q\ll 1$, by recovering the constant-psi solution that
describes forced reconnection for times $t\gg S^{1/3}$. The corresponding
Bromwich contour is shown in Figure 2.

For large frequencies, the time-asymptotic response is dominated by the
residues at the forcing frequency (see Fig. 2), and if the condition
$\omega_{0}S^{1/3}k^{-2/3}\gg 1$ is satisfied, then a different expansion of
the stream function for their evaluation is required. In this case, the
solution is given by the asymptotic expansion of the inner layer solution for
$Q\gg 1$ with terms in $z/Q^{3/2}\ll 1$ retained to higher order. This regime
represents Alfvén resonance in the time-asymptotic limit.

A comparison of the theoretical predictions and numerical simulations is given
in Figure 3 and 4, testing our model’s ability to capture both forced
reconnection and Alfvén resonance. We show that for large frequencies the
resonant poles of the boundary driver were sufficiently far from zero
(accordingly, the resonant points are far from the neutral point) and a pair
of anti-symmetric resonances, of width $\delta\sim S^{-1/3}$, are clearly
defined in that case. As these poles are brought towards $\omega_{0}=0$ (the
resonant points approach the neutral point $x=0$) in the small frequency
limit, the widths of the resonant layers overlap across $x=0$ by coupling to
the reconnection mode. The rate of reconnection and its signatures are
strongest when the poles overlap entirely, that is for $\omega_{0}=0$.

Finally, we demonstrated that our solution recovers known scalings of the
resonant layer $\delta$ as a function of the Lundquist number $S$ under
specific limits. Additionally, we demonstrated new scalings of the
reconnection rate and reconnected flux as functions of the driving frequency
of the boundary condition. In Figure 5 we show that, in the constant-psi
regime, the maximum reconnection rate and reconnected flux are largest for
frequencies less than $\omega_{0}\lesssim S^{-3/5}$ and is reduced at larger
frequencies until reconnection is completely suppressed for $\omega_{0}\gg
S^{-1/3}$.

We conclude that forced reconnection is indeed a particular case of Alfvén
resonant absorption and there is a transition from Alfvén resonance to forced
reconnection when the forcing frequency $\omega_{0}\ll S^{-1/3}$. However, the
reconnection rate and reconnected flux are strongly reduced for $\omega_{0}\gg
S^{-3/5}$. When the forcing frequency is $\omega_{0}\gg S^{-1/3}$, Alfvén
resonance decouples from the reconnecting mode and the formation of resonant
layers away from the neutral line effectively shields reconnection.

Here we have studied the interaction between Alfvén resonances and forced
reconnection by considering an idealized symmetric forcing away from the
neutral line. However, our approach is general and provides the first step for
further investigations of resonances in diverse environments and within
different plasma models. In particular, more realistic perturbations can be
considered, such as asymmetric and impulsive. This type of forcing is
especially of interest for reconnection onset at Earth’s magnetopause. In this
regard, recent observations from the MMS spacecraft within Earth’s
magnetosphere [27] have identified reconnection events taking place at
electron scales, termed electron-only reconnection [28, 29]. By building on
prior investigations into forced reconnection within the EMHD framework [30],
our approach can be extended to investigate the effect of a time-dependent
forcing on electron-only reconnection, with applications to reconnection onset
within Earth’s magnetosphere. Additionally, our unified approach holds
potential implications for studying the transition from wave-based to
reconnection-based coronal heating scenarios [31, 14], which will be explored
in future work.

###### Acknowledgements.

This research was supported by NSF grant 2108320. One of us (AT) was also
supported by NSF CAREER grant 2141564 and another (FW) by DOE grant DE-
FG02-04ER54742. We also acknowledge the Texas Advanced Computing Center (TACC)
at The University of Texas at Austin for providing HPC resources that have
contributed to the research results reported within this paper.

##  Appendix A Numerical Simulations

The numerical code used to validate our analytical solutions is a 2D
incompressible MHD code. The normalized equations that we integrated are:

$\frac{\partial\psi}{\partial t}=x\frac{d\phi}{dy}+\frac{1}{S}\nabla^{2}\psi$
(55) $\frac{\partial\Omega}{\partial t}=-x\frac{d}{dy}\nabla^{2}\psi$ (56)

where $\Omega=-\nabla^{2}\phi$ is the vorticity. We impose a boundary
condition at $x=\pm 1$ with small amplitude given by:

$\psi(x=\pm
1,y,t)=\Xi_{0}\cos(k_{0}y)\cos(\omega_{0}t)\left[1-\left(1+\frac{t}{\tau}\right)e^{-t/\tau}\right],$
(57)

$\phi(x=\pm
1,y,t)=\pm\Xi_{0}\frac{\sin(k_{0}y)}{k_{0}}\left\\{\left(\frac{t}{\tau}\right)^{2}e^{-t/\tau}\cos(\omega_{0}t)-\omega_{0}\sin(\omega_{0}t)\left[1-\left(1+\frac{t}{\tau}\right)e^{-t/\tau}\right]\right\\},$
(58)

with $\tau$ being a small growth parameter, following [3]. The code is
periodic in the $y$ direction, where the Fast Fourier Transform is used to
calculate spatial derivatives. The $x$ direction is non-periodic and a finite
difference scheme of the sixth order is used to calculate spatial derivatives
[32]. We use an explicit Runge-Kutta of the third order to advance in time. We
impose the 2/3 rule for dealiasing with respect to the periodic variable, and
a pseudo-spectral filter in $x$ [32]. Simulations are performed by using a
grid $L_{x}\times L_{y}=[-1,1]\times[0,2\pi]$ with mesh points $N_{x}\times
N_{y}=512\times 32$ for $S=10^{4}$, $N_{x}\times N_{y}=1024\times 32$ for
$S=10^{5}$, $N_{x}\times N_{y}=2048\times 16$ for $S=10^{6}$, and $N_{x}\times
N_{y}=4096\times 16$ for $S=10^{7}$, with $N_{x}$ chosen to resolve their
respective boundary layers. In all simulations we have fixed $\tau=0.001$,
$\Xi_{0}=0.001$, and wavenumber $k_{0}=1$ while varying $S$ and the frequency
$\omega_{0}$.

## References

  * [1] Chanchal Uberoi and Ellen G Zweibel.  Alfvén resonances and forced reconnection.  Journal of plasma physics, 62(3):345–350, 1999.
  * [2] TS Hahm and RM Kulsrud.  Forced magnetic reconnection.  The Physics of fluids, 28(8):2412–2418, 1985.
  * [3] Andrew Cole and Richard Fitzpatrick.  Forced magnetic reconnection in the inviscid taylor problem.  Physics of plasmas, 11(7):3525–3529, 2004.
  * [4] Xiaogang Wang and A Bhattacharjee.  Forced reconnection and current sheet formation in taylor’s model.  Physics of Fluids B: Plasma Physics, 4(7):1795–1799, 1992.
  * [5] Richard Fitzpatrick.  A numerical study of forced magnetic reconnection in the viscous taylor problem.  Physics of Plasmas, 10(6):2304–2312, 2003.
  * [6] Luca Comisso, Daniela Grasso, and François L Waelbroeck.  Extended theory of the taylor problem in the plasmoid-unstable regime.  Physics of Plasmas, 22(4), 2015.
  * [7] J Tataronis and W Grossmann.  Decay of mhd waves by phase mixing: I. the sheet-pinch in plane geometry.  Zeitschrift für Physik A Hadrons and nuclei, 261(3):203–216, 1973.
  * [8] JM Kappraff and JA Tataronis.  Resistive effects on alfvén wave heating.  Journal of Plasma Physics, 18(2):209–226, 1977.
  * [9] Y Mok and G Einaudi.  Resistive decay of alfvén waves in a non-uniform plasma.  Journal of plasma physics, 33(2):199–208, 1985.
  * [10] G Bertin, G Einaudi, and F Pegoraro.  Alfvén modes in inhomogeneous plasmas.  Comments Plasma Phys. Controll. Fus.;(United Kingdom), 10(4), 1986.
  * [11] Chris C Hegna, James D Callen, and Robert J LaHaye.  Dynamics of seed magnetic island formation due to geometrically coupled perturbations.  Physics of Plasmas, 6(1):130–136, 1999.
  * [12] Liu Chen and Akira Hasegawa.  Plasma heating by spatial resonance of alfvén wave.  The Physics of Fluids, 17(7):1399–1403, 1974.
  * [13] Joseph M Davila.  Heating of the solar corona by the resonant absorption of alfvén waves.  Astrophysical Journal, Part 1 (ISSN 0004-637X), vol. 317, June 1, 1987, p. 514-521. NASA-supported research., 317:514–521, 1987.
  * [14] Francesco Malara and M Velli.  Wave-based heating mechanisms for the solar corona.  In IAU Colloq. 144: Solar Coronal Structures, pages 443–451, 1994.
  * [15] Akira Hasegawa and Chanchal Uberoi.  The Alfvén Wave: Prepared for the Office of Fusion Energy, Office of Energy Research, US Department of Energy, volume 2.  Technical Information center, US Department of energy, 1982.
  * [16] C Uberoi.  Resonant absorption of alfvén waves near a neutral point.  Journal of plasma physics, 52(2):215–221, 1994.
  * [17] Wolfram Research, Inc.  Mathematica, Version 14.0.  Champaign, IL, 2024.
  * [18] Harold P Furth, John Killeen, and Marshall N Rosenbluth.  Finite-resistivity instabilities of a sheet pinch.  The physics of Fluids, 6(4):459–484, 1963.
  * [19] Francesco Pegoraro, F Porcelli, and TJ Schep.  Internal kink modes in the ion-kinetic regime.  Physics of Fluids B: Plasma Physics, 1(2):364–374, 1989.
  * [20] Francesco Porcelli.  Viscous resistive magnetic reconnection.  The Physics of fluids, 30(6):1734–1742, 1987.
  * [21] Francesco Pegoraro and TJ Schep.  Theory of resistive modes in the ballooning representation.  Plasma physics and controlled fusion, 28(4):647, 1986.
  * [22] Milton Abramowitz, Irene A Stegun, and Robert H Romer.  Handbook of mathematical functions with formulas, graphs, and mathematical tables, 1988.
  * [23] Marshall N. Rosenbluth, R. Y. Dagazian, and P. H. Rutherford.  Nonlinear properties of the internal m = 1 kink instability in the cylindrical tokamak.  Physics of Fluids, 16(11):1894–1902, Nov 1973.
  * [24] B. Coppi, R. Galvao, R. Pellat, M. Rosenbluth, and P. Rutherford.  Resistive internal kink modes.  Sov. J. Plasma Phys., 2(6):533–535, November 1976.
  * [25] Q Luan and X Wang.  On transition from alfvén resonance to forced magnetic reconnection.  Physics of Plasmas, 21(7), 2014.
  * [26] R Fitzpatrick and TC Hender.  The interaction of resonant magnetic perturbations with rotating plasmas.  Physics of Fluids B: Plasma Physics, 3(3):644–673, 1991.
  * [27] JL Burch, RB Torbert, TD Phan, L-J Chen, TE Moore, RE Ergun, JP Eastwood, DJ Gershman, PA Cassak, MR Argall, et al.  Electron-scale measurements of magnetic reconnection in space.  Science, 352(6290):aaf2939, 2016.
  * [28] TD Phan, Jonathan P Eastwood, MA Shay, JF Drake, BU Ö Sonnerup, Masaki Fujimoto, PA Cassak, M Øieroset, JL Burch, RB Torbert, et al.  Electron magnetic reconnection without ion coupling in earth’s turbulent magnetosheath.  Nature, 557(7704):202–206, 2018.
  * [29] SY Huang, QY Xiong, LF Song, J Nan, ZG Yuan, K Jiang, XH Deng, and L Yu.  Electron-only reconnection in an ion-scale current sheet at the magnetopause.  The Astrophysical Journal, 922(1):54, 2021.
  * [30] K Avinash, SV Bulanov, T Esirkepov, P Kaw, Francesco Pegoraro, PV Sasorov, and A Sen.  Forced magnetic field line reconnection in electron magnetohydrodynamics.  Physics of Plasmas, 5(8):2849–2860, 1998.
  * [31] M Velli, F Pucci, F Rappazzo, and A Tenerani.  Models of coronal heating, turbulence and fast reconnection.  Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 373(2042):20140262, 2015.
  * [32] Sanjiva K Lele.  Compact finite difference schemes with spectral-like resolution.  Journal of computational physics, 103(1):16–42, 1992.


11institutetext: 1University College London, London, United Kingdom
2RayCloud Inc., Hangzhou, China
11email<EMAIL_ADDRESS><EMAIL_ADDRESS>
# Deep Learning over Multi-field Categorical Data

– A Case Study on User Response Prediction

Weinan Zhang 11 Tianming Du 1122 Jun Wang 11

###### Abstract

Predicting user responses, such as click-through rate and conversion rate, are
critical in many web applications including web search, personalised
recommendation, and online advertising. Different from continuous raw features
that we usually found in the image and audio domains, the input features in
web space are always of multi-field and are mostly discrete and categorical
while their dependencies are little known. Major user response prediction
models have to either limit themselves to linear models or require manually
building up high-order combination features. The former loses the ability of
exploring feature interactions, while the latter results in a heavy
computation in the large feature space. To tackle the issue, we propose two
novel models using deep neural networks (DNNs) to automatically learn
effective patterns from categorical feature interactions and make predictions
of users’ ad clicks. To get our DNNs efficiently work, we propose to leverage
three feature transformation methods, i.e., factorisation machines (FMs),
restricted Boltzmann machines (RBMs) and denoising auto-encoders (DAEs). This
paper presents the structure of our models and their efficient training
algorithms. The large-scale experiments with real-world data demonstrate that
our methods work better than major state-of-the-art models.

##  1 Introduction

User response (e.g., click-through or conversion) prediction plays a critical
part in many web applications including web search, recommender systems,
sponsored search, and display advertising. In online advertising, for
instance, the ability of targeting individual users is the key advantage
compared to traditional offline advertising. All these targeting techniques,
essentially, rely on the system function of predicting whether a specific user
will think the potential ad is “relevant”, i.e., the probability that the user
in a certain context will click a given ad [6]. Sponsored search, contextual
advertising, and the recently emerged real-time bidding (RTB) display
advertising all heavily rely on the ability of learned models to predict ad
click-through rates (CTR) [32, 41]. The applied CTR estimation models today
are mostly linear, ranging from logistic regression [32] and naive Bayes [14]
to FTRL logistic regression [28] and Bayesian probit regression [12], all of
which are based on a huge number of sparse features with one-hot encoding [1].
Linear models have advantages of easy implementation, efficient learning but
relative low performance because of the failure of learning the non-trivial
patterns to catch the interactions between the assumed (conditionally)
independent raw features [12]. Non-linear models, on the other hand, are able
to utilise different feature combinations and thus could potentially improve
estimation performance. For example, factorisation machines (FMs) [29] map the
user and item binary features into a low dimensional continuous space. And the
feature interaction is automatically explored via vector inner product.
Gradient boosting trees [38] automatically learn feature combinations while
growing each decision/regression tree. However, these models cannot make use
of all possible combinations of different features [20]. In addition, many
models require feature engineering that manually designs what the inputs
should be. Another problem of the mainstream ad CTR estimation models is that
most prediction models have shallow structures and have limited expression to
model the underlying patterns from complex and massive data [15]. As a result,
their data modelling and generalisation ability is still restricted.

Deep learning [25] has become successful in computer vision [22], speech
recognition [13], and natural language processing (NLP) [19, 33] during recent
five years. As visual, aural, and textual signals are known to be spatially
and/or temporally correlated, the newly introduced unsupervised training on
deep structures [18] would be able to explore such _local_ dependency and
establish a _dense_ representation of the feature space, making neural network
models effective in learning high-order features directly from the raw feature
input. With such learning ability, deep learning would be a good candidate to
estimate online user response rate such as ad CTR. However, most input
features in CTR estimation are of multi-field and are discrete categorical
features, e.g., the user location city (London, Paris), device type (PC,
Mobile), ad category (Sports, Electronics) etc., and their local dependencies
(thus the sparsity in the feature space) are unknown. Therefore, it is of
great interest to see how deep learning improves the CTR estimation via
learning feature representation on such large-scale multi-field discrete
categorical features. To our best knowledge, there is no previous literature
of ad CTR estimation using deep learning methods thus far111Although the
leverage of deep learning models on ad CTR estimation has been claimed in
industry (e.g., [42]), there is no detail of the models or implementation.. In
addition, training deep neural networks (DNNs) on a large input feature space
requires tuning a huge number of parameters, which is computationally
expensive. For instance, unlike image and audio cases, we have about 1 million
binary input features and 100 hidden units in the first layer; then it
requires 100 million links to build the first layer neural network.

In this paper, we take ad CTR estimation as a working example to study deep
learning over a large multi-field categorical feature space by using embedding
methods in both supervised and unsupervised fashions. We introduce two types
of deep learning models, called Factorisation Machine supported Neural Network
(FNN) and Sampling-based Neural Network (SNN). Specifically, FNN with a
supervised-learning embedding layer using factorisation machines [31] is
proposed to efficiently reduce the dimension from sparse features to dense
continuous features. The second model SNN is a deep neural network powered by
a sampling-based restricted Boltzmann machine (SNN-RBM) or a sampling-based
denoising auto-encoder (SNN-DAE) with a proposed negative sampling method.
Based on the embedding layer, we build multiple layers neural nets with full
connections to explore non-trivial data patterns. Our experiments on multiple
real-world advertisers’ ad click data have demonstrated the consistent
improvement of CTR estimation from our proposed models over the state-of-the-
art ones.

##  2 Related Work

Click-through rate, defined as the probability of the ad click from a specific
user on a displayed ad, is essential in online advertising [39]. In order to
maximise revenue and user satisfaction, online advertising platforms must
predict the expected user behaviour for each displayed ad and maximise the
expectation that users will click. The majority of current models use logistic
regression based on a set of sparse binary features converted from the
original categorical features via one-hot encoding [26, 32]. Heavy engineering
efforts are needed to design features such as locations, top unigrams,
combination features, etc. [15].

Embedding very large feature vector into low-dimensional vector spaces is
useful for prediction task as it reduces the data and model complexity and
improves both the effectiveness and the efficiency of the training and
prediction. Various methods of embedding architectures have been proposed [37,
23]. Factorisation machine (FM) [31], originally proposed for collaborative
filtering recommendation, is regarded as one of the most successful embedding
models. FM naturally has the capability of estimating interactions between any
two features via mapping them into vectors in a low-rank latent space.

Deep Learning [2] is a branch of artificial intelligence research that
attempts to develop the techniques that will allow computers to handle complex
tasks such as recognition and prediction at high performance. Deep neural
networks (DNNs) are able to extract the hidden structures and intrinsic
patterns at different levels of abstractions from training data. DNNs have
been successfully applied in computer vision [40], speech recognition [8] and
natural language processing (NLP) [7, 19, 33]. Furthermore, with the help of
unsupervised pre-training, we can get good feature representation which guides
the learning towards basins of attraction of minima that support better
generalisation from the training data [10]. Usually, these deep models have
two stages in learning [18]: the first stage performs model initialisation via
unsupervised learning (i.e., the restricted Boltzmann machine or stacked
denoising auto-encoders) to make the model catch the input data distribution;
the second stage involves a fine tuning of the initialised model via
supervised learning with back-propagation. The novelty of our deep learning
models lies in the first layer initialisation, where the input raw features
are high dimensional and sparse binary features converted from the raw
categorical features, which makes it hard to train traditional DNNs in large
scale. Compared with the word-embedding techniques used in NLP [19, 33], our
models deal with more general multi-field categorical features without any
assumed data structures such as word alignment and letter-n-gram etc.

##  3 DNNs for CTR Estimation given Categorical Features

Figure 1: A 4-layer FNN model structure.

In this section, we discuss the two proposed DNN architectures in detail,
namely Factorisation-machine supported Neural Networks (FNN) and Sampling-
based Neural Networks (SNN). The input categorical features are field-wise
one-hot encoded. For each field, e.g., city, there are multiple units, each of
which represents a specific value of this field, e.g., city=London, and there
is only one positive (1) unit, while all others are negative (0). The encoded
features, denoted as $\boldsymbol{x}$, are the input of many CTR estimation
models [32, 26] as well as our DNN models, as depicted at the bottom layer of
Figure 1.

###  3.1 Factorisation-machine supported Neural Networks (FNN)

Our first model FNN is based on the factorisation machine as the bottom layer.
The network structure is shown in Figure 1. With a top-down description, the
output unit is a real number $\hat{y}\in(0,1)$ as predicted CTR, i.e., the
probability of a specific user clicking a given ad in a certain context:

$\displaystyle\hat{y}=\text{sigmoid}(\boldsymbol{W}_{3}\boldsymbol{l}_{2}+b_{3}),$
(1)

where $\text{sigmoid}(x)=1/(1+e^{-x})$ is the logistic activation function,
$\boldsymbol{W}_{3}\in\mathbb{R}^{1\times L}$, $b_{3}\in\mathbb{R}$ and
$\boldsymbol{l}_{2}\in\mathbb{R}^{L}$ as input for this layer. The calculation
of $\boldsymbol{l}_{2}$ is

$\displaystyle\boldsymbol{l}_{2}=\text{tanh}(\boldsymbol{W}_{2}\boldsymbol{l}_{1}+\boldsymbol{b}_{2}),$
(2)

where $\tanh(x)=(1-e^{-2x})/(1+e^{-2x})$,
$\boldsymbol{W}_{2}\in\mathbb{R}^{L\times M}$,
$\boldsymbol{b}_{2}\in\mathbb{R}^{L}$ and
$\boldsymbol{l}_{1}\in\mathbb{R}^{M}$. We choose $\tanh(\cdot)$ as it has
optimal empirical learning performance than other activation functions, as
will be discussed in Section 4.3. Similarly,

$\displaystyle\boldsymbol{l}_{1}=\text{tanh}(\boldsymbol{W}_{1}\boldsymbol{z}+\boldsymbol{b}_{1}),$
(3)

where $\boldsymbol{W}_{1}\in\mathbb{R}^{M\times J}$,
$\boldsymbol{b}_{1}\in\mathbb{R}^{M}$ and $\boldsymbol{z}\in\mathbb{R}^{J}$.

$\displaystyle\boldsymbol{z}=(w_{0},\boldsymbol{z}_{1},\boldsymbol{z}_{2},...\boldsymbol{z}_{i},...,\boldsymbol{z}_{n}),$
(4)

where $w_{0}\in\mathbb{R}$ is a global scalar parameter and $n$ is the number
of fields in total. $\boldsymbol{z}_{i}\in\mathbb{R}^{K+1}$ is a parameter
vectors for the $i$-th field in factorisation machines:

$\boldsymbol{z}_{i}=\boldsymbol{W}_{0}^{i}\cdot\boldsymbol{x}[\text{start}_{i}:\text{end}_{i}]=(w_{i},v_{i}^{1},v_{i}^{2},\ldots,v_{i}^{K}),$
(5)

where $\text{start}_{i}$ and $\text{end}_{i}$ are starting and ending feature
indexes of the $i$-th field,
$\boldsymbol{W}_{0}^{i}\in\mathbb{R}^{(K+1)\times(\text{end}_{i}-\text{start}_{i}+1)}$
and $\boldsymbol{x}$ is the input vector as described at beginning. All
weights $\boldsymbol{W}_{0}^{i}$ are initialised with the bias term $w_{i}$
and vector $\boldsymbol{v}_{i}$ respectively (e.g.,
$\boldsymbol{W}_{0}^{i}[0]$ is initialised by $w_{i}$,
$\boldsymbol{W}_{0}^{i}[1]$ is initialised by $v_{i}^{1}$,
$\boldsymbol{W}_{0}^{i}[2]$ is initialised by $v_{i}^{2}$, etc.). In this way,
$\boldsymbol{z}$ vector of the first layer is initialised as shown in Figure 1
via training a factorisation machine (FM) [31]:

$y_{\text{FM}}(\boldsymbol{x}):=\text{sigmoid}\Big{(}w_{0}+\sum_{i=1}^{N}w_{i}x_{i}+\sum_{i=1}^{N}\sum_{j=i+1}^{N}\langle\boldsymbol{v}_{i},\boldsymbol{v}_{j}\rangle
x_{i}x_{j}\Big{)},$ (6)

where each feature $i$ is assigned with a bias weight $w_{i}$ and a
$K$-dimensional vector $\boldsymbol{v}_{i}$ and the feature interaction is
modelled as their vectors’ inner product
$\langle\boldsymbol{v}_{i},\boldsymbol{v}_{j}\rangle$. In this way, the above
neural nets can learn more efficiently from factorisation machine
representation so that the computational complexity problem of the high-
dimensional binary inputs has been naturally bypassed. Different hidden layers
can be regarded as different internal functions capturing different forms of
representations of the data instance. For this reason, this model has more
abilities of catching intrinsic data patterns and leads to better performance.

The idea using FM in the bottom layer is ignited by Convolutional Neural
Networks (CNNs) [11], which exploit spatially local correlation by enforcing a
local connectivity pattern between neurons of adjacent layers. Similarly, the
inputs of hidden layer 1 are connected to the input units of a specific field.
Also, the bottom layer is not fully connected as FM performs a field-wise
training for one-hot sparse encoded input, allowing local sparsity,
illustrated as the dash lines in Figure 1. FM learns good structural data
representation in the latent space, helpful for any further model to build on.
A subtle difference, though, appears between the product rule of FM and the
sum rule of DNN for combination. However, according to [21], if the
observational discriminatory information is highly ambiguous (which is true in
our case for ad click behaviour), the posterior weights (from DNN) will not
deviate dramatically from the prior (FM).

Furthermore, the weights in hidden layers (except the FM layer) are
initialised by layer-wise RBM pre-training [3] using contrastive divergence
[17], which effectively preserves the information in input dataset as detailed
in [18, 16]. The initial weights for FMs are trained by stochastic gradient
descent (SGD), as detailed in [31]. Note that we only need to update weights
which connect to the positive input units, which largely reduces the
computational complexity. After pre-training of the FM and upper layers,
supervised fine-tuning (back propagation) is applied to minimise loss function
of cross entropy:

$\displaystyle L(y,\hat{y})=-y\log\hat{y}-(1-y)\log(1-\hat{y}),$ (7)

where $\hat{y}$ is the predicted CTR in Eq. (1) and $y$ is the binary click
ground-truth label. Using the chain rule of back propagation, the FNN weights
including FM weights can be efficiently updated. For example, we update FM
layer weights via

$\displaystyle\frac{\partial L(y,\hat{y})}{\partial\boldsymbol{W}_{0}^{i}}$
$\displaystyle=\frac{\partial
L(y,\hat{y})}{\partial\boldsymbol{z}_{i}}\frac{\partial\boldsymbol{z}_{i}}{\partial\boldsymbol{W}_{0}^{i}}=\frac{\partial
L(y,\hat{y})}{\partial\boldsymbol{z}_{i}}\boldsymbol{x}[\text{start}_{i}:\text{end}_{i}]$
(8) $\displaystyle\boldsymbol{W}_{0}^{i}$
$\displaystyle\leftarrow\boldsymbol{W}_{0}^{i}-\eta\cdot\frac{\partial
L(y,\hat{y})}{\partial\boldsymbol{z}_{i}}\boldsymbol{x}[\text{start}_{i}:\text{end}_{i}].$
(9)

Due to the fact that the majority entries of
$\boldsymbol{x}[\text{start}_{i}:\text{end}_{i}]$ are 0, we can accelerate
fine-tuning by updating weights linking to positive units only.

Figure 2: A 4-layer SNN architecture and two first-layer pre-training methods.

###  3.2 Sampling-based Neural Networks (SNN)

The structure of the second model SNN is shown in Figure 2(a). The difference
between SNN and FNN lies in the structure and training method in the bottom
layer. SNN’s bottom layer is fully connected with sigmoid activation function:

$\boldsymbol{z}=\text{sigmoid}(\boldsymbol{W}_{0}\boldsymbol{x}+\boldsymbol{b}_{0}).$
(10)

To initialise the weights of the bottom layer, we tried both restricted
Boltzmann machine (RBM) [16] and denoising auto-encoder (DAE) [4] in the pre-
training stage. In order to deal with the computational problem of training
large sparse one-hot encoding data, we propose a sampling-based RBM (Figure
2(b), denoted as SNN-RBM) and a sampling-based DAE in (Figure 2(c), denoted as
SNN-DAE) to efficiently calculate the initial weights of the bottom layer.

Instead of modelling the whole feature set for each training instance set, for
each feature field, e.g., city, there is only one positive value feature for
each training instance, e.g., city=London, we sample $m$ negative units, e.g.,
city=Paris when $m=1$, randomly with value 0. Black units in Figure 2(b) and
2(c) are unsampled and thus ignored when pre-training the data instance. With
the sampled units, we can train an RBM via contrastive divergence [17] and a
DAE via SGD with unsupervised approaches to largely reduce the data dimension
with high recovery performance. The real-value dense vector is used as the
input of the further layers in SNN.

In this way, computational complexity can be dramatically reduced and, in
turn, initial weights can be calculated quickly and back-propagation is then
performed to fine-tune SNN model.

###  3.3 Regularisation

To prevent overfitting, the widely used L2 regularisation term is added to the
loss function. For example, the L2 regularisation for FNN in Figure 1 is

$\displaystyle\Omega(\boldsymbol{w})=||\boldsymbol{W}_{0}||_{2}^{2}+\sum_{l=1}^{3}\Big{(}||\boldsymbol{W}_{l}||_{2}^{2}+||\boldsymbol{b}_{l}||_{2}^{2}\Big{)}.$
(11)

On the other hand, _dropout_ [35] is a technique which becomes a popular and
effective regularisation technique for deep learning during the recent years.
We also implement this regularisation and compare them in our experiment.

##  4 Experiment

###  4.1 Experiment Setup

Data. We evaluate our models based on iPinYou dataset [27], a public real-
world display ad dataset with each ad display information and corresponding
user click feedback. The data logs are organised by different advertisers and
in a row-per-record format. There are 19.50M data instances with 14.79K
positive label (click) in total. The features for each data instance are all
categorical. Feature examples in the ad log data are user agent, partially
masked IP, region, city, ad exchange, domain, URL, ad slot ID, ad slot
visibility, ad slot size, ad slot format, creative ID, user tags, etc. After
one-hot encoding, the number of binary features is 937.67K in the whole
dataset. We feed each compared model with these binary-feature data instances
and the user click (1) and non-click (0) feedback as the ground-truth labels.
In our experiments, we use training data from advertiser 1458, 2259, 2261,
2997, 3386 and the whole dataset, respectively.

Models. We compare the performance of the following CTR estimation models:

  * LR:

Logistic Regression [32] is a linear model with simple implementation and fast
training speed, which is widely used in online advertising estimation.

  * FM:

Factorisation Machine [31] is a non-linear model able to estimate feature
interactions even in problems with huge sparsity.

  * FNN:

Factorisation-machine supported Neural Network is our proposed model as
described in Section 3.1.

  * SNN:

Sampling-based Neural Network is also our proposed model with sampling-based
RBM and DAE pre-training methods for the first layer in Section 3.2, denoted
as SNN-RBM and SNN-DAE respectively.

Our experiment code222The source code with demo data:
https://github.com/wnzhang/deep-ctr of both FNN and SNN is implemented with
Theano333Theano: http://deeplearning.net/software/theano/.

Metric. To measure the CTR estimation performance of each model, we employ the
area under ROC curve (AUC)444Besides AUC, root mean square error (RMSE) is
also tested. However, positive/negative examples are largly unbalanced in ad
click scenario, and the empirically best regression model usually provides the
predicted CTR close to 0, which results in very small RMSE values and thus the
improvement is not well captured.. The AUC [12] metric is a widely used
measure for evaluating the CTR performance.

###  4.2 Performance Comparison

Table 1: Overall CTR estimation AUC performance. |  LR |  FM |  FNN | SNN-DAE | SNN-RBM
---|---|---|---|---|---
1458 | 70.42% | 70.21% | 70.52% | 70.46% | 70.49%
2259 | 69.66% | 69.73% | 69.74% | 68.08% | 68.34%
2261 | 62.03% | 60.97% | 62.99% | 63.72% | 63.72%
2997 | 60.77% | 60.87% | 61.41% | 61.58% | 61.45%
3386 | 80.30% | 79.05% | 80.56% | 79.62% | 80.07%
all | 68.81% | 68.18% | 70.70% | 69.15% | 69.15%

Table 1 shows the results that compare LR, FM, FNN and SNN with RBM and DAE on
5 different advertisers and the whole dataset. We observe that FM is not
significantly better than LR, which means 2-order combination features might
not be good enough to catch the underlying data patterns. The AUC performance
of the proposed FNN and SNN is better than the performance of LR and FM on all
tested datasets. Based on the latent structure learned by FM, FNN further
learns effective patterns between these latent features and provides a
consistent improvement over FM. The performance of SNN-DAE and SNN-RBM is
generally consistent, i.e., the relative order of the results of the SNN are
almost the same.

###  4.3 Hyperparameter Tuning

Due to the fact that deep neural networks involve many implementation details
and need to tune a fairly large number of hyper-parameters, following details
show how we implement our models and tune hyperparameters in the models.

We use stochastic gradient descent to learn most of our parameters for all
proposed models. Regarding selecting the number of training epochs, we use
early stopping [30], i.e., the training stops when the validation error
increases. We try different learning rate from 1, 0.1, 0.01, 0.001 to 0.0001
and choose the one with optimal performance on the validation dataset.

For negative unit sampling of SNN-RBM and SNN-DAE, we try the negative sample
number $m=1,2$ and $4$ per field as described in Section 3.2, and find $m=2$
produces the best results in most situations. For the activation functions in
both models on the hidden layers (as Eqs. (3) and (2)), we try linear
function, sigmoid function and tanh function, and find the result of tanh
function is optimal. This might be because the hyperbolic tangent often
converges faster than the sigmoid function.

###  4.4 Architecture Selection

In our models, we investigate architectures with 3, 4 and 5 hidden layers by
fixing all layer sizes and find the architecture with 3 hidden layers (i.e., 5
layers in total) is the best in terms of AUC performance. However, the range
of choosing their layer sizes is exponential in the number of hidden layers.
Suppose there is a deep neural network with $L$ hidden layers and each of the
hidden layers is trained with a range of hidden units from 100 to 500 with
increments of 100, thus there are $5^{L}$ models in total to compare.

(a) Architecture

(b) SNN-RBM on 2997

(c) SNN-RBM on 3386

Figure 3: AUC Performance with different architectures.

Instead of trying all combinations of hidden units, in our experiment we use
another strategy by starting tuning the different hidden layer sizes with the
same number of hidden units in all three hidden layers555Some advanced
Bayesian methods for hyperparameter tuning [34] are not considered in this
paper and may be investigated in the future work. since the architecture with
equal-size hidden layers is empirically better than the architecture with
increasing width or decreasing width in [24]. For this reason, we start tuning
layer sizes with equal hidden layer sizes. In fact, apart from increasing,
constant, decreasing layer sizes, there is a more effective structure, which
is the diamond shape of neural networks, as shown in Figure 3(a). We compare
our diamond shape network with other three shapes of networks and tune the
total number of total hidden units on two different datasets shown in Figures
3(b) and 3(c). The diamond shape architecture outperforms others in almost all
layer size settings. The reason why this diamond shape works might be because
this special shape of neural network has certain constraint to the capacity of
the neural network, which provides better generalisation on test sets. On the
other hand, the performance of diamond architecture picks at the total hidden
unit size of 600, i.e., the combination of (200, 300, 100). This depends on
the training data observation numbers. Too many hidden units against a limited
dataset could cause overfitting.

###  4.5 Regularisation Comparison

Neural network training algorithms are very sensitive to the overfitting
problem since deep networks have multiple non-linear layers, which makes them
very expressive models that can learn very complicated functions. For DNN
models, we compared L2 regularisation (Eq. (11)) and dropout [35] for
preventing complex co-adaptations on the training data. The dropout rate
implemented in this experiment refers to the probability of each unit being
active.

Figure 4(a) shows the compared AUC performance of SNN-RBM regularised by L2
norm and dropout. It is obvious that dropout outperforms L2 in all compared
settings. The reason why dropout is more effective is that when feeding each
training case, each hidden unit is stochastically excluded from the network
with a probability of dropout rate, i.e., each training case can be regarded
as a new model and these models are averaged as a special case of bagging [5],
which effectively improves the generalisation ability of DNN models.

(a) Dropout vs. L2

(b) FNN on 2997 dataset

(c) SNN on 2997 dataset

Figure 4: AUC performance w.r.t difference regularisation settings.

###  4.6 Analysis of Parameters

As a summary of Sections 4.4 and 4.5, for both FNN and SNN, there are two
important parameters which should be tuned to make the model more effective:
(i) the parameters of layer size decide the architecture of the neural network
and (ii) the parameter of dropout rate changes generalisation ability on all
datasets compared to neural networks just with L2 regularisation.

Figures 4(b) and 4(c) show how the AUC performance changes with the increasing
of dropout in both FNN and SNN. We can find that there is an upward trend of
performance in both models at the beginning and then drop sharply with
continuous decreasing of dropout rate. The distinction between two models is
the different sensitivities of the dropout. From Figure 4(c), we can see the
model SNN is sensitive to the dropout rate. This might be caused by the
connectivities in the bottom layer. The bottom layer of the SNN is fully
connected with the input vector while the bottom layer for FNN is partially
connected and thus the FNN is more robust when some hidden units are dropped
out. Furthermore, the sigmoid activation function tend to more effective than
the linear activation function in terms of dropout. Therefore, the dropout
rates at the best performance of FNN and SNN are quite different. For FNN the
optimal dropout rate is around 0.8 while for SNN is about 0.99.

##  5 Conclusion

In this paper, we investigated the potential of training deep neural networks
(DNNs) to predict users’ ad click response based on multi-field categorical
features. To deal with the computational complexity problem of high-
dimensional discrete categorical features, we proposed two DNN models: field-
wise feature embedding with supervised factorisation machine pre-training, and
fully connected DNN with field-wise sampling-based RBM and DAE unsupervised
pre-training. These architectures and pre-training algorithms make our DNNs
trained very efficiently. Comprehensive experiments on a public real-world
dataset verifies that the proposed DNN models successfully learn the
underlying data patterns and provide superior CTR estimation performance than
other compared models. The proposed models are very general and could enable a
wide range of future works. For example, the model performance can be improved
by momentum methods in that it suffices for handling the curvature problems in
DNN training objectives without using complex second-order methods [36]. In
addition, the partial connection in the bottom layer could be extended to
higher hidden layers as partial connectivities have many advantages such as
lower complexity, higher generalisation ability and more similar to human
brain [9].

## References

  * [1] Beck, J.E., Woolf, B.P.: High-level student modeling with machine learning. In: Intelligent tutoring systems. pp. 584–593. Springer (2000)
  * [2] Bengio, Y.: Learning deep architectures for ai. Foundations and trends® in Machine Learning 2(1), 1–127 (2009)
  * [3] Bengio, Y., Lamblin, P., Popovici, D., Larochelle, H., et al.: Greedy layer-wise training of deep networks. NIPS 19, 153 (2007)
  * [4] Bengio, Y., Yao, L., Alain, G., Vincent, P.: Generalized denoising auto-encoders as generative models. In: NIPS. pp. 899–907 (2013)
  * [5] Breiman, L.: Bagging predictors. Machine learning 24(2), 123–140 (1996)
  * [6] Broder, A.Z.: Computational advertising. In: SODA. vol. 8, pp. 992–992 (2008)
  * [7] Collobert, R., Weston, J., Bottou, L., Karlen, M., Kavukcuoglu, K., Kuksa, P.: Natural language processing (almost) from scratch. JMLR 12, 2493–2537 (2011)
  * [8] Deng, L., Abdel-Hamid, O., Yu, D.: A deep convolutional neural network using heterogeneous pooling for trading acoustic invariance with phonetic confusion. In: ICASSP. pp. 6669–6673. IEEE (2013)
  * [9] Elizondo, D., Fiesler, E.: A survey of partially connected neural networks. International journal of neural systems 8(05n06), 535–558 (1997)
  * [10] Erhan, D., Bengio, Y., Courville, A., Manzagol, P.A., Vincent, P., Bengio, S.: Why does unsupervised pre-training help deep learning? JMLR 11 (2010)
  * [11] Fukushima, K.: Neocognitron: A self-organizing neural network model for a mechanism of pattern recognition unaffected by shift in position. Biological cybernetics 36(4), 193–202 (1980)
  * [12] Graepel, T., Candela, J.Q., Borchert, T., Herbrich, R.: Web-scale bayesian click-through rate prediction for sponsored search advertising in microsoft’s bing search engine. In: ICML. pp. 13–20 (2010)
  * [13] Graves, A., Mohamed, A.r., Hinton, G.: Speech recognition with deep recurrent neural networks. In: ICASSP. pp. 6645–6649. IEEE (2013)
  * [14] Hand, D.J., Yu, K.: Idiot’s bayes not so stupid after all? International statistical review 69(3), 385–398 (2001)
  * [15] He, X., Pan, J., Jin, O., Xu, T., Liu, B., Xu, T., Shi, Y., Atallah, A., Herbrich, R., Bowers, S., et al.: Practical lessons from predicting clicks on ads at facebook. In: ADKDD. pp. 1–9. ACM (2014)
  * [16] Hinton, G.: A practical guide to training restricted boltzmann machines. Momentum 9(1), 926 (2010)
  * [17] Hinton, G.E.: Training products of experts by minimizing contrastive divergence. Neural computation 14(8), 1771–1800 (2002)
  * [18] Hinton, G.E., Salakhutdinov, R.R.: Reducing the dimensionality of data with neural networks. Science 313(5786), 504–507 (2006)
  * [19] Huang, P.S., He, X., Gao, J., Deng, L., Acero, A., Heck, L.: Learning deep structured semantic models for web search using clickthrough data. In: CIKM. pp. 2333–2338 (2013)
  * [20] Juan, Y.C., Zhuang, Y., Chin, W.S.: 3 idiots approach for display advertising challenge. In: Internet and Network Economics, pp. 254–265. Springer (2011)
  * [21] Kittler, J., Hatef, M., Duin, R.P., Matas, J.: On combining classifiers. PAMI 20(3), 226–239 (1998)
  * [22] Krizhevsky, A., Sutskever, I., Hinton, G.E.: Imagenet classification with deep convolutional neural networks. In: NIPS (2012)
  * [23] Kurashima, T., Iwata, T., Takaya, N., Sawada, H.: Probabilistic latent network visualization: inferring and embedding diffusion networks. In: KDD. pp. 1236–1245. ACM (2014)
  * [24] Larochelle, H., Bengio, Y., Louradour, J., Lamblin, P.: Exploring strategies for training deep neural networks. JMLR 10, 1–40 (2009)
  * [25] LeCun, Y., Bengio, Y., Hinton, G.: Deep learning. Nature 521(7553) (2015)
  * [26] Lee, K.c., Orten, B., Dasdan, A., Li, W.: Estimating conversion rate in display advertising from past performance data. In: KDD. pp. 768–776. ACM (2012)
  * [27] Liao, H., Peng, L., Liu, Z., Shen, X.: ipinyou global rtb bidding algorithm competition dataset. In: ADKDD. pp. 1–6. ACM (2014)
  * [28] McMahan, H.B., Holt, G., Sculley, D., Young, M., Ebner, D., Grady, J., Nie, L., Phillips, T., Davydov, E., Golovin, D., et al.: Ad click prediction: a view from the trenches. In: KDD. pp. 1222–1230. ACM (2013)
  * [29] Oentaryo, R.J., Lim, E.P., Low, D.J.W., Lo, D., Finegold, M.: Predicting response in mobile advertising with hierarchical importance-aware factorization machine. In: WSDM (2014)
  * [30] Prechelt, L.: Automatic early stopping using cross validation: quantifying the criteria. Neural Networks 11(4), 761–767 (1998)
  * [31] Rendle, S.: Factorization machines with libfm. ACM TIST 3(3), 57 (2012)
  * [32] Richardson, M., Dominowska, E., Ragno, R.: Predicting clicks: estimating the click-through rate for new ads. In: WWW. pp. 521–530. ACM (2007)
  * [33] Shen, Y., He, X., Gao, J., Deng, L., Mesnil, G.: A latent semantic model with convolutional-pooling structure for information retrieval. In: CIKM (2014)
  * [34] Snoek, J., Larochelle, H., Adams, R.P.: Practical bayesian optimization of machine learning algorithms. In: NIPS. pp. 2951–2959 (2012)
  * [35] Srivastava, N., Hinton, G., Krizhevsky, A., Sutskever, I., Salakhutdinov, R.: Dropout: A simple way to prevent neural networks from overfitting. JMLR 15(1), 1929–1958 (2014)
  * [36] Sutskever, I., Martens, J., Dahl, G., Hinton, G.: On the importance of initialization and momentum in deep learning. In: ICML. pp. 1139–1147 (2013)
  * [37] Tang, J., Qu, M., Wang, M., Zhang, M., Yan, J., Mei, Q.: Line: Large-scale information network embedding. In: WWW. pp. 1067–1077 (2015)
  * [38] Trofimov, I., Kornetova, A., Topinskiy, V.: Using boosted trees for click-through rate prediction for sponsored search. In: WINE. p. 2. ACM (2012)
  * [39] Wang, X., Li, W., Cui, Y., Zhang, R., Mao, J.: Click-through rate estimation for rare events in online advertising. Online Multimedia Advertising: Techniques and Technologies pp. 1–12 (2010)
  * [40] Zeiler, M.D., Taylor, G.W., Fergus, R.: Adaptive deconvolutional networks for mid and high level feature learning. In: ICCV. pp. 2018–2025. IEEE (2011)
  * [41] Zhang, W., Yuan, S., Wang, J.: Optimal real-time bidding for display advertising. In: KDD. pp. 1077–1086. ACM (2014)
  * [42] Zou, Y., Jin, X., Li, Y., Guo, Z., Wang, E., Xiao, B.: Mariana: Tencent deep learning platform and its applications. VLDB 7(13), 1772–1777 (2014)