# HyperTime: Implicit Neural Representation
for Time Series

Elizabeth Fons
J. P. Morgan AI Research
&Alejandro Sztrajman
University College London
&Yousef El-laham
J. P. Morgan AI Research
&Alexandros Iosifidis
Aarhus University
&Svitlana Vyetrenko
J. P. Morgan AI Research

###### Abstract

Implicit neural representations (INRs) have recently emerged as a powerful
tool that provides an accurate and resolution-independent encoding of data.
Their robustness as general approximators has been shown in a wide variety of
data sources, with applications on image, sound, and 3D scene representation.
However, little attention has been given to leveraging these architectures for
the representation and analysis of time series data. In this paper, we analyze
the representation of time series using INRs, comparing different activation
functions in terms of reconstruction accuracy and training convergence speed.
We show how these networks can be leveraged for the imputation of time series,
with applications on both univariate and multivariate data. Finally, we
propose a hypernetwork architecture that leverages INRs to learn a compressed
latent representation of an entire time series dataset. We introduce an FFT-
based loss to guide training so that all frequencies are preserved in the time
series. We show that this network can be used to encode time series as INRs,
and their embeddings can be interpolated to generate new time series from
existing ones. We evaluate our generative method by using it for data
augmentation, and show that it is competitive against current state-of-the-art
approaches for augmentation of time series.

##  1 Introduction

Modeling time series data has been a key topic of research for many years,
constituting a crucial component of applications in a wide variety of areas
such as climate modeling, medicine, biology, retail and finance [22].
Traditional methods for time series modeling have relied on parametric models
informed by expert knowledge. However, the development of modern machine
learning methods has provided purely data-driven techniques to learn temporal
relationships. In particular, neural network-based methods have gained
popularity in recent times, with applications on a wide range of tasks, such
as time series classification [18], clustering [27, 2], segmentation [32, 50],
anomaly detection [11, 47, 17], upsampling [30, 6], imputation [24, 25, 9],
forecasting [22, 44] and synthesis [1, 48, 21]. In particular, the generation
of time series data for augmentation has remained as an open problem, and is
currently gaining interest due to the large number of potential applications
such as in medical and financial datasets, where data cannot be shared, either
for privacy reasons or for proprietary restrictions [20, 21, 3, 12].

In recent years, implicit neural representations (INRs) have gained popularity
as an accurate and flexible method to parameterize signals, such as from
image, video, audio and 3D scene data [37, 29]. Conventional methods for data
encoding often rely on discrete representations, such as data grids, which are
limited by their spatial resolution and present inherent discretization
artifacts. In contrast, implicit neural representations encode data in terms
of continuous functional relationships between signals, and thus are uncoupled
to spatial resolution. In practical terms, INRs provide a new data
representation framework that is resolution-independent, with many potential
applications on time series data, where irregularly sampled and missing data
are common occurrences [14]. However, there are currently no works exploring
the suitability of INRs on time series representation and analysis.

In this work, we propose an implicit neural representation for univariate and
multivariate time series data. We compare the performance of different
activation functions in terms of reconstruction accuracy and training
convergence, and we formulate and compare different strategies for data
imputation in time series, relying on INRs (Section 4.2). Finally, we combine
these representations with a hypernetwork architecture, in order to learn a
prior over the space of time series. The training of our hypernetwork takes
into account the accurate reconstruction of both the time series signals and
their respective power spectra. This motivates us to propose a Fourier-based
loss that proves to be crucial in guiding the learning process. The advantage
of employing such a Fourier-based loss is that it allows our hypernetwork to
preserve all frequencies in the time series representation. In Section 4.3, we
leverage the latent embeddings learned by the hypernetwork for the synthesis
of new time series by interpolation, and show that our method performs
competitively against recent state-of-the-art methods for time series
augmentation.

##  2 Related Work

#### Implicit Neural Representations

INRs (or coordinate-based neural networks) have recently gained popularity in
computer vision applications. The usual implementation of INRs consists of a
fully-connected neural network (MLP) that maps coordinates (_e.g_. xyz-
coordinates) to the corresponding values of the data, essentially encoding
their functional relationship in the network. One of the main advantages of
this approach for data representation, is that the information is encoded in a
continuous/grid-free representation, that provides a built-in non-linear
interpolation of the data. This avoids the usual artifacts that arise from
discretization, and has been shown to combine flexible and accurate data
representation, with high memory efficiency [37, 43]. Whilst INRs have been
shown to work on data from diverse sources, such as video, images and audio
[37, 10, 33], their recent popularity has been motivated by multiple
applications in the representation of 3D scene data, such as 3D geometry [31,
28, 36, 38] and object appearance [29, 41].

In early architectures, INRs showed a lack of accuracy in the encoding of
high-frequency details of signals. Mildenhall  _et al_. [29] proposed
positional encodings to address this issue, and Tancik  _et al_. [43] further
explored them, showing that by using Fourier-based features in the input
layer, the network is able to learn the full spectrum of frequencies from
data. Concurrently, Sitzmann  _et al_. [37] tackled the encoding of high-
frequency data by proposing the use of sinusoidal activation functions (SIREN:
Sinusoidal Representation Networks), and Benbarka  _et al_. [7] showed the
equivalence between Fourier features and single-layer SIRENs. Our INR
architecture for time series data (Section 3.1) is based on the SIREN
architecture by Sitzmann  _et al_. In Section 4.1 we compare the performance
of different activation layers, in terms of reconstruction accuracy and
training convergence speed, for both univariate and multivariate time series.

#### Hypernetworks

A hypernetwork is a neural network architecture designed to predict the weight
values of a secondary neural network, denominated a HypoNetwork [36]. The
concept of hypernetwork was formalized by Ha  _et al_. [16], drawing
inspiration from methods in evolutionary computing [40]. Moreover, while
convolutional encoders have been likened to the function of the human visual
system [39], the analogy cannot be extended to convolutional decoders, and
some authors have argued that hypernetworks much more closely match the
behavior of the prefrontal cortex [34]. Hypernetworks have been praised for
their expressivity, compression due to weight sharing, and for their fast
inference times[39]. They have been leveraged for multiple applications,
including few-shot learning [35, 52], continual learning [45] and architecture
search [51, 8]. Moreover, in the last two years some works have started to
leverage hypernetworks for the training of INRs, enabling the learning of
latent encodings of data, while also maintaining the flexible and accurate
reconstruction of signals provided by INRs. This approach has been implemented
with different hypernetwork architectures, to learn priors over image data
[37, 39], 3D scene geometry [23, 38, 36] and material appearance [41]. Tancik
_et al_. [42] leverage hypernetworks to speed-up the training of INRs by
providing learned initializations of the network weights. Sitzmann  _et al_.
[37] combine a set encoder with a hypernetwork decoder to learn a prior over
INRs representing image data, and apply it for image in-painting. Our
hypernetwork architecture from Section 3 is similar to Sitzmann  _et al_.’s,
however we learn a prior over the space of time series and leverage it for new
data synthesis through interpolation of the learned embeddings. Furthermore,
our architecture implements a Fourier-based loss, which we show is crucial for
the accurate reconstruction of time series datasets (Section 4.3).

#### Time Series Generation

Realistic time series generation has been previously studied in the literature
by using the generative adversarial networks (GANs). In TimeGAN architecture
[49], realistic generation of temporal patterns was achieved by jointly
optimizing with both supervised and adversarial objectives to learn an
embedding space. QuantGAN [46] consists of a generator and discriminator
functions represented by temporal convolutional networks, which allows it to
synthesize long-range dependencies such as the presence of volatility clusters
that are characteristic of financial time series. TimeVAE [13] was recently
proposed as a variational autoencoder alternative to GAN-based time series
generation. GANs and VAEs are typically used for creating statistical replicas
of the training data, and not the distributionally new scenarios needed for
data augmentation. More recently, Alaa  _et al_. [1] presented Fourier Flows,
a normalizing flow model for time series data that leverages the frequency
domain representation, currently considered together with TimeGAN as state-of-
the-art for time series augmentation.

Data augmentation is well established in computer vision tasks due to the
simplicity of label-preserving geometric image transformation techniques, but
it is still not widely used for time series with some early work being
discussed in the literature [19]. For example, simple augmentation techniques
applied to financial price time series such as adding noise or time warping
were shown to improve the quality of next day price prediction model. [15]

##  3 Formulation

In this Section we describe the network architectures that we use to encode
time series data (Subsection 3.1), and the hypernetwork architecture
(HyperTime) leveraged for prior learning and new data generation (Subsection
3.2).

###  3.1 Time Series Representation

In Figure 1 we present a diagram of the INR used for univariate time series.
The network is composed of fully-connected layers of dimensions $1\times
60\times 60\times 60\times 1$, with sine activations (SIREN [37]):

$\phi_{i}(\mathbf{x}_{i})=\text{sin}(\omega_{0}\mathbf{W}_{i}\mathbf{x}_{i}+\mathbf{b}_{i})$
(1)

where $\phi_{i}$ corresponds to the $i^{th}$ layer of the network. A general
factor $\omega_{0}$ multiplying the network weights determines the order of
magnitude of the frequencies that will be used to encode the signal. Input and
output of the INR are uni-dimensional, and correspond to the time coordinate
$t$ and the time series evaluation $f(t)$. Training of the network is done in
a supervised manner, with MSE loss, and takes less than $10s$ in a GeForce GTX
1080 Ti GPU. After training, the network encodes a continuous representation
of the functional relationship $f(t)$ for a single time series.

Figure 1: Diagram of the implicit neural representation (INR) for univariate
time series. Neurons with a black border use sine activations. The network is
composed of fully-connected layers with dimensions $1\times 60\times 60\times
60\times 1$.

The architecture from Figure 1 can be modified to encode multivariate time
series, by simply increasing the number of neurons of the output layer to
match the number of channels of the signal. Due to weight-sharing, this adds a
potential for data compression of the time series. In addition, the
simultaneous encoding of multiple correlated channels can be leveraged for
channel imputation, as we will show in Section 4.2.

###  3.2 Time Series Generation with HyperTime

In Figure 2 we display a diagram of the HyperTime architecture, which allows
us to leverage INRs to learn priors over the space of time series. The Set
Encoder (green network), composed of SIREN layers [37] with dimensions
$2\times 128\times 128\times 40$, takes as input a pair of values,
corresponding to the time-coordinate $t$ and the time series signal $f(t)$.
Each pair of input values is thus encoded into a full $40$-values embedding
and fed to the HyperNet decoder (red network), composed of fully-connected
layers with ReLU activations (MLP), with dimensions $40\times 128\times 7500$.
The output of the HyperNet is a one-dimensional $7500$-values embedding that
contains the network weights of an INR which encodes the time series data from
the input. The INR architecture used within HyperTime is the same described in
the previous section, and illustrated in Figure 1. Following previous works
[36], in order to avoid ambiguities we refer to these predicted INRs as
HypoNets.

Figure 2: Diagram of HyperTime network architecture. Each pair of time-
coordinate $t$ and time series $f(t)$ is encoded as a 40-values embedding $Z$
by the Set Encoder. The HyperNet decoder learns to predict HypoNet weights
from the embeddings. During training, the output of the HyperNet is used to
build a HypoNet and evaluate it on in the input time-coordinates. The loss is
computed as a difference between $f(t)$ and the output of the HypoNet
$\hat{f}(t)$.

During the training of HyperTime, we use the weights predicted by the HyperNet
decoder to instantiate a HypoNet and evaluate it on the input time-coordinate
$t$, to produce the predicted time series value $\hat{f}(t)$. The entire chain
of operations is implemented within the same differentiable pipeline, and
hence the training loss can be computed as the difference between the ground
truth time series signal $f(t)$ and the value predicted by the HypoNet
$\hat{f}(t)$.

After the training of HyperTime, the Set Encoder is able to generate latent
embeddings $Z$ for entire time series. In Section 4.3, we show that these
embeddings can be interpolated to synthesize new time series signals from
known ones, which can be leveraged for data augmentation (see additional
material for a pseudo-code of the procedure).

#### Loss

The training of HyperTime is done by optimizing the following loss, which
contains an MSE reconstruction term $\mathcal{L}_{\text{rec}}$ and two
regularization terms $\mathcal{L}_{\text{weights}}$ and
$\mathcal{L}_{\text{latent}}$, for the network weights and the latent
embeddings respectively:

$\mathcal{L}=\underbrace{\frac{1}{N}\sum_{i=1}^{N}\left\|f(t_{i})-\hat{f}(t_{i})\right\|^{2}}_{\mathcal{L}_{\text{rec}}}+\lambda_{1}\underbrace{\frac{1}{W}\sum_{j=1}^{W}w^{2}_{j}}_{\mathcal{L}_{\text{weights}}}+\lambda_{2}\underbrace{\frac{1}{Z}\sum_{k=1}^{Z}z^{2}_{k}}_{\mathcal{L}_{\text{latent}}}+\lambda_{3}\mathcal{L}_{\text{FFT}}$
(2)

In addition, we introduce a Fourier-based loss $\mathcal{L}_{\text{FFT}}$ that
focuses on the accurate reconstruction of the power spectrum of the ground
truth signal (see Supplement for more details):

$\mathcal{L}_{\text{FFT}}=\frac{1}{N}\sum_{i=1}^{N}\left\|\text{FFT}[f(t)]_{i}-\text{FFT}[\hat{f}(t)]_{i}\right\|.$
(3)

In Section 4.3, we show that $\mathcal{L}_{\text{FFT}}$ is crucial for the
accurate reconstruction of the time series signals.

##  4 Experiments

###  4.1 Reconstruction

We start by showing that encoding time series using SIRENS leads to a better
reconstruction error than using implicit networks with other activations. We
use univariate and multivariate time series datasets from the UCR archive
[4].111The datasets can be downloaded from the project’s website:
www.timeseriesclassification.com [5] We selected datasets with different
characteristics, either short length time series or long, or in the case of
the multivariate datasets, with many features (in some cases, more features
than time series length). We sample 300 time series (or the maximum number
available) from each dataset, train a single SIREN for each time series and
calculate the reconstruction error. For comparison we train implicit networks
using ReLU, Tanh and Sigmoid activations. As a sample case, we show in Figure
3 the losses (left) and reconstruction plots (right) for one of the univariate
datasets (NonInvasiveFetalECGThorax1). Here we observe that sine activations
converge much faster, and to lower error values, than other activation
functions (for the full set of loss and reconstruction plots, along with a
description of each dataset, see additional materials). A summary of results
can be found in Table 1, where we observe that the MSE error is at least an
order of magnitude lower for sine activations, with respect to other
activation layers.

Losses Reconstructions

MSE (log scale) Time-series value

Iterations Time Steps

Figure 3: Comparison of implicit networks using different activation functions, MSE loss (left) and reconstructed time-series (right). Table 1: Comparison using MSE of implicit networks using different activation functions on different univariate and multivariate time series from the UCR dataset. | Sine | ReLU | Tanh | Sigmoid
---|---|---|---|---
Univariate |  |  |  |
Crop | 5.1e-06 | 5.4e-03 | 2.8e-02 | 5.1e-01
NonInvasiveFetalECGThorax1 | 2.3e-05 | 2.8e-02 | 5.7e-02 | 8.1e-02
PhalangesOutlinesCorrect | 7.5e-06 | 1.9e-02 | 1.4e-01 | 3.3e-01
FordA | 9.2e-06 | 1.4e-01 | 1.5e-01 | 1.5e-01
Multivariate |  |  |  |
Cricket | 1.6e-04 | 4.2e-03 | 5.1e-03 | 1.6e-02
DuckDuckGeese | 9.1e-05 | 8.0e-04 | 8.7e-04 | 9.1e-04
MotorImagery | 1.7e-03 | 1.1e-02 | 1.1e-02 | 1.8e-02
PhonemeSpectra | 1.1e-06 | 6.0e-03 | 1.6e-02 | 1.8e-02

###  4.2 Imputation

Real world time series data tend to suffer from missing values, in some cases
achieving missing rates of up to $90\%$, making the data difficult to use
[14]. Motivated by this, we first use a SIREN trained with the available data
points to infer the missing values. Besides fitting a SIREN to the time
series, we include a simple regularization term based on a total variation
(TV) prior [26]. The TV prior samples random points $t_{i}\in[-1,1]$ from the
available points of the time series (those that do not have missing values)
and the regularization consists of the L1 norm on the gradient

$p_{TV}=\frac{1}{N}\sum_{i=1}^{N}|\nabla\Phi_{\theta}(t_{i})|.$

We compare this approach with common time series imputation methods, mean,
that fills missing values with the average, kNN, cubic spline and linear. We
are not solely interested in the reconstruction error, but also in the
accurate and smooth reproduction of the Fourier spectrum composition of the
original time series. Therefore, we evaluate both the MSE of the time series
reconstruction, and the MAE between the Fourier spectra of ground truth and
reconstructed, that we name Fourier Error (FFTE).

Table 2 shows the comparison of SIREN with and without prior with the
classical methods using different ratios of missing values. We can see that
the reconstruction error from baseline methods tends to be low, although when
we compare the spectrum through FFTE, we can see that the lowest errors are
achieved by SIREN, except when we have a fraction of missing values of $90\%$,
where in small time series such as Crop or Phalanges, that have less than 100
time steps, both reconstruction error and FFTE are poor. With regard to adding
a prior over the gradient, we can see that for very short time series such as
Crop and Phalanges (which are 46 and 80 time steps of length, respectively),
the prior improves the reconstruction error when compared to SIREN without a
prior, and this is also the case for FFTE. Figure 4 shows the comparison of
SIREN, SIREN plus TV prior and Linear for a randomly selected time series from
the NonInvasiveFetalECGThorax1 (NonInv) dataset with $70\%$ of data points
missing (the plot is zoomed in a segment of 250 timesteps to highlight the
imputation characteristics). We can see that the linear imputation is not
smooth, so it’s not a good representation of the series, while SIREN is smooth
but tends to have higher deviations from the series. SIREN plus the TV prior
is a good compromise between maintaining the smoothness of the series while
also deviating less from the ground truth.

Table 2: Comparison of implicit neural representations (INRs) with classical
imputation methods using different fractions of missing values. Each method
shows the MSE reconstruction errors and the Fourier reconstruction errors.

|  | SIREN | SIREN_TV | Mean | kNN | Cubic Spline | Linear
---|---|---|---|---|---|---|---
|  | MSE | FFTE | MSE | FFTE | MSE | FFTE | MSE | FFTE | MSE | FFTE | MSE | FFTE
Crop | 0.0 | 5.5e-07 | 0.00 | 1.7e-06 | 0.01 | - | - | - | - | - | - | - | -
0.1 | 6.5e-03 | 0.43 | 1.2e-03 | 0.18 | 1.7e-02 | 3.01 | 1.7e-02 | 3.03 | 6.7e-07 | 2.25 | 2.2e-06 | 2.25
0.5 | 1.2e-01 | 1.80 | 1.3e-02 | 0.63 | 1.7e-02 | 3.02 | 1.6e-02 | 2.97 | 8.3e-05 | 2.25 | 7.6e-05 | 2.25
0.7 | 2.9e-01 | 2.52 | 3.8e-02 | 0.95 | 1.7e-02 | 3.03 | 1.7e-02 | 3.03 | 9.5e-03 | 2.25 | 8.6e-04 | 2.24
0.9 | 5.4e-01 | 2.74 | 3.5e-01 | 2.02 | 1.7e-02 | 3.02 | 1.7e-02 | 3.00 | 1.1e+02 | 2.49 | 1.6e-02 | 2.24
NonInv | 0.0 | 1.5e-06 | 0.02 | 3.5e-06 | 0.04 | - | - | - | - | - | - | - | -
0.1 | 2.1e-06 | 0.03 | 4.7e-06 | 0.04 | 1.7e-02 | 3.00 | 1.6e-02 | 2.97 | 3.5e-07 | 2.25 | 2.2e-06 | 2.25
0.5 | 9.2e-05 | 0.15 | 7.6e-05 | 0.14 | 1.7e-02 | 3.00 | 1.7e-02 | 3.02 | 1.1e-04 | 2.25 | 1.1e-04 | 2.25
0.7 | 1.2e-03 | 0.46 | 1.3e-03 | 0.45 | 1.7e-02 | 2.99 | 1.7e-02 | 2.99 | 1.7e-02 | 2.25 | 6.9e-04 | 2.24
0.9 | 1.5e-02 | 1.40 | 1.6e-02 | 1.55 | 1.6e-02 | 2.98 | 1.6e-02 | 2.97 | 2.8e-01 | 2.27 | 1.8e-02 | 2.24
Phalanges | 0.0 | 4.3e-07 | 0.00 | 1.7e-06 | 0.01 | - | - | - | - | - | - | - | -
0.1 | 9.0e-04 | 0.22 | 7.6e-04 | 0.20 | 1.7e-02 | 3.00 | 1.7e-02 | 3.05 | 3.4e-07 | 2.25 | 1.7e-06 | 2.25
0.5 | 2.8e-02 | 1.12 | 1.1e-02 | 0.74 | 1.7e-02 | 3.00 | 1.7e-02 | 3.04 | 3.5e-05 | 2.25 | 1.1e-04 | 2.25
0.7 | 1.1e-01 | 1.94 | 4.3e-02 | 1.27 | 1.6e-02 | 2.99 | 1.7e-02 | 2.99 | 2.3e-03 | 2.25 | 6.2e-04 | 2.24
0.9 | 3.3e-01 | 2.75 | 2.2e-01 | 2.26 | 1.6e-02 | 2.99 | 1.7e-02 | 2.99 | 7.1e-02 | 2.26 | 1.4e-02 | 2.24
FordA | 0.0 | 2.1e-06 | 0.02 | 4.1e-06 | 0.03 | - | - | - | - | - | - | - | -
0.1 | 7.2e-06 | 0.04 | 3.3e-05 | 0.09 | 1.7e-02 | 3.00 | 1.7e-02 | 3.02 | 5.3e-07 | 2.25 | 1.7e-06 | 2.25
0.5 | 2.0e-03 | 0.54 | 5.1e-03 | 0.97 | 1.7e-02 | 3.02 | 1.7e-02 | 3.01 | 9.4e-05 | 2.25 | 8.4e-05 | 2.25
0.7 | 1.9e-02 | 1.57 | 2.9e-02 | 2.06 | 1.6e-02 | 2.96 | 1.6e-02 | 2.97 | 3.4e-03 | 2.25 | 6.4e-04 | 2.24
0.9 | 1.3e-01 | 3.52 | 1.5e-01 | 3.78 | 1.7e-02 | 3.01 | 1.7e-02 | 2.99 | 1.2e-01 | 2.26 | 1.6e-02 | 2.24

NonInv - Linear NonInv - SIREN NonInv - SIREN + TV

Time series value

Time step Time step Time step

Figure 4: Comparison of imputation methods, zoomed in 250 time steps for the
NonInv dataset, with $70\%$ of missing values.

###  4.3 Time Series Generation

To evaluate the utility of learning a prior over the space of implicit
functions, we use the set encoder network and the hypernetwork to generate new
time series. We do so by projecting time series into the latent vector of the
HyperTime network and interpolating the latent vector. This is similar to
training an autoencoder and interpolating the latent space, but the output of
the decoder of HyperTime are the weights of the SIREN networks. We follow the
experimental set up proposed in [1] for the evaluation, were the performance
of the synthetic data is evaluated using a predictive score (MAE) that
corresponds to the prediction accuracy of an off-the-shelf neural network
trained on the synthetic data and tested on the real data. Additionally, to
measure the quality of the synthetic data, we use the precision and recall
averaged over all time steps, which are then combined into a single F-score.
We use the same datasets as before, and we add two datasets that were used in
Fourier Flows [1] and TimeGAN [49], Google stocks data and UCI Energy data. We
compare our HyperTime model with generating data using PCA, with Fourier Flows
and TimeGAN, two state-of-the-art methods for time series generation. Table 3
shows the performance scores for all models and datasets. Additionally, we
visualize the generated samples using t-SNE plots 5 where we can see that the
generated data from HyperTime exhibits the same patterns as the original data.
In the case of Fourier Flows, in the UCR datasets we see that NonInv and
Phalanges do not show a good agreement.

The synthesis of time series via principal component analysis is performed in
a similar fashion as our HyperTime generation pipeline. We apply PCA to
generate a decomposition of time series into a basis of $40$ principal
components. The coefficients of these components constitute a latent
representation for each time series of the dataset, and we can interpolate
between embeddings of known time series to synthesize new ones. The main
limitation of this procedure, besides its linearity, is that it can only be
applied to datasets of equally sampled time series.

Table 3: Performance scores for data generated with HyperTime and for all baselines. | Crop | NonInv | Phalanges | Energy | Stock
---|---|---|---|---|---
PCA |  |  |  |  |
MAE | 0.050 | 0.019 | 0.050 | 0.007 | 0.110
F1 Score | 0.999 | 0.999 | 0.999 | 0.998 | 0.999
HyperTime (Ours) |  |  |  |  |
MAE | 0.040 | 0.005 | 0.026 | 0.058 | 0.013
F1 Score | 0.999 | 0.996 | 0.998 | 0.999 | 0.995
TimeGAN |  |  |  |  |
MAE | 0.048 | 0.028 | 0.108 | 0.056 | 0.173
F1 Score | 0.831 | 0.914 | 0.960 | 0.479 | 0.938
Fourier Flows |  |  |  |  |
MAE | 0.040 | 0.018 | 0.056 | 0.029 | 0.008
F1 Score | 0.991 | 0.990 | 0.992 | 0.945 | 0.992
Figure 5: t-SNE visualization on univariate datasets (in rows: Stocks, Energy,
Crop, NonInv and Phalanges), using different time series generation methods
(in columns: HyperTime, PCA, Fourier Flows and TimeGAN). Blue corresponds to
original data and orange to synthetic data.

Figure 6: _Left:_ t-SNE visualization of ground truth and generated data on
two univariate datasets (NonInv and FordA), using HyperNet with and without
the Fourier-based loss $\mathcal{L}_{\text{FFT}}$ (Eq. 3). _Right:_ Standard
deviation of the power spectra for the time series of the same two datasets.
FordA shows a considerably larger number of variations in the distributions of
the power spectra, which explains the difficulty of HyperTime to learn
patterns from the data.

Finally, we analyze the importance of the Fourier-based loss
$\mathcal{L}_{\text{FFT}}$ from equation 3 on the training of HyperTime. In
Figure 6-left we display t-SNE visualizations of time series synthesized by
HyperTime with and without the use of the FFT loss during training, for two
datasets (NonInv and FordA). In both cases, the addition of the
$\mathcal{L}_{\text{FFT}}$ loss results in an improved matching between ground
truth and generated data. However, in the case of FordA, the addition of this
loss becomes crucial to guide the learning process. This is also reflected in
the numerical evaluations from Table 4, which shows steep improvements in
performance for the FordA dataset.

A likely explanation for the difficulty of the network to learn meaningful
patterns from the data of this dataset is provided by the right plot in Figure
6. Here we show the standard deviation of the power spectrum for both
datasets, as a function of the frequency. The difference in the distributions
indicates that FordA is composed of spectra that present larger variability,
while NonInv’s spectra are considerably more clustered. Further research on
the characteristics of the datasets that benefit the more from the
$\mathcal{L}_{\text{FFT}}$ loss should be further investigated, especially
focusing on non-stationary time series.

Table 4: Performance scores for data generated with HyperTime, with and without the Fourier-based loss $\mathcal{L}_{\text{FFT}}$, for two datasets (NonInvasiveFetalECGThorax1, FordA). | NonInv | FordA
---|---|---
HyperTime + FFT loss |  |
MAE | 0.0053 | 0.0076
F1 Score | 0.9962 | 0.9987
HyperTime (no FFT) |  |
MAE | 0.0058 | 0.1647
F1 Score | 0.9960 | 0.0167

##  5 Conclusions

In this paper we explored the use of implicit neural representations for the
encoding and analysis of both univariate and multivariate time series data. We
compared multiple activation functions, and showed that periodic activation
layers outperform traditional activations in terms of reconstruction accuracy
and training speed. Additionally, we showed that INRs can be leveraged for
data imputation, resulting in good reconstructions of both the original data
and its power spectrum, when compared with classical imputation methods.
Finally, we presented HyperTime, a hypernetwork architecture to generate
synthetic data which enforces not only learning an accurate reconstruction
over the learned space of time series, but also preserving the shapes of the
power distributions. For this purpose, we introduced a new Fourier-based loss
on the training, and showed that for some datasets it becomes crucial to
enable the network to find meaningful patterns in the data. We leveraged the
latent representations learned by the hypernetwork for the generation of new
time-series data, and compared favorably against current state-of-the-art
methods for time-series augmentation. Besides reconstruction, imputation and
synthesis, we believe that both INRs and HyperTime open the door for a large
number of potential applications in time-series, such as upsampling,
forecasting and privacy preservation.

#### Disclaimer

This paper was prepared for informational purposes in part by the Artificial
Intelligence Research group of JPMorgan Chase & Coȧnd its affiliates (“JP
Morgan”), and is not a product of the Research Department of JP Morgan. JP
Morgan makes no representation and warranty whatsoever and disclaims all
liability, for the completeness, accuracy or reliability of the information
contained herein. This document is not intended as investment research or
investment advice, or a recommendation, offer or solicitation for the purchase
or sale of any security, financial instrument, financial product or service,
or to be used in any way for evaluating the merits of participating in any
transaction, and shall not constitute a solicitation under any jurisdiction or
to any person, if such solicitation under such jurisdiction or to such person
would be unlawful.

## References

  * [1] Ahmed Alaa, Alex James Chan, and Mihaela van der Schaar.  Generative time-series modeling with fourier flows.  In International Conference on Learning Representations, 2021.
  * [2] Ali Alqahtani, Mohammed Ali, Xianghua Xie, and Mark W. Jones.  Deep time-series clustering: A review.  Electronics, 10(23), 2021.
  * [3] Samuel A. Assefa, Danial Dervovic, Mahmoud Mahfouz, Robert E. Tillman, Prashant Reddy, and Manuela Veloso.  Generating synthetic data in finance: Opportunities, challenges and pitfalls.  In Proceedings of the First ACM International Conference on AI in Finance, ICAIF ’20, New York, NY, USA, 2020. Association for Computing Machinery.
  * [4] A. Bagnall, J. Lines, A. Bostrom, J. Large, and E. Keogh.  The great time series classification bake off: a review and experimental evaluation of recent algorithmic advances.  Data Mining and Knowledge Discovery, 31:606–660, 2017.
  * [5] Anthony Bagnall, Jason Lines, William Vickers, and Eamonn Keogh.  The uea $\&$ ucr time series classification repository.  www.timeseriesclassification.com.  Accessed: 2022-05-10.
  * [6] Dimitrios Bellos, Mark Basham, Tony P. Pridmore, and Andrew P. French.  A convolutional neural network for fast upsampling of undersampled tomograms in x-ray ct time-series using a representative highly sampled tomogram.  Journal of Synchrotron Radiation, 26:839 – 853, 2019.
  * [7] Nuri Benbarka, Timon Höfer, Hamd ul Moqeet Riaz, and Andreas Zell.  Seeing implicit neural representations as fourier series.  2022 IEEE/CVF Winter Conference on Applications of Computer Vision (WACV), pages 2283–2292, 2022.
  * [8] Andrew Brock, Theo Lim, J.M. Ritchie, and Nick Weston.  SMASH: One-shot model architecture search through hypernetworks.  In International Conference on Learning Representations, 2018.
  * [9] Wei Cao, Dong Wang, Jian Li, Hao Zhou, Lei Li, and Yitan Li.  Brits: Bidirectional recurrent imputation for time series.  In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
  * [10] Yinbo Chen, Sifei Liu, and Xiaolong Wang.  Learning continuous image representation with local implicit image function.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 8628–8638, 2021.
  * [11] Kukjin Choi, Jihun Yi, Changhwa Park, and Sungroh Yoon.  Deep learning for anomaly detection in time-series data: Review, analysis, and guidelines.  IEEE Access, 9:120043–120065, 2021.
  * [12] Andrea Coletta, Matteo Prata, Michele Conti, Emanuele Mercanti, Novella Bartolini, Aymeric Moulin, Svitlana Vyetrenko, and Tucker Balch.  Towards realistic market simulations: A generative adversarial networks approach.  In Proceedings of the Second ACM International Conference on AI in Finance, ICAIF ’21, New York, NY, USA, 2021. Association for Computing Machinery.
  * [13] Abhyuday Desai, Cynthia Freeman, Zuhui Wang, and Ian Beaver.  Timevae: A variational auto-encoder for multivariate time series generation.  arXiv preprint arXiv:2111.08095, 2021.
  * [14] Chenguang Fang and Chen Wang.  Time series data imputation: A survey on deep learning approaches.  ArXiv, abs/2011.11347, 2020.
  * [15] Elizabeth Fons, Paula Dawson, Xiao jun Zeng, John Keane, and Alexandros Iosifidis.  Adaptive weighting scheme for automatic time-series data augmentation, 2021.
  * [16] David Ha, Andrew M. Dai, and Quoc V. Le.  Hypernetworks.  In 5th International Conference on Learning Representations, ICLR 2017, Toulon, France, April 24-26, 2017, Conference Track Proceedings. OpenReview.net, 2017.
  * [17] Kyle Hundman, Valentino Constantinou, Christopher Laporte, Ian Colwell, and Tom Soderstrom.  Detecting spacecraft anomalies using lstms and nonparametric dynamic thresholding.  In Proceedings of the 24th ACM SIGKDD International Conference on Knowledge Discovery $\&$; Data Mining, KDD ’18, page 387–395, New York, NY, USA, 2018. Association for Computing Machinery.
  * [18] Hassan Ismail Fawaz, Benjamin Lucas, Germain Forestier, Charlotte Pelletier, Daniel F. Schmidt, Jonathan Weber, Geoffrey I. Webb, Lhassane Idoumghar, Pierre-Alain Muller, and François Petitjean.  Inceptiontime: Finding alexnet for time series classification.  Data Mining and Knowledge Discovery, 2020.
  * [19] Brian Kenji Iwana and Seiichi Uchida.  An empirical survey of data augmentation for time series classification with neural networks.  PLOS ONE, 16(7):e0254841, Jul 2021.
  * [20] James Jordon, Daniel Jarrett, Evgeny Saveliev, Jinsung Yoon, Paul Elbers, Patrick Thoral, Ari Ercole, Cheng Zhang, Danielle Belgrave, and Mihaela van der Schaar.  Hide-and-seek privacy challenge: Synthetic data generation vs. patient re-identification.  In Hugo Jair Escalante and Katja Hofmann, editors, Proceedings of the NeurIPS 2020 Competition and Demonstration Track, volume 133 of Proceedings of Machine Learning Research, pages 206–215. PMLR, 06–12 Dec 2021\.
  * [21] James Jordon, Jinsung Yoon, and Mihaela van der Schaar.  Pate-gan: Generating synthetic data with differential privacy guarantees.  In ICLR, 2019.
  * [22] Bryan Lim and Stefan Zohren.  Time-series forecasting with deep learning: a survey.  Phylosophical Transactions of the Royal Society A, 2021.
  * [23] Gidi Littwin and Lior Wolf.  Deep meta functionals for shape representation.  In IEEE/CVF International Conference on Computer Vision (ICCV), pages 1824–1833, 10 2019.
  * [24] Yan Liu.  Recurrent neural networks for multivariate time series with missing values.  Scientific Reports, 8(1):6085, 2018.
  * [25] Yonghong Luo, Xiangrui Cai, Ying ZHANG, Jun Xu, and Yuan xiaojie.  Multivariate time series imputation with generative adversarial networks.  In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
  * [26] Marius Lysaker.  Iterative image restoration combining total variation minimization and a second-order functional.  International Journal of Computer Vision, 66:5–18, 01 2006.
  * [27] Qianli Ma, Jiawei Zheng, Sen Li, and Gary W Cottrell.  Learning representations for time series clustering.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
  * [28] Lars Mescheder, Michael Oechsle, Michael Niemeyer, Sebastian Nowozin, and Andreas Geiger.  Occupancy networks: Learning 3d reconstruction in function space.  In Proceedings IEEE Conf. on Computer Vision and Pattern Recognition (CVPR), 2019.
  * [29] Ben Mildenhall, Pratul P. Srinivasan, Matthew Tancik, Jonathan T. Barron, Ravi Ramamoorthi, and Ren Ng.  Nerf: Representing scenes as neural radiance fields for view synthesis.  In ECCV, 2020.
  * [30] Cheolhwan Oh, Seungmin Han, and Jongpil Jeong.  Time-series data augmentation based on interpolation.  Procedia Computer Science, 175:64–71, 2020.  The 17th International Conference on Mobile Systems and Pervasive Computing (MobiSPC),The 15th International Conference on Future Networks and Communications (FNC),The 10th International Conference on Sustainable Energy Information Technology.
  * [31] Jeong Joon Park, Peter Florence, Julian Straub, Richard Newcombe, and Steven Lovegrove.  Deepsdf: Learning continuous signed distance functions for shape representation.  In The IEEE Conference on Computer Vision and Pattern Recognition (CVPR), June 2019.
  * [32] Mathias Perslev, Michael Jensen, Sune Darkner, Poul Jø rgen Jennum, and Christian Igel.  U-time: A fully convolutional network for time series segmentation applied to sleep staging.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
  * [33] Tamar Rott Shaham, Michael Gharbi, Richard Zhang, Eli Shechtman, and Tomer Michaeli.  Spatially-adaptive pixelwise networks for fast image translation.  In Computer Vision and Pattern Recognition (CVPR), 2021.
  * [34] Jacob Russin Russin, Randall O’Reilly, and Yoshua Bengio Bengio.  Deep learning needs a prefrontal cortex.  In Bridging AI and Cognitive Science ICLR 2020 Workshop, 2020.
  * [35] Andrei A. Rusu, Dushyant Rao, Jakub Sygnowski, Oriol Vinyals, Razvan Pascanu, Simon Osindero, and Raia Hadsell.  Meta-learning with latent embedding optimization.  In International Conference on Learning Representations, 2019.
  * [36] Vincent Sitzmann, Eric R. Chan, Richard Tucker, Noah Snavely, and Gordon Wetzstein.  Metasdf: Meta-learning signed distance functions.  In Proc. NeurIPS, 2020.
  * [37] Vincent Sitzmann, Julien N.P. Martel, Alexander W. Bergman, David B. Lindell, and Gordon Wetzstein.  Implicit neural representations with periodic activation functions.  In Proc. NeurIPS, 2020.
  * [38] Vincent Sitzmann, Michael Zollhöfer, and Gordon Wetzstein.  Scene representation networks: Continuous 3d-structure-aware neural scene representations.  In Advances in Neural Information Processing Systems, 2019.
  * [39] Ivan Skorokhodov, Savva Ignatyev, and Mohamed Elhoseiny.  Adversarial generation of continuous images.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), pages 10753–10764, June 2021.
  * [40] Kenneth O. Stanley, David B. D’Ambrosio, and Jason Gauci.  A Hypercube-Based Encoding for Evolving Large-Scale Neural Networks.  Artificial Life, 15(2):185–212, 04 2009.
  * [41] Alejandro Sztrajman, Gilles Rainer, Tobias Ritschel, and Tim Weyrich.  Neural brdf representation and importance sampling.  Computer Graphics Forum, 40(6):332–346, 2021.
  * [42] Matthew Tancik, Ben Mildenhall, Terrance Wang, Divi Schmidt, Pratul P. Srinivasan, Jonathan T. Barron, and Ren Ng.  Learned initializations for optimizing coordinate-based neural representations.  In CVPR, 2021.
  * [43] Matthew Tancik, Pratul P. Srinivasan, Ben Mildenhall, Sara Fridovich-Keil, Nithin Raghavan, Utkarsh Singhal, Ravi Ramamoorthi, Jonathan T. Barron, and Ren Ng.  Fourier features let networks learn high frequency functions in low dimensional domains.  NeurIPS, 2020.
  * [44] José F. Torres, Dalil Hadjout, Abderrazak Sebaa, Francisco Martínez-Álvarez, and Alicia Troncoso Lora.  Deep learning for time series forecasting: A survey.  Big data, 2021.
  * [45] Johannes von Oswald, Christian Henning, Benjamin F. Grewe, and João Sacramento.  Continual learning with hypernetworks.  In International Conference on Learning Representations, 2020.
  * [46] Magnus Wiese, Robert Knobloch, Ralf Korn, and Peter Kretschmer.  Quant gans: deep generation of financial time series.  Quantitative Finance, page 1–22, Apr 2020.
  * [47] Haowen Xu, Wenxiao Chen, Nengwen Zhao, Zeyan Li, Jiahao Bu, Zhihan Li, Ying Liu, Youjian Zhao, Dan Pei, Yang Feng, Jie Chen, Zhaogang Wang, and Honglin Qiao.  Unsupervised anomaly detection via variational auto-encoder for seasonal kpis in web applications.  In Proceedings of the 2018 World Wide Web Conference, WWW ’18, page 187–196, 2018.
  * [48] Jinsung Yoon, Daniel Jarrett, and Mihaela van der Schaar.  Time-series generative adversarial networks.  In NeurIPS, 2019.
  * [49] Jinsung Yoon, Daniel Jarrett, and Mihaela van der Schaar.  Time-series generative adversarial networks.  In NeurIPS, 2019.
  * [50] Li Zeng, Baifan Zhou, Mohammad Al-Rifai, and Evgeny Kharlamov.  Segtime: Precise time series segmentation without sliding window, 2022\.
  * [51] Chris Zhang, Mengye Ren, and Raquel Urtasun.  Graph hypernetworks for neural architecture search.  In International Conference on Learning Representations, 2019.
  * [52] Dominic Zhao, Seijin Kobayashi, João Sacramento, and Johannes von Oswald.  Meta-learning via hypernetworks.  In 4th Workshop on Meta-Learning at NeurIPS 2020 (MetaLearn 2020). NeurIPS, 2020.


The NA64 Collaboration

#  EUROPEAN LABORATORY FOR PARTICLE PHYSICS
CERN-EP-2021-017

Constraints on New Physics in the Electron $g-2$ from a Search for Invisible
Decays
of a Scalar, Pseudoscalar, Vector, and Axial Vector

Yu. M. Andreev  Institute for Nuclear Research, 117312 Moscow, Russia D.
Banerjee  CERN, European Organization for Nuclear Research, CH-1211 Geneva 23,
Switzerland J. Bernhard  CERN, European Organization for Nuclear Research,
CH-1211 Geneva 23, Switzerland V. E. Burtsev  Joint Institute for Nuclear
Research, 141980 Dubna, Russia A. G. Chumakov  Tomsk Polytechnic University,
634050 Tomsk, Russia Tomsk State Pedagogical University, 634061 Tomsk, Russia
D. Cooke  UCL Departement of Physics and Astronomy, University College London,
Gower St. London WC1E 6BT, United Kingdom P. Crivelli  ETH Zürich, Institute
for Particle Physics and Astrophysics, CH-8093 Zürich, Switzerland E. Depero
ETH Zürich, Institute for Particle Physics and Astrophysics, CH-8093 Zürich,
Switzerland A. V. Dermenev  Institute for Nuclear Research, 117312 Moscow,
Russia S. V. Donskov  State Scientific Center of the Russian Federation
Institute for High Energy Physics of National Research Center ’Kurchatov
Institute’ (IHEP), 142281 Protvino, Russia R. R. Dusaev  Tomsk Polytechnic
University, 634050 Tomsk, Russia T. Enik  Joint Institute for Nuclear
Research, 141980 Dubna, Russia N. Charitonidis  CERN, European Organization
for Nuclear Research, CH-1211 Geneva 23, Switzerland A. Feshchenko  Joint
Institute for Nuclear Research, 141980 Dubna, Russia V. N. Frolov  Joint
Institute for Nuclear Research, 141980 Dubna, Russia A. Gardikiotis  Physics
Department, University of Patras, 265 04 Patras, Greece S. G. Gerassimov
Technische Universität München, Physik Department, 85748 Garching, Germany
P.N. Lebedev Physical Institute, Moscow, Russia, 119 991 Moscow, Russia S. N.
Gninenko111Corresponding author<EMAIL_ADDRESS>Institute for Nuclear
Research, 117312 Moscow, Russia M. Hösgen  Universität Bonn, Helmholtz-
Institut für Strahlen-und Kernphysik, 53115 Bonn, Germany V. A. Kachanov
State Scientific Center of the Russian Federation Institute for High Energy
Physics of National Research Center ’Kurchatov Institute’ (IHEP), 142281
Protvino, Russia A. E. Karneyeu  Institute for Nuclear Research, 117312
Moscow, Russia G. Kekelidze  Joint Institute for Nuclear Research, 141980
Dubna, Russia B. Ketzer  Universität Bonn, Helmholtz-Institut für Strahlen-und
Kernphysik, 53115 Bonn, Germany D. V. Kirpichnikov  Institute for Nuclear
Research, 117312 Moscow, Russia M. M. Kirsanov  Institute for Nuclear
Research, 117312 Moscow, Russia V. N. Kolosov  State Scientific Center of the
Russian Federation Institute for High Energy Physics of National Research
Center ’Kurchatov Institute’ (IHEP), 142281 Protvino, Russia I. V. Konorov
Technische Universität München, Physik Department, 85748 Garching, Germany
P.N. Lebedev Physical Institute, Moscow, Russia, 119 991 Moscow, Russia S. G.
Kovalenko  Departamento de Ciencias Físicas, Universidad Andres Bello, Sazié
2212, Piso 7, Santiago, Chile V. A. Kramarenko  Joint Institute for Nuclear
Research, 141980 Dubna, Russia Skobeltsyn Institute of Nuclear Physics,
Lomonosov Moscow State University, 119991 Moscow, Russia L. V. Kravchuk
Institute for Nuclear Research, 117312 Moscow, Russia N. V. Krasnikov  Joint
Institute for Nuclear Research, 141980 Dubna, Russia Institute for Nuclear
Research, 117312 Moscow, Russia S. V. Kuleshov  Departamento de Ciencias
Físicas, Universidad Andres Bello, Sazié 2212, Piso 7, Santiago, Chile SAPHIR
Millennium Institute of ANID, Chile V. E. Lyubovitskij  Tomsk Polytechnic
University, 634050 Tomsk, Russia Tomsk State Pedagogical University, 634061
Tomsk, Russia Universidad Técnica Federico Santa María, 2390123 Valparaíso,
Chile V. Lysan  Joint Institute for Nuclear Research, 141980 Dubna, Russia V.
A. Matveev  Joint Institute for Nuclear Research, 141980 Dubna, Russia Yu. V.
Mikhailov  State Scientific Center of the Russian Federation Institute for
High Energy Physics of National Research Center ’Kurchatov Institute’ (IHEP),
142281 Protvino, Russia L. Molina Bueno  ETH Zürich, Institute for Particle
Physics and Astrophysics, CH-8093 Zürich, Switzerland D. V. Peshekhonov  Joint
Institute for Nuclear Research, 141980 Dubna, Russia V. A. Polyakov  State
Scientific Center of the Russian Federation Institute for High Energy Physics
of National Research Center ’Kurchatov Institute’ (IHEP), 142281 Protvino,
Russia B. Radics  ETH Zürich, Institute for Particle Physics and Astrophysics,
CH-8093 Zürich, Switzerland R. Rojas  Universidad Técnica Federico Santa
María, 2390123 Valparaíso, Chile A. Rubbia  ETH Zürich, Institute for Particle
Physics and Astrophysics, CH-8093 Zürich, Switzerland V. D. Samoylenko  State
Scientific Center of the Russian Federation Institute for High Energy Physics
of National Research Center ’Kurchatov Institute’ (IHEP), 142281 Protvino,
Russia H. Sieber  ETH Zürich, Institute for Particle Physics and Astrophysics,
CH-8093 Zürich, Switzerland D. Shchukin  P.N. Lebedev Physical Institute,
Moscow, Russia, 119 991 Moscow, Russia V. O. Tikhomirov  P.N. Lebedev Physical
Institute, Moscow, Russia, 119 991 Moscow, Russia I. Tlisova  Institute for
Nuclear Research, 117312 Moscow, Russia A. N. Toropin  Institute for Nuclear
Research, 117312 Moscow, Russia A. Yu. Trifonov  Tomsk Polytechnic University,
634050 Tomsk, Russia Tomsk State Pedagogical University, 634061 Tomsk, Russia
B. I. Vasilishin  Tomsk Polytechnic University, 634050 Tomsk, Russia P. V.
Volkov  Joint Institute for Nuclear Research, 141980 Dubna, Russia Skobeltsyn
Institute of Nuclear Physics, Lomonosov Moscow State University, 119991
Moscow, Russia V. Yu. Volkov  Skobeltsyn Institute of Nuclear Physics,
Lomonosov Moscow State University, 119991 Moscow, Russia

###### Abstract

We performed a search for a new generic $X$ boson, which could be a scalar
($S$), pseudoscalar ($P$), vector ($V$) or an axial vector ($A$) particle
produced in the 100 GeV electron scattering off nuclei, $e^{-}Z\to e^{-}ZX$,
followed by its invisible decay in the NA64 experiment at CERN. No evidence
for such process was found in the full NA64 data set of $2.84\times 10^{11}$
electrons on target. We place new bounds on the $S,P,V,A$ coupling strengths
to electrons, and set constraints on their contributions to the electron
anomalous magnetic moment $a_{e}$, $|\Delta a_{X}|\lesssim 10^{-15}-10^{-13}$
for the $X$ mass region $m_{X}\lesssim 1$ GeV. These results are an order of
magnitude more sensitive compared to the current accuracy on $a_{e}$ from the
electron $g-2$ experiments and recent high-precision determination of the fine
structure constant.

Searching for new physics (NP) with mass below the electroweak scale ($\ll
100$ GeV) at the high-intensity and high-precision frontiers has received
significant attention in recent years jr ; Essig:2013lka ; report1 ; report2 ;
pbc-bsm ; esu ; berlin ; gaia . Motivations for searches of feebly-coupled
particles in the low-mass range come from the evidence for NP in the neutrino
and dark matter sectors, and are well supported by theoretical arguments, see,
e.g. Refs.jr ; berlin ; gaia ; bf ; pf ; prv ; mp ; arkani . Existing
anomalies observed in particle experiments also contribute to the field. Well-
known examples are the current muon $g-2$ anomaly - the $\simeq 3.6\sigma$
discrepancy between the predicted and observed value of the muon anomalous
magnetic moment g-2_mu , or the X17 anomaly - an excess of $e^{+}e^{-}$ events
in the 8Be and 4He nuclei transitions atomki1 ; atomki2 , which might be
explained by NP models at low-mass scale, see, e.g. Refs.gk ; feng . These
anomalies are being scrutinized in the upcoming experiments at Fermilab and
JPARC fnal ; jpark , and with NA64 at CERN na64be1 ; na64be2 ; na64be3 ,
respectively.

Recently, a new puzzle indicating the possible presence of NP in the electron
$g-2$ has emerged. The precise measurements performed at Laboratoire Kastler
Brossel (LKB) with 87Rb rubidium atoms report a new value for the fine-
structure constant $\alpha^{-1}$ = 137.035999206(11) with a relative accuracy
of 81 parts per trillion lkb . This result improves the accuracy on $\alpha$
by 2.5 over the previous measurements performed at Berkeley with 137Cs atoms
berkeley  but, surprisingly, it reveals a $5.4\sigma$ difference from this
latest result. Using these measurements of the fine-structure constant, the
Standard Model (SM) prediction of the anomalous magnetic moment of the
electron, $a_{e}=(g-2)_{e}/2$ aoyama1 ; aoyama2 , is $1.6~{}\sigma$ lower and
$-2.4~{}\sigma$ higher than the direct experimental measurement of
$a_{e}^{exp}$ gabri :

$\displaystyle\Delta a_{e}=a^{exp}_{e}-a^{LKB}_{e}=(4.8\pm 3.0)\times
10^{-13}$ (1) $\displaystyle\Delta a_{e}=a^{exp}_{e}-a^{B}_{e}=(-8.8\pm
3.6)\times 10^{-13}$ (2)

for the LKB and Berkeley measurements, respectively. The errors on $\Delta
a_{e}$ are dominated mostly by the uncertainty in $a^{exp}_{e}$. As the SM
predicts a certain value of the $a_{e}$ aoyama1 ; aoyama2  the measurements of
this parameter in different processes should be consistent with each other.
With new measurements and improved SM calculations, one hopes to clarify
whether the deviations of Eqs.(1,2) are a result of yet unknown experimental
errors, or it is a sign of new physics in the electron $g-2$ giud . This
motivates recent significant efforts towards possible explanation of the
deviation, in particular the discrepancy of Eq.(2), with a NP effect, see,
e.g., Refs.nvk -gkm .

In this Letter, we study the question of whether a new light $X$ boson could
contribute to the electron $g-2$. We consider models with a generic $X$ in
sub-GeV mass range, which could be a scalar ($S$), pseudoscalar ($P$), vector
($V$), or an axial vector ($A$) particle feebly coupled to electrons. It is
assumed that the $X$ decays predominantly invisibly, $\Gamma(X\rightarrow
invisible)/\Gamma_{tot}\simeq 1$, e.g. into dark sector particles, thus
escaping stringent constraints placed today on the visible decay modes of the
$X$ into SM particles from collider, fixed-target, and atomic experiments pdg
. The most stringent limits on the invisible $X$ in the sub-GeV mass range are
obtained, so far, for the $V$ case of dark photons coupled to electrons
through the mixing with the ordinary photons by the NA64 na64prl19  and
$BABAR$ babarg-2  experiments, leaving a large area of the parameter space for
the generic $X$ still unexplored. Various aspects of such invisible $X$ weakly
coupled to leptons including possible phenomenological implications can be
found in Refs.jr ; Essig:2013lka ; report1 ; report2 ; pbc-bsm ; esu ; berlin
; gaia ; gkm ; Gninenko:2018ter ; Kirpichnikov:2020tcf .

The $e-X$-interaction with the coupling strength $g_{X}$ defined as
$g_{X}=\varepsilon_{X}e$ (here $\varepsilon_{X}$ is a parameter and $e$ is the
charge of the electron) is given for the $S,P,V,A$ cases by phenomenological
Lagrangians:

$\displaystyle\mathcal{L}_{S}$ $\displaystyle=$ $\displaystyle
g_{S}\overline{e}eS$ $\displaystyle\mathcal{L}_{P}$ $\displaystyle=$
$\displaystyle ig_{P}\overline{e}\gamma_{5}eP$ $\displaystyle\mathcal{L}_{V}$
$\displaystyle=$ $\displaystyle g_{V}\overline{e}\gamma_{\mu}eV_{\mu}$
$\displaystyle\mathcal{L}_{A}$ $\displaystyle=$ $\displaystyle
g_{A}\overline{e}\gamma_{\mu}\gamma_{5}eA_{\mu}$ (3)

Figure 1: One-loop contribution of the $S$ and $P$ (left panel) and the $V$
and $A$ (right panel) to $\Delta a_{e}$.

The corresponding one-loop contributions to the $(g-2)_{e}$ factor induced by
diagrams shown in Fig. 1 are given by:

$\displaystyle\Delta a_{S}$ $\displaystyle=$
$\displaystyle\frac{g_{S}^{2}}{4\pi^{2}}\bigl{(}\frac{m_{e}}{m_{X}}\bigr{)}^{2}\bigl{[}{\rm
ln}\frac{m_{X}}{m_{e}}-\frac{7}{12}\bigr{]}$ (4) $\displaystyle\Delta a_{P}$
$\displaystyle=$
$\displaystyle\frac{g_{P}^{2}}{4\pi^{2}}\bigl{(}\frac{m_{e}}{m_{X}}\bigr{)}^{2}\bigl{[}-{\rm
ln}\frac{m_{X}}{m_{e}}+\frac{11}{12}\bigr{]}$ (5) $\displaystyle\Delta a_{V}$
$\displaystyle=$
$\displaystyle\frac{g_{V}^{2}}{4\pi^{2}}\bigl{(}\frac{m_{e}}{m_{X}}\bigr{)}^{2}\frac{1}{3}$
(6) $\displaystyle\Delta a_{A}$ $\displaystyle=$
$\displaystyle\frac{g_{A}^{2}}{4\pi^{2}}\bigl{(}\frac{m_{e}}{m_{X}}\bigr{)}^{2}\bigl{(}-\frac{5}{3}\bigr{)}$
(7)

assuming that $m_{X}\gg m_{e}$. One can see that presumably a scalar and a
vector can explain the positive deviation of Eq.((1)), while only a
pseudoscalar and an axial vector could explain the negative value of Eq.((2)).
The required couplings $g_{X}$ to explain deviations of Eqs.(1,2) are in the
range $10^{-3}\lesssim|g_{X}|\lesssim 10^{-4}$ which is accessible to the NA64
search, thus making it interesting.

The method of the search, discussed in this work and proposed in Refs.
Gninenko:2013rka ; Andreas:2013lya , is based on the detection of the missing
energy, carried away by the hard bremsstrahlung $X$ produced in the process
$e^{-}Z\to e^{-}ZX;~{}X\to invisible$ of high-energy electrons scattering in
an active beam dump. The NA64 experiment employed a 100 GeV pure electron
beam, using the H4 beam-line of the CERN’s North Area. The beam was slowly
extracted towards NA64 in 4.8 s spills, and had an intensity up to $\simeq
10^{7}$ electrons per spill. The e- beam was defined by the scintillator
(${\it S}$) and veto (${\it V}_{1}$) counters. A magnetic spectrometer
consisting of two successive dipole magnets with the integral magnetic
strength of $\simeq$7 T$\cdot$m and a low-material-budget tracker consisting
of a set of Micromegas (MM), Straw-Tube (ST) and Gaseous Electron Multiplier
(GEM) chambers allowed to measure the incoming $e^{-}$ momenta with the
precision $\delta p/p\simeq 1\%$ Banerjee:2015eno . The synchrotron radiation
(SR) emitted in the magnets was used for the electron identification and their
efficient tagging with a SR detector (SRD)na64srd , which was an array of a
Pb-Sc sandwich calorimeter of a fine segmentation. By using the SRD the
intrinsic hadron contamination of the beam of the order of $\sim 1\%$ was
further suppressed to a negligible level. The downstream part of the detector
was equipped with an electromagnetic (e-m) calorimeter (ECAL), a matrix of
$6\times 6$ Shashlik-type modules assembled from Pb and Sc plates serving as
an active beam-dump target for measurement of the electron energy $E_{{\rm
ECAL}}$. Each ECAL module has $\simeq 40$ radiation lengths ($X_{0}$) with the
first 4$X_{0}$ serving as a preshower detector (PS). Further downstream the
detector was equipped with a high-efficiency veto counter (${\it V}_{2}$), and
a hermetic hadronic calorimeter (HCAL) of $\simeq 30$ nuclear interaction
lengths in total. The HCAL was used as an efficient veto against hadronic
secondaries and also to detect muons produced in $e^{-}$ interactions in the
target.

The search described in this paper uses the data samples of
$n_{EOT}=2.84\times 10^{11}$ electrons on target (EOT), collected in the years
2016, 2017 and 2018 (runs I,II, and III, respectively) at the beam intensities
mostly in the range $\simeq(5-9)\times 10^{6}$ e- per spill with the hardware
trigger na64prl17 ; na64prd18 ; na64prl19

$Tr(X)=\Pi{\it S}_{i}\cdot{\it V}_{1}\cdot{\rm PS}(>E^{th}_{{\rm
PS}})\cdot\overline{{\rm ECAL}}(<E^{th}_{{\rm ECAL}}),$ (8)

accepting events with in-time hits in beam-defining counters $S_{i}$ and
clusters in the PS and ECAL with the energy exceeding the thresholds
$E^{th}_{{\rm PS}}\simeq 0.3$ GeV and $E^{th}_{{\rm ECAL}}\lesssim 80$ GeV,
respectively. The missing energy events have the signature

$S(X)=Tr(X)\cdot{\rm Track}(P_{e})\cdot{\it V}_{2}(<E^{th}_{{\it
V}_{2}})\cdot{\rm HCAL}(<E^{th}_{{\rm HCAL}})$ (9)

with the incoming track momentum $P_{e}\simeq 100\pm 3$ GeV, and $V_{2}$ and
HCAL zero-energy deposition, defined as energy below the thresholds
$E^{th}_{{\it V}_{2}}\simeq 1$ MIP (minimum ionizing particle) and
$E^{th}_{{\rm HCAL}}\simeq 1$ GeV, respectively. Data from these three runs,
were processed with selection criteria similar to the one used in Refs.
na64prd18 ; na64prl19  and finally analysed as described below.

A detailed Geant4 Agostinelli:2002hh ; geant  based Monte Carlo (MC)
simulation was used to study detector performance and signal acceptance, to
simulate backgrounds and selection cuts. For the calculations of the signal
yield we used the fully Geant4 compatible package DMG4 dmg4 . Using this
package the production of $X$ in the process $e^{-}Z\to e^{-}ZX;~{}X\to
invisible$ has been simulated for each type of interactions listed in Eq.(3)
with cross-sections obtained from exact tree-level (ETL) calculations, see,
e.g., Refs. gkkk ; gkkketl ; dkk . The produced signal samples were processed
by the same reconstruction program as the real data and passed through the
same selection criteria. The total number ${n}_{X}$ of the produced $X$ per
single electron on target (EOT) was calculated as

${n}_{X}(g_{X},~{}m_{X},~{}E_{0})=\frac{\rho N_{A}}{A_{{\rm
Pb}}}\sum_{i}n(E_{0},E_{e},s)\sigma_{X}(E_{e})\Delta s_{i}$ (10)

where $\rho$ is density of the target, $N_{A}$ is the Avogadro’s number,
$A_{{\rm Pb}}$ is the Pb atomic mass, $n(E_{0},E_{e},s)$ is the number of
$e^{\pm}$ in the e-m shower at the depth $s$ (in radiation lengths) with the
energy $E_{e}$ within the target of total thickness $T$, and $\sigma(E_{e})$
is the cross section of the $X$ production in the kinematically allowed region
up to $E_{X}\simeq E_{e}$ by an electron with the energy $E_{e}$ in the
reaction $e^{-}Z\to e^{-}ZX;~{}X\to invisible$. The latter depends in
particular on the coupling and mass $g_{X},~{}m_{X}$, and the beam energy
$E_{0}$. The $X$ energy distribution $\frac{dn_{X}}{dE_{X}}$ was calculated
for each case by taking into account the corresponding differential cross-
section $\frac{d\sigma(E_{e},E_{X})}{dE_{X}}$, as described in Ref.gkkketl .
An example of the simulated X (or missing) energy spectrum in the target
calculated by using the detailed simulation of e-m shower development by
Geant4 is shown for the $P$ and $V$ cases in Fig. 2 for the mass $m_{X}=20$
MeV.

Figure 2: The emission spectra of the 20 MeV P(solid line) and V(dashed line)
particles produced from the interactions of the 100 GeV electron beam in the
ECAL target obtained from the ETL calculations. The spectra are normalized to
the same number of EOT.

.

The expected number of $X$ events in our detector from the reaction $e^{-}Z\to
e^{-}ZX;~{}X\to invisible$ was determined for each $X$ interaction type also
by comparison to the rare process of dimuon production, $e^{-}Z\to
e^{-}Z\gamma;\gamma\to\mu^{+}\mu^{-}$, which has a well-known reaction rate.
These events originate from the QED reaction in the ECAL, dominated by the
hard bremsstrahlung photon conversion into dimuon pairs on a target nucleus
and accompanied by small energy deposition in the HCAL, thus mimicking the
$X\to invisible$ decay events below the two-MIP threshold. The reaction was
also used as a benchmark process allowing us to verify the reliability of the
MC simulation, correct the signal acceptance, cross-check systematic
uncertainties and background estimate na64prd18 ; na64prl19 . Good agreement
was found between the observations and simulations. Using rare dimuon events
as a crosscheck for normalization to the signal modes cancels many systematic
uncertainties by keeping selection cuts identical whenever possible.

In order to avoid biases in the determination of the selection criteria for
signal events, a blind analysis similar to the one described in Ref.na64prl19
was performed. The signal box ($E_{{\rm ECAL}}<50~{}{\text{G}eV};E_{{\rm
HCAL}}<1~{}{\text{G}eV}$) was defined based on the energy spectrum
calculations for $X$s emitted by $e^{\pm}$ from the e-m shower generated by
the primary $e^{-}$s in the ECAL gkkk ; gkkketl  and the HCAL zero-energy
threshold determined mostly by the noise of the read-out electronics. Finally,
to maximize the acceptance for signal events and to minimize backgrounds we
used the following selection criteria: (i) The incoming electron track
momentum should be within $100\pm 3$ GeV; (ii) The SRD energy should be within
the range of the SR energy emitted by $e^{-}$s in the magnets and in time with
the trigger; (iii) The shower shape in the ECAL should be consistent with the
one expected for the signal shower gkkk ; (iv) There should be only a single
track activity in the tracker chambers upstream of the dump in order to reject
interactions in the beam line materials, and no activity in ${\it V}_{2}$.

The dominant background for $e^{-}Z\to e^{-}ZX;~{}X\to invisible$ arises from
the interactions of the $e^{-}$ beam in the downstream part of the detector
resulting in hadron electro-production in the beam line materials. In rare
cases, these reactions are accompanied by the emission of large-angle (high
$p_{T}$) hadronic secondaries faking the signal due to the insufficient
downstream detector coverage. Charged secondaries were rejected by requiring
no additional tracks or hits in the downstream ST chambers, which have the
largest transverse acceptance in our setup. We also requested no extra in-time
hits upstream of the magnets and at most one extra in-time hit downstream of
the magnets in the MM chambers. The remaining background from the large-angle
neutral hadronic secondaries was evaluated mainly from data by the
extrapolation of events from the sideband ($E_{{\rm
ECAL}}>50~{}{\text{G}eV};E_{{\rm HCAL}}<1~{}{\text{G}eV}$) into the signal
region and assessing the systematic errors by varying the fit functions
selected as described in Ref. na64prd18 . The shape of the extrapolation
functions was evaluated from the study of a larger data sample of events from
hadronic $e^{-}$ interactions in the dump, which was also cross-checked with
simulations. Another background from punch-through of leading (with energy
$\gtrsim 0.5~{}E_{0}$) neutral hadrons $(n,K^{0}_{L})$ produced in the $e^{-}$
interactions in the target, was studied by using events from the region
($E_{{\rm ECAL}}<50~{}{\text{G}eV};E_{{\rm HCAL}}>1~{}{\text{G}eV}$), which
were pure neutral hadronic secondaries produced in the ECAL. Its level was
estimated from the data by using the longitudinal segmentation of the HCAL and
the punch-through probability estimated conservatively and was found to be
negligible. Several other background sources that may fake the signal, such as
loss of dimuons due to statistical fluctuations of the signal or muon decays,
and decays in flight of mistakenly SRD tagged beam $\pi$, $K$ were simulated
with the full statistics of the data and also were found to be negligible.
After determining all the selection criteria and background levels, we
unblinded the signal region and found 0 events consistent with $0.53\pm 0.17$
events from the conservative background estimations na64prl19  allowing us to
obtain the $m_{X}$-dependent upper limits on the $e-X$ coupling strengths.

Figure 3: The 90% C.L. upper limits on the coupling parameter
$\varepsilon_{X}$ in the ($m_{X},\varepsilon_{X}$) plane obtained by NA64 and
presented in comparison with the bounds derived from the results of the LKB
lkb  and Berkeley (B) berkeley  experiments. The limits are shown by lines
labeled with the $X$ type of the same color.

The overall signal efficiency $\epsilon_{X}$ defined as the product of signal
efficiencies accounting for the geometrical acceptance, the track, SRD, ${\it
V}_{2}$ and HCAL reconstruction, and the DAQ dead time was found to be
slightly dependent on $m_{X},E_{X}$ values na64prl19 . The signal-event
reconstruction efficiency $\epsilon_{{\rm ECAL}}$ was estimated as a function
of energy deposited in the ECAL for different $X$ masses. Compared to the
ordinary e-m shower, the $\epsilon_{{\rm ECAL}}$ value for a shower from $X$
event has to be corrected due to difference in the e-m showers development at
the early stage in the ECAL PS gkkk . Depending on the energy threshold in the
PS ($E^{th}_{{\rm PS}}$) used in trigger (8) this correction was
$\lesssim(5\pm 3)\%$ dominated by the errors due to the $E^{th}_{{\rm PS}}$
variation during the run. The ${\it V}_{2}$ and HCAL efficiency defined by the
leak of the signal shower energy from the ECAL to these detectors, was studied
for different $X$ masses with simulations that were validated with
measurements at the $e^{-}$ beam. The uncertainty in the efficiencies
dominated mostly by the pileup effect was estimated to be $\lesssim 4\%$. The
trigger efficiency was found to be $0.95\pm 0.02$. The $X$ signal-event
acceptance was estimated by taking into account the efficiency of selection
cuts for the signal shower shape in the ECAL gkkk . The dominant uncertainty
in the signal yield $\simeq 10\%$ was conservatively accounted for the
difference between the predicted and measured dimuon yield na64prd18 . The
total signal efficiency $\epsilon_{X}$ was in the range 0.5 - 0.7 depending on
the beam intensity and the $X$ mass.

To set the limits we analysed runs I-III simultaneously using the technique
based on the RooStats package root  allowing multibin limit setting na64prd18
. For each of $X=S,~{}P,~{}V,~{}A$ cases, we tried to optimize the size of the
signal box by comparing sensitivities defined as an average expected limit
calculated using the profile likelihood method. The calculations were done by
taking into account the background estimate, efficiencies, and their
corrections with uncertainties used as nuisance parameters Gross:2007zz . For
this optimization, the most important inputs came from the background
extrapolation into the signal region from the data samples of runs I-III with
their errors estimated from the extrapolation procedure. The optimal signal
box size was found to be weakly dependent on the $e-X$ type of interaction and
$X$ mass varying with a few GeV, and was finally set to $E_{{\rm
ECAL}}\lesssim 50$ GeV for all four cases of Eq.(3) and the whole mass range.

Figure 4: Shown are the NA64 90% C.L. exclusion region in the ($m_{X},|\Delta
a_{X}|$) plane for the $S,P,V$ and $A$ contributions to $a_{e}$ together with
the bands of Eqs.(1,2), representing the results of the LKB lkb  (black
dashed) and Berkeley berkeley  (blue solid) experiments. The legend is the
same as for Fig. 3.

The total number of signal events was the sum of expected events from the all
three runs in the signal box:

$N_{X}=\sum_{i=1}^{3}N_{X}^{i}=\sum_{i=1}^{3}n_{EOT}^{i}\epsilon_{X}^{i}n_{X}^{i}(g_{X},m_{X},\Delta
E_{e})$ (11)

where $\epsilon_{X}^{i}$ and $n_{X}^{i}(\epsilon,m_{X},\Delta E_{X})$ is the
signal efficiency and the signal yield per EOT in the energy range $\Delta
E_{e}$, respectively. These values were calculated from simulations and
processing of signal events through the reconstruction program with the same
selection cuts and efficiency corrections as for the data sample from run $i$.

The combined 90% C.L. exclusion limits on the coupling parameter
$\varepsilon_{X}$ as a function of the $X$ mass, calculated by using the
modified frequentist approach na64prl19 ; junk ; limit ; Read:2002hq  are
shown in Fig. 3. By using Eqs.(1), (2) and (4) - (7), it is also possible to
translate the measurements of Refs.lkb ; berkeley  into constraints on the
coupling $\varepsilon_{X}$ which are shown in Fig. 3 for comparison. The
limits were calculated by taking into account the sign of the contributions
$\Delta a_{X}$ in Eqs.(4) - (7) assuming that the $S$ and $V$ contribute to
the deviation of Eq.(1) , while only the $P$ and $A$ can resolve the
discrepancy of Eq.(2). Our bounds are more stringent than those derived from
the results of high-precision measurements of Refs.gabri ; lkb ; berkeley .
Using Eqs.(4) - (7) and obtained limits on the $X$ coupling strength we can
derive constraints on the $X$ contribution $\Delta a_{X}$ to $a_{e}$. This
results in stringent bounds in the range $|\Delta a_{X}|\lesssim
10^{-15}-10^{-13}$ for S, P, V and A with sub-GeV masses, which are shown in
the $(m_{X};|\Delta a_{X}|)$ plane in Fig. 4 together with the experimental
bands of the $\Delta a_{X}$ values defined by Eqs.(1, 2). For the low mass
region $m_{X}\lesssim 10$ MeV the limits were obtained by taking into account
corrections from the exact calculations. These results are an order of
magnitude more sensitive compared to the current accuracy on $a_{e}$ from the
electron $g-2$ experiments and recent high-precision determination of the fine
structure constant, thus demonstrating the strength of the NA64 approach on
probing new physics in the electron $g-2$.

We gratefully acknowledge the support of the CERN management and staff and the
technical staff of the participating institutions for their vital
contributions. This work was supported by the Helmholtz-Institut für Strahlen-
und Kern-physik (HISKP), University of Bonn, the Carl Zeiss Foundation
0653-2.8/581/2, and Verbundprojekt-05A17VTA-CRESST-XENON (Germany), Joint
Institute for Nuclear Research (JINR) (Dubna), the Ministry of Science and
Higher Education (MSHE) and RAS (Russia), ETH Zurich and SNSF Grant No.
197346, 186181 and 186158 (Switzerland), FONDECYT Grants No.1191103 and No.
1190845, SAPHIR Millennium Institute of ANID and ANID PIA/APOYO AFB180002
(Chile).

## References

  * (1) J. Jaeckel and A. Ringwald, Annu. Rev. Nucl. Part. Sci. 60, 405 (2010).
  * (2) R. Essig et al., arXiv:1311.0029.
  * (3) J. Alexander et al., arXiv:1608.08632.
  * (4) M. Battaglieri et al., arXiv:1707.04591.
  * (5) J. Beacham et al., J. Phys. G 47, 010501 (2020); arXiv:1901.09966.
  * (6) R. K. Ellis at el. (European Strategy for Particle Physics Preparatory Group), arXiv:1910.11775.
  * (7) A. Berlin, N. Blinov, G. Krnjaic, P. Schuster, and N. Toro, Phys. Rev. D 99, 075001 (2019).
  * (8) G. Lanfranchi, M. Pospelov, and P. Schuster, arXiv: 2011.02157.
  * (9) C. Boehm and P. Fayet, Nucl. Phys. B 683, 219 (2004).
  * (10) P. Fayet, Phys. Rev. D 75, 115017 (2007).
  * (11) M. Pospelov, A. Ritz, and M. B. Voloshin, Phys. Lett. B 662, 53 (2008).
  * (12) M. Pospelov, Phys. Rev. D 80, 095002 (2009).
  * (13) N. Arkani-Hamed, D. P. Finkbeiner, T. R. Slatyer, and N. Weiner, Phys. Rev. D 79, 015014 (2009).
  * (14) G. W. Bennett et al. (Muon g-2), Phys. Rev. D 73, 072003 (2006).
  * (15) A. J. Krasznahorkay, et al., Phys. Rev. Lett. 116, 042501 (2016).
  * (16) A. J. Krasznahorkay, et al., arXiv:1910.10459.
  * (17) S. N. Gninenko and N. V. Krasnikov, Phys. Lett. B 513, 119 (2001).
  * (18) J. L. Feng, B. Fornal, I. Galon, S. Gardner, J. Smolinsky, T. M. P. Tait, and P. Tanedo, Phys. Rev. Lett. 117, 071803 (2016).
  * (19) J. Grange et al. (Muon g-2), arXiv:1501.06858.
  * (20) T. Mibe (J-PARC g-2), Chin. Phys. C 34, 745 (2010).
  * (21) D. Banerjee et al. (NA64 Collaboration), Phys. Rev. Lett. 120, 231802 (2018).
  * (22) D. Banerjee et al. (NA64 Collaboration), Phys. Rev. D 101, 071101 (2020).
  * (23) E. Depero et al. (NA64 Collaboration), Eur. Phys. J. C 80, 1159 (2020).
  * (24) L. Morel, Zh. Yao, P. Cladé, and S. Guellati-Khélifa, Nature 588, 61 (2020).
  * (25) R. H. Parker, C. Yu, W. Zhong, B. Estey, and H. Müller, Science 360, 191 (2018).
  * (26) T. Aoyama, M. Hayakawa, T. Kinoshita, and M. Nio, Phys. Rev. Lett. 109, 111807 (2012).
  * (27) T. Aoyama, T. Kinoshita, and M. Nio, Atoms 7, 28 (2019).
  * (28) D. Hanneke, S. Fogwell, and G. Gabrielse, Phys. Rev. Lett. 100, 120801 (2008).
  * (29) G. F. Giudice, P. Paradisi, and M. Passera, J. High Energy Phys. 11 (2012) 113.
  * (30) N. V. Krasnikov, Mod. Phys. Lett. A 35, 2050116 (2020).
  * (31) H. S. Lee, Phys. Rev. D 90, 091702(R) (2014).
  * (32) H. Davoudiasl and W. J. Marciano, Phys. Rev. D 98, 075011 (2018).
  * (33) W. J. Marciano, A. Masiero, P. Paradisi, and M. Passera, Phys. Rev. D 94, 115033 (2016).
  * (34) A. Crivellin, M. Hoferichter, and P. Schmidt-Wellenburg, Phys. Rev. D 98, 113002 (2018).
  * (35) E. J. Chun, T. Mondal, J. High Energy Phys. 11 (2020) 077.
  * (36) J. Liu, C. E. M. Wagner, and X.-P. Wang, J. High Energy Phys. 03, (2019) 008.
  * (37) M. Bauer, M. Neuber, S. Renner, M. Schnubel, and A. Thamm, Phys. Rev. Lett. 124, 211803 (2020).
  * (38) S. Gardner, X. Yan, Phys. Rev. D 102, 075016 (2020).
  * (39) L. Darme, F. Giacchino, E. Nardi, and M. Raggi, arXiv: 2012.07894.
  * (40) L. D. Rose, S. Khalil, and S. Moretti, arXiv: 2012.06911 [hep-ph].
  * (41) M. Endo, W. Yin, J. High Energy Phys. 08 (2019) 122.
  * (42) I. Dorsner, S. Fajfer, and S. Saad, Phys. Rev. D 102, 075007 (2020).
  * (43) K.-F. Chen, C.-W. Chiang, and K. Yagyu, J. High Energy Phys. 09 (2020) 119.
  * (44) C. Boehm, X. Chu, J.-L. Kuo, J. Pradler, arXiv: 2010.02954.
  * (45) S. N. Gninenko, N. V. Krasnikov, and V. A. Matveev, Phys. Part. Nucl. 51, 829 (2020).
  * (46) P. A. Zyla et al. (Particle Data Group), Prog. Theor. Exp. Phys. 2020, 083C01 (2020).
  * (47) D. Banerjee et al. (NA64 Collaboration), Phys. Rev. Lett. 123, 121801 (2019).
  * (48) J. P. Lees et al. (BABAR Collaboration), Phys. Rev. Lett. 119, 131804 (2017).
  * (49) S. N. Gninenko, D. V. Kirpichnikov, and N. V. Krasnikov, Phys. Rev. D 100, 035003 (2019).
  * (50) D. V. Kirpichnikov, V. E. Lyubovitskij, and A. S. Zhevlakov, Phys. Rev. D 102, 095024 (2020).
  * (51) S. N. Gninenko, Phys. Rev. D 89, 075008 (2014).
  * (52) S. Andreas et al., arXiv:1312.3309.
  * (53) D. Banerjee, P. Crivelli, and A. Rubbia, Adv. High Energy Phys. 2015, 105730 (2015).
  * (54) E. Depero et al., Nucl. Instrum. Methods Phys. Res., Sect. A 866, 196 (2017).
  * (55) D. Banerjee et al. (NA64 Collaboration), Phys. Rev. Lett. 118, 011802 (2017).
  * (56) D. Banerjee et al. (NA64 Collaboration), Phys. Rev. D 97, 072002 (2018).
  * (57) S. Agostinelli et al. [GEANT4 Collaboration], Nucl. Instrum. Methods Phys. Res., Sect. A 506, 250 (2003).
  * (58) J. Allison et al., IEEE Trans. Nucl. Sci. 53, 270 (2006).
  * (59) A. Celentano, M. Bondi, R. R. Dusaev, D. V. Kirpichnikov, M. M. Kirsanov, N. V. Krasnikov, L. Marsicano, and D. Shchukin, arXiv: 2101.12192.
  * (60) S. N. Gninenko, N. V. Krasnikov, M. M. Kirsanov, and D. V. Kirpichnikov, Phys. Rev. D 94, 095025 (2016).
  * (61) S. N. Gninenko, D. V. Kirpichnikov, M. M. Kirsanov, and N. V. Krasnikov, Phys. Lett. B 782, 406 (2018).
  * (62) R. R. Dusaev, D. V. Kirpichnikov, and M. M. Kirsanov, Phys. Rev. D 102, 055018 (2020).
  * (63) I. Antcheva et al., Comput. Phys. Commun. 180, 2499 (2009).
  * (64) E. Gross, “LHC statistics for pedestrians,” CERN Report No. CERN-2008-001, 2008, p.71.
  * (65) T. Junk, Nucl. Instrum. Methods Phys. Res., Sect. A 434, 435 (1999).
  * (66) G. Cowan, K. Cranmer, E. Gross, and O. Vitells, Eur. Phys. J. C 71, 1 (2011).
  * (67) A. L. Read, J. Phys. G 28, 2693 (2002).