# Summation method in optimal control problem with delay

Pavel Barkhayev  Pavel Barkhayev,
Department of Mathematical Sciences Norwegian University of Science and
Technology, 7491, Trondheim, Norway and
Institute for Low Temperature Physics and Engineering, 47 Nauki ave, 61103,
Kharkiv, Ukraine,
<EMAIL_ADDRESS>and  Yurii Lyubarskii  Yurii Lyubarskii,
St. Petersburg State University, St. Petersburg, Russia, and
Department of Mathematical Sciences, Norwegian University of Science and
Technology, NO-7491 Trondheim, Norway,
<EMAIL_ADDRESS>
The work of P.Barkhayev was supported by Norwegian Research Council project
No. 275113, "COMAN"
The work of Yu.Lyubarskii was supported by a grant of the Goverment of the
Russian Federation for the state support of scientific research, carried under
the supervision of leading scientists, agreement 075-15-2021-602

##  1\. Introduction

We consider the control problem for differential equation with distributed
delay

(1.1) $\dot{x}(t)=x(t-1)+\int_{-1}^{0}x(t+\tau)\varphi(\tau)d\tau+u(t),\ t>0$

with the initial conditions

(1.2) $x(t)=x_{0}(t),\ t\in(-1,0);\ x(0)=\xi,\quad\xi\in\mathbb{C},\ x_{0}\in
L^{2}(-1,0).$

The distribution of delay is defined by the function $\varphi\in L^{2}(-1,0)$
and the control function $u$ is to be chosen to steer the initial condition
$(\xi,x_{0}(t))$ to the desired state.

Let $\mathcal{M}=\mathbb{C}\times L^{2}(-1,0)$ be the space of all initial
data $\mathcal{M}=\\{\vec{x}=\bigl{(}\begin{smallmatrix}\xi\\\
x(t)\end{smallmatrix}\bigr{)},\xi\in\mathbb{C},\ x_{0}\in L^{2}(-1,0)\\}$.

The equation (1.1) always has a solution $x(t)=x(t;\vec{x},u)$ and we deal
with null controlability problem: given $T>1$ and $\vec{x}\in\mathcal{M}$,
find the control function $u$ which steers the initial data to zero, i.e.
$x(t;\vec{x},u)=0$, $t\in(T-1,T)$.

It is known (see e.g. [1]) that such control function always exists, and can
be given by a bounded operator. Actually there exist many such functions. We
consider the problem of finding the optimal control. In other we look for the
problem

(1.3) $\|u\|_{L^{2}(0,T)}\to\min,\quad u\in U(x_{0},\xi;T),$

here $U(x_{0},\xi,T)$ stays for all control functions which steers the initial
data to zero at the moment $T$.

Such kind of problems were investigated in various settings by many authors.
We refer the reader to ([8], [7],[10], [5] ) to mention a few. We also refer
the reader to [7], [5] for comprehensive survey and description of the current
stay of art.

We combine the semigroup techniques from [6] with the projection techniques,
see e.g. [2] and [9]. These techniques lead one to the study of the
eigenfunction expansions for the corresponding non-selfadjoint operator. Since
the corresponding systems of eigenfunctions do not form a Riesz basis, the
expansions diverge, generally speaking. In [2] the authors prove the point-
wise convergence on the inner intervals of the segment $[T-1,T]$, yet the
problem of construction the optimal control on the whole interval remains
open. In this article we consider a simplest model of time delay equation and
suggest the reconstruction procedure which converges in the whole
$L^{2}(0,T)$. We believe the same procedure can be applied in much more
general settings.

We adjust for this case the techniques in [4], this techniques cannot be
applied directly because the generating function does not satisfy the
Muckenhoupt condition. Also it is worth to mention that when constructing the
optimal control for a single eigenvector we arrived to already familiar
problem which appeared in study of the hereditary completeness of sets of
exponential functions: given a set of exponents $\Lambda\in\mathbb{C}$ such
that, the corresponding exponential functions $\\{\exp(i\lambda t)\\}$ are
incomplete in $L^{2}(0,T)$ we need to construct a function which is in the
span of $\\{\exp(i\lambda t)\\}$ and is orthogonal to all but one of these
exponential functions. This problem has been considered in [3] in a very
general setting, however the explicit construction of such functions is a
complicated problem. In the present article we restrict ourselves to the
simplest model case when $\varphi=0$.

##  2\. The homogeneous problem

The corresponding homogeneous problem has the form

$\dot{x}(t)=x(t-1)+\int_{-1}^{0}x(t+\tau)\varphi(\tau)d\tau,\ t>0$

with the initial conditions

$\bigl{(}\begin{smallmatrix}x(0)\\\
x(\cdot)\end{smallmatrix}\bigr{)}=\vec{x}_{0}\in\mathcal{M}$

This problem admits the unique solution $x(t)=x(t,\vec{x}_{0})$ and also
generates the natural semigroup

$S_{t}=e^{At}:\mathcal{M}\to\mathcal{M},\
S_{t}\vec{x}_{0}=\bigl{(}\begin{smallmatrix}x(t)\\\
x(t+\cdot)\end{smallmatrix}\bigr{)}.$

The infinitesimal operator of the group is defined by the relation

(2.1)
$A:\vec{x}_{0}\mapsto\begin{pmatrix}x_{0}(-1)+\int_{-1}^{0}x_{0}(\tau)\varphi(\tau)d\tau\\\
\dot{x}_{0}(\tau)\end{pmatrix}.$

Its domain

$\mathcal{D}_{A}=\bigl{\\{}(\begin{smallmatrix}\xi\\\
x(t)\end{smallmatrix}\bigr{)},\ x\in W^{2}_{1}(-1,0),\ \xi=x(0)\bigr{\\}}.$

Let

$D(z)=-iz+e^{-iz}+\int_{-1}^{0}e^{i\tau z}\varphi(\tau)d\tau$

be the characteristic function of the homogeneous problem. Its indicator
diagram is $[0,i]$. We denote by $\Lambda=\\{\lambda_{n}\\}$ the set of all
zeros of $D$. It can be easily seen that

$\lambda_{n}=-\pi/2+2\pi n+i\log|n|+o(1),n\in\mathbb{Z}.$

Here and in what follows we refer the reader to [11] for facts about entire
functions. We denote by For simplicity we assume
$\Lambda\subset\mathbb{C}_{+}$ and that $D$ has no multiple zeroes, this may
be not true for a finite number of zeros only.

We also have

(2.2) $|D(z)|\asymp\begin{cases}|z|,\ \text{for}\ z=x+iy\in\mathbb{C}_{+},\
y<\log x,\ \text{dist}(z,\Lambda)>\epsilon\\\ e^{y},\ \text{for}\
z=x+iy\in\mathbb{C}_{+},\ y>\log x,\
\text{dist}(z,\Lambda)>\epsilon.\end{cases}$

The propositions below are straightforward (again we refer to [11] for the
proof of completeness and minimality).

######  Proposition 2.1.

The spectra of $A$ is simple and coincides with $i\Lambda$. The corresponding
eigenvectors are of the form
$\vec{e}_{\lambda}=\bigl{(}\begin{smallmatrix}1\\\
\exp(i\lambda\tau)\end{smallmatrix}\bigr{)},\lambda\in\Lambda$.

Let
$\tilde{\mathcal{E}}_{\Lambda}=\\{\vec{e}_{\lambda}\\}_{\lambda\in\Lambda}$ be
the set of all eigenfunctions.

######  Proposition 2.2.

The system $\tilde{\mathcal{E}}_{\Lambda}$ is complete and minimal in
$\mathcal{M}$.

By
$\mathcal{X}_{\Lambda}=\\{\vec{x}_{\lambda}\\}_{\lambda\in\Lambda}\subset\mathcal{M}$,
$\vec{x}_{\lambda}=\bigl{(}\begin{smallmatrix}\xi_{\lambda}\\\
x_{\lambda}(\tau)\end{smallmatrix}\bigr{)}$. we denote the system biorthogonal
to $\tilde{\mathcal{E}}_{\Lambda}$.

We will use the relations

$e^{A\tau}\vec{e}_{\lambda}=e^{i\lambda\tau}\vec{e}_{\lambda},\
e^{A^{*}\tau}\vec{x}_{\lambda}=e^{-i\bar{\lambda}\tau}\vec{x}_{\lambda},\
\lambda\in\Lambda.$

######  Proposition 2.3.

The following relations take place

(2.3) $\overline{\xi_{\lambda}}=\frac{-i}{D^{\prime}(\lambda)},$ (2.4)
$\int_{-1}^{0}e^{izt}\overline{x_{\lambda}(t)}dt=\frac{1}{D^{\prime}(\lambda)}\frac{e^{-iz}+\int_{-1}^{0}e^{iz\tau}\varphi(\tau)d\tau+i\lambda}{z-\lambda}.$

Since $D(\lambda)=0$, one can rewrite (2.4 ) as

(2.5)
$\int_{-1}^{0}e^{izt}\overline{x_{\lambda}(t)}dt=\frac{1}{D^{\prime}(\lambda)}\bigl{[}\frac{e^{-izt}-e^{-i\lambda
t}}{z-\lambda}+\int_{-1}^{0}\frac{e^{-izt}-e^{-i\lambda
t}}{z-\lambda}\varphi(t)dt\bigr{]}.$

In what follows we use both (2.4) and (2.5).

Even after being normalised, the system $\tilde{\mathcal{E}}_{\Lambda}$ does
not form a Riesz basis in $\mathcal{M}$. However for each
$\vec{x}\in\mathcal{M}$ we can consider it formal Fourier series

(2.6)
$\vec{x}\sim\sum_{\lambda\in\Lambda}\langle\vec{x},\vec{x}_{\lambda}\rangle\vec{e}_{\lambda}.$

##  3\. Structure of the optimal solution

Denote the operator $\mathcal{B}:\mathbb{C}\to\mathcal{M}$ as
$\mathcal{B}a=\bigl{(}\begin{smallmatrix}a\\\ 0\end{smallmatrix}\bigr{)}$,
$a\in\mathbb{C}$. The problem (1.1), (1.2) can now be written as

(3.1) $\begin{cases}\dot{\vec{x}}(t)=A\vec{x}(t)+\mathcal{B}u(t),\ t>0,\\\
\vec{x}(0)=\vec{x}_{0}.\end{cases}$

Its solution has the form

$\vec{x}(t)=e^{At}\vec{x}_{0}+\int_{0}^{t}e^{A(t-\tau)}\mathcal{B}u(\tau)d\tau.$

Fix now $T>1$. The admissibility condition $u\in U(T;\vec{x}_{0})$ reads as

$e^{AT}\vec{x}_{0}=-\int_{0}^{T}e^{A(T-\tau)}\mathcal{B}u(\tau)d\tau.$

The set of exponential functions ${\mathcal{E}}_{\Lambda}=\\{e^{i\lambda
t}\\}$ is incomplete in $L^{2}(0,T)$. Let ${\mathcal{E}}_{\Lambda}(T)$ be the
closure of its span in $L^{2}(0,T)$ and
${\mathcal{E}}_{\Lambda}(T)^{\perp}=L^{2}(0,T)\ominus{\mathcal{E}}_{\Lambda}(T)$.

We follow the result from [2], see also [9], Ch.7.4.

######  Theorem 3.1.

Let $u_{\lambda}\in L^{2}(0,T)$ solve the problem (1.3) for
$\vec{x}=\vec{e}_{\lambda}$ and
$v_{\lambda}(t)=-u_{\lambda}(T-t)$. Then

$v_{\lambda}\in{\mathcal{E}}_{\Lambda}(T)$

and also

$\int_{0}^{T}e^{i\mu t}\overline{v_{\lambda}(t)}dt=\delta_{\lambda,\mu},\
\mu\in\Lambda.$

######  Lemma 3.2.

Let $v\in L^{2}(0,T)$ be such that

(3.2) $0=\int_{0}^{T}e^{At}\mathcal{B}v(t)dt.$

Then $v\in{\mathcal{E}}_{\Lambda}(T)^{\perp}$.

Proof Given $\lambda\in\Lambda$ relation (3.2) yields

$0=\int_{0}^{T}\langle e^{At}\mathcal{B}v(t),\vec{x}_{\lambda}\rangle
dt=\int_{0}^{T}\langle\mathcal{B}v(t),e^{A^{*}t}\vec{x}_{\lambda}\rangle
dt=\int_{0}^{T}\langle\mathcal{B}v(t),e^{-i\bar{\lambda}t}\vec{x}_{\lambda}\rangle
dt.$

It remains to mention that $\mathcal{B}v(t)=\bigl{(}\begin{smallmatrix}v(t)\\\
0\end{smallmatrix}\bigr{)}$,
$\vec{x}_{\lambda}=\bigl{(}\begin{smallmatrix}\xi_{\lambda}\\\
x_{\lambda}(t)\end{smallmatrix}\bigr{)}$, therefore
$\langle\mathcal{B}v(t),e^{-i\bar{\lambda}t}\vec{x}_{\lambda}\rangle=\overline{\xi_{\lambda}}v(t)e^{-i\bar{\lambda}t}$
and $\xi_{\lambda}\neq 0$. ∎

######  Lemma 3.3.

Let $\lambda\in\Lambda$ and $v\in L^{2}(0,T)$ be such that

$e^{AT}\vec{e}_{\lambda}=\int_{0}^{T}e^{At}\mathcal{B}v(t)dt.$

Then

$\int_{0}^{T}v(t)e^{-i\bar{\mu}t}dt=\bar{y}_{\lambda}^{-1}e^{i\lambda
t}\delta_{\lambda,\mu},\ \mu\in\Lambda.$

Proof is similar to that of the previous lemma. We write

$e^{i\lambda t}\delta_{\lambda,\mu}=\langle
e^{AT}\vec{e}_{\lambda},\vec{x}_{\mu}\rangle=\int_{0}^{T}\langle
e^{At}\mathcal{B}v(t),\vec{x}_{\mu}\rangle dt$

and then apply similar arguments. ∎

######  Corollary 3.4.

The optimal control linearly depends on initial data: if $u_{1}$ and $u_{2}$
are solutions of the problem (1.3) for $\vec{x}=\vec{x}_{1},\vec{x}_{2}$
respectively, then $u_{1}+u_{2}$ solves this problem for
$\vec{x}=\vec{x}_{1}+\vec{x}_{2}$.

Now in order to solve the general problem we need to reconstruct an arbitrary
initial data $\vec{x}\in\mathcal{M}$ from its expansion (2.6).

##  4\. Summation method

###  4.1. Formulation of the result

In this section we describe linear summation method for the series (2.6), i.e
we describe the matrices $(w_{n}(\lambda))_{n\in\mathbb{N},\lambda\in\Lambda}$
with the following properties:

(4.1) $\displaystyle w_{n}(\lambda)\to 1\ \text{as}\ n\to\infty,\ \text{for
each}\ \lambda\in\Lambda;$ (4.2) $\displaystyle\\#\\{n;w_{n}(\lambda)\neq
0\\}<\infty,\ \text{for each}\ n\in\mathbb{N};$ (4.3) $\displaystyle
S_{n}\vec{x}\to\vec{x}\ \text{as}\ n\to\infty,\ \text{here}\
S_{n}\vec{x}=\sum_{\lambda\in\Lambda}w_{n}(\lambda)\langle\vec{x},\vec{x}_{\lambda}\rangle\vec{e}_{\lambda}.$

Given the sequences $\\{l_{n}\\}$, $\\{R_{n}\\}$, satisfying the conditions

(4.4) $\displaystyle l_{n},\ R_{n}\to\infty\ \text{as}\ n\to\infty,$ (4.5)
$\displaystyle l_{n}^{2}/n\to 0,\ n/R_{n}\to 0\ \text{as}\ n\to\infty,$ (4.6)
$\displaystyle e^{-\pi l_{n}/2}R_{n}\to 0\ \text{as}\ n\to\infty.\ $

######  Remark 4.1.

It suffices, for example, to take $l_{n}=|n|^{1/4}$, $R_{n}=n^{4}$.

Consider the function

$W_{n}(z)=e^{-l_{n}\pi-il_{n}\log\frac{z-n}{z+n}},\ z\in\mathbb{C}_{+}.$

######  Theorem 4.2.

Let, for $\lambda\in\Lambda$

(4.7) $w_{n}(\lambda)=\begin{cases}W_{n}(\lambda),&|\lambda|<R_{n},\\\
0,&\text{otherwise}\end{cases}$

and the operator $S_{n}$ be defined by

(4.8) $\
S_{n}\vec{x}=\sum_{\lambda\in\Lambda}w_{n}(\lambda)\langle\vec{x},\vec{x}_{\lambda}\rangle\vec{e}_{\lambda}.$

Then $S_{n}\vec{x}\to\vec{x}$, as $n\to\infty.$

###  4.2. Preliminary estimates

######  Lemma 4.3.

(4.9) $W_{n}(z)\rightrightarrows 1,\ {\text{as}}\ n\to\infty\ {\text{on each
compact set}}\ K\subset\mathbb{C}_{+}$

######  Lemma 4.4.

(4.10) $|W_{n}(z)|<e^{-\pi l_{n}/2},\ z\in\mathbb{C}_{+},|z|=R_{n}.$

These lemmas are straightforward.

On the real line we have

$W_{n}(t)=\begin{cases}e^{il_{n}\log\bigl{|}\frac{t-n}{t+n}\bigr{|}},&|t|<n,\\\
e^{il_{n}\log\bigl{|}\frac{t-n}{t+n}\bigr{|}d-l_{n}\pi},&|t|>n\end{cases}$

######  Lemma 4.5.

$\sup_{n}\bigl{|}\int_{1<|t|<n}W_{n}(t)\frac{dt}{t}\bigr{|}<\infty$

######  Lemma 4.6.

$\sup_{n}\bigl{|}\int_{|t|>n}W_{n}(t)\frac{dt}{t}\bigr{|}<\infty$

Proof of Lemma 4.5

We have

$I_{n}:=\int_{1<|t|<n}W_{n}(t)\frac{dt}{t}=\bigl{(}\int_{1<|t|<n/l_{n}}+\int_{n/l_{n}<|t|<n}\bigr{)}=I_{n,1}+I_{n,2}$

The estimate of $I_{n,1}$ is straightforward. We have for $|t|<n/l_{n}$

$l_{n}\log\bigl{|}\frac{t-n}{t+n}|\prec\frac{l_{n}}{n}|t|<<1,$

Respectively

$I_{n,1}=\int_{1<|t|<n/l_{n}}\bigl{[}1+O\bigl{(}\frac{l_{n}}{n}|t|\bigr{)}\bigr{]}\frac{dt}{t}=O(1).$

When estimating $I_{2,n}$ we consider each interval $\pm[n/l_{n},n]$
separately. Consider the integral

$I_{n,2}^{\prime}=\int_{n/l_{n}<t<n}e^{il_{n}\log\bigl{|}\frac{t-n}{t+n}\bigr{|}}\frac{dt}{t}$

By changing of variables

$u=\frac{n-t}{n+t};\ \tau=1-u;\ v=\log(1-\tau),$

and extracting the main term on each step, we finally arrive to the integral
$\int_{1/l_{n}}^{\infty}e^{il_{n}v}\frac{dv}{v},$ which is bounded uniformly
in $n$. ∎

Proof of Lemma 4.6 goes in the same vien. We write

$\int_{|t|>n}W_{n}(t)\frac{dt}{t}=\bigl{(}\int_{|t|>nl_{n}}+\int_{n<|t|<nl_{n}}\bigr{)}W_{n}(t)\frac{dt}{t}$

and use the similar reasonings. ∎

###  4.3. Proof of Theorem 4.2

a. Let operator $S_{n}$ be defined by (4.3). It follows from (4.9) that
$S_{n}\vec{x}\to\vec{x}$ as $n\to\infty$ on the dense set of finite linear
combinations of $\vec{e}_{\lambda}$. So it suffices to prove the uniform
boundedness of all operators $S_{n}$. Consider instead the adjoint operator

$S^{*}\vec{y}=\sum_{\lambda}\overline{w_{n}(\lambda)}\langle\vec{y},\vec{e}_{\lambda}\rangle\vec{x}_{\lambda}.$

Actually it is more convenient to deal with the conjugate objects:

$\overline{S^{*}\vec{y}}=\sum_{\lambda}{w_{n}(\lambda)}\langle\vec{e}_{\lambda},\vec{y}\rangle\overline{\vec{x}_{\lambda}}=:\vec{s}_{n}.$

Let $\vec{s}_{n}=\bigl{(}\begin{smallmatrix}\sigma_{n}\\\
s_{n}(t)\end{smallmatrix}\bigr{)}$ and
$\vec{y}=\bigl{(}\begin{smallmatrix}\eta\\\ y(t)\end{smallmatrix}\bigr{)}$.
Since $\vec{x}_{\lambda}=\bigl{(}\begin{smallmatrix}\xi_{\lambda}\\\
x_{\lambda}(t)\end{smallmatrix}\bigr{)}$ we have

(4.11)
$\sigma_{n}=\sum_{\lambda}w_{n}(\lambda)\bar{\xi}_{\lambda}+\sum_{\lambda}w_{n}(\lambda)\bar{\xi}_{\lambda}\int_{-1}^{0}e^{i\lambda
t}\overline{y(t)}dt=\sigma_{n}^{(1)}+\sigma_{n}{(2)},$ (4.12)
$s_{n}(t)=\sum_{\lambda}w_{n}(\lambda)\overline{x_{\lambda}(t)}+\sum_{\lambda}w_{n}(\lambda)\int_{-1}^{0}e^{i\lambda\tau}\overline{y(\tau)}d\tau\
\overline{x_{\lambda}(t)}=s_{n}^{(1)}(t)+s_{n}^{(2)}(t).$

Each of the four summands in the right of (4.11), (4.12) should be estimates
separately.

Consider the contour $\Gamma_{n}=L_{n}\cup C_{n}$. Here
$L_{n}=[-R_{n},R_{n}]$, $C_{n}=\\{R_{n}e^{i\theta};\theta\in[0,\pi]\\}$. The
arcs $C_{n}$ uniformly in $n$ satisfy the Carleson condition, hence, for any
function $h$ from the Hardy space $\mathcal{H}(\mathbb{C}_{+})$,

$\Bigl{(}\int_{C_{n}}|h(\zeta)|^{2}|d\zeta|\Bigr{)}^{1/2}\leq{\rm
Const}\|h\|.$

We also assume that ${\rm dist}(\Lambda,\Gamma_{n})>0$ uniformly in $n$.

b. Estimate of $\sigma_{n}^{(1)}$. We use (2.3). By the residues

(4.13)
$2i\pi\sigma_{n}^{(1)}=\int_{\Gamma_{n}}\frac{W_{n}(\zeta)}{D(\zeta)}d\zeta=\bigl{(}\int_{C_{n}}+\int_{L_{n}}\bigr{)}\frac{W_{n}(\zeta)}{D(\zeta)}d\zeta=I_{n}+J_{n}.$

That $I_{n}$ is bounded (and even tends to zero as $n\to 0$) is very clear,
this follows from (4.6), (4.10), and (2.2). In order to see the boundedness of
$J_{n}$ it suffices to mention that, according Lemma 4.6 the integral

$\int_{L_{n}}\frac{W_{n}(t)}{it+1}dt$

is bounded uniformly in $n$. Therefore it suffices to establish the
boundedness of

$\int_{L_{n}}W_{n}(t)\bigl{\\{}\frac{1}{D(t)}-\frac{1}{it+1}\bigr{\\}}dt.$

The later is evident since the expression in parenthesis uniformly $O(t^{-2})$
as $t\to\infty$.

c. Estimate of $\sigma_{n}^{(2)}$ Denote

(4.14) $Y(z)=\int_{-1}^{0}e^{iz\tau}\overline{y(\tau)}d\tau.$

As before we have

$2i\pi\sigma_{n}^{(2)}=\int_{\Gamma_{n}}W_{n}(\zeta)\frac{Y(\zeta)}{D(\zeta)}d\zeta=\Bigl{(}\int_{C_{n}}+\int_{L_{n}}\Bigr{)}W_{n}(\zeta)\frac{Y(\zeta)}{D(\zeta)}d\zeta=I_{n}+J_{n}.$

The estimate of $J_{n}$ is now straightforward because $Y\in
L^{2}(\mathbb{R})$ uniformly in $y$ in the unit ball of $L^{2}(-1,0)$ and
$|D(t)|\asymp|t|$ as $t\to\pm\infty$.

In order to estimate $I_{n}$ we observe that the function
$\Psi(\zeta)=e^{iz}Y(\zeta)$ belongs to the Hardy space
$H^{2}(\mathbb{C}_{+})$ and has the same $L^{2}(\mathbb{R})$ norm as $Y$.
Also, it follows from (2.2) that

$\Bigl{|}\frac{Y(\zeta)}{D(\zeta)}\Bigr{|}\prec|\Psi(\zeta)|,\ \zeta\in
C_{n}.$

Now we use that all $C_{n}$’s are uniformly Carleson curves in
$\mathbb{C}_{+}$ and also (4.10). We obtain

$[J_{n}|\prec e^{-l_{n}\pi/2}\int_{C_{n}}|\Psi(\zeta)||d\zeta|\leq
e^{-l_{n}\pi/2}|C_{n}|^{1/2}\bigl{(}\int_{C_{n}}|\Psi(\zeta)|^{2}|d\zeta|\bigr{)}^{1/2}\prec
R_{n}^{1/2}e^{-l_{n}\pi/2}\|\Psi\|.$

It remains to refer to (4.6).

d. Estimate of $\|s_{n}^{(1)}\|$. We apply the Fourier transform and use
representation (2.5):

$({s_{n}^{(1)}})\hat{}(x)=\sum_{\lambda}w_{n}(\lambda)\frac{1}{D^{\prime}(\lambda)}\Bigl{(}\frac{e^{-izt}-e^{-i\lambda
t}}{z-\lambda}+\int_{-1}^{0}\frac{e^{-izt}-e^{-i\lambda
t}}{z-\lambda}\varphi(t)dt\Bigr{)}$

Therefore it suffices to obtain uniform in $t$ estimates of
$\|\Phi_{n,t}\|_{L^{2}(\mathbb{R})}$ where

$\Phi_{n,t}(x)=\sum_{\Lambda}w_{n}(\lambda)\frac{1}{D^{\prime}(\lambda)}\frac{e^{-ixt}-e^{-i\lambda
t}}{x-\lambda}$

Since $\Lambda\subset\mathbb{C}_{+}$ this function belongs to the Hardy space
$\mathcal{H}^{2}(\mathbb{C}_{-})$ and

$\|\Phi_{n,t}\|_{L^{2}(\mathbb{R})}=\sup\Bigl{\\{}\bigl{|}\int_{-\infty}^{\infty}\Phi_{n,t}(x)h(x)dx\Big{|}\
;\ h\in\mathcal{H}^{2}(\mathbb{C}_{+}),\ \|h\|=1\Bigr{\\}}.$

We obtain two terms to be estimated:

$\sum_{\Lambda}\frac{w_{n}(\lambda)}{D^{\prime}(\lambda)}e^{i\lambda
t}\int_{-\infty}^{\infty}\frac{h(x)}{x-\lambda}dx=2i\pi\sum_{\Lambda}\frac{w_{n}(\lambda)}{D^{\prime}(\lambda)}e^{i\lambda
t}h(\lambda)=:a_{1,n},$

and

(4.15)
$\sum_{\Lambda}\frac{w_{n}(\lambda)}{D^{\prime}(\lambda)}\int_{-\infty}^{\infty}\frac{e^{ixt}h(x)}{x-\lambda}dx=:a_{2,n}.$

Once again we have

$a_{1,n}=\Bigl{(}\int_{L_{n}}+\int_{C_{n}}\Bigr{)}W_{n}(\zeta)\frac{e^{i\zeta
t}}{D(\zeta)}h(\zeta)d\zeta.$

The estimate of the first integral follows blantly from. the Schwartz
inequality. In order to estimate the second one we observe that
$\bigl{|}e^{it\zeta}/D(\zeta)\bigr{|}\prec 1$, $\zeta\in C_{n}$ uniformly in
$t$ and $n$ and then again use the fact that all $C_{n}$’s are uniformly
Carlesonian together with (4.6):

$\Bigl{|}\int_{C_{n}}W_{n}(\zeta)\frac{e^{i\zeta
t}}{D(\zeta)}h(\zeta)d\zeta\Bigr{|}\prec\int_{C_{n}}|W_{n}(\zeta)||h(\zeta)||d\zeta|\prec
e^{-\pi l_{n}/2}R_{n}^{1/2}\|h\|\prec\|h\|.$

In order to estimate $a_{n,2}$ we represent $h$ as

$h(x)=\int_{0}^{\infty}e^{i\xi x}\omega(\xi)d\xi$

for some $\omega\in L^{2}(0,\infty)$.

For $t\in[-1,0]$ we have

(4.16) $h(x)=\int_{t}^{0}e^{i\xi x}\omega(\xi-t)d\xi+\int_{0}^{\infty}e^{i\xi
x}\omega(\xi-t)d\xi=h_{t}^{-}(x)+h_{t}^{+}(x).$

Then $h_{t}^{\pm}\in\mathcal{H}(\mathbb{C}_{\pm})$, $\|h_{t}^{+}\|\leq\|h\|$
and

$a_{2,n}=\sum_{\lambda}w_{n}(\lambda)\frac{h_{t}^{+}(\lambda)}{D^{\prime}(\lambda)}.$

The rest of the estimate is similar to that of $a_{n,1}$.

e. Estimate of $\|s_{n}^{(2)}\|$. Again we apply the Fourier transform but now
use representation (2.4):

$(\widehat{s_{n}^{(2)}})(x)=\sum_{\lambda}w_{n}(\lambda)\frac{Y(\lambda)}{D^{\prime}(\lambda)}\
\frac{e^{-ix}+\int_{-1}^{0}e^{ix\tau}\varphi(\tau)d\tau+i\lambda}{x-\lambda},$

here $Y(\lambda)$ is defined in (4.14).

Once again we use the duality reasons

$\|\widehat{s_{n}^{(2)}}\|=\sup\bigl{\\{}\bigl{|}\int_{-\infty}^{\infty}\widehat{s_{n}^{(2)}}(x)h(x)dx\bigr{|}\
,\ h\in\mathcal{H}^{2}(\mathbb{C}_{+}),\ \|h\|\leq 1\bigr{\\}}.$

Let the functions $h_{t}$ be defined as in (4.16). Consider the function

$f(x)=h_{-1}+\int_{-1}^{0}h_{t}(x)\varphi(t)dt.$

We have $f\in\mathcal{H}^{2}(\mathbb{C}^{+})$, $\|f\|\prec 1$ and

$\int_{-\infty}^{\infty}\widehat{s_{n}^{(2)}}(x)h(x)dx=\sum_{\lambda}w_{n}(\lambda)\frac{1+i\lambda}{D^{\prime}(\lambda)}Y(\lambda)f(\lambda).$

This expression can be estimated in a similar way as when estimating
$\|s_{n,1}\|$. ∎

We have proved the following statement.

######  Theorem 4.7.

Given the control problem (1.1), (1.2), let $u_{0}(t)$ be the optimal solution
defined by (1.3). Let, further, $\Lambda$ and $\tilde{\mathcal{E}}_{\Lambda}$
be the corresponding eigenvalues and eigenvectors and also

$\vec{x}_{0}\sim\sum_{\lambda\in\Lambda}\langle\vec{x}_{0},\vec{x}_{\lambda}\rangle\vec{e}_{\lambda}$

be the formal Fourier expansion of the initial data $\vec{x}_{0}$.

Then

$u_{0}(t)=\lim_{n\to\infty}\sum_{\lambda}w_{n}(\lambda)\langle\vec{x}_{0},\vec{x}_{\lambda}\rangle
u_{\lambda}(t),$

here $w_{n}(\lambda)$ are defined in (4.7) and the functions $u_{\lambda}\in
L^{2}(0,T)$ are defined by the relations

$u_{\lambda}(t)=-\overline{v_{\lambda}(T-t)};\
v_{\lambda}\in{\mathcal{E}}(\Lambda),\ \int_{0}^{T}e^{i\mu
t}\overline{v_{\lambda}(t)}dt=\delta_{\lambda,\mu},\ \mu\in\Lambda.$

The limit exists in $L^{2}(0,T)$ sense.

##  5\. A model case

a. In this section we consider the simplest model case

(5.1) $\displaystyle\dot{x}(t)=x(t-1)+u(t),\ t>0;$ (5.2) $\displaystyle
x(t)=x_{0}(t),\ t\in(-1,0),\ x(0)=\xi_{0}.$

For this case we give an explicit construction of the corresponding functions
$v_{\lambda}$ for $T\in(1,2)$. We denote $\delta=T-1$. The characteristic
function has the form

(5.3) $D(z)=-iz+e^{-iz}.$

b. Description of ${\mathcal{E}}(\Lambda)^{\perp}$. Let $g\in L^{2}(0,T)$,
$g\in{\mathcal{E}}(\Lambda)^{\perp}$. Then the entire function

$G(z)=\int_{0}^{T}e^{izt}\overline{g(t)}dt$

vanishes on $\Lambda$ and the standard reasonings on comparing its growth with
that of $D(z)$ yield

(5.4) $G(z)=e^{iz}D(z)\Omega(z),$

where $\Omega$ runs through the set of all functions of the form

(5.5) $\Omega(z)=\int_{0}^{\delta}e^{izt}\overline{\omega(t)}dt,\
x\Omega(x)\in L^{2}(\mathbb{R})$

The later yields

(5.6) $\omega\in W^{2}_{1}(0,\delta),\ \omega(0)=\omega(\delta)=0,$

and, by (5.3) - (5.4), we obtain

$g(t)=\begin{cases}\omega^{\prime}(t-1),&t\in(1,1+\delta),\\\
\omega(t),&t\in(0,\delta),\\\ 0,&\rm{otherwise}\end{cases}$

This is the description of ${\mathcal{E}}(\Lambda)^{\perp}$.

c. Description of the biorthogonal system. Fix $\lambda\in\Lambda$ and let a
function $v_{\lambda}\in L^{2}(0,T)$ satisfy

$\int_{0}^{T}e^{i\mu t}\overline{v_{\lambda}(t)}dt=\delta_{\lambda,\mu},\
\mu\in\Lambda.$

Denote

$V_{\lambda}(z)=\int_{0}^{T}e^{izt}\overline{v_{\lambda}(t)}dt.$

As before we have

$V_{\lambda}(z)=\frac{1}{e^{i\lambda}D^{\prime}(\lambda)}\frac{e^{iz}D(z)}{z-\lambda}P_{\lambda}(z),$

where the function

$P_{\lambda}(z)=\int_{0}^{\delta}e^{izt}\overline{p_{\lambda}(t)}dt,\
p_{\lambda}\in L^{2}(0,\delta)$

Satisfies

(5.7) $P_{\lambda}(\lambda)=\int_{0}^{\delta}e^{i\lambda
t}\overline{p_{\lambda}(t)}dt=1.$

We remind that $D(\lambda)=0$. Therefore

$\frac{e^{ix}D(x)}{x-\lambda}=-ie^{ix}-i\frac{e^{i(x-\lambda)}-1}{x-\lambda}=-i\left[e^{ix}+\int_{0}^{1}e^{ixt}e^{-i\lambda
t}dt\right]$

and, denoting $a_{\lambda}=1/(e^{i\lambda}D^{\prime}(\lambda)$ we obtain

$v_{\lambda}(t)=-ia_{\lambda}\left[\overline{p_{\lambda}(t-1)}+\int_{t-1}^{t}e^{-i\lambda(t-u)}\overline{p_{\lambda}(u)}du\right]$

When $p_{\lambda}$ runs through $L^{2}(0,\delta)$ the function $v_{\lambda}$
runs through all functions in $L^{2}(0,T)$ orthogonal to
${\mathcal{E}}_{\Lambda}\setminus\\{e^{i\lambda t}\\}$. The condition
$v_{\lambda}\in{\mathcal{E}}_{\lambda}^{\perp}$ leads one to condition on
$p_{\lambda}$.

d. We assume that the functions $\omega$, $p_{\lambda}$ are defined on the
whole $\mathbb{R}$ and vanish outside of $(0,\delta)$. The condition
$v_{\lambda}\in{\mathcal{E}}_{\Lambda}$ takes the form

(5.8)
$0=\int_{0}^{T}[\omega(t)+\omega^{\prime}(t-1)]\left[\overline{p_{\lambda}(t-1)}+e^{-i\lambda
t}\int_{t-1}^{t}e^{-i\lambda\tau}\overline{p_{\lambda}(\tau)}d\tau\right]dt,$

for all $\omega\in W^{2}_{1}(0,\delta)$, $\omega(0)=\omega(\delta)=0$.

When opening the parenthesis we obtain four summands $A_{i}$, $i=1,2,3,4$.
Consider each of them separately.

$A_{1}=\int_{0}^{T}\omega(t)\overline{p_{\lambda}(t-1)}dt=0,$

since the supports of $\omega(t)$ and $p_{\lambda}(t-1)$ are disjoint.

Denote $\varepsilon_{\lambda}(t)=e^{-i\lambda t}\chi_{(0,\infty)}(t)$. Then

$A_{2}=\int_{0}^{T}\omega(t)\int_{t-1}^{t}\overline{p_{\lambda}(\tau)}e^{-i\lambda(t-\tau)}d\tau=\int_{0}^{\delta}\omega(t)(\overline{p_{\lambda}}*\varepsilon_{\lambda})(t)dt;$
$A_{3}=\int_{0}^{T}\omega^{\prime}(t-1)\overline{p_{\lambda}(t-1)}dt=-\int_{0}^{T}\omega(t-1)\overline{p_{\lambda}^{\prime}(t-1)}dt.$

Here we assume that $p_{\lambda}^{\prime}$ exists. With this assumption we
will arrive to equation with respect to $p_{\lambda}^{\prime}$ which has
smooth solution. This suffices since we already know that there is only one
$p_{\lambda}$, satisfying all conditions.

When considering $A_{4}$ we use (5.7). We have

$\displaystyle A_{4}=\int_{0}^{T}\omega^{\prime}(t-1)e^{-i\lambda
t}\int_{t-1}^{t}\overline{p_{\lambda}(\tau)}e^{i\lambda\tau}d\tau
dt=\int_{1}^{T}\omega^{\prime}(t-1)e^{-i\lambda
t}\int_{t-1}^{\delta}\overline{p_{\lambda}(\tau)}e^{i\lambda\tau}d\tau dt=$
$\displaystyle e^{-i\lambda}\int_{0}^{\delta}\omega^{\prime}(t)e^{-i\lambda
t}\left(1-\int_{0}^{t}\overline{p_{\lambda}(\tau)}e^{i\lambda\tau}d\tau\right)dt=e^{-i\lambda}\left[i\lambda\int_{0}^{\delta}\omega(t)e^{-i\lambda
t}dt-\int_{0}^{\delta}\omega(t)(\overline{p_{\lambda}}*\varepsilon_{\lambda})^{\prime}(t)dt\right].$

We collecting these relations together and remark that $\omega$ can be any
function in $W^{2}_{1}(0,\delta)$, satisfying (5.6). We then arrive to the
equation

(5.9)
$-\overline{p_{\lambda}}^{\prime}(t)-e^{-i\lambda}(\overline{p_{\lambda}}*\varepsilon_{\lambda})^{\prime}(t)+(\overline{p_{\lambda}}*\varepsilon_{\lambda})(t)-e^{-i\lambda}\varepsilon_{\lambda}^{\prime}(t)=0,\
0<t<\delta.$

The left-hand side of this equation is well-defined for all $t>0$. Denote it
by $\alpha(t)$. Let $\mathcal{P}_{\lambda}(s)$ and $e^{-\delta
s}\mathcal{A}(s)$ be the Laplace transform of $\overline{p_{\lambda}}$ and
$\alpha$ respectively. Taking once again into account that $D(\lambda)=0$, we
transform (5.9) into equation for $\mathcal{P}_{\lambda}(s)$ :

$\mathcal{P}_{\lambda}(s)(1-s^{2})=(\lambda^{2}-i\lambda\overline{p_{\lambda}(0)})+s\overline{p_{\lambda}(0)}+e^{-\delta
s}\mathcal{A}(s).$

The last summand in the right-hand side does not influence the values of
$\overline{p_{\lambda}}$ on $(0,\delta)$. Therefore

$\overline{p_{\lambda}}(t)=[\lambda^{2}-i\lambda\overline{p_{\lambda}}(0)]\sinh
t+\overline{p_{\lambda}}(0)\cosh t,\ t\in(0,\delta).$

The value $\overline{p_{\lambda}}(0)$ should now be chosen to meet (5.7). We
omit the corresponding calculation.

## References

  * [1] H. T. Banks, M. Q. Jacobs, C. E. Langenhop, Characterization of the controlled states in W(1)2 of linear hereditary systems, SIAM J. Control 13 (1975), 611–649.
  * [2] H. T. Banks, A. Manitius, Projection series for retarded functional differential equations with applications to optimal control problems J. Differential Equations 18 (1975), no. 2, 296–332.
  * [3] A. Baranov, Yu. Belov, A. Borichev, Hereditary completeness for systems of exponentials and reproducing kernels. Adv. Math. 235 (2013), 525–554.
  * [4] Yu. Belov, Yu. Lyubarskii, On summation of nonharmonic Fourier series, Constr. Approx. 43 (2016), no. 2, 291–309.
  * [5] A. Boccia, R. B.V̇inter, The maximum principle for optimal control problems with time delays. SIAM J. Control Optim. 55 (2017), no. 5, 2905–2935.
  * [6] J. A. Burns, T. L.; Herdman, H. W. Stech, Linear functional-differential equations as semigroups on product spaces, SIAM J. Math. Anal. 14 (1983), no. 1, 98–116.
  * [7] M. C. Delfour, The linear-quadratic optimal control problem with delays in state and control variables: a state space approach, SIAM J. Control Optim. 24 (1986), no. 5, 835–883.
  * [8] A. Halanay, Optimal controls for systems with time lag;SIAM J. Control, 6, 1968, 215–234.
  * [9] J. K. Hale, S. M. Verduyn Lunel Introduction to functional-differential equations, Springer-Verlag, New York, 1993. x+447 pp
  * [10] C. H. Lee, S. P. Yung, S. P. Sufficient conditions for optimal control problems with time delay, J. Optim. Theory Appl. 88 (1996), no. 1, 157–176.
  * [11] B. Ya. Levin, Lectures on entire functions, American Mathematical Society, Providence, RI, 1996.


# Adversarial methods to reduce simulation bias in neutrino interaction event
filtering at Liquid Argon Time Projection Chambers

M. Babicz<EMAIL_ADDRESS>H. Niewodniczanski Institute of Nuclear
Physics PAN, ul. Radzikowskiego 152, Krakow, Poland  CERN, Espl. des
Particules 1, Meyrin, Switzerland  S. Alonso-Monsalve
<EMAIL_ADDRESS>The author’s work presented in the manuscript
was done when the author was a member of the CERN neutrino group. ETH Zurich,
Rämistrasse 101, Zürich, Switzerland  CERN, Espl. des Particules 1, Meyrin,
Switzerland  S. Dolan<EMAIL_ADDRESS>CERN, Espl. des
Particules 1, Meyrin, Switzerland  K. Terao<EMAIL_ADDRESS>SLAC
National Accelerator Laboratory, Menlo Park, CA, 94025, USA

###### Abstract

For current and future neutrino oscillation experiments using large Liquid
Argon Time Projection Chambers (LAr-TPCs), a key challenge is identifying
neutrino interactions from the pervading cosmic-ray background. Rejection of
such background is often possible using traditional cut-based selections, but
this typically requires the prior use of computationally expensive
reconstruction algorithms. This work demonstrates an alternative approach of
using a 3D Submanifold Sparse Convolutional Network trained on low-level
information from the scintillation light signal of interactions inside LAr-
TPCs. This technique is applied to example simulations from ICARUS, the far
detector of the Short Baseline Neutrino (SBN) program at Fermilab. The results
of the network, show that cosmic background is reduced by up to 76.3% whilst
neutrino interaction selection efficiency remains over 98.9%. We further
present a way to mitigate potential biases from imperfect input simulations by
applying Domain Adversarial Neural Networks (DANNs), for which modified
simulated samples are introduced to imitate real data and a small portion of
them are used for adverserial training. A series of mock-data studies are
performed and demonstrate the effectiveness of using DANNs to mitigate biases,
showing neutrino interaction selection efficiency performances significantly
better than that achieved without the adversarial training.

##  I Introduction

The current and next generation of neutrino oscillation experiments offer a
tantalising opportunity to explore physics beyond the Standard Model. However,
as detectors grow larger and neutrino beams more powerful, a pre-filtering of
relevant neutrino interaction data becomes increasingly important for
experiments to be computationally viable. This is particularly crucial for
liquid Argon time projection chambers (LAr-TPCs), which are commonly used
neutrino detectors (see e.g., Rubbia:1977zz ; Willis:1974gi ; ArgoneuT:2016wjb
; Antonello:2015lea ; Abi:2020wmh ). LAr-TPCs offer precise spatial and
calorimetric measurements based on the electron drift signal from ionisation
and the scintillation photons from the excitation of Argon atoms caused by
interacting particles. However, detecting neutrino interactions with this
technology becomes challenging due to the significant background of incoming
cosmic rays.

For example, at the ICARUS detector of the SBN experiment, cosmic rays are
expected to outnumber neutrino interactions within the Booster Neutrino Beam’s
(BNB) spill gate by more than three to one Antonello:2015lea . Even located
1.5 km underground, the DUNE far detectors will experience a comparable rate
of cosmic rays and neutrinos Abi:2020wmh .

In LAr-TPCs, the ionisation electrons, stimulated by propagating charged
particles, are drifted by an applied electric field to be collected by TPC
anode wires, whilst the emitted LAr scintillation light is recorded by
photodetectors, often using Photomultiplier Tubes (PMT). As such, the TPC
records charged particle trajectories as images with a high spatial resolution
($\sim$mm/pixel) and the photodetectors provide event timing information with
nanosecond resolution. The scintillation light signal thereby provides an
easily accessible means to classify events requiring little or no processing,
which may help distinguish cosmic rays from neutrino interactions before
running any reconstruction algorithms.

The rejection of cosmic-ray backgrounds in LAr-TPCs typically starts at the
online stage. A trigger to record data is issued only if a fast signal from
the photodetectors is observed in the beam spill window, the time window
during which neutrino signal is expected. However, it does not prevent the
selection of background events caused by an accidental coincidence of a beam
spill window with incident cosmic rays. More sophisticated filtering may be
achieved through multi-dimensional analyses, as discrimination power can be
found by analysing which PMTs received light, at what time, and for how long
relative to all other PMTs inside the LAr-TPC.

Machine-learning methods, capable of analysing such high dimensional
information, are therefore excellently suited to classifying events using the
available PMT information. In particular, when detector data is represented as
images, the use of Convolutional Neural Networks (CNNs) firstCNN ; LeCun-et-
al-1998-gradient  are especially effective. The use of CNNs for event
classification is well established across the field of neutrino physics
Aurisano:2016jvx ; DUNE:2020gpm ; MicroBooNE:2018kka . Although CNNs are
broadly used in the neutrino community, images of neutrino interactions are
typically very sparse such that most of the pixels have empty values which can
render some standard methods ineffective. A straightforward solution to this
issue is to use submanifold sparse convolutional networks (SSCNs)
(graham2017submanifold, ), a variation of standard CNNs that use a new sparse
convolutional operator to efficiently handle sparse inputs, with already a
remarkable number of successful applications in neutrino physics
Adams-2019-deep ; Domine-2020-scalable . All the neural network architectures
shown in this paper belong to the class of SSCNs. Nevertheless, other
architectures, such as graph neural networks SperdutiFirstGNN ; zhou2018graph
, were also considered but were found to produce significantly worse results
in comparison to SSCNs for our application (although further optimisation of
these alternatives may have been able to improve performance).

CNNs111We refer to CNNs and SSCNs indistinguishably in the rest of this
section. can be optimized to discriminate signal images against backgrounds
through a supervised training process. This is often done using simulated
images (e.g., simulated neutrino and cosmic images) where the true labels are
available. However, when this model is applied to images from the real
detector, its performance is typically worse than what is observed on
simulated images because of discrepancies between two data domains (i.e.,
physics of the real world v.s. simulation) due to imperfect simulation. To
address this problem, the CNN classifier may rely on domain adaptation (DA)
techniques ben2010theory ; redko2019advances  so that the classifier learned
from the training domain (i.e., simulated data) can also be applied to the
testing domain (i.e., eventual experimental data). This DA can be achieved
through the application of Domain Adversarial Neural Networks (DANNs) dann ,
in which the detector data is used in an unsupervised (or semi-supervised)
manner to prevent the CNN exploiting features that differ between data and
simulation. DANNs were first used in neutrino physics by the MINER$\nu$A
experiment, where the bias of a deep-learning-based neutrino vertex
identification method was mitigated using these techniques perdue2018reducing
. In this manuscript, we present the first application of DANN for a CNN as an
event classifier for a LAr-TPC to discriminate neutrino signal against cosmic
backgrounds.

To test the effectiveness of CNNs and DANNs at distinguishing cosmic-ray
backgrounds from neutrino interactions using only information from the
scintillation signal, we consider the ICARUS detector of the SBN program
Antonello:2015lea  as a case study. ICARUS is currently the world’s largest
LAr-TPC employed in neutrino physics and operates close to the surface, and so
is subject to a particularly challenging cosmic-ray background rejection. The
details of the ICARUS detector and simulation are summarised in Sec. II. The
CNN approach to event filtering is detailed and demonstrated in Sec III. The
application of DANNs to reduce the CNN sensitivity to input simulation
dependence, and a method of using mock-data studies to test their
effectiveness, is then described and applied in Sec IV. Finally, the results
and the main conclusions of this work are presented in Sec V.

##  II Event filtering at the ICARUS detector

The ICARUS detector Amerio:2004ze  is a 760-ton LAr-TPC, serving as the far
detector of the SBN program Antonello:2015lea , positioned 600 m away from the
Booster Neutrino Beam (BNB) at FNAL. The detector consists of two identical
adjacent modules, each housing two TPCs separated by a common cathode used to
generate the electric field that directs the Argon ionisation signal to the
anode. The prompt (order of nanosecond) LAr scintillation light signal from
charged particles propagating within ICARUS is readout by 360 8 inches PMTs
Babicz:2018svg  arranged on the walls of the TPCs, placed as shown in Fig. 1.
The PMT system provides the means to trigger the readout of signals within the
1.6 $\mu$s beam spill windows whilst also enabling fast spatial localisation
of neutrino beam associated events. The placement, performance and timing
resolution of the PMTs are expected to allow the localisation of the
associated charged particle track with accuracy better than 1 m Babicz:2018svg
; Ali-Mohammadzadeh:2020fbd ; Babicz:2020rpf . The signals of all PMTs are
continuously readout, where pairs of adjacent PMTs are typically used together
within the ICARUS trigger system (the pairing scheme is shown in Fig. 1). For
each beam spill window, the ICARUS trigger system can assess which PMTs have a
signal exceeding a predefined threshold, at what time that signal was recorded
with respect to the start of the beam window and how many times the PMT
recorded an opening above the threshold. Note that a new opening is counted
every 0.16 $\mu$s that the PMT signal remains above a pre-set threshold and so
the number of openings acts as a discretised measurement of the time over the
threshold which itself is highly correlated with the signal amplitude.

Figure 1: A schematic view of a one wall segment of ICARUS TPCs showing how
the PMTs are distributed. The PMTs are arranged in groups of 15 that are
connected to the same digitiser board. The hollow ellipses show the pairing
system of adjacent PMTs. Figure 2: Example images of one ICARUS cryostat used
as an input to the CNN. Each dot represents a PMT pair position (taken as the
pair’s barycentre) which are distributed across the walls of the TPCs. The
colour of the dots represents the number of openings (left) or the opening
time (right).

The BNB delivers neutrinos in a bunch at the rate up to 5 Hz. The accelerator
complex informs the arrival time of each bunch, and hence the timing of the
beam spill window, during which the detector may issue a trigger if there is a
signal in photodetectors above the threshold. This results in about 5% of beam
spills being recorded. Despite a large reduction factor, the recorded events
are still dominated by cosmic-ray backgrounds due to an accidental coincidence
of the beam window and optical signals produced by cosmic rays. Our goal is to
further reduce the cosmic-ray backgrounds via CNN-based fast event filter that
only requires optical data.

Given the possibility of better separation of neutrino interactions from
cosmic rays after full reconstruction of the TPC signal, the primary aim of
the low-level event filter is to reduce the vast majority of cosmic rays while
not compromising neutrino selection efficiency. In this way, the amount of
data that needs to be processed for higher-level analyses is greatly reduced
while avoiding the risk of losing the neutrino signal. In order to allow the
possibility of online event filtering, only the information available to the
ICARUS event trigger is used. This means that information from each PMT pair
is concatenated where the earliest opening time from each pair is stored
alongside the total number of openings across both PMTs within the pair.

###  II.1 Simulation

The cosmic-ray particles impinging the ICARUS detector were generated with
CORSIKA event generator Heck:1998vt . These particles are then propagated
through the ICARUS detector and the surrounding material using GEANT4
Agostinelli:2002hh  implementation in LArSoft Snider:2017wjd . Scintillation
photons are then propagated to the PMTs using a parameterised model based on
pre-calculated tables (also derived from GEANT4). A further parameterised PMT
readout model, constrained from test-beam data, is then used to simulate the
digitised data-like signal from the detector.

The incoming flux of neutrinos is modelled using a GEANT4-based simulation of
the BNB beamline Antonello:2015lea ; AguilarArevalo:2008yp  whilst their
interactions with nuclei (and electrons) within ICARUS are described using
GENIE version 3 Andreopoulos:2009rq . The particle propagation and detector
response are simulated identically to the case of cosmic rays.

For this study, we simulate 396,200 PMT readout windows (events) containing
cosmic rays and 120,000 containing a single neutrino interaction. For this
work, only one of ICARUS’ two cryostats is considered. Whilst it is possible
to have a PMT readout window containing cosmic rays and a neutrino interaction
or multiple neutrino interactions, this is not particularly common, and such
details are beyond the scope of this study. The number of events used in this
work was reduced to those that passed the ICARUS trigger conditions, resulting
in 114,589 neutrino and 46,115 cosmic events.

##  III CNN event filter

The goal of the CNN is to classify whether events are from neutrino or cosmic-
ray interactions. To train the CNN, the simulated PMT data is presented as 3D
images, where the position of each PMT pair, alongside its opening time and a
number of openings, are stored. The image pixel size is chosen to be 40 cm,
which is the maximum distance such that two PMT pairs do not appear within the
same pixel. An example of the image provided as input, divided into two sub-
images representing the two PMT pair observable for better visualisation, is
shown in Fig. 2. Each event contains one image, which expresses both the
opening time and number of openings for each PMT pair, of which 80% are used
for training, 10% for validation, and the remaining 10% for testing. Each
image in the training sample is labelled as a cosmic-ray or neutrino event.

Figure 3: The sparse convolutional network architecture used for this
analysis. It was developed using the MinkowskiEngine package choy20194d  to
handle sparse inputs more efficiently.

The main feature of CNNs is that they learn a series of filters (using
convolutions), applied in sequence to extract increasingly powerful and
abstract features that allow the CNN to learn a mapping between input images
and target labels. Once the CNN is trained, it can be applied to new images to
make accurate predictions on unseen examples during the training. The designed
CNN architecture is depicted in Fig. 3. As introduced in Sec. I, it is based
on 3D Submanifold Sparse Convolutions graham2017submanifold  to deal with the
sparse images used in this work. The CNN is trained222The rest of this
paragraph applies to all the neural networks analysed in this article. for 50
epochs333Epoch: one forward pass and one backward pass of all the training
examples. In other words, an epoch is one pass over the entire dataset. For
this paper, we show the test results for the epoch that performs best on the
validation set. using Python 3.6.9 and PyTorch 2.1.0 Abadi-et-
al-2016-tensorflow , as well as the Minkowski Engine package, version 0.5.4
choy20194d , on an NVIDIA Tesla V100 GPUs. Stochastic Gradient Descent (SGD)
is used as the optimiser, with a mini-batch size of 32 events, a learning rate
of 0.1 (divided by 10 when the error plateaus, as suggested in He-et-
al-2015-deep ), a weight decay of 0.0001, and a momentum of 0.9444See Ref.
Goodfellow-et-al-2016-deep  for a description of optimisers and associated
terminology.. The model weights used for the analysis correspond to those at
the epoch that maximises the overall accuracy on the validation set. There are
no signs of overtraining.

###  III.1 Performance

Once trained, the output of the CNN is a continuous score for each event
between 0 (neutrino-like) and 1 (cosmic-like). The distribution of CNN scores
for each true event type in the test sample is shown in Fig. 4. If a selection
of neutrino events is made by cutting at a CNN score of 0.5, a 98.9% selection
efficiency is maintained whilst 76.3% of cosmic-ray backgrounds are rejected.
The charged-current selection efficiency was found to be flat (i.e., unbiased
by kinematics) in various tested observables. An example of outgoing lepton
angle is shown in Fig. 4.

##  IV Reducing model-dependence with DANN-based training

Figure 4:  The trained CNN’s classification of events in the test sample
(top). Note that the relative normalisation of cosmic and neutrino events is
fixed to approximate what would be expected in data. The distribution of
outgoing lepton angle (blue) with respect to the incoming neutrino from GENIE
is shown alongside the neutrino selection efficiency of the CNN (red)
following a cut at a classification score of 0.5 (bottom).  Figure 5: Domain-
adversarial neural network architecture. The feature extractor (blue) and the
label predictor (green) form the standard neural network classifier shown in
Fig. 3. The domain classifier (purple) provides the domain adaptation part
since it is connected to the feature extractor through a gradient reversal
layer, allowing the alignment of feature distributions across the source and
target domains. Figure adapted from dann .

Whilst the CNN presented in Sec. III.1 shows excellent performance, the
results assume perfect modelling of the neutrino and cosmic-ray events, the
particle propagation and the detector response. If the CNN is trained with
events that do not suitably represent what is in the real data, then the test
sample’s performance will not be reliable. Modern deep neural networks consist
of millions (and sometimes billions) of parameters and have a strong
representative capability with which they may exploit every detailed feature
present in the simulation, including those that may not be present in data as
well as others that may not follow the true physics model behind real data.
Thus, it is not easy to ever be sure that the pertinent aspects of the events
are well modelled. To alleviate this issue, adversarial training methods can
be employed to prevent neural networks from exploiting features that are only
present in one of two domains. As a result, the performance can be made
consistent in both domains. In this analysis, we show that it is possible to
mitigate challenges associated with domain discrepancies through the
application of DANNs.

In DANNs, the neural network model is trained on examples from two domains:
(a) the source domain, which consists of labelled simulated data; and (b) the
target domain, which consists of unlabelled true experimental data. The goal
is to learn a discriminator from the labelled source domain examples and use
the unlabelled target domain examples to ensure the discriminator relies on
only domain-invariant features to perform the predictions. Regarding the
implementation of the neural network, the classifier architecture remains
identical, and it can be seen as the combination of a feature extractor (i.e.,
the bulk of the CNN, in our case) and a label predictor (i.e., the sparse
linear layer(s) at the end). However, this alternative neural network has a
second path, which connects the output of the feature extractor through a
gradient reversal layer with a few linear layers that form a domain
classifier. The gradient reversal layer performs an identity transformation
during the forward propagation process but multiplies the gradient by a
negative constant during the back-propagation, guarantying that the parameters
learnt by the feature extractor are made similar for the source and target
distributions. In other words, with this approach, the features learnt by this
model are simultaneously discriminative - thanks to the label predictor - and
domain-invariant - thanks to the domain classifier. This behaviour is shown in
Fig. 5. Furthermore, if some events from the target distribution are labelled
(e.g., experimental data cosmic rays produced without a neutrino beam), those
events might be used for the feature extractor learning too, making the domain
adaptation semi-supervised, in contrast to the unsupervised case where all the
events from the target distribution are unlabelled.

In order to test the effectiveness of DANNs as a method of reducing simulation
dependence, we perform a series of mock-data studies. For these studies,
statistically independent simulations of events (from neutrinos and cosmic
rays) are produced before being modified to simulate possible mis-modelling
bias. Since the coarse PMT information used in this analysis is likely not
sensitive to the exact details of the neutrino interaction or cosmic-ray
production, we focus primarily on applying distortions to the simulated
detector response. The details of the mock data are as follows:

Figure 6:  Results of the classification score of the nominal and
adversarially trained CNN/DANN model applied to the original and: local noise
mock data using a 5% noise spread (left) and the global noise mock data using
a 10% noise spread (right). None of the cosmic-ray-mock-data events are
labelled by event type.

#### “Global noise” data:

in this mock data, noise, which is uncorrelated with the event content, is
randomly added to each PMT with some pre-specified ‘global’ probability that
is common to all PMTs. The global noise probabilities considered are 2%, 5%
and 10%. The timing of the noise is considered as uniform distribution, and if
noise is simulated to arrive before a PMT is opened by a simulated signal, the
opening time of the signal is overwritten by that of the noise.

#### “Local noise” data:

similarly to the global noise data, this mock data set considers the addition
of random noise to each PMT but where the probability of producing noise is
different for every PMT. Noise probabilities per event for each PMT were
generated randomly using a uniform distribution between 0 and either 2%, 5% or
10%.

For each mock-data study, the DANN is trained as described in Sec III but with
an addition of 22,778 cosmic-ray-mock-data events and 94,306 neutrino-mock-
data events (40% train, 10% validation, 50% test555Following the suggestion in
perdue2018reducing , we do not use all the available mock data events for
training. Moreover, the large test set size provides enough statistics for the
analysis.), which are labelled by the domain (i.e., mock data or original
simulation), all the simulated events are labelled by event type (i.e., cosmic
or neutrino), and we consider three scenarios for the event type of mock-data
events: (1) none of the events are labelled by event type (unsupervised domain
adaptation), (2) 10% of the cosmic-ray-mock-data events are labelled (semi-
supervised domain adaptation), and (3) 50% of the cosmic-ray-mock-data events
are labelled (semi-supervised domain adaptation). This method could equally be
applied to real data instead of mock data, using run periods with no neutrino
beam to label the cosmic-ray events.

Both the originally trained CNN (as described in Sec. III) and the newly
trained DANNs are used to attempt to classify events from the original sample
and from the new mock data sample. An example of the classification scores for
each model applied to the original and mock data sets is shown for two mock-
data studies in Fig. 6. A summary of the neutrino selection efficiency and the
background rejection performance for the nominal simulation as well as for
each mock data set is shown in Tab. 1. The presented numbers are provided for
a selection cut set to 0.5 of the network classification score in all cases.
It should be noted that there is additional freedom to optimise performance as
desired by tuning the cut value applied. For example, it was found that for
the case of mock data with global noise of 10%, the results using an
unsupervised domain adaptation improve with respect to the original model by
$\sim$14% for neutrino selection efficiency and $\sim$1% for cosmic-background
rejection by setting a cut at 0.25. However, with the cut at 0.5, the neutrino
selection efficiency improves dramatically (by $\sim$22%), but the background
rejection performance decreases (by $\sim$8%). An alternative assessment of
the CNN and DANN performance where the cut is varied to keep the background
rejection factor constant is presented in Appendix A.

Dataset | Model type | $\mathcal{E}_{\nu}$[%] | $\mathcal{R}_{cos}$[%] | Domain adaptation
---|---|---|---|---
Unsupervised | Semisup. 10% | Semisup. 50%
$\mathcal{E}_{\nu}$[%] |  $\mathcal{R}_{cos}$[%] |  $\mathcal{E}_{\nu}$[%] |  $\mathcal{R}_{cos}$[%] |  $\mathcal{E}_{\nu}$[%] |  $\mathcal{R}_{cos}$[%]
Nominal | CNN | 98.9 | 76.3 | - | - | -
DANNG2 | - | 97.2 | 77.2 | 99.0 | 74.8 | 97.9 | 77.2
DANNG5 | - | 98.3 | 76.8 | 98.0 | 76.6 | 97.9 | 76.9
DANNG10 | - | 98.7 | 75.8 | 98.1 | 76.4 | 98.5 | 76.4
DANNL2 | - | 98.2 | 76.1 | 98.1 | 75.8 | 98.7 | 76.1
DANNL5 | - | 98.8 | 75.9 | 98.7 | 76.0 | 98.8 | 76.2
DANNL10 | - | 98.7 | 75.2 | 98.1 | 76.8 | 98.0 | 76.8
Mock data |  Global noise | 2% | CNN | 91.7 | 74.8 | - | - | -
DANNG2 | - | 92.7 | 76.5 | 92.7 | 75.1 | 87.2 | 78.7
5% | CNN | 81.0 | 75.8 | - | - | -
DANNG5 | - | 89.7 | 72.1 | 89.9 | 76.6 | 84.1 | 78.7
10% | CNN | 66.4 | 79.0 | - | - | -
DANNG10 | - | 88.8 | 71.0 | 87.2 | 69.8 | 88.1 | 75.9
Local noise | 2% | CNN | 95.5 | 75.1 | - | - | -
DANNL2 | - | 97.6 | 74.2 | 98.2 | 74.8 | 96.5 | 76.6
5% | CNN | 90.2 | 74.8 | - | - | -
DANNL5 | - | 89.9 | 75.0 | 90.2 | 75.9 | 90.3 | 78.9
10% | CNN | 81.9 | 75.7 | - | - | -
DANNL10 | - | 90.2 | 73.3 | 88.7 | 77.9 | 88.4 | 78.9

Table 1: Efficiency ($\mathcal{E}_{\nu}$) and proportion of rejected cosmic-
ray background events ($\mathcal{R}_{cos}$) using the original and
adversarially trained CNN/DANN to classify events in the nominal simulations
and mock-data studies. Table 2 details the model-type nomenclature.

Network name | Description
---|---
CNN | Original neural network trained on the nominal simulation.
DANNG2 | Adversarial network trained on nominal simulation + mock data (global noise = 2%).
DANNG5 | Adversarial network trained on nominal simulation + mock data (global noise = 5%).
DANNG10 | Adversarial network trained on nominal simulation + mock data (global noise = 10%).
DANNL2 | Adversarial network trained on nominal simulation + mock data (local noise = 2%).
DANNL5 | Adversarial network trained on nominal simulation + mock data (local noise = 5%).
DANNL10 | Adversarial network trained on nominal simulation + mock data (local noise = 10%).

Table 2: Networks legend.

These results show that, without the adversarial training, the original CNN
can reject a sizeable portion of neutrino interactions in the mock data.
However, once the adversarial training is used, the network is able to
mitigate the bias and maintain a very high neutrino selection efficiency (the
main goal of the filter) whilst continuing to achieve a significant rejection
of cosmic-ray backgrounds. This occurs for the unsupervised and the semi-
supervised domain adaptations, and it is more visible for larger noises (i.e.,
10% of global noise and 10% of local noise). It can equally be observed that
the use of a DANN does not degrade the performance of the baseline CNN when
applied to the nominal simulation. This demonstrates that even if the
simulation very well describes the data, the use of a DANN over a CNN is not
expected to degrade the performance.

Concerning the unsupervised domain adaptation compared to the semi-supervised
cases, we find a small but non-negligible improvement in the rejection of
mock-data cosmic-ray backgrounds for some mock data studies, while suffering a
slight reduction in mock-data neutrino selection efficiency. This behaviour is
expected since, for the semi-supervised cases, the models have more labelled
cosmic events to learn to reject from (especially for the case where 50% of
the mock-data cosmic rays are labelled). It is possible that labelling a
larger portion of the training sample may allow improved performance from the
semi-supervision (recall that the majority of events in the training sample
are neutrino interactions, which cannot be labelled in real data).

##  V Conclusion

The studies presented in this manuscript demonstrate that easily accessible
information from LAr-TPC experiment’s light detection systems, which requires
very little processing, may be used to effectively separate neutrino from
cosmic-ray induced signals within a neutrino beam spill. The use of a
specially adapted CNN ensures that the majority of cosmic-ray interactions can
be filtered out without the rejection of almost any neutrino induced
interactions.

Whilst the use of a CNN trained on simulated event samples is susceptible to
bias due to mis-modelling, potentially causing the inadvertent rejection of
neutrino events, it is demonstrated that adversarial training via a DANN can
mitigate the loss of efficiency at the cost of some reduced background
rejection. It is further shown that in some cases, the background rejection
performance may be improved through semi-supervised domain adaptation of the
DANN using labelled real cosmic ray events.

Overall the techniques presented in this manuscript demonstrate a method for
providing a significant rejection of cosmic-ray events without the need for
computationally expensive reconstruction algorithms. These methods are shown
to be effective when applied to simulations from the ICARUS experiment, but
are easily adaptable and could likely achieve similar success if applied to
other LAr-TPC experiments.

## References

  * [1] C. Rubbia.  The Liquid Argon Time Projection Chamber: A New Concept for Neutrino Detectors.  5 1977.
  * [2] W. J. Willis and V. Radeka.  Liquid Argon Ionization Chambers as Total Absorption Detectors.  Nucl. Instrum. Meth., 120:221–236, 1974.
  * [3] R. Acciarri et al.  First Observation of Low Energy Electron Neutrinos in a Liquid Argon Time Projection Chamber.  Phys. Rev. D, 95(7):072005, 2017.
  * [4] M. Antonello et al.  A Proposal for a Three Detector Short-Baseline Neutrino Oscillation Program in the Fermilab Booster Neutrino Beam.  3 2015.
  * [5] B. Abi et al.  Deep Underground Neutrino Experiment (DUNE), Far Detector Technical Design Report, Volume I Introduction to DUNE.  JINST, 15(08):T08008, 2020.
  * [6] Y. LeCun et al.  Backpropagation applied to handwritten zip code recognition.  Neural Computation, 1(4):541–551, Dec 1989.
  * [7] Y. LeCun, L. Bottou, Y. Bengio, and P. Haffner.  Gradient-based learning applied to document recognition.  Proceedings of the IEEE, 86(11):2278–2324, 1998.
  * [8] A. Aurisano et al.  A Convolutional Neural Network Neutrino Event Classifier.  JINST, 11(09):P09001, 2016.
  * [9] B. Abi et al.  Neutrino interaction classification with a convolutional neural network in the DUNE far detector.  Phys. Rev. D, 102(9):092003, 2020.
  * [10] C. Adams et al.  Deep neural network for pixel-level electromagnetic particle identification in the MicroBooNE liquid argon time projection chamber.  Phys. Rev. D, 99(9):092001, 2019.
  * [11] B. Graham and L. van der Maaten.  Submanifold sparse convolutional networks.  2017\.
  * [12] C. Adams et al.  Deep neural network for pixel-level electromagnetic particle identification in the MicroBooNE liquid argon time projection chamber.  Phys. Rev. D, 99:092001, May 2019.
  * [13] L. Dominé and K. Terao.  Scalable deep convolutional neural networks for sparse, locally dense liquid argon time projection chamber data.  Phys. Rev. D, 102:012005, Jul 2020.
  * [14] A. Sperduti and A. Starita.  Supervised neural networks for the classification of structures.  IEEE Transactions on Neural Networks, 8(3):714–735, 1997.
  * [15] Jie Zhou, Ganqu Cui, Zhengyan Zhang, Cheng Yang, Zhiyuan Liu, Lifeng Wang, Changcheng Li, and Maosong Sun.  Graph neural networks: A review of methods and applications, 2018.
  * [16] S. Ben-David et al.  A theory of learning from different domains.  Machine learning, 79(1):151–175, 2010.
  * [17] I. Redko et al.  Advances in domain adaptation theory.  Elsevier, 2019.
  * [18] Y. Ganin et al.  Domain-adversarial training of neural networks.  J. Mach. Learn. Res., 17(1):2096–2030, January 2016.
  * [19] GN Perdue, A Ghosh, M Wospakrik, F Akbar, DA Andrade, M Ascencio, L Bellantoni, A Bercellie, M Betancourt, GFR Caceres Vera, et al.  Reducing model bias in a deep learning classifier using domain adversarial neural networks in the minerva experiment.  Journal of Instrumentation, 13(11):P11020, 2018.
  * [20] S. Amerio et al.  Design, construction and tests of the ICARUS T600 detector.  Nucl. Instrum. Meth. A, 527:329–410, 2004.
  * [21] M. Babicz et al.  Test and characterization of 400 Hamamatsu R5912-MOD photomultiplier tubes for the ICARUS T600 detector.  JINST, 13(10):P10030, 2018.
  * [22] B. Ali-Mohammadzadeh et al.  Design and implementation of the new scintillation light detection system of ICARUS T600.  JINST, 15(10):T10007, 2020.
  * [23] M. Babicz et al.  A particle detector that exploits Liquid Argon scintillation light.  Nucl. Instrum. Meth. A, 958:162421, 2020.
  * [24] D. Heck, J. Knapp, J. N. Capdevielle, G. Schatz, and T. Thouw.  CORSIKA: A Monte Carlo code to simulate extensive air showers.  2 1998.
  * [25] S. Agostinelli et al.  GEANT4–a simulation toolkit.  Nucl. Instrum. Meth. A, 506:250–303, 2003.
  * [26] E. L. Snider and G. Petrillo.  LArSoft: Toolkit for Simulation, Reconstruction and Analysis of Liquid Argon TPC Neutrino Detectors.  J. Phys. Conf. Ser., 898(4):042057, 2017.
  * [27] A. A. Aguilar-Arevalo et al.  The Neutrino Flux prediction at MiniBooNE.  Phys. Rev. D, 79:072002, 2009.
  * [28] C. Andreopoulos et al.  The GENIE Neutrino Monte Carlo Generator.  Nucl. Instrum. Meth. A, 614:87–104, 2010.
  * [29] C. Choy, J. Gwak, and S. Savarese.  4d spatio-temporal convnets: Minkowski convolutional neural networks.  2019\.
  * [30] M. Abadi et al.  Tensorflow: A system for large-scale machine learning.  In OSDI, volume 16, pages 265–283, 2016.
  * [31] K. He, X. Zhang, S. Ren, and J. Sun.  Deep residual learning for image recognition.  CoRR, 2015.
  * [32] I. Goodfellow, Y. Bengio, and A. Courville.  Deep Learning.  MIT Press, 2016.  http://www.deeplearningbook.org.

## Acknowledgements

The authors would like to thank the ICARUS collaboration for supporting this
work. Particular thanks is given to François Drieslma for providing detailed
comments on a draft version of this manuscript. The work of M. Babicz was
supported by the National Science Center, Poland, research project No.
2019/33/N/ST2/02874.

##  Appendix A Complementary results

Dataset | Model type | $\mathcal{E}_{\nu}$[%] | Domain adaptation
---|---|---|---
Unsupervised | Semisup. 10% | Semisup. 50%
$\mathcal{E}_{\nu}$[%] |  $\mathcal{E}_{\nu}$[%] |  $\mathcal{E}_{\nu}$[%]
Nominal | CNN | 99.2 | - | - | -
DANNG2 | - | 98.3 | 98.9 | 98.8
DANNG5 | - | 99.0 | 98.7 | 98.9
DANNG10 | - | 99.0 | 98.6 | 98.7
DANNL2 | - | 98.6 | 98.3 | 99.1
DANNL5 | - | 99.1 | 98.9 | 99.0
DANNL10 | - | 98.7 | 98.8 | 98.8
Mock data |  Global noise | 2% | CNN | 91.5 | - | - | -
DANNG2 | - | 93.6 | 92.7 | 89.6
5% | CNN | 81.7 | - | - | -
DANNG5 | - | 88.0 | 90.8 | 87.3
10% | CNN | 74.7 | - | - | -
DANNG10 | - | 86.0 | 84.5 | 88.6
Local noise | 2% | CNN | 95.4 | - | - | -
DANNL2 | - | 97.1 | 97.9 | 97.1
5% | CNN | 89.9 | - | - | -
DANNL5 | - | 89.9 | 90.6 | 93.7
10% | CNN | 82.6 | - | - | -
DANNL10 | - | 89.1 | 90.5 | 91.3

Table 3: Efficiency ($\mathcal{E}_{\nu}$) using the original and adversarially
trained CNN/DANN to classify events in the nominal simulations and mock-data
studies. In contrast to what is shown in Table 1 and for a better comparison
of the efficiencies, the classification score cuts were tuned so that the
proportion of rejected cosmic-ray background events for each model is always
75%. Table 2 details the model-type nomenclature.

As discussed in Sec. IV, the performance of the CNN and DANNs’ application to
the mock-data sets may be demonstrated in alternative ways to as presented in
Tab. 1, which shows the efficiency and background rejection achieved for a
fixed cut in the networks’ output classification score distribution. As such,
Tab. 3 instead changes the cut such that the background rejection remains
fixed (at 75%) so that the efficiencies can be more directly compared. The
conclusions remain unchanged from those presented in Sec. IV: the performance
improvement offered by the adversarial training of the DANNs is substantial
with respect to naively applying CNN trained only on the input simulation. The
improvement can be seen to be stronger for more extreme fake data studies.
Labelling some proportion of the cosmic ray events in the mock-data samples to
provide a semi-supervised of the DANNs can offer a small additional
improvement in some cases.


* (18) G. S. Bali, S. Collins, A. Cox and A. Schäfer, _Masses and decay constants of the $D_{s0}^{*}(2317)$ and $D_{s1}(2460)$ from $N_{f}=2$ lattice QCD close to the physical point_, _Phys. Rev._ D96 (2017) 074501 [1706.01247].
  * (19) UKQCD collaboration, M. Foster and C. Michael, _Quark mass dependence of hadron masses from lattice QCD_ , _Phys. Rev._ D59 (1999) 074503 [hep-lat/9810021].
  * (20) G. Martinelli and C. T. Sachrajda, _A Lattice Study of Nucleon Structure_ , _Nucl. Phys._ B316 (1989) 355.
  * (21) Y.-B. Yang, A. Alexandru, T. Draper, M. Gong and K.-F. Liu, _Stochastic method with low mode substitution for nucleon isovector matrix elements_ , _Phys. Rev._ D93 (2016) 034503 [1509.04616].
  * (22) S. Gusken, _A Study of smearing techniques for hadron correlation functions_ , _Nucl. Phys. Proc. Suppl._ 17 (1990) 361.
  * (23) G. S. Bali, B. Lang, B. U. Musch and A. Schäfer, _Novel quark smearing for hadrons with high momenta in lattice QCD_ , _Phys. Rev._ D93 (2016) 094515 [1602.05525].
  * (24) QCDSF/UKQCD collaboration, D. Brömmel et al., _The Pion form-factor from lattice QCD with two dynamical flavours_ , _Eur. Phys. J._ C51 (2007) 335 [hep-lat/0608021].
  * (25) C. Thron, S. J. Dong, K. F. Liu and H. P. Ying, _Padé - $Z(2)$ estimator of determinants_, _Phys. Rev._ D57 (1998) 1642 [hep-lat/9707001].
  * (26) SESAM collaboration, G. S. Bali, H. Neff, T. Düssel, T. Lippert and K. Schilling, _Observation of string breaking in QCD_ , _Phys. Rev._ D71 (2005) 114513 [hep-lat/0505012].
  * (27) G. S. Bali, S. Collins and A. Schäfer, _Effective noise reduction techniques for disconnected loops in Lattice QCD_ , _Comput. Phys. Commun._ 181 (2010) 1570 [0910.3970].
  * (28) M. Göckeler et al., _Perturbative and Nonperturbative Renormalization in Lattice QCD_ , _Phys. Rev._ D82 (2010) 114511 [1003.5756].
  * (29) S. Capitani, M. Göckeler, R. Horsley, H. Perlt, P. E. L. Rakow, G. Schierholz et al., _Renormalization and off-shell improvement in lattice perturbation theory_ , _Nucl. Phys._ B593 (2001) 183 [hep-lat/0007004].
  * (30) S. Sint and P. Weisz, _Further results on O(a) improved lattice QCD to one loop order of perturbation theory_ , _Nucl. Phys._ B502 (1997) 251 [hep-lat/9704001].
  * (31) Y. Taniguchi and A. Ukawa, _Perturbative calculation of improvement coefficients to $O(g^{2}a)$ for bilinear quark operators in lattice QCD_, _Phys. Rev._ D58 (1998) 114503 [hep-lat/9806015].
  * (32) P. Fritzsch, J. Heitger and N. Tantalo, _Non-perturbative improvement of quark mass renormalization in two-flavour lattice QCD_ , _JHEP_ 08 (2010) 074 [1004.3978].
  * (33) R. A. Briceño, J. V. Guerrero, M. T. Hansen and C. J. Monahan, _Finite-volume effects due to spatially non-local operators_ , 1805.01034.
  * (34) K. Cichy, K. Jansen and P. Korcyl, _Non-perturbative renormalization in coordinate space for $N_{f}=2$ maximally twisted mass fermions with tree-level Symanzik improved gauge action_, _Nucl. Phys._ B865 (2012) 268 [1207.0628].
  * (35) M. Burkardt, _Impact parameter dependent parton distributions and off forward parton distributions for $\zeta\to 0$_, _Phys. Rev._ D62 (2000) 071503 [hep-ph/0005108].
  * (36) S. Dürr, _Validity of ChPT – is $M_{\pi}$=135 MeV small enough?_, _PoS_ LATTICE2014 (2015) 006 [1412.6434].
  * (37) S. Aoki et al., _Review of lattice results concerning low-energy particle physics_ , _Eur. Phys. J._ C74 (2014) 2890 [1310.8555].
  * (38) J. F. Donoghue, C. Ramirez and G. Valencia, _The Spectrum of QCD and Chiral Lagrangians of the Strong and Weak Interactions_ , _Phys. Rev._ D39 (1989) 1947.
  * (39) G. Ecker, J. Gasser, A. Pich and E. de Rafael, _The Role of Resonances in Chiral Perturbation Theory_ , _Nucl. Phys._ B321 (1989) 311.
  * (40) SciDAC, LHPC, UKQCD collaboration, R. G. Edwards and B. Joo, _The Chroma software system for lattice QCD_ , _Nucl. Phys. Proc. Suppl._ 140 (2005) 832 [hep-lat/0409003].
  * (41) M. Lüscher and S. Schaefer, _Lattice QCD with open boundary conditions and twisted-mass reweighting_ , _Comput. Phys. Commun._ 184 (2013) 519 [1206.2809].
  * (42) Y. Nakamura and H. Stüben, _BQCD - Berlin quantum chromodynamics program_ , _PoS_ LATTICE2010 (2010) 040 [1011.0199].
  * (43) Y. Nakamura, A. Nobile, D. Pleiter, H. Simma, T. Streuer, T. Wettig et al., _Lattice QCD Applications on QPACE_ , 1103.1363.
  * (44) M. Matteo Frigo and S. G. Johnson, _The design and implementation of FFTW3_ , _Proc. IEEE_ 93 (2005) 216–231.
  * (45) D. Binosi and L. Theussl, _JaxoDraw: A Graphical user interface for drawing Feynman diagrams_ , _Comput. Phys. Commun._ 161 (2004) 76 [hep-ph/0309015].
  * (46) D. Binosi, J. Collins, C. Kaufhold and L. Theussl, _JaxoDraw: A Graphical user interface for drawing Feynman diagrams. Version 2.0 release notes_ , _Comput. Phys. Commun._ 180 (2009) 1709 [0811.4113].