Fast production of large area patterns with nanometre resolution is crucial for the established semiconductor industry and for enabling industrial scale production of next generation quantum devices. Metastable atom lithography with binary holography masks has been suggested as a higher resolution/low cost alternative to the current state of the art: extreme ultra violet (EUV) lithography. However, it was recently shown that the interaction of the metastable atoms with the mask material (SiN) leads to a strong perturbation of the wavefront, not included in existing mask generation theory, which is based on classical scalar waves. This means that the inverse problem  (creating a mask based on a desired pattern) cannot be solved analytically even in 1D. Here we present a machine learning approach to mask generation targeted for metastable atoms. Our algorithm uses a combination of genetic optimisation and deep learning to obtain the mask. A novel deep neural architecture is trained to produce an initial approximation of the mask. This approximation is then used to generate the initial population of the genetic optimisation algorithm that is able to converge to arbitrary precision. We demonstrate generation of arbitrary 1D patterns for system dimensions within the Fraunhofer approximation limit.

§ INTRODUCTION

All semiconductor device fabrication is currently based on pattern generation using mask based photolithography. The standing aim is to create patterns with smaller and smaller features at higher and higher information densities. The state of the art is EUV photolithography, which, with photons (light-waves) of a wavelength of 13.5 nm, should be able to produce patterns with a resolution of 6.75 nm according to the Abbe criterion. Smaller features could potentially be achieved using immersion and/or over-exposure or under-development, as is done with the 193-nm lightsource, EUVs predecessor. However, the problem is that in EUV lithography, due to the high energy of the photons  (91.8 eV), the pattern generation process is mediated by photo-generated secondary electrons, which can travel for several nm before inducing a reaction. Current experiments and theory indicate that the secondary electron blur radius for EUV is around 3 nm, which would limit the feature size that can be achieved to around 6 nm [1]. This means that devices based on small quantum dots and individual atoms and molecules cannot be produced with EUV. Moving to wavelengths even
shorter than 13.5 nm, would just exacerbate the secondary electron issue - an alternative to photons is needed.

Already more than three decades ago lithography with atoms (matter-waves) was suggested as an alternative to photolithography. The reason for this is that for a given wavelength, the energy of an atom is much less than the energy of a photon.  A helium atom with the same wavelength as an EUV photon has a kinetic energy of only $E = h^2/(2m\lambda^2) \approx 0.011 \,\rm{eV}$, where $h$ denotes Planck’s constant and $m$ is the mass of the helium atom. This is almost a factor 10.000 less than the EUV photon. Some of the first atom lithography experiments created the patterns by depositing atoms directly on a substrate.  However, this method strongly limits the application areas, because it put restrictions on the materials which the pattern can be made of. This can be avoided by using metastable  atoms. When a metastable atom hits the substrate, it decays, and the decay energy is transferred to the substrate [2, 3, 4]. In 1995, Berggren et al. [2] demonstrated pattern generation in a thiol-based resist using a beam of metastable argon atoms manipulated by a light-field mask.  The energy released when a metastable atom decays is about 10 eV for argon and 20 eV for helium [3, 5]. A range of experiments have been done generating patterns by manipulating  atomic beams either by light or electrostatic fields [2, 6, 7, 8, 9]. However, these methods cannot be used to generate complex patterns with high resolution.  For high resolution pattern generation experiments have  been done focusing atom beams with lenses [10, 11, 12, 13, 14, 15, 16, 17, 18, 19], which were then used for serial writing of arbitrary patterns. However, for mass-scale production serial writing is not a suitable method.

The particular challenge for pattern generations with metastable atoms is that they decay when they impinge on a surface, so it is not possible to use masks made on substrates as is done in photolithography. The pattern must be generated by open areas in the substrate, where the beam can go through. This means any pattern with a closed path (i.e. a circle) would lead to a segment of the mask falling out.

In 1996, Fujita et al. [20] came up with an idea to circumvent this problem. They demonstrated pattern generation with metastable neon atoms, using a solid mask consisting of a distribution of uniformly sized holes (30 nm in diameter), etched into a silicon nitride membrane. The hole distribution was calculated using the theory of grid-based binary holography developed by Lohmann and Paris [21], and later by Onoe and Kaneko [22] for scalar waves: This theory imposes the limitation  that the openings are all of the same size and positioned on a regular grid structure. The hole distribution is an approximated Fourier transform of the final, desired pattern. In recent publications binary philography was further explored, investigating issues such as  how many holes are needed to generate a particular pattern and what maximum resolution can be achieved with typical matter-wave wavelength scales [23, 9]. However, earlier this year it was shown that the wavefront of metastable helium atoms is perturbed by dispersion forces when it passes through a hole in a silicon nitride membrane [24], an effect not considered in the previous publications. Dispersion forces, in this case, the Casimir–Polder forces, are caused by the quantum-mechanical ground-state fluctuations of the electromagnetic field in the absence of charges [25]. Due to the field fluctuations, the helium atom will be polarised for a short amount of time. The resulting induced dipole moment then interacts with the dielectric membrane via dipole-dipole interactions. These forces decay dramatically with the distance, $r^{-3}$-power law [26], but play a significant role on the nanometre length scale and has been observed in several matter-wave diffraction experiments [12, 27, 28, 29, 30].

For Fujita et al., the mask holes were so big that the dispersion force interaction was negligible. Hence, they could use the standard binary holography theory to generate the mask design for their patterns. However, for smaller mask holes, which are required for high resolution, the effect of the dispersion force interaction becomes significant. In Ref. [24], it is shown that metastable helium atoms cannot penetrate holes in 5 nm thick silicon nitride membranes with less than 2 nm diameter.

The aim of this paper is to establish a theoretical framework for mask-based matter-wave lithography, initially in one dimension, so that for a given desired pattern, atom wavelength and dispersion force interaction a mask can be generated.  Mathematically this is a so called inverse problem [31]. We now explain why machine learning is the only realistic approach to solving this problem even in the one dimensional case:

As discussed above, existing binary holography theory for scalar waves, uses a grid of holes with equal size as the base for the masks and assumes system dimensions and wavelength that fulfil the well known optics condition: Fraunhofer approximation [32]. In the Fraunhofer approximation the inverse problem is reduced to a Fourier transform of the desired pattern, which yields the required mask. The mask needs to be sampled according to the Nyquist–Shannon sampling theorem to account for the finite size of the  holes and mask. In one dimension, this can be done analytically, see Sec. <ref>. However, this solution  cannot be adapted to matter waves, because the dispersion force interaction between the particles and the diffraction object induces a complex phase distribution of the matter wave [19].

An alternative approach to Fourier transformation, starts from convolution of the mask function with the so-called point-spread functions [33] representing a basis of elementary openings. This approach however rises several issues: (i) the complex structure of the required integral kernel (analogously to the Airy disc for matter waves [29]) which (ii) does not scale linearly with the dimension of the holes [24] and (iii) would require additional optimisation algorithms to consider a variation of the hole sizes.

In this paper, we use machine learning as a new approach to mask generation for metastable atoms, in particular we include the dispersion force interaction and opens up for the possibility that holes can be of different sizes and not positioned on a regular grid structure. Our machine learning architecture is based on a combination of deep learning and genetic algorithms. The problem is solved in two steps: first, a deep neural network [34] is trained on large amounts of generated data. Second, the neural network resolves the desired diffraction pattern, producing a mask. This mask is then mutated to produce an initial population that is fed to a genetic algorithm [35]. The genetic algorithm then is able to quickly converge to a mask that produces the desired diffraction pattern with arbitrary precision (see Fig. <ref>).

Deep learning has been used extensively to solve inverse problems due to neural networks being universal function approximators [36]. For example, it has been applied to light inverse scattering problems [37]. The high-level combination used in this paper (deep learning plus genetics algorithm) has been used in photonic device design [38], inverse molecular design [39], and robot manipulators [40], among others. However, to our knowledge, it has never been used to solve an inverse diffraction problem, where a diffraction mask is recovered from a diffraction pattern.

§ RESULTS

§.§ Mask generation framework

The general idea of our mask generation framework follows: first, a deep convolutional neural network is used to approximate a general solution to the inverse problem. This neural network is trained with a large data set of randomly generated examples that map a mask to a diffraction pattern (see Sec. <ref> for details of how the examples are generated). Masks are represented as binary sequences of finite length where 1 represents openings in the mask and 0 blocked space.

Once the approximate solution is obtained, this solution (a binary mask) is randomly mutated. The random mutations of this approximate solution form the initial population of a genetic algorithm. The genetic algorithm — known to be an efficient optimiser for binary sequences [41] — further refines the mask provided by the neural network until a convergence constraint is satisfied.

General mask generation framework used in this paper. First, the desired intensity pattern is fed to . then computes an approximated mask that is expected to produce this pattern. From this mask, an initial population is generated and fed to a genetic algorithm that is run until the accuracy condition is met.

§.§.§ : a neural network for approximate inversion of the Kirchhoff diffraction formula

is a deep convolutional neural network that takes as features (inputs) the desired intensity pattern and its fast Fourier transform and has to output the mask pattern, a binary intensity map. We choose the real-input fast Fourier transform as one of the network's features because the Fourier transform of the diffraction pattern plays a significant role in far-field diffraction [32, 9]. Thus, the diffraction patterns have periodic components that can be efficiently compressed by using its Fourier decomposition, and for this reason, allowing for better learning by our neural network.

We split the network graph in three towers, each specialised in one of the three different features (the intensity map, and the real and intensity parts of the fast Fourier transform). Within each tower, we use skip connections - a mechanism that has been used effectively in the U-Net architecture, an inverse architecture for light scattering problems [42]. Skip connections simply consist of taking the input to a neural network layer (or succession of layers) and summing it after the neural network layer has been applied to it [43].

Finally, we use a fully connected layer to learn the predicted mask (see fig. <ref> for the full architecture). Note that unlike many popular computer vision architectures, we do not use batch normalisation due to the reason that after extensive experimentation, we did not see any benefit for our task.

architecture. A three-tower architecture . The two towers in the left are 1-D convolutional layers with skip blocks that ingest the angle and the modulus of the right fast Fourier transform of the input. The tower in the right is also composed by one skip block with two convolutions. Fully connected layers and concatenations are used to combine the information coming from the three blocks.

We train the neural network with the focal loss. A loss commonly used to train segmentation masks in computer vision [44]. The focal loss is derived from the binary cross entropy loss $l_{\rm BC}$:

\begin{equation}
l_{\rm BC}(y,\hat{y})=\left[y\log\hat{y}+(1-y)\log(1-\hat{y})\right] \,.
\end{equation}

By summing over the all mask points with the total number $N_{\rm points}$, we reach the total loss per mask $L$:

\begin{equation}\label{eq:BCE_loss}
L(y,\hat{y})=-\frac{1}{N_{\rm pixels}}\sum_{\rm mask}l_{\rm BC}(y,\hat{y}) \,,
\end{equation}

with $y$ is the true value of each pixel in the mask sequence and $\hat{y}$ is the value predicted by the neural network.

The focal loss is a variation of the binary cross-entropy loss and reads

\begin{equation}\label{eq:focal_loss}
l_{\rm FO}^i(y,\hat{y})=-\alpha_i\left[1-\mathrm e^{l_{\rm BC}(y,\hat{y})}\right]^\gamma l_{\rm BC}(y,\hat{y}) \,.
\end{equation}

This loss function is designed to achieve two purposes: (I) weight the classes (in our case 1, and 0) according to their rarity through the weights $\alpha_i$, so that they contribute equally to the loss, and (II) reduce the loss of easily-classified examples, so that the network can focus on the harder parts of the mask.

can be defined as a deterministic non-linear function that depends on its parameters $w$. The optimisation problem that solves as training is to find the set of parameters $w$ that minimise the loss function over the training set $\Psi$,

\begin{equation}
w^*=\min_w\left[\sum_{\Psi} L(y,\hat{y})\right]\,.
\end{equation}

is trained on a dataset $\mathcal{D}$ of 300k samples, randomly generated as described in Sec. <ref>. As it's standard in deep learning, this dataset is split in three smaller subsets: a training set $\mathcal{D}_{\rm train}$, a test set $\mathcal{D}_{\rm test}$, and a validation set $\mathcal{D}_{\rm val}$. The split used is 81/9/10 (train/validation/test). One tenth of the data is reserved for test and one tenth of the remaining training data is reserved to validate hyper-parameters. The training of is split in two parts: hyperparameter tuning (which optimises for hyperparameters of the neural network such as learning rate, loss function type etc), and the training of the network itself (i.e. the training of its weights and biases).

We use a variation of Bayesian search to train the hyperparameters of the neural network, and evaluate the validation dataset using Mean Squared Error over the result of the integrals of the produced mask. We choose to evaluate using Mean Squared Error because that is the metric that really tells whether the mask produces the required pattern. In the discretised space, the inverse may not be unique and small differences in the binary sequence that represents the mask might have a strong effect on the focal loss, while actually affecting the mean squared error loss of the integrated mask to a much lesser extent.

For the training of the neural network, we use the best optimisation algorithm (typically Adam [45] or RMSProp [46]) that is returned by the hyperparameter search. Training is stopped using early stopping with patience 10 and delta 0.001. Full details of this procedure and the space over which hyperparameter search is performed can be found in the supplementary information. The hyperparameters used to produce the results presented here are optimiser: Adam, learning rate: 0.00016, batch size: 225, $\alpha$: 0.439, $\gamma$: 5.952.
Figure <ref> shows the focal loss and the mean squared error after integrating the masks evaluated over the validation dataset.

Training (red) and validation (green) focal loss (eq. <ref>) compared with mean squared error computed over the validation set after performing the integral of the solution (blue). Note how the focal loss over the validation seems to exhibit a clear over-fitting effect. However, when calculating the mean squared error of intensity produced by the mask, the network is actually seen to improve gradually as the epochs advance. This is likely due to two effects (i) the intensity integral does not have a unique inverse when discretised, (ii) after a certain number of epochs, the masks produced become very stable as the 0.5 thresholds convert the continuous output of the network into a discrete binary mask.

§.§.§ Genetic algorithm

Genetic algorithms are stochastic classical evolutionary algorithms, that is algorithms that dynamically change in order to optimise a fitness function $F$. Genetic algorithms are inspired by Darwin's theory of evolution that describes how the genes in the population evolve according to their capacity of reproducing and mutating. In genetic algorithms, the population is not an actual population of living beings but instead a set of solutions for a particular optimisation problem (known as chromosomes). Furthermore, mutation is not necessarily the result of natural stochastic processes but occurs according to different rules that can be decided by the programmer.  For each chromosome in the population a fitness value can be calculated using a pre-defined formula set at the choice of the researcher. Finally, the reproductive process of life is simulated by combining different chromosomes (parents) into next generations of offspring. The reproductive process can also be defined by the researcher [47].

Genetic algorithms are very robust optimisation algorithms that are especially suited to work with sequences of categorical variables (such as the four types of bases present in DNA: adenine (A), cytosine (C), guanine (G), and thymine (T)). The masks that we want to obtain in this paper are a perfect example of that: they are fully defined by a sequence of binary values.

The genetic algorithm that we use here is applied to the solution by , the fitness function is simply the inverse of the mean absolute error of the generated diffraction pattern plus a small numerical constant $\alpha$.

\begin{equation}
F=\frac{1}{\alpha+\sum_{r_i}  \lvert\tilde{P}(r_i)-\psi(r_i;{\bm{x}},{\bm{d}})\psi^\star(r_i;{\bm{x}},{\bm{d}})\lvert}.
\end{equation}

$\tilde{P}(r_i)$ here is the absolute square value of the  wave function discretised over a grid of radial coordinates $r_i$ measuring the distance to the centre of the diffraction pattern. $\psi(r_i;{\bm{x}},{\bm{d}})$ is the wave function that is made to depend on the the positions ${\bm{x}}=\left(x_1,x_2,\dots\right)$ and thicknesses ${\bm{d}}=(d_1,d_2,\dots)$ of the grating (mask) openings. A detailed description of how the wave function is computed can be found in Sec. <ref>. To run our genetic algorithm we use a pygad: a well-known genetic algorithm solver for python [48]. Within pygad we use the following hyper-parameters (obtained through grid search). (i) uniform crossover, in which offspring is formed by random recombination of parents' chromosomes. (ii) fit parent persistence: the fittest chromosomes are carried on to the next generation. (iii) initial population size of 50 chromosomes. (iv):  The number of solutions to be selected as parents is set to 7.

figureMean difference in percentage between using to generate the initial population for the genetic algorithm and using a random approach. Results averaged over 100 never-seen randomly sampled test masks. Each mask is approximated by the generated algorithm 50 times. Error bars (shaded area) are 4$\sigma$. Note how outperforms the randomly-generated population significantly, with the percentage difference growing larger as the genetic algorithm runs more generations.

figureAverage fitnesses $F$ for an arbitrary mask after several generations of the genetic algorithm. The plot shows the fitness of the solution when the initial population is given by (blue) compared with an initial population sampled randomly (red). Here shown an example where helps converging to a much more precise solution. Error bars (shaded area) are 4$\sigma$, equivalent to 99.99% certainty.

The initial population is not generated by pygad but by our own bespoke script. The script takes as an input the number of mutations allowed in each chromosome and generates a population by randomly mutating the output of . 's solution is also kept as part of the initial population. The results presented here have 15 mutations in the chromosome. Mutations are assigned sampling randomly from an discrete uniform distribution - which means that a gen randomly set to mutate can maintain its initial value if it is sampled again from he distribution.

Figure <ref> shows the average difference in percentage when comparing using to generate the initial solution population and using a random initial population. Note that within $4\sigma$ error bars outperforms a naive genetic approach. More importantly perhaps, becomes better compared to using a randomly initiated genetic algorithm the more generations that the genetic algorithm is run, showing that the edge provided by the neural network increases with computational time.

Figure <ref> shows an example for a particular mask in which performs significantly better than random. Results for all 100 randomly generated masks can be found in Appendix <ref>. As shown above, when all masks are averaged, is found to outperform random significantly (by more $4\sigma$). Figure <ref> compares the solutions returned by and the genetic algorithm for a particular inverse problem. Note how provides a reasonable first approximation to mask design and how further Genetic iteration incorporates small changes (mutations) that make the solution closer to the ground truth.

Example of a mask returned by (blue) and by enhanced by the Genetic Algorithm after 2000 generations compared (red) to the ground truth (green). Note how the genetic algorithm introduces certain openings in the mask that were missed by the neural network, thus increasing the general fitness of the solution. Note as well how (blue) provides a reasonable first approximation to the mask design.

§.§ Forward propagation

Here we describe how the training patterns are generated by forward propagation of matter-waves through randomly generated masks.  The diffraction at a dielectric interface can be described via classical waves due to the duality of waves and particle [49], which can be derived via Kirchhoff's diffraction formula [32]. This formula determines the propagation of a classical wave with amplitude $a_0$

\begin{eqnarray}
\psi_{\rm{cl}}({\bm{r}}_{\rm D}) = \frac{a_0 k_0}{2\pi \mi}\int \md^2 r_{\rm{i}}\, T({\bm{r}}_{\rm i}) \frac{\me^{\mi k_0\left(r_{\rm S,i}+r_{\rm i,D}\right)}}{r_{\rm S,i}r_{\rm i,D}} \frac{\cos\vartheta+\cos\vartheta'}{2}\,, \label{eq:Kirchhoff}
\end{eqnarray}

through an interface, in analogy to Huygens–Fresnel principle [50]. The interface is usually described by a transmission function $T$, which varies (for electromagnetic waves) between $1$ for transmission and $0$ for absorption. The source of the wave is located at $\rr_{\rm S}$ (subscript S for source), and the grating  is placed at the intermediate points $\rr_{\rm i}$ (subscript i for intermediate), leading to the relative coordinates $\rr_{\rm S,i}={\bm{r}}_{\rm S} - {\bm{r}}_{\rm i}$ to the grating and continues the propagation to the detector $\rr_{\rm D}$ (subscript D for detector) $\rr_{\rm i,D} =\rr_{\rm i}-\rr_{\rm D}$. The wave vector is related via the de-Broglie wavelength $k_0 = 2\pi/\lambda_{\rm dB} = p/\hbar$ with the particle's momentum ${\bm{p}}= p {\bm{e}}_p$. [49] The geometric correction angles $\vartheta$ and $\vartheta'$ are the angles between the aperture's normal ${\bm{n}}$ and the wave's propagation directions, $\rr_{\rm S,i}$ and $\rr_{\rm i,D}$, respectively.

By considering plane waves passing the obstacle, the propagation lengths are dominated by the distances between the source and the interface $r_{\rm S,i}\approx L_1$ and between the interface and the detector $r_{\rm i,D} \approx L_2$. This approach simplifies Kirchhoff's diffraction formula (<ref>) to the Fourier transform of the transmission function and is known as Fraunhofer approximation

\begin{eqnarray}
\psi_{\rm{cl}}(r) = \frac{a_0 k_0}{2\pi\mi}\frac{\me^{\mi k_0 \left[L_1 + L_2 + r^2/(2L_2)\right]}}{L_1L_2}\int\md A\, T \me^{\mi \frac{k_0}{L_2}{\bm{r}}\cdot{\bm{s}}}\,,
\end{eqnarray}

where $r$ denotes the position at the screen. Furthermore, the transmission function is invariant along the direction of the slit and, thus, the wave at the screen is just given by the one-dimensional Fourier transform of the transmission function

\begin{eqnarray}
\psi_{\rm{cl}}(r) = \frac{a_0 k_0}{2\pi\mi}\frac{\me^{\mi k_0 \left[L_1 + L_2 + r^2/(2L_2)\right]}}{L_1L_2}\int\md x\, T(x) \me^{\mi \frac{k_0}{L_2}rx}\,. \label{eq:1d}
\end{eqnarray}

Assuming the mask to be an array of slit openings the total transmission function can be determined by the superposition of the single slits

\begin{equation}
T(x) = \sum_n t(x-x_n;d_n)\,,\label{eq:transmission}
\end{equation}

with the centre coordinate of each opening $x_n$. We allow each opening to vary in its width $d_n$. By inserting the transmission function (<ref>) into the wave equation (<ref>), the field at the screen decomposes into the superposition of single slits with a phase correction due to the corresponding spatial shift

\begin{eqnarray}
\psi(r) = \frac{a_0 k_0}{2\pi\mi}\frac{\me^{\mi k_0 \left[L_1 + L_2 + r^2/(2L_2)\right]}}{L_1L_2}\sum_n\me^{-\mi \frac{k_0}{L_2}rx_n}\int\md x\, t(x;d_n) \me^{\mi \frac{k_0}{L_2}rx}\,. \label{eq:wave}
\end{eqnarray}

The interference pattern is given by the absolute square value of the wave function

\begin{equation}
P(r) = \psi(r)\psi^\star(r) \,.
\end{equation}

This yields that a target pattern $\tilde{P}(r)$ can be approximated by minimising error function

\begin{eqnarray}
E({\bm{x}},{\bm{d}}) = \int \md r \left(\tilde{P}(r)-\psi(r;{\bm{x}},{\bm{d}})\psi^\star(r;{\bm{x}},{\bm{d}})\right)^2\,,\label{eq:error}
\end{eqnarray}

with respect to the positions ${\bm{x}}=\left(x_1,x_2,\dots\right)$ and thicknesses ${\bm{d}}=(d_1,d_2,\dots)$ of the grating openings. Ordinarily, one would apply a least-square fitting algorithm to reduce the error (<ref>) concerning the positions and grating openings.

§.§.§ Propagation of electromagnetic waves

In the absence of dispersion forces, which is achieved for electromagnetic waves, the wave propagation (<ref>) can be calculated explicitly leading to

\begin{equation}
\psi(r) = \frac{a_0 }{\pi\mi}\frac{\me^{\mi k_0 \left(L_1 + L_2 + r^2/(2L_2)\right)}}{L_1 r}\sum_n\me^{-\mi \frac{k_0}{L_2}rx_n}\sin\left(\frac{k_0 r d_n}{2L_2}\right) \,, \label{eq:emw}
\end{equation}

and, thus, the interference pattern reads

\begin{equation}
P(r) =\frac{a_0^2}{\pi^2 L_1^2 r^2}\sum_{n,m} \me^{\mi \frac{k_0}{L_2}r\left(x_m-x_n\right)}\sin\left(\frac{k_0 r d_n}{2L_2}\right) \sin\left(\frac{k_0 r d_m}{2L_2}\right) \,.
\end{equation}

To analyse the relation between the interference pattern $P(r)$ and the transmission function parameters ${\bm{x}}$ and ${\bm{d}}$ concerning the resolution and contrast of the interference pattern, a Fourier analysis is required, which yields [51]

\begin{eqnarray}
\lefteqn{P(\kappa) = -\frac{a_0^2}{2 L_1^2} \sum_{n,m} \left\lbrace \left(\kappa-\frac{c_1}{2\pi}\right)\left[2\Theta\left(\kappa-\frac{c_1}{2\pi}\right)-1\right]+\left(\kappa-\frac{c_2}{2\pi}\right)\left[2\Theta\left(\kappa-\frac{c_2}{2\pi}\right)-1\right]\right.}\nonumber\\
\right\rbrace\,,
\end{eqnarray}

with the Heaviside function $\Theta$ and the frequency shifts

\begin{eqnarray}
c_1= \frac{k_0}{2L_2}\left(2x_m-2x_n-d_n+d_m\right) \,,&&
c_2= -\frac{k_0}{2L_2}\left(-2x_m-2x_n-d_n+d_m\right) \,,\\
c_3= \frac{k_0}{2L_2}\left(2x_m-2x_n+d_n+d_m\right) \,,&&
c_4= -\frac{k_0}{2L_2}\left(-2x_m-2x_n+d_n+d_m\right) \,.
\end{eqnarray}

This analysis illustrates that the relevant contributions is located in the range $\kappa =\mathcal{O}\left(10^{-5} /\lambda_{\rm dB}\right)$ for typical parameters used in matter-wave lithography [spatial resolution of the mask $\Delta x =\mathcal{O}(\rm{nm})$ and distance between the mask and the object plane $L_2=\mathcal{O}(100\,\rm{\mu m})$]. By increasing the mask's extension to be in the same order as $L_2$, the entire spectral resolution will be in the range $10^{-5}/\lambda_{\rm dB} \le \kappa  \le 1/\lambda_{\rm dB}$.

However, due to the presence of dispersion forces, the application of such method would result in several numerical issues caused by the complexity of the dependence on these quantities, which will be illuminated in the next section.

§.§.§ Dispersion  force interactions

The dispersion force interaction between an atom and a dielectric membrane can be approximated by [24]

\begin{equation}
U_{\rm CP,app}(\rr) = -\frac{9 C_3}{\pi}\int \frac{\md^3s}{\left[\left({\bm{s}}-\rr\right)\cdot \left({\bm{s}}-\rr\right)\right]^3}\,,\label{eq:Ham}
\end{equation}

with the $C_3$-coefficient denoting the interaction strength between the atom and a plane built of the same material as the membrane, the position of the atom  $\rr$, and the integration volume $\md^3s$ is bounded by the membrane's surface.

Comparison of the interference patterns obtain by passing a light wave (electromagnetic wave) (orange) and a metastable helium atom matter wave (blue) through the same SiN double-slit mask of 5 nm thickness with 8 nm openings separated by an 8 nm wall. The inset shows the mask used for both cases.

The reduction of the slit opening can be found analogously to Ref. [24] (table 3) due to the applied assumption of a single-wall interface. However, the phase shift needs to be adapted and reads for a membrane of thickness $t$ with a slit opening $d$

\begin{eqnarray}
\varphi(x) \approx - \frac{m\lambda_{\rm dB}}{2\pi\hbar^2} \int U(x,z) \md z= -\frac{12 C_3 m\lambda_{\rm dB} t d \left(d^2+12 x^2\right)}{\hbar^2 \pi \left(d+2x\right)^3\left(d-2x\right)^3}\,.
\end{eqnarray}

To this end, the transmission function for a single slit adapts due to the dispersion force interaction to

\begin{eqnarray}
t(x) =\Theta\left(\frac{d-2\Delta R}{2}-x\right)\Theta\left(\frac{d-2\Delta R}{2}+x\right) \mathrm e^{\mi \varphi(x)}\,. \label{eq:transmw}
\end{eqnarray}

The impact of the dispersion forces (Casimir–Polder)  is illustrated in Fig. <ref>, which shows the difference in the interference pattern from the classical double slit experiment using a scalar wave and a matter-wave with dispersion force interaction. A much higher population of the higher diffraction orders can be observed.

§.§.§ Generation of training data

In this proof of concept paper we have chosen to consider a  1D system with dimensions so that the Fraunhofer approximation applies. We use a typical wavelength for metastable helium $\lambda_{\rm dB} = 0.1 \,\rm{nm}$, assuming a helium point-source located $L_1 = 1 \,\rm{m}$ away from the mask, and a substrate patterning plane at $L_2 = 300\,\rm{\mu m}$ behind the mask. The mask is considered to be made of SiN with a thickness of $5\,{\rm{nm}}$. The maximum extension of the mask is restricted to $200\,\rm{nm}$ to stay within the Fraunhofer approximation (satisfying a Fresnel number below 1). Furthermore the mask is separated into 50 sections of $4\,\rm{nm}$ width each, corresponding to the minimum opening width. This width was selected based on Ref. [24], which shows that the effective entrance width is reduced by about 2  nm due to the dispersion force interaction.

We randomly generate masks represented via arrays with entries 0  for close or 1 for open sections. Neighbouring open sections will be combined to a single opening with larger thickness. Hence, the randomly generated array will be transferred into a set of width and position for each opening. Thus, the allowed openings varies between 1 and 50 units. To this end, we calculate the single slit diffraction [integrant in Eq. (<ref>)]

\begin{equation}\label{eq:data_generation}
\int \mathrm d x t(x;d_n) \mathrm e^{\mi \frac{k_0}{L_2}r x} \,,
\end{equation}

for all possible thicknesses $d_n$ for matter waves (<ref>) and electromagnetic waves (<ref>) via standard numerical integration techniques. These results are tabled and a second program calculates the superposition of multislits to generate the training data.

§.§ An application example: the double slit pattern

To demonstrate the performance of , we  invert a well-known diffraction pattern: the pattern resulting from a double slit mask, see Fig. <ref>. Note that a double slit is very different from the 's training set, which consist of randomly generated masks with no particular preset structure.

successfully inverts the pattern to a good degree of accuracy after 10000 generations of the Genetic Algorithm. Interestingly, it does not obtain the exact same mask as the ground truth but the mask that it produces matches very closely to the desired pattern (see Figs. <ref> and <ref>).

figureMask returned by after attempting to invert the pattern for a metastable helium beam corresponding to a double slit mask with slit width 12 nm, spacing 28 nm (measured between the end of the first slit and the beginning of the second), and an assumed membrane material of 5 nm SiN. Note how the mask is slightly different. However, the pattern produced by it is very similar t the desired one.

figurePattern returned by the mask produced by after attempting to solve the double slit inversion problem. Note the near-perfect match with the ground truth despite the significant difference in the returned mask (Fig. <ref>).

§ DISCUSSION

Despite the known capabilities of deep learning as a powerful inverse strategy, we found that the best strategy to obtain accurate result was to couple it with an additional optimisation step (in our case, a Genetic Algorithm). 's difficulty in fully solving the inverse could be due to (i) insufficient data, (ii) insufficient network complexity, or (iii) the very hard inverse problem that it is trying to solve.

We performed several experiments aimed to determine the cause of 's difficulty with some examples. Interestingly, increasing the complexity of the network (by adding layers, and increasing the number of neurons in the layers), did not in general produce a better generalisation. With regards to point (ii), we did see a clear improvement with increasing the size of our data-set (which we can arbitrarily scale up by solving for Eq. (<ref>)). The experiments presented here use a data set of 300k examples. By using the data generation procedure described in Sec. <ref>, one can extend the data set arbitrarily, potentially increasing the generalisation capabilities of .

Finally, we saw some examples for which initialising the mask using a random sequence or using made no significant difference. By observing these examples, it seems to be the case that this happens more often for masks that are containing many small openings. These masks are more likely to be automatically randomly produced by the Genetic Algorithm. Therefore, it is possible that for some cases the random initialisation is as good as, or even better than 's guess.

§ CONCLUSION

We present a machine learning architecture for solving the inverse problem for matter wave lithography: . The architecture of consists of a convolutional deep neural network followed by a Genetic Algorithm.  This is a key technological step for industrial scale matter-wave lithography, needed to achieve the ultimate goal of fast patterning of individual atoms and molecules over large areas. The dispersion force interaction for matter waves presents a great computational challenge compared to scalar waves due to the complex integral and its nonlinear solution.

In this proof of concept paper we have restricted ourselves to the 1D case in the Fraunhofer approximation regime. The focus of future work will be to modify to solve the general problem in 2D.

§ ACKNOWLEDGEMENTS

The computations were performed on the Norwegian Research and Education
Cloud (NREC), using resources provided by the University of
Bergen. http://www.nrec.no/
Special thanks to Tore Burheim and Tor Lædre
This work has received funding from the European Union’s Horizon 2020 research and innovation
programme. B.H. gratefully acknowledges support from the  H2020-FETOPEN-2018-2019-2020-01 under Grant Agreement No. 863127 Nanolace (www.nanolace.eu). J.F. gratefully acknowledges support from H2020-MSCA-IF-2020, Grant Agreement No. 101031712.

§ HYPERPARAMETER TUNING

We optimise our hyperparameters using the Tree-structured Parzen Estimator algorithm [52] - a bayesian algorithm commonly found in hyperparameter optimisation packages. To do so, we use optuna [53], an open source hyperparameter search library for python. To find our hyperparameters we do two searches: one in which we focus in finding the best batch size for our problem - which the hyper-parameter optimiser finds to be 225 and another one in which we focus on optimising the parameters of the focal loss within a more constrained space. The complete list of trials and their reported Mean Squared Error can be found in Table <ref>.
Results of the hyperparameter tuning rounds in

optimizer            lr     batch_size         alpha         gamma         MSE

optimizer            lr     batch_size         alpha         gamma         MSE

11          Adam      0.000163            225      0.439672      5.951962      1.399777

43       RMSprop      0.000288            225      0.448443      6.431573      1.404423

47       RMSprop      0.000123             87      0.439803      6.291578      1.421799

48       RMSprop      0.000283            200      0.448358      6.074838      1.433143

39       RMSprop      0.000099             66      0.466304      7.131291      1.450584

4           Adam      0.000098            225      0.584123      5.253094      1.456233

22          Adam      0.000173            225       0.40145      5.910984      1.456791

23          Adam      0.000109            225      0.450265      5.525971      1.465228

12          Adam      0.000187            225      0.412308       6.12187      1.465743

32       RMSprop      0.000089            114      0.179152      7.315897      1.466564

40       RMSprop      0.000082            570      0.471317      3.003355      1.469744

46       RMSprop      0.000059            153      0.295671      7.276354      1.472631

16          Adam      0.000133            225      0.521078       5.03229      1.477234

6           Adam      0.000097            225      0.467451       5.42483      1.478946

35       RMSprop      0.000179             66      0.540199      7.791814      1.485383

2        RMSprop      0.000052            225      0.487136      7.756889      1.485913

0        RMSprop      0.000075            225      0.505513      6.664503      1.488749

15          Adam      0.000035            225      0.436437      5.701631      1.490682

37       RMSprop      0.000151             73      0.560933      6.776026      1.495198

28       RMSprop      0.000035            265      0.364326      4.117896      1.497929

19          Adam      0.000144            225      0.582637      6.224196      1.499658

21          Adam      0.000261            225      0.422763      6.114556      1.502561

20          Adam      0.000025            225      0.401408      6.951192      1.506756

7           Adam      0.000069            225      0.496068      7.374401       1.51577

13          Adam      0.000031            225      0.452886      5.889216      1.530536

49       RMSprop      0.000265            195       0.52036       6.28119      1.544203

18          Adam      0.000016            225       0.43794       5.56616      1.549521

38       RMSprop      0.000434            605      0.494833      6.954664      1.549902

29          Adam      0.000039            400      0.360762      7.630682      1.560381

45       RMSprop      0.000329            332      0.413058      4.932252      1.580241

1        RMSprop      0.000285            225      0.555782      7.176087      1.583419

26          Adam      0.000047            987       0.32565       4.25768      1.593861

33          Adam      0.000032            155      0.120995      5.757826      1.606285

36       RMSprop       0.00018            142       0.59511      7.965855      1.611265

27       RMSprop       0.00003            391      0.234176      1.175896      1.626354

34       RMSprop      0.000013            789       0.38773       4.12605      1.627522

31          Adam      0.000076            900      0.167727      5.418842      1.630291

10          Adam      0.000013            225      0.593005      6.110873      1.638141

41       RMSprop      0.000013            241      0.270598      5.188218      1.647733

25          Adam       0.00029            784      0.230935      2.162501      1.686709

8        RMSprop       0.00037            225      0.478398      5.746003      1.688195

14          Adam      0.000184            225      0.533729      6.477521       1.69213

42       RMSprop      0.000112            435      0.101017      6.930886       1.71168

30       RMSprop      0.000706            285      0.407261      6.188819      1.725915

24          Adam      0.000514            225      0.400875      5.872189      1.773344

5        RMSprop      0.000757            225      0.415991      5.231926      1.793177

44       RMSprop      0.000935            226      0.460966      6.326114      1.799258

17          Adam       0.00037            225      0.553661      6.577556      1.904586

3           Adam      0.000923            225      0.498035       5.16874      2.041645

9           Adam      0.000786            225       0.57071      5.326067      2.277016

§ EXTENDED RESULTS

In Figure <ref>, we show a randomly selected example from the 100 examples that are averaged in Fig. <ref>. In this section we show all other examples. There is no inherent physical difference between the figures shown here and  Fig. <ref> - they simply correspond to different randomly-generated masks. The figures are vectorised, so they can be zoomed in digitally and explored at full resolution. It is evident from the figures presented here that most of them follow the same trend as Fig. <ref>.
in figure_netvsrandom_0,

in figure_netvsrandom_8,

in figure_netvsrandom_80,

in figure_netvsrandom_88,

[1]
R. A. Lawson and A.P. Robinson.
Materials and processes for next generation lithography.
Elsevier, Amsterdam, 1th edition, 2016.

[2]
Karl K. Berggren, Andreas Bard, James L. Wilbur, John D. Gillaspy, Andreas G.
Helg, Jabez J. McClelland, Steven L. Rolston, William D. Phillips, Mara
Prentiss, and George M. Whitesides.
Microlithography by using neutral metastable atoms and self-assembled
Science, 269(5228):1255–1257, 1995.

[3]
Kenneth Baldwin.
Metastable helium: Atom optics with nano-grenades.
Contemporary Physics, 46(2):105–120, 2005.

[4]
B Ueberholz, S Kuhr, D Frese, V Gomer, and D Meschede.
Cold collisions in a high-gradient magneto-optical trap.
Journal of Physics B: Atomic, Molecular and Optical Physics,
35(23):4899–4914, nov 2002.

[5]
L. Colli.
Ultraviolet photons in the decay of metastable argon atoms.
Phys. Rev., 95:892–894, Aug 1954.

[6]
Jamie R. Gardner, Erik M. Anciaux, and Mark G. Raizen.
Communication: Neutral atom imaging using a pulsed electromagnetic
The Journal of Chemical Physics, 146(8):081102, 2017.

[7]
C. S. Adams, M. Sigel, and J. Mlynek.
Atom optics.
Physics Reports, 240(3):143–210, 1994.

[8]
H. Hinderthür, A. Pautz, F. Ruschewitz, K. Sengstock, and W. Ertmer.
Atom interferometry with polarizing beam splitters.
Phys. Rev. A, 57:4730–4735, Jun 1998.

[9]
Torstein Nesse, Ingve Simonsen, and Bodil Holst.
Nanometer-resolution mask lithography with matter waves: Near-field
binary holography.
Phys. Rev. Applied, 11:024009, Feb 2019.

[10]
S. D. Eder, A. K. Ravn, B. Samelin, G. Bracco, A. Salvador Palau, T. Reisinger,
E. B. Knudsen, K. Lefmann, and B. Holst.
Zero-order filter for diffractive focusing of de broglie matter
Phys. Rev. A, 95:023618, Feb 2017.

[11]
M. KOCH, S. REHBEIN, G. SCHMAHL, T. REISINGER, G. BRACCO, W. E. ERNST, and
B. HOLST.
Imaging with neutral atoms—a new matter-wave microscope.
Journal of Microscopy, 229(1):1–5, 2008.

[12]
R. B. Doak, R. E. Grisenti, S. Rehbein, G. Schmahl, J. P. Toennies, and Ch.
Towards realization of an atomic de broglie microscope: Helium atom
focusing using fresnel zone plates.
Phys. Rev. Lett., 83:4229–4232, Nov 1999.

[13]
O. Carnal, M. Sigel, T. Sleator, H. Takuma, and J. Mlynek.
Imaging and focusing of atoms by a fresnel zone plate.
Phys. Rev. Lett., 67:3231–3234, Dec 1991.

[14]
S. D. Eder, X. Guo, T. Kaltenbacher, M. M. Greve, M. Kalläne, L. Kipp, and
B. Holst.
Focusing of a neutral helium beam with a photon-sieve structure.
Phys. Rev. A, 91:043608, Apr 2015.

[15]
S D Eder, T Reisinger, M M Greve, G Bracco, and B Holst.
Focusing of a neutral helium beam below one micron.
New Journal of Physics, 14(7):073014, jul 2012.

[16]
Forest S. Patton, Daniel P. Deponte, Greg S. Elliott, and Stephen D. Kevan.
Speckle patterns with atomic and molecular de broglie waves.
Phys. Rev. Lett., 97:013202, Jul 2006.

[17]
M Barr, A Fahy, J Martens, A P Jardine, D J Ward, J Ellis, W Allison, and P C
Unlocking new contrast in a scanning helium microscope.
Nature Communications, 7(1):10189, 2016.

[18]
Thomas Reisinger and Bodil Holst.
Neutral atom and molecule focusing using a fresnel zone plate.
Journal of Vacuum Science & Technology B: Microelectronics and
Nanometer Structures Processing, Measurement, and Phenomena,
26(6):2374–2379, 2008.

[19]
Johannes Fiedler, Wijnand Broer, and Stefan Scheel.
Reconstruction of casimir—polder interactions from
matter-wave interference experiments.
Journal of Physics B: Atomic, Molecular and Optical Physics,
50(15):155501, jul 2017.

[20]
J Fujita, M Morinaga, T Kishimoto, M Yasuda, S Matsui, and F Shimizu.
Manipulation of an atomic beam by a computer-generated hologram .
Nature, 380(6576):691–694, 1996.

[21]
A. W. Lohmann and D. P. Paris.
Binary fraunhofer holograms, generated by computer.
Appl. Opt., 6(10):1739–1748, Oct 1967.

[22]
M. Onoe and M. Kaneko.
Computer generated pure binary hologram.
Electronics Communications of Japan, 62:118–126, January 1979.

[23]
Torstein Nesse, Jean-Philippe Banon, Bodil Holst, and Ingve Simonsen.
Optimal design of grid-based binary holograms for matter-wave
Phys. Rev. Applied, 8:024011, Aug 2017.

[24]
Johannes Fiedler and Bodil Holst.
An atom passing through a hole in a dielectric membrane: impact of
dispersion forces on mask-based matter-wave lithography.
Journal of Physics B: Atomic, Molecular and Optical Physics,
55(2):025401, jan 2022.

[25]
Stefan Scheel and Stefan Yoshi Buhmann.
Macroscopic QED - concepts and applications.
Acta Physica Slovaca, 58(5):675–809, feb 2008.

[26]
Stefan Yoshi Buhmann, Stefan Scheel, Simen Å. Ellingsen, Klaus Hornberger,
and Andreas Jacob.
Casimir-polder interaction of fullerene molecules with surfaces.
Phys. Rev. A, 85:042513, Apr 2012.

[27]
Markus Arndt, Olaf Nairz, Julian Vos-Andreae, Claudia Keller, Gerbrand van der
Zouw, and Anton Zeilinger.
Wave–particle duality of C60 molecules.
Nature, 401(6754):680–682, oct 1999.

[28]
Christian Brand, Johannes Fiedler, Thomas Juffmann, Michele Sclafani, Christian
Knobloch, Stefan Scheel, Yigal Lilach, Ori Cheshnovsky, and Markus Arndt.
A Green's function approach to modeling molecular diffraction in the
limit of ultra-thin gratings.
Annalen der Physik, 527(9-10):580–591, oct 2015.

[29]
Joshua L Hemmerich, Robert Bennett, Thomas Reisinger, Stefan Nimmrichter,
Johannes Fiedler, Horst Hahn, Herbert Gleiter, and Stefan Yoshi Buhmann.
Impact of Casimir-Polder interaction on Poisson-spot diffraction at
a dielectric sphere.
Physical Review A, 94(2):1–11, 2016.

[30]
Nicolas Gack, Christian Reitz, Joshua Leo Hemmerich, Max Könne, Robert
Bennett, Johannes Fiedler, Herbert Gleiter, Stefan Yoshi Buhmann, Horst Hahn,
and Thomas Reisinger.
Signature of short-range van der waals forces observed in poisson
spot diffraction with indium atoms.
Phys. Rev. Lett., 125:050401, Jul 2020.

[31]
Richard C. Aster, Brian Borchers, and Clifford H. Thurber, editors.
Parameter Estimation and Inverse Problems (Second Edition).
Academic Press, Boston, second edition edition, 2013.

[32]
M. Born and E. Wolf.
Principles of Optics: Electromagnetic Theory of Propagation,
Interference and Diffraction of Light (7th Edition).
Cambridge University Press, Cambridge, 7th edition, 1999.

[33]
Joseph J.M. Braat, Sven van Haver, Augustus J.E.M. Janssen, and Peter
Chapter 6 assessment of optical systems by means of point-spread
volume 51 of Progress in Optics, pages 349–468. Elsevier,

[34]
Ian Goodfellow, Yoshua Bengio, and Aaron Courville.
Deep learning.
MIT press, 2016.

[35]
Antonio Gonzalez and Raúl Perez.
A learning system of fuzzy control rules based on genetic algorithms.
Genetic algorithms and soft computing, Physica-Verlag, pages
202–225, 1996.

[36]
Kurt Hornik, Maxwell Stinchcombe, and Halbert White.
Multilayer feedforward networks are universal approximators.
Neural networks, 2(5):359–366, 1989.

[37]
Xudong Chen, Zhun Wei, Maokun Li, and Paolo Rocca.
A review of deep learning approaches for inverse scattering problems
(invited review).
Progress In Electromagnetics Research, 167:67–81, 2020.

[38]
Yangming Ren, Linxuan Zhang, Weiqiang Wang, Xinyu Wang, Yufang Lei, Yulong Xue,
Xiaochen Sun, and Wenfu Zhang.
Genetic-algorithm-based deep neural networks for highly efficient
photonic device design.
Photonics Research, 9(6):B247–B252, 2021.

[39]
AkshatKumar Nigam, Robert Pollice, and Alán Aspuru-Guzik.
Janus: parallel tempered genetic algorithm guided by deep neural
networks for inverse molecular design.
arXiv preprint arXiv:2106.04011, 2021.

[40]
RaşIt KöKer.
A genetic algorithm approach to a neural-network-based inverse
kinematics solution of robotic manipulators based on error minimization.
Information Sciences, 222:528–543, 2013.

[41]
M Zhao, Zhe Chen, and Frede Blaabjerg.
Optimisation of electrical system for offshore wind farms via genetic
IET Renewable Power Generation, 3(2):205–216, 2009.

[42]
Zhun Wei and Xudong Chen.
Deep-learning schemes for full-wave nonlinear inverse scattering
IEEE Transactions on Geoscience and Remote Sensing,
57(4):1849–1860, 2018.

[43]
A Emin Orhan and Xaq Pitkow.
Skip connections eliminate singularities.
arXiv preprint arXiv:1701.09175, 2017.

[44]
Shruti Jadon.
A survey of loss functions for semantic segmentation.
In 2020 IEEE Conference on Computational Intelligence in
Bioinformatics and Computational Biology (CIBCB), pages 1–7. IEEE, 2020.

[45]
Diederik P Kingma and Jimmy Ba.
Adam: A method for stochastic optimization.
arXiv preprint arXiv:1412.6980, 2014.

[46]
Geoffrey Hinton, Nitish Srivastava, and Kevin Swersky.
Neural networks for machine learning lecture 6a overview of
mini-batch gradient descent.
Cited on, 14(8):2, 2012.

[47]
Agoston E Eiben, James E Smith, et al.
Introduction to evolutionary computing, volume 53.
Springer, 2003.

[48]
Ahmed Fawzy Gad.
Pygad: An intuitive genetic algorithm python library.
arXiv preprint arXiv:2106.06158, 2021.

[49]
A. Messiah.
Quantum Mechanics.
Dover Books on Physics. Dover Publications, Mineola, New York, 2014.

[50]
L. D. Landau and E. M. Lifshitz.
The Classical Theory of Fields. Vol. 2.
Pergamon Press, Oxford, 4th edition, 1975.

[51]
D.W. Kammler.
A First Course in Fourier Analysis.
Prentice Hall, 2000.

[52]
James Bergstra, Rémi Bardenet, Yoshua Bengio, and Balázs Kégl.
Algorithms for hyper-parameter optimization.
Advances in neural information processing systems, 24, 2011.

[53]
Takuya Akiba, Shotaro Sano, Toshihiko Yanase, Takeru Ohta, and Masanori Koyama.
Optuna: A next-generation hyperparameter optimization framework.
In Proceedings of the 25th ACM SIGKDD international conference
on knowledge discovery & data mining, pages 2623–2631, 2019.


# Region Invariant Normalizing Flows for Mobility Transfer

Vinayak Gupta  IIT Delhi<EMAIL_ADDRESS>and  Srikanta Bedathur
IIT Delhi<EMAIL_ADDRESS>
(2021)

###### Abstract.

There exists a high variability in mobility data volumes across different
regions, which deteriorates the performance of spatial recommender systems
that rely on region-specific data. In this paper, we propose a novel transfer
learning framework called Reformd, for continuous-time location prediction for
regions with sparse checkin data. Specifically, we model user-specific
checkin-sequences in a region using a marked temporal point process (MTPP)
with normalizing flows to learn the inter-checkin time and geo-distributions.
Later, we transfer the model parameters of spatial and temporal flows trained
on a data-rich origin region for the next check-in and time prediction in a
target region with scarce checkin data. We capture the evolving region-
specific checkin dynamics for MTPP and spatial-temporal flows by maximizing
the joint likelihood of next checkin with three channels (1) checkin-category
prediction, (2) checkin-time prediction, and (3) travel distance prediction.
Extensive experiments on different user mobility datasets across the U.S. and
Japan show that our model significantly outperforms state-of-the-art methods
for modeling continuous-time sequences. Moreover, we also show that Reformd
can be easily adapted for product recommendations _i.e._ , sequences without
any spatial component.

Normalizing Flows, POI Recommendation, Transfer Learning

††journalyear: 2021††copyright: acmlicensed††conference: Proceedings of the
30th ACM International Conference on Information and Knowledge Management;
November 1–5, 2021; Virtual Event, Australia.††booktitle: Proceedings of the
30th ACM Int’l Conf. on Information and Knowledge Management (CIKM ’21),
November 1–5, 2021, Virtual Event, Australia††price: 15.00††isbn:
978-1-4503-8446-9/21/11††doi: 10.1145/3459637.3482169††ccs: Information
systems Location based services

##  1\. Introduction

Recent research has shown that accurate advertisements on Points-of-Interest
(POI) networks, such as Foursquare and Instagram, can achieve up to 25 times
the return-on-investment (Metev, 2021). Consequently, predicting the time-
evolving mobility of users, _i.e.__where_ and _when_ , is of utmost importance
to power systems relying on spatial data. Current approaches (Likhyani et al.,
2017; Cho et al., 2011; Yuan et al., 2013) overlook the temporal aspect of a
recommender system as it involves modeling continuous-time checkin sequences –
which is non-trivial with standard neural architectures (Du et al., 2016; Zuo
et al., 2020; Mei and Eisner, 2017). The problem is further aggravated by the
variation in volumes of mobility data across regions due to the growing
awareness for personal data privacy (Barkhuus and Dey, 2003; Wei et al.,
2016). Therefore, there exists an underlying region-based data-scarcity, which
further exacerbates the training of large neural models.

In recent years, Marked Temporal Point Processes(MTPP) have outperformed other
neural architectures for characterizing asynchronous events localized in
continuous time and are even used in a wide range of applications, including
healthcare (Rizoiu et al., 2018), finance (Zhang et al., 2020; Bacry et al.,
2015), and social networks (Gupta et al., 2021; Mei and Eisner, 2017; Zuo et
al., 2020). Recent works that deploy MTPP for predicting user mobility
patterns are either: (i) limited to predicting the time of user-location
interactions rather than actual locations (Yang et al., 2018), (ii) restricted
to _one_ dataset without a foreseeable way to easily utilize external
information (Likhyani et al., 2020), or (iii) disregard the opportunity to
reuse trained parameters from external datasets by jointly embedding the
checkin and time distributions (Chen et al., 2020). Thus, none of these
approaches can be used for designing mobility prediction models for limited
data regions.

In this paper, we present Reformd (Reusable Flows for Mobility Data), a novel
transfer learning framework that learns spatial and temporal distribution of
checkins using normalizing flows(NFs) on a checkin-rich _origin_ region and
transfers them for efficient prediction in a checkin-scarce _target_ region.
Specifically, we consider the series of checkins made by a user as her checkin
sequence and model these sequences for all users from a region using a neural
MTPP and learn the inter-checkin time interval and spatial-distance
distributions as two independent NFs (McLachlan and Peel, 2004; Rezende and
Mohamed, 2015). To make the learned spatial and temporal NFs _invariant_ of
the underlying region, we restrict our model to learn the distribution of
inter-checkin time intervals and spatial distance. These features are
unaffected by the network characteristics that vary across regions – POI
categories and user affinities towards these POIs. Therefore, these NFs can be
easily extended for prediction in other mobility regions. The ability of NFs
to provide faster sampling and closed-form training for continuous-time event
sequences (Shchur et al., 2020) make them a perfect medium to transfer
mobility information. Moreover, for transferring across regions, we cluster
the checkin sequences of each region, with each cluster containing checkin-
sequences with _similar_ spatial and temporal checkin patterns and only
transfer the parameters across these clusters.

In summary, the key contributions we make in this paper via Reformd are three-
fold:  (i) We propose Reformd, a transfer-learning model for predicting
mobility dynamics in checkin-scarce datasets by incorporating mobility
parameters trained on a checkin-rich region.  (ii) We present a novel NF-based
transfer over the MTPP that not only enables a faster sampling of time and
distance features of next checkin, but also achieves high performance even
with limited fine-tuning on the target region.  (iii) Finally, we empirically
show that Reformd outperforms the state-of-the-art models by up to 20% and 23%
for checkin-category and time prediction and can easily be extended to product
recommendation datasets.

##  2\. Related Work

The key related works fall into the following categories:

Mobility Prediction: Recent sequential POI prediction models consider the
checkin trajectory for each user as a sequence of events and utilize an RNN
based learning (Cheng et al., 2013; Manotumruksa et al., 2018; Yuan et al.,
2013) with some variants that incorporate the spatial features as well (Liu et
al., 2016; Feng et al., 2018). Another approach (Figueiredo et al., 2016) is a
generic model for predicting user trajectories as well as next product
recommendation. Recent approaches for checkin time prediction are limited to a
single dataset (Chen et al., 2020; Yang et al., 2018; Likhyani et al., 2020).
They also model event-times as random variables rather than sequential flows
and thus cannot be used for transfer across regions.

Temporal Point Process: In recent years TPPs have emerged as a powerful tool
to model asynchronous events localized in continuous time (Daley and Vere-
Jones, 2007; Hawkes, 1971), which have a wide variety of applications, e.g
information diffusion, disease modeling, finance, etc. Driven by these
motivations, in recent years, there has been a surge of works on TPPs (Rizoiu
et al., 2018; Gupta et al., 2021; Farajtabar et al., 2017). Modeling the event
sequences via a neural network led to further developments including neural
Hawkes process (Mei and Eisner, 2017) and several other neural models of TPPs
(Xiao et al., 2017b, a; Omi et al., 2019), but cannot incorporate
heterogeneous features as in a spatial networks. The approach most similar to
our model is (Shchur et al., 2020) that learns the inter-event time intervals
using NFs, but ignores the spatial dynamics and is limited to a single data
source.

##  3\. Problem Setup

We consider the mobility records for two regions with non-overlapping
locations and users, origin and target as $\mathcal{R}^{\text{org}}$ and
$\mathcal{R}^{\text{tgt}}$ respectively. For any region, we represent a user
trajectory as a sequence of checkins represented by
$\mathcal{S}_{k}=\\{e_{i}=(c_{i},t_{i},d_{i})|i\in[k],t_{i}<t_{i+1},d_{i}<d_{i+1}\\}$,
where $t_{i}\in\mathbb{R}^{+}$ is the checkin time, $d_{i}\in\mathbb{R}^{+}$
is the total distance traveled, and $c_{i}\in\mathcal{C}$ is a discrete
category of the $i$-th checkin with $\mathcal{C}$ as the set of all
categories, and $\mathcal{S}_{k}$ denotes the first $k$ checkins. We represent
the inter-checkin times and distances as, $\Delta_{t,k}=t_{k}-t_{k-1}$ and
$\Delta_{d,k}=d_{k}-d_{k-1}$ respectively and model their distribution using
NFs. Our goal is to capture these region invariant dynamics in origin region
for mobility prediction in target region, _i.e._ given the checkin sequence
for target region, $\mathcal{S}^{\text{tgt}}_{K}$ and the MTPP trained on
origin, we aim to predict the time and category of the next checkin,
$e^{\text{tgt}}_{K+1}$.

##  4\. Model Description

We divide the working of Reformd into two parts: (i) the neural MTPP to
capture mobility dynamics specific to a region, and (ii) transfer of NFs
trained on the origin region to the target region. The overall schematic of
Reformd is given in Figure 1.

###  4.1. Region-Specific MTPP

We model the checkin sequences using an MTPP that we build on a recurrent
neural network (RNN). The RNN is used to obtain time-conditioned vector
representation of sequences as in (Du et al., 2016; Mei and Eisner, 2017; Omi
et al., 2019). Later, via these embeddings we estimate the mark distribution
and inter-event time and space densities sing a three-stage architecture:

Input stage: In this stage we represent the incoming checkin at index $k$,
$e_{k}$ using a suitable vector embedding, $\boldsymbol{v}_{k}$ as:

(1)
$\boldsymbol{v}_{k}=\boldsymbol{w}_{c}c_{k}+\boldsymbol{w}_{t}\Delta_{t,k}+\boldsymbol{w}_{d}\Delta_{d,k}+\boldsymbol{b}_{v},$

where $\boldsymbol{w}_{\bullet},\boldsymbol{b}_{\bullet}$ are trainable
parameters and $\boldsymbol{v}_{k}$ denotes the vector embedding for checkin
$e_{k}$ respectively.

Update stage: In this stage, we update the hidden state representation of the
RNN to include the current checkin $e_{k}$ as:

(2)
$\boldsymbol{s}_{k}=\tanh(\boldsymbol{G}_{s}\boldsymbol{s}_{k-1}+\boldsymbol{G}_{v}\boldsymbol{v}_{k}+\boldsymbol{g}_{t}\Delta_{t,k}+\boldsymbol{g}_{d}\Delta_{d,k}+\boldsymbol{b}_{s}),$

where
$\boldsymbol{G}_{\bullet},\boldsymbol{g}_{\bullet},\boldsymbol{b}_{\bullet}$
are trainable parameters and $\boldsymbol{s}_{k}$ denotes the RNN hidden
state, i.e. a cumulative embedding for all previous checkins till the current
time $t_{k}$.

Figure 1.  Architecture of Reformd with flow-based transfer between origin
region (red) and target region (blue).

Output stage: Given the trajectory embedding $\boldsymbol{s}_{k}$, we predict
the next checkin time and the checkin category. Unlike (Du et al., 2016; Mei
and Eisner, 2017) that learn the time distribution using the RNN hidden state,
we model the density of arrival times using a LogNormal (Shchur et al., 2020)
flow denoted as $p_{t}(\Delta_{t,k+1})$ conditioned on $s_{k}$ as:

(3)
$p_{t}(\Delta_{t,k+1}|\boldsymbol{s}_{k})=\texttt{LogNormal}\big{(}\mu_{t}(\boldsymbol{s}_{k}),\sigma^{2}_{t}(\boldsymbol{s}_{k})\big{)},$

with
$[\mu_{t}(\boldsymbol{s}_{k}),\sigma^{2}_{t}(\boldsymbol{s}_{k})]=[\boldsymbol{W}_{\mu}\boldsymbol{s}_{k}+\boldsymbol{\mu}_{t},\boldsymbol{W}_{\sigma^{2}}\boldsymbol{s}_{k}+\boldsymbol{\sigma^{2}}_{t}]$
denote the mean and variance of the time distribution. Such a formulation
reduces model complexity, facilitates faster training and sampling in a
closed-form (Shchur et al., 2020).

To predict the time of the next checkin, we sample the probable time
difference between the current and the next checkin as
$\Delta_{t,k+1}\thicksim\mathcal{L}\mathcal{N}_{t}\big{(}\mu_{t}(\boldsymbol{s}_{k}),\sigma^{2}_{t}(\boldsymbol{s}_{k})\big{)}$,
where $\mathcal{L}\mathcal{N}_{t}$ denotes the learned log-normal parameters.
The actual time of the next checkin is the sum of the sampled time difference
and the current checkin time, $\widehat{t_{k+1}}=t_{k}+\Delta_{t,k+1}$.
Similar to the temporal flow, we also model the inter-checkin density of
spatial distances using a log-normal denoted as
$p_{d}(\Delta_{d,k+1}|\boldsymbol{s}_{k})$. We interpret this distribution as
the spatial flow for a region.

The inter-location spatial distance plays a crucial role in determining the
next POI (Cho et al., 2011; Manotumruksa et al., 2018). Unlike time, the
distances between two checkin locations are unchanged throughout the data.
Previous approaches (Du et al., 2016; Mei and Eisner, 2017) ignore these
spatial features and rely solely on the past checkin-categories. Moreover, in
a sequential setting the distance that the user will travel for her next
checkin is not known. Our MTPPs, being generative models, and spatial flows
overcome this drawback as we can sample the probable travel distance for the
next checkin from the spatial flow as
$\Delta_{d,k+1}\thicksim\mathcal{L}\mathcal{N}_{d}\big{(}\mu_{d}(\boldsymbol{s}_{k}),\sigma^{2}_{d}(\boldsymbol{s}_{k})\big{)}$.
Then, for predicting the next checkin, we use the sampled distance
$\Delta_{d,k+1}$ and RNN hidden state $\boldsymbol{s}_{k}$ via an attention
weighted embedding (Bahdanau et al., 2015).

(4)
$\boldsymbol{s}^{*}_{k}=\boldsymbol{s}_{k}+\alpha\cdot\boldsymbol{w}_{f}\Delta_{d,k+1},$

where $\alpha,\boldsymbol{w}_{f}$ denote the attention weight, a trainable
parameter and $\boldsymbol{s}^{*}_{k}$ denotes the updated hidden state. We
then predict the next checkin category as:

(5)
$\mathbb{P}(c_{k+1}=c|\boldsymbol{s}^{*}_{k})=\frac{\exp(\boldsymbol{V}_{s,c}\boldsymbol{s}^{*}_{k}+\boldsymbol{b}_{s,c})}{\sum_{\forall
c^{\prime}\in\mathcal{C}}\exp(\boldsymbol{V}_{s,c^{\prime}}\boldsymbol{s}^{*}_{k}+\boldsymbol{b}_{s,c^{\prime}})},$

where $\boldsymbol{V}_{s,\bullet},\boldsymbol{b}_{s,\bullet}$ are trainable
parameters and $\bullet$ denotes the entry corresponding to a category.
$\mathbb{P}(c_{k+1}=c|\boldsymbol{s}^{*}_{k})$ denotes the probability of next
checkin being of category $c$ with $c\in\mathcal{C}$.

Optimization: Given the set of all sequences $\mathcal{S}$ for a region
$\mathcal{R}$, we maximize the joint likelihood for the next checkin, the log-
normal density distribution of spatial, and temporal normalizing flows.

(6)
$\mathscr{L}=\sum_{\forall\mathcal{S}}\sum_{k=1}^{|\mathcal{S}|}\log\big{(}\mathbb{P}(c_{k+1}|\boldsymbol{s}^{*}_{k})\cdot
p_{t}(\Delta_{t,k+1}|\boldsymbol{s}_{k})\cdot
p_{d}(\Delta_{d,k+1}|\boldsymbol{s}_{k})\big{)}.$

where $\mathscr{L}$ denotes the joint likelihood, which we represent as the
sum of the likelihoods for all user sequences. We learn the parameters of
Reformd using Adam (Kingma and Ba, 2015) optimizer.

###  4.2. Flow-based Transfer

For transferring the mobility parameters across the regions, we follow the
standard transfer learning procedure (Pan and Yang, 2010; Szegedy et al.,
2015) of training exclusively on the origin region and then fine-tuning for
the target region. However, the affinity of a user towards a POI evolves with
time (Cho et al., 2011; Likhyani et al., 2017). For example, a POI with
frequent user-checkins during the summer season might not be an attractive
option in winters. We include these insights by training multiple independent
normalizing flows each for spatial and temporal densities. Specifically, we
cluster the checkin sequences in the origin region into $\mathcal{M}$ equal
clusters based on the median of the occurrence times for all the checkins.
Later, for each cluster of sequences, we train spatial and temporal flows
independently. In this setting our net likelihood changes as to include the
sum of all $\mathcal{M}$ likelihoods,
$\mathscr{L}=\sum_{i=1}^{\mathcal{M}}\mathscr{L}_{i}$, where $\mathscr{L}_{i}$
is the joint likelihood for trajectories in clusters $\mathcal{M}$.

As in origin region, we divide the user trajectories in the target region as
well into $\mathcal{M}$ clusters and for trajectories in target-cluster
$m^{\text{tgt}}_{i}$ we attentively factor the spatial and temporal flows
corresponding to origin-cluster $m^{\text{org}}_{i}$. Mathematically, for the
temporal flows in the target region our density of arrival times changes to:

(7)
$[\mu_{t}(\boldsymbol{s}_{k}),\sigma^{2}_{t}(\boldsymbol{s}_{k})]^{\text{tgt}}=[\boldsymbol{W}_{\mu}\boldsymbol{s}_{k}+\boldsymbol{b}_{\mu}+\phi_{t}\boldsymbol{\mu}^{\text{org}}_{t},\boldsymbol{W}_{\sigma}\boldsymbol{s}_{k}+\boldsymbol{b}_{\sigma}+\phi_{t}\boldsymbol{\sigma}^{\text{org}}_{t}],$

where
$\boldsymbol{s}_{k},\phi_{t},\boldsymbol{\mu}^{\text{org}}_{t},\boldsymbol{\sigma}^{\text{org}}_{t}$
denotes the hidden state representation for the target region, attention
parameter for temporal flow and the learned flow parameters of mean and
variance for cluster $m_{i}$ in origin region. Similarly, our spatial flows
for target region include the origin flow parameters with an attention
parameter $\phi_{d}$. For faster convergence, we share $\phi_{t}$ and
$\phi_{d}$ across all $\mathcal{M}$. Other model components are same as in
origin region and we maximize the joint likelihood for target region as in
Equation 6.

We highlight that our choice to divide the sequences based on median of user
trajectories rather than the individual checkin locations is driven by the
following technical point: in the latter case, the _net flow_ –be it spatial
or temporal– would be the sum of lognormal flows for each set. Such a
formulation is undesirable since the result is neither _closed_ nor does it
remain a lognormal (Daniel, 2004), thus requiring involved techniques to
approximate them (Barouch et al., 1986), which we would like to explore in
future work. However with the current formulation, we can learn the parameters
of different flows independently.

##  5\. Evaluation

In this section, we conduct an empirical evaluation of Reformd. Specifically,
we address the following research questions.

  * RQ1

Can Reformd outperform state-of-the-art baselines for time and checkin
prediction?

  * RQ2

What is the advantage of transferring via normalizing flows?

  * RQ3

Can we extend Reformd for non-spatial datasets?

For evaluating mobility prediction, we consider six POI datasets from the U.S.
and Japan. All our models are implemented in Tensorflow on an NVIDIA Tesla
V100 GPU and are made public at https://github.com/data-iitd/reformd.

Table 1. Statistics of datasets used in our experiments. The origin region
columns are followed by target regions.

Property | NY | MI | NE | VI | TY | AI | CH | SA
---|---|---|---|---|---|---|---|---
#Users or #Sequences ($|\mathcal{S}|$) | 25.6k | 6.7k | 4.1k | 6.5k | 32.1k | 10.9k | 7.5k | 11.4k
#Categories ($|\mathcal{C}|$) | 403 | 364 | 311 | 357 | 376 | 319 | 286 | 289
Mean Length ($\mu_{|\mathcal{S}|}$) | 57.17 | 66.21 | 48.56 | 56.33 | 61.72 | 56.60 | 63.60 | 53.08

###  5.1. Experimental Settings

Dataset Description: We use POI data from Foursquare (Yang et al., 2019) in
United States(US) and Japan(JP) and for each country we construct 4 datasets:
one with large check-in data and three with limited data. The statistics of
all datasets is given in Table 1 with each acronym denoting the following
region: (i) NY: New York(US), (ii) MI: Michigan(US), (iii) NV: Nevada(US),
(iv) VI: Virginia(US), (v) TY: Tokyo(JP), (vi) CH: Chiba(JP), (vii) SA:
Saitama(JP) and (viii) AI: Aichi(JP). We consider NY and TY as the origin
regions and MI, NV, VI and CH, SA, AI as the corresponding target regions. For
each region, we consider the time of checkin and category as event time and
mark and normalize the times based on the minimum and maximum event times. We
set the embedding and RNN hidden dimension to $64$ and $\mathcal{M}=3$ for all
our experiments. Other values for the model parameter had negligible
differences.

Evaluation Protocol: We split each stream of say $N$ checkins
$\mathcal{S}_{N}$ into training and test set, where the training set (test
set) consists of first 80% (last 20%) checkin. We evaluate models using
standard metrics (Du et al., 2016) of (i) mean absolute error (MAE) of
predicted and actual checkin times,
$\frac{1}{|\mathcal{S}|}\sum_{e_{i}\in\mathcal{S}}[|t_{i}-\widehat{t}_{i}|]$
and (ii) mark (checkin category) prediction accuracy (MPA), i.e,
$\frac{1}{|\mathcal{S}|}\sum_{e_{i}\in\mathcal{S}}\\#(c_{i}=\widehat{c}_{i})$.
Here $\widehat{t_{i}}$ and $\widehat{c_{i}}$ are the predicted time and
category of the $i$-th checkin. Moreover, the clustering of sequences into
different sets is done based solely on the training data and using these
thresholds we assign clusters to sequences in the test data.

Table 2. Performance of all the methods in terms of mark prediction accuracy
(MPA) and mean absolute error (MAE) across all datasets. Results marked
$\dagger$ are statistically significant (i.e. two-sided $t$-test with $p\leq
0.1$) over the best baseline.

| Mark Prediction Accuracy (MPA) | Mean Absolute Error (MAE)
---|---|---
$\mathcal{R}^{\text{org}}\rightarrow\mathcal{R}^{\text{tgt}}$ | NY $\rightarrow$ MI | NY $\rightarrow$ NE | NY $\rightarrow$ VI | TY $\rightarrow$ AI | TY $\rightarrow$ CH | TY $\rightarrow$ SA | NY $\rightarrow$ MI | NY $\rightarrow$ NE | NY $\rightarrow$ VI | TY $\rightarrow$ AI | TY $\rightarrow$ CH | TY $\rightarrow$ SA
NHP (Mei and Eisner, 2017) | 0.1745 | 0.1672 | 0.1348 | 0.2162 | 0.4073 | 0.3820 | 0.0920 | 0.1710 | 0.1482 | 0.1146 | 0.1217 | 0.1288
RMTPP (Du et al., 2016) | 0.1761 | 0.1684 | 0.1377 | 0.2293 | 0.4250 | 0.4036 | 0.0817 | 0.1581 | 0.1360 | 0.1058 | 0.1162 | 0.1205
SAHP (Zhang et al., 2020) | 0.1587 | 0.1529 | 0.1303 | 0.1968 | 0.3864 | 0.3943 | 0.1132 | 0.1958 | 0.1705 | 0.1671 | 0.1696 | 0.1574
THP (Zuo et al., 2020) | 0.1793 | 0.1545 | 0.1493 | 0.2361 | 0.4229 | 0.4057 | 0.0983 | 0.1735 | 0.1652 | 0.1445 | 0.1426 | 0.1468
Reformd | 0.2159$\dagger$ | 0.1868$\dagger$ | 0.1631$\dagger$ | 0.2588$\dagger$ | 0.4474 | 0.4208 | 0.0672$\dagger$ | 0.1317$\dagger$ | 0.1089$\dagger$ | 0.0858$\dagger$ | 0.0897$\dagger$ | 0.0973$\dagger$
$\Delta$ (%) | 20.41 | 10.92 | 9.24 | 9.61 | 5.27 | 3.72 | 17.74 | 16.69 | 19.92 | 18.90 | 22.80 | 19.25

###  5.2. Baselines

We compare the prediction performance of Reformd with the following state-of-
the-art methods:

  * RMTPP (Du et al., 2016): A recurrent neural network that models time-differences to learn a representation of the past events.

  * NHP (Mei and Eisner, 2017): Models an MTPP using continuous-time LSTMs for capturing the temporal evolution of sequences.

  * SAHP (Zhang et al., 2020): A self-attention model to learn the temporal dynamics using an aggregation of historical events.

  * THP (Zuo et al., 2020): Extends the transformer model (Vaswani et al., 2017) to include the conditional intensity of event arrival and the inter-mark influences.

We omit comparison with other continuous-time models (Omi et al., 2019; Shchur
et al., 2020; Xiao et al., 2017a, b; Hawkes, 1971) as they have already been
outperformed by these approaches.

Table 3. Prediction performance of all the methods for product recommendation
in Amazon datasets. Results marked $\dagger$ are statistically significant as
in Table 2.

| Mark Prediction Accuracy | Mean Absolute Error
---|---|---
$\mathcal{R}^{\text{org}}\rightarrow\mathcal{R}^{\text{tgt}}$ | DM $\rightarrow$ AP | DM $\rightarrow$ BY | DM $\rightarrow$ AP | DM $\rightarrow$ BY
NHP (Mei and Eisner, 2017) | 0.8773 | 0.5711 | 0.0903 | 0.1795
RMTPP (Du et al., 2016) | 0.8975 | 0.5530 | 0.0884 | 0.1758
SAHP (Zhang et al., 2020) | 0.8931 | 0.5517 | 0.1439 | 0.2214
THP (Zuo et al., 2020) | 0.9084 | 0.5879 | 0.1253 | 0.2035
Reformd | 0.9129 | 0.6035 | 0.0756$\dagger$ | 0.1564$\dagger$
$\Delta$ (%) | 0.49 | 2.65 | 14.47 | 11.03

(a) Virginia

(b) Aichi

Figure 2.  Real life true and predicted inter-arrival times $\Delta_{t,k}$ of
different events $e_{k}$ for (a) Virginia and (b) Aichi.

###  5.3. Prediction Performance (RQ1)

We report the prediction performance of different methods across our target
datasets in Table 2 and make the following observations:

  * $\bullet$

Reformd consistently yields the best performance on all the datasets. In
particular, it improves over the strongest baselines by 10% and 19% for
category and time prediction respectively. These results indicate the
importance of spatial and temporal flow-based transfer from external data for
prediction in limited-data regions.

  * $\bullet$

RMTPP (Du et al., 2016) is the second-best performer in terms of MAE of time
prediction almost for all the datasets. For some datasets THP (Zuo et al.,
2020) outperforms RMTPP for mark category prediction. However Reformd
significantly outperforms it across all metrics.

Qualitative Analysis: We also perform a qualitative analysis to demonstrate
how Reformd is able to model the checkin time distribution. For this, we plot
the actual inter-checkin time differences and the difference time predicted by
Reformd in Figure 2 for Virginia and Aichi datasets. From the results we note
that the predicted inter-arrival times closely match with the true inter-
arrival times and Reformd is even able to capture large time differences
(peaks). For brevity, we omit the results for other datasets.

Runtime: For all datasets, the times for training on the origin and later on
target regions is within 3 hours, thus are within the range for practical
deployment.

###  5.4. Transfer Advantage (RQ2)

Reformd outperforms other baselines and also brings exhibits a key feature of
transfer learning, i.e. quick parameter learning (Szegedy et al., 2015; Pan
and Yang, 2010). We highlight this by plotting the time prediction error (MAE)
corresponding to the epochs trained on the target region for Reformd and the
best time prediction model, i.e. RMTPP. Figure 3 summarizes the results for
Virginia and Aichi. We note that Reformd exhibits faster convergence than
RMTPP for both datasets. More specifically, the flow-based transfer procedure
of Reformd can outperform most baselines even with a fine-tuning of a few
epochs. The results also highlight the stable learning procedure of Reformd.

###  5.5. Product Recommendation (RQ3)

We further evaluate the performance of Reformd in product recommendation, i.e.
without spatial coordinates. Consequently, we use purchase records for three
item categories from Amazon (Ni et al., 2019), namely Digital Music(DM),
Appliances(AP) and Beauty(BY). For each item we use the user reviews as the
events in a sequence with the time of the written review as the event time and
the rating (1 to 5) as the corresponding mark. As in this case, we do not have
a spatial density function $p_{d}(\Delta_{d,k+1})$, we change the fusion
equation 4 to include the predicted time of next purchase as:

(8)
$\boldsymbol{s}^{*}_{k}=\boldsymbol{s}_{k}+\alpha\cdot\boldsymbol{w}_{f}\Delta_{t,k+1},$

We consider Digital Music($|\mathcal{S}|=12k$) as origin and
Appliances($|\mathcal{S}|=7k$) and Beauty($|\mathcal{S}|=6k$) as target. From
the results in Table 3, we note that even in the absence of spatial flows,
Reformd outperforms other baselines across all metrics.

(a) Virginia

(b) Aichi

Figure 3.  Training curves of Reformd and RMTPP for time prediction with best
MAE for (a) Virginia and (b) Aichi.

##  6\. Conclusion and Future Work

In this paper, we present Reformd, a novel method for transferring mobility
knowledge across regions by sharing the spatial and temporal NFs for
continuous-time checkin prediction. As a future work, we plan to incorporate
meta-learning based transfer (Finn et al., 2017).

## References

  * (1)
  * Bacry et al. (2015) Emmanuel Bacry, Iacopo Mastromatteo, and Jean-François Muzy. 2015\.  Hawkes processes in finance. In _Market Microstructure and Liquidity_.
  * Bahdanau et al. (2015) Dzmitry Bahdanau, Kyunghyun Cho, and Yoshua Bengio. 2015\.  Neural machine translation by jointly learning to align and translate. In _ICLR_.
  * Barkhuus and Dey (2003) Louise Barkhuus and Anind K Dey. 2003.  Location-Based Services for Mobile Telephony: a Study of Users Privacy Concerns. In _Interact_.
  * Barouch et al. (1986) E Barouch, GM Kaufman, and ML Glasser. 1986.  On sums of lognormal random variables.  _Studies in Applied Mathematics_ 75, 1 (1986), 37–55.
  * Chen et al. (2020) Yile Chen, Cheng Long, Gao Cong, and Chenliang Li. 2020\.  Context-aware Deep Model for Joint Mobility and Time Prediction. In _WSDM_.
  * Cheng et al. (2013) Chen Cheng, Haiqin Yang, Michael R Lyu, and Irwin King. 2013\.  Where you like to go next: Successive point-of-interest recommendation. In _IJCAI_.
  * Cho et al. (2011) Eunjoon Cho, Seth A. Myers, and Jure Leskovec. 2011\.  Friendship and mobility: User movement in location-based social networks. In _KDD_.
  * Daley and Vere-Jones (2007) Daryl J Daley and David Vere-Jones. 2007.  _An introduction to the theory of point processes: volume II: general theory and structure_.  Springer Science & Business Media.
  * Daniel (2004) Dufresne Daniel. 2004\.  The log-normal approximation in financial and other computations. In _Advances in Applied Probability_.
  * Du et al. (2016) Nan Du, Hanjun Dai, Rakshit Trivedi, Utkarsh Upadhyay, Manuel Gomez-Rodriguez, and Le Song. 2016\.  Recurrent marked temporal point processes: Embedding event history to vector. In _KDD_.
  * Farajtabar et al. (2017) Mehrdad Farajtabar, Jiachen Yang, Xiaojing Ye, Huan Xu, Rakshit Trivedi, Elias Khalil, Shuang Li, Le Song, and Hongyuan Zha. 2017\.  Fake news mitigation via point process based intervention. In _ICML_.
  * Feng et al. (2018) Jie Feng, Yong Li, Chao Zhang, Funing Sun, Fanchao Meng, Ang Guo, and Depeng Jin. 2018.  DeepMove: Predicting Human Mobility with Attentional Recurrent Networks. In _WWW_.
  * Figueiredo et al. (2016) Flavio Figueiredo, Bruno Ribeiro, Jussara M Almeida, and Christos Faloutsos. 2016.  TribeFlow: mining & predicting user trajectories. In _WWW_.
  * Finn et al. (2017) Chelsea Finn, Pieter Abbeel, and Sergey Levine. 2017\.  Model-agnostic meta-learning for fast adaptation of deep networks. In _ICML_.
  * Gupta et al. (2021) Vinayak Gupta, Srikanta Bedathur, Sourangshu Bhattacharya, and Abir De. 2021.  Learning Temporal Point Processes with Intermittent Observations. In _AISTATS_.
  * Hawkes (1971) Alan G Hawkes. 1971\.  Spectra of some self-exciting and mutually exciting point processes. In _Biometrika_.
  * Kingma and Ba (2015) Diederik P. Kingma and Jimmy Ba. 2015.  Adam: A Method for Stochastic Optimization. In _ICLR_.
  * Likhyani et al. (2017) Ankita Likhyani, Srikanta Bedathur, and Deepak P. 2017\.  LoCaTe: Influence Quantification for Location Promotion in Location-based Social Networks. In _IJCAI_.
  * Likhyani et al. (2020) Ankita Likhyani, Vinayak Gupta, PK Srijith, P Deepak, and Srikanta Bedathur. 2020.  Modeling Implicit Communities from Geo-tagged Event Traces using Spatio-Temporal Point Processes. In _WISE_.
  * Liu et al. (2016) Qiang Liu, Shu Wu, Liang Wang, and Tieniu Tan. 2016\.  Predicting the Next Location: A Recurrent Model with Spatial and Temporal Contexts. In _AAAI_.
  * Manotumruksa et al. (2018) Jarana Manotumruksa, Craig Macdonald, and Iadh Ounis. 2018\.  A Contextual Attention Recurrent Architecture for Context-Aware Venue Recommendation. In _SIGIR_.
  * McLachlan and Peel (2004) Geoffrey McLachlan and David Peel. 2004.  _Finite mixture models_.  John Wiley & Sons.
  * Mei and Eisner (2017) Hongyuan Mei and Jason M Eisner. 2017.  The neural hawkes process: A neurally self-modulating multivariate point process. In _NeurIPS_.
  * Metev (2021) Denis Metev. 2021\.  41+ Must Know Foursquare Statistics in 2020.  Available At: https://review42.com/resources/foursquare-statistics/.
  * Ni et al. (2019) Jianmo Ni, Jiacheng Li, and Julian McAuley. 2019.  Justifying recommendations using distantly-labeled reviews and fine-grained aspects. In _EMNLP-IJCNLP_.
  * Omi et al. (2019) Takahiro Omi, Naonori Ueda, and Kazuyuki Aihara. 2019\.  Fully Neural Network based Model for General Temporal Point Processes. In _NeurIPS_.
  * Pan and Yang (2010) Sinno Jialin Pan and Qiang Yang. 2010.  A Survey on Transfer Learning. In _IEEE TKDE_.
  * Rezende and Mohamed (2015) Danilo Rezende and Shakir Mohamed. 2015.  Variational inference with normalizing flows. In _ICML_.
  * Rizoiu et al. (2018) Marian-Andrei Rizoiu, Swapnil Mishra, Quyu Kong, Mark Carman, and Lexing Xie. 2018\.  SIR-Hawkes: on the Relationship Between Epidemic Models and Hawkes Point Processes. In _WWW_.
  * Shchur et al. (2020) Oleksandr Shchur, Marin Biloš, and Stephan Günnemann. 2020\.  Intensity-Free Learning of Temporal Point Processes. In _ICLR_.
  * Szegedy et al. (2015) C. Szegedy, W. Liu, Y. Jia, P. Sermanet, S. Reed, D. Anguelov, D. Erhan, V. Vanhoucke, and A. Rabinovich. 2015\.  Going deeper with convolutions. In _CVPR_.
  * Vaswani et al. (2017) Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin. 2017.  Attention is all you need. In _NeurIPS_.
  * Wei et al. (2016) Ying Wei, Yu Zheng, and Qiang Yang. 2016.  Transfer knowledge between cities. In _KDD_.
  * Xiao et al. (2017a) Shuai Xiao, Mehrdad Farajtabar, Xiaojing Ye, Junchi Yan, Le Song, and Hongyuan Zha. 2017a.  Wasserstein Learning of Deep Generative Point Process Models. In _NeurIPS_.
  * Xiao et al. (2017b) Shuai Xiao, Junchi Yan, Xiaokang Yang, Hongyuan Zha, and Stephen M Chu. 2017b.  Modeling the Intensity Function of Point Process Via Recurrent Neural Networks. In _AAAI_.
  * Yang et al. (2019) Dingqi Yang, Bingqing Qu, Jie Yang, and Philippe Cudre-Mauroux. 2019\.  Revisiting User Mobility and Social Relationships in LBSNs: A Hypergraph Embedding Approach. In _WWW_.
  * Yang et al. (2018) Guolei Yang, Ying Cai, and Chandan K. Reddy. 2018.  Recurrent Spatio-Temporal Point Process for Check-in Time Prediction. In _CIKM_.
  * Yuan et al. (2013) Quan Yuan, Gao Cong, Zongyang Ma, Aixin Sun, and Nadia Magnenat-Thalmann. 2013.  Time-aware point-of-interest recommendation. In _SIGIR_.
  * Zhang et al. (2020) Qiang Zhang, Aldo Lipani, Omer Kirnap, and Emine Yilmaz. 2020\.  Self-attentive Hawkes processes. In _ICML_.
  * Zuo et al. (2020) Simiao Zuo, Haoming Jiang, Zichong Li, Tuo Zhao, and Hongyuan Zha. 2020.  Transformer Hawkes Process. In _ICML_.