# Agriculture-Vision Challenge 2022 – The Runner-Up Solution for Agricultural
Pattern Recognition via Transformer-based Models

Zhicheng Yang
PAII Inc.
Palo Alto, CA, USA
<EMAIL_ADDRESS>Jui-Hsin Lai
PAII Inc.
Palo Alto, CA, USA
<EMAIL_ADDRESS>Jun Zhou
Ping An Bank Co., Ltd
Shenzhen, Guangdong, China
<EMAIL_ADDRESS>Hang Zhou
PAII Inc.
Palo Alto, CA, USA
<EMAIL_ADDRESS>Chen Du
PAII Inc.
Palo Alto, CA, USA
<EMAIL_ADDRESS>Zhongcheng Lai
Ping An Bank Co., Ltd
Shenzhen, Guangdong, China
<EMAIL_ADDRESS>
###### Abstract

The Agriculture-Vision Challenge in CVPR is one of the most famous and
competitive challenges for global researchers to break the boundary between
computer vision and agriculture sectors, aiming at agricultural pattern
recognition from aerial images. In this paper, we propose our solution to the
third Agriculture-Vision Challenge in CVPR 2022. We leverage a data pre-
processing scheme and several Transformer-based models as well as data
augmentation techniques to achieve a mIoU of 0.582, accomplishing the 2nd
place in this challenge.

##  1 Introduction

Computer vision applications in agricultural domain has become one of hot
topics nowadays, especially using remote sensing satellite images and aerial
images. With the rapid development of deep learning methods, numerous research
studies have proposed pioneer and practical solutions to various computer
vision problems in agriculture [2, 9, 4, 3]. Aside from fruitful research
achievements, various algorithm challenges have been held at top-tier
conferences for global researchers in recent years, in order to explore more
effective algorithms to solve the specific problems. The _Agriculture-Vision
Challenge_ in CVPR since 2020 is one of most famous and competitive challenges
in this interdisciplinarity study. It aims at applying computer vision
algorithms to agricultural pattern recognition from high-resolution aerial
images. This year, CVPR 2022 holds the 3rd Agriculture-Vision Challenge, and
we form our team “PAII-RS” to participate in this contest.

##  2 Materials and Methods

In this section, we elaborate on the given datasets, the pre-processing
method, the proposed deep learning-based framework, and the test-time
augmentation (TTA) strategy.

###  2.1 Description of Dataset

Table 1: Information of the given and resampled datasets for training and validation categories. Class Index | Class Name | Original Amount (Train/Val) | Resampled Amount (Train/Val)
---|---|---|---
0 | Background | 56944 / 18334 | 75121 / 13642
1 | Double Plant | 6234 / 2322 | 10961 / 2294
2 | Drydown | 16806 / 5800 | 19320 / 3383
3 | Endrow | 4481 / 1755 | 8544 / 1858
4 | Nutrient Deficiency | 13308 / 3883 | 14859 / 2610
5 | Planter Skip | 2599 / 1197 | 5361 / 1015
6 | Water | 2155 / 987 | 4132 / 721
7 | Waterway | 3899 / 696 | 6024 / 1109
8 | Weed Cluster | 11111 / 2834 | 14423 / 2773
Figure 1: Framework of our deep learning-based model.

The challenge this year provides the entire Agriculture-Vision dataset
released in [1]. It contains 94,986 aerial farmland images collected
throughout 2019 across the U.S. Each image has a size of 512$\times$512 pixels
and has 4 channels (RGB and NIR). A total of 9 label classes are manually
labeled for every image. Table 1 shows the given amount of images in each
class. Note that many images have multiple labels, and even have overlapped
labels (one pixel has multiple labels).

Although the amount of the given training data is considerable, we still
generate more data following the data augmentation scheme of the winner
solution last year. They conducted an image mosaic scheme to enable the model
to have multi-scale views during the training. To fit the model input size, we
create two new datasets using mosaicked images with down-sampling 2X (2 times)
and down-sampling 3X as shown in Fig. 1. The down-sampling dataset has the
same image size of 512$\times$512 pixels that the recognition model can share
the same network architecture among 1X, 2X, and 3X imagery.

###  2.2 Data Pre-Processing

We observe that the image counts in each category are uneven. For example, the
image count of the background class is 25 times larger than the water class.
To tackle the unbalance issue, we try to sample more images in the few-shot
classes. The re-sampled image counts are listed in Table 1.

###  2.3 Framework

Fig. 1 shows our deep learning-based framework. _SegFormer_ is a Transformer-
based efficient segmentation model [7]. It designs a hierarchical Transformer
encoder with multi-level feature outputs. Unlike other cumbersome decoders,
SegFormer’s decoder adopts MLP layers to aggregate multi-scale feature outputs
from different layers. One of the key advantages of SegFormer is that its
model size is relatively small but the performance keeps outstanding.
Therefore, SegFormer is suitable for this challenge due to the model size
parameter limit of 150M.

SegFormer provides six versions with various settings of Transformer encoders,
leading to different model sizes. These six models are named from B0 to B5,
with the increased model size. To follow the policy, we select Mix Transformer
(MiT) B3 and Mix Transformer B2 as our training models. Their model size
information can be found in Table 7 “Mix Transformer Encoder” in [7]. After
obtaining the individual inference result from each model, the model ensemble
is performed to predict the final segmentation results.

###  2.4 Test-Time Augmentation

Since our models are trained with 1X, 2X, and 3X down-sampling imagery, we
conduct the same processing on the test dataset. In addition to the scale
augmentation, we include image rotation and flip.

##  3 Results

###  3.1 Evaluation Metric

The required evaluation metric is the average Intersection over Union metric
(mIoU), which is defined as Eq. 1 to measure the performance.

$\text{mIoU}=\frac{1}{c}\sum\frac{\text{ Area }\left(P_{c}\cap
T_{c}\right)}{\text{ Area }\left(P_{c}\cup T_{c}\right)}$ (1)

where $c$ is the number of label classes (8 foreground classes + 1 background
class for this challenge); $P_{c}$ and $T_{c}$ are the predicted label mask
and ground truth label mask of the class $c$, respectively.

Table 2: Performance comparisons among various models. The bold font of numeric results indicates the best performance on the test set. BG: Background; DP: Double Plant; D: Drydown; E: Endrow; ND: Nutrient Deficiency; PS: Planter Skip; W: Water; WW: Waterway; WC: Weed Cluster. The number in the parentheses following the class name refers to the class index. Models | mIoU | BG(0) | DP(1) | D(2) | E(3) | ND(4) | PS(5) | W(6) | WW(7) | WC(8)
---|---|---|---|---|---|---|---|---|---|---
(Other methods, on the val set) |  |  |  |  |  |  |  |  |  |
Agriculture-Vision baseline(RGBN) [1] | 0.434 | 0.743 | 0.285 | 0.574 | 0.217 | 0.389 | 0.336 | 0.736 | 0.344 | 0.283
MiT-B3(RGBN)[5] | 0.454 | 0.768 | 0.371 | 0.609 | 0.245 | 0.424 | 0.413 | 0.692 | 0.269 | 0.299
MiT-B5(RGB)[6] | 0.464 | 0.755 | 0.370 | 0.585 | 0.227 | 0.313 | 0.414 | 0.802 | 0.401 | 0.304
MiT-B5(RGBN)[6] | 0.490 | 0.762 | 0.373 | 0.618 | 0.246 | 0.428 | 0.420 | 0.813 | 0.437 | 0.318
(Our implementation, on the test set) |  |  |  |  |  |  |  |  |  |
HRNet-W48+OCR[8](RGB baseline) | 0.413 | 0.717 | 0.316 | 0.567 | 0.233 | 0.269 | 0.283 | 0.718 | 0.289 | 0.326
MiT-B3[7](RGB baseline) | 0.448 | 0.720 | 0.395 | 0.557 | 0.325 | 0.364 | 0.330 | 0.687 | 0.293 | 0.358
MiT-B2[7](RGBN+Our method) | 0.554 | 0.778 | 0.483 | 0.632 | 0.476 | 0.570 | 0.403 | 0.768 | 0.410 | 0.466
MiT-B3[7](RGBN+Our method) | 0.563 | 0.773 | 0.471 | 0.640 | 0.452 | 0.569 | 0.442 | 0.782 | 0.463 | 0.475
(Our implementation, on the test set) |  |  |  |  |  |  |  |  |  |
Model Ensemble(RGBN+Our method) | 0.582 | 0.777 | 0.485 | 0.646 | 0.481 | 0.573 | 0.471 | 0.779 | 0.547 | 0.479

###  3.2 Experiment Results

Table 2 presents our results, the baseline provided by the host Agriculture-
Vision organizers, and the results of other methods. Note that other baselines
evaluate their performance on the validation set due to the unavailable test
set. As we can see, while our single model baselines are competitive with
other baselines, our proposed method effectively improves the single model
performance. Even though some single models have peak performance in some
classes (0.778 for “Background” and 0.782 for “Water”), our model ensemble
enjoys the merits of multiple single models’ strength to achieve the mIoU of
0.582. It also shows that our ensemble results significantly outperform other
baselines and our implementation of various single models.

##  4 Conclusions

In this paper, we propose our solution to the 3rd Agriculture-Vision Challenge
in CVPR 2022. For data usage, we perform data pre-processing and test data
augmentation schemes. Several SegFormer models are leveraged. We finally
accomplish a mIoU of 0.582, achieving the 2nd place in this challenge.

Future Directions. The potential applications of our proposed algorithm
include crop type identification in precision agriculture, agricultural asset
estimation and agricultural insurance product design in the Environmental,
Social, and Governance (ESG) domain. These future directions can illuminate
the revitalization of rural areas and facilitate the service of inclusive
finance in an eco-friendly way.

## References

  * [1] Mang Tik Chiu, Xingqian Xu, Yunchao Wei, Zilong Huang, Alexander G Schwing, Robert Brunner, Hrant Khachatrian, Hovnatan Karapetyan, Ivan Dozier, Greg Rose, et al.  Agriculture-vision: A large aerial image database for agricultural pattern analysis.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 2828–2838, 2020.
  * [2] Fulin Huang, Zhicheng Yang, Hang Zhou, Chen Du, Andy JY Wong, Yuchuan Gou, Mei Han, and Jui-Hsin Lai.  Unsupervised superpixel-driven parcel segmentation of remote sensing images using graph convolutional network.  In Companion Proceedings of the Web Conference 2022, 2022.
  * [3] Andreas Kamilaris and Francesc X Prenafeta-Boldú.  Deep learning in agriculture: A survey.  Computers and electronics in agriculture, 147:70–90, 2018.
  * [4] Nan Qiao, Yi Zhao, Ruei-Sung Lin, Bo Gong, Zhongxiang Wu, Mei Han, and Jiashu Liu.  Generative-discriminative crop type identification using satellite images.  In 2019 IEEE Global Conference on Signal and Information Processing (GlobalSIP), pages 1–5. IEEE, 2019.
  * [5] Yao Shen, Lei Wang, and Yue Jin.  Aaformer: A multi-modal transformer network for aerial agricultural images.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 1705–1711, 2022.
  * [6] Antonio Tavera, Edoardo Arnaudo, Carlo Masone, and Barbara Caputo.  Augmentation invariance and adaptive sampling in semantic segmentation of agricultural aerial images.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 1656–1665, 2022.
  * [7] Enze Xie, Wenhai Wang, Zhiding Yu, Anima Anandkumar, Jose M Alvarez, and Ping Luo.  Segformer: Simple and efficient design for semantic segmentation with transformers.  Advances in Neural Information Processing Systems, 34, 2021.
  * [8] Yuhui Yuan, Xilin Chen, and Jingdong Wang.  Object-contextual representations for semantic segmentation.  In European conference on computer vision, pages 173–190. Springer, 2020.
  * [9] Liheng Zhong, Lina Hu, and Hang Zhou.  Deep learning based multi-temporal crop classification.  Remote sensing of environment, 221:430–443, 2019.


# Nonlocal Correlation Mediated by Weyl Orbits

Zhe Hou  International Center for Quantum Materials, School of Physics, Peking
University, Beijing 100871, China  Qing-Feng Sun<EMAIL_ADDRESS>International Center for Quantum Materials, School of Physics, Peking
University, Beijing 100871, China  Collaborative Innovation Center of Quantum
Matter, Beijing 100871, China  CAS Center for Excellence in Topological
Quantum Computation, University of Chinese Academy of Sciences, Beijing
100190, China

###### Abstract

Nonlocality is an interesting topic in quantum physics and is usually mediated
by some unique quantum states. Here we investigate a Weyl semimetal slab and
find an exotic nonlocal correlation effect when placing two potential wells
merely on the top and bottom surfaces. This correlation arises from the
peculiar Weyl orbit in Weyl semimetals and is a consequence of the bulk-
boundary correspondence in topological band theory. A giant nonlocal transport
signal and a body breakdown by Weyl fermions are further uncovered, which can
serve as signatures for verifying this nonlocal correlation effect
experimentally. Our results extend a new member in the nonlocality family and
have potential applications for designing new electric devices with fancy
functions.

_Introduction._ Nonlocality has always been an intriguing topic in quantum
physics. The most famous one is the quantum entanglement Einstein ; Horodecki
of two or more particles concerning their spin, momentum or position in the
entangled state. The quantum entanglement violates the Bell inequality Bell ;
Gisin , comes along with a spooky action at a distance Popescu , and has been
utilized for quantum dense encoding Bennett1  and information teleportation
Bennett2 . In solid state materials, nonlocality also exists in some zero mode
excitations whose properties mimic the fundamental particles in quantum field
theory, such as the Majorana zero modes in a $p$-wave superconductor chain
Kitaev ; Fu ; Lutchyn  and the Jackiw-Rebbi zero mode in a quantum spin Hall
insulator constriction Klinovajia ; Alicea ; Wu . Conceptually the zero modes
can be regarded as half an intact particle, and similar to the positive and
negative monopoles in a magnet, they have an intrinsic nonlocal correlation
and can only create or annihilate in pairs.

Basically in a nonlocal phenomenon, applying a perturbation at one site of the
system could arouse a remote response at the other site, which realizes a
“communication” between particles over long distances. Yet nonlocality of
fermion in three dimensional bulk materials are rarely reported. The main
reason stems from the electrostatic screening to the perturbation and the
overall distribution of the bulk wave-functions which make the spatial
response decay exponentially or proportionally to $r^{-3}$ with the distance
$r$ Fetter .

Recently a series of unconventional quantum oscillations in topological Dirac
semimetals $\rm Cd_{3}As_{2}$ were reported in experiments Moll ; Zhang1 ;
Zhang2 . In those systems, there exists a closed orbit in momentum space
connecting the surface and bulk states, namely the Weyl orbit. Under the
magnetic field, Weyl fermions can move along the Weyl orbit and tunnel between
the top and bottom surfaces through the bulk states. The bulk and surface
connection by the Weyl orbit implies that there may exist a deeper and
intrinsic nonlocality in topological semimetals.

Here we report a three dimensional nonlocal correlation effect in topological
Weyl semimetals (WSMs) Wan ; Burkov ; Xu1 ; Lv ; Huang . We show that by
putting two attractive potential wells (PWs) merely on the top and bottom
surfaces, a giant inter-surface correlation across the body would happen. The
PW in the WSM works as a particle beam deflector which injects the Weyl
fermions on the surface into the bulk or extract them from the bulk to the
surface. With both PWs, the surface Weyl fermions would first be pumped into
the bulk and then appear on the other surface, which realizes a remote
nonlocal correlation between the PWs. This phenomenon arises from the unique
band structure, i.e. the Weyl orbit, and is an intrinsic property in WSMs
because no externally global electric or magnetic fields are applied. For Weyl
fermions taking a round trip between the top/bottom surface and contributing
to a constructive interference, a body breakdown by a newly formed resonant
bridge state would happen, which can be a strong signature for the
experimental verification on the remote correlation effect.

_Results_ We consider a WSM slab which is infinitely large in $x$ and $z$
directions but has finite width ($N_{y}$ layers) in $y$ direction, as shown in
Fig. 1(a). The Hamiltonian of the WSM has a two band form Hou1 ; Yang ; Wang :
$\hat{H}=\sum_{\bf k}c_{\bf k}^{\dagger}H({\bf k})c_{\bf k}$ with $c_{\bf
k}=(c_{\bf k\uparrow},c_{\bf k\downarrow})^{T}$ and $H({\bf
k})=t_{z}\sigma_{z}(2+\gamma-\cos{k_{x}a}-\cos{k_{y}a}-\cos{k_{z}a})\sigma_{z}+t_{x}\sigma_{x}\sin{k_{x}a}+t_{y}\sigma_{y}\sin{k_{y}a}$.
Here $\bf k$ is the wave-vector, $a$ is the lattice constant and
$t_{i}(i=x,y,z)$ is the hopping energy of nearest sites in the $i$-direction.
The Pauli matrix $\sigma_{i}$ acts on the spin space. In this model, there
exist two Weyl nodes locating at ${\bf K}_{\pm}=(0,0,\pm k_{0})$ with
$\cos{k_{0}a}=\gamma$ in the Brillouin zone. The Weyl node plays the role of
magnetic monopoles with positive or negative chirality in the momentum space,
and due to the bulk-boundary correspondence, the chiral surface states are
formed in the WSM slab sketched in Fig. 1(a) and propagate in $\pm x$
direction in the bottom (top) surface, as shown with the yellow and green
arrows. With equal energies, the surface states cannot get closed themselves
in momentum space but have to transit into bulk states near the Weyl nodes,
which constitute a closed loop, i.e. the Weyl orbit in Fig. 1(b).

Figure 1: (Color online) Schematic diagram of the nonlocal inter-surface
correlation mediated by the Weyl orbit. (a) The real space trajectory of the
Weyl fermions in the WSM slab with an elliptical PW (the brown circle) on the
top and bottom surfaces. The slab is infinitely large in $x$ and $z$
directions but has a thickness $N_{y}a$ in $y$ direction. (b) The Weyl orbit
in the momentum space for the WSM slab in (a). Here the up and down line
segments denote the bottom and top surface states, respectively. The black
half circles denote the trivial bulk states which are connected with the
surface states near the Weyl nodes. The black arrows on the Weyl orbit
indicate the directions of the group velocity.

Next we show how the applied PWs induce a nonlocal inter-surface response
within the Weyl orbit regime. The dark brown circles in Fig. 1(a) indicate the
PWs which are on the same position of the top/bottom surfaces and yield an
elliptical form: $U({\bf
r})=\omega_{x}(x-x_{C})^{2}+\omega_{z}(z-z_{C})^{2}-\omega_{x}r_{0}^{2}$, if
$\omega_{x}(x-x_{C})^{2}+\omega_{z}(z-z_{C})^{2}\leq\omega_{x}r_{0}^{2}$ and
$y=0,N_{y}a$. For other positions the potential is assumed to be zero. Here
$(x_{C},z_{C}),\omega_{x,z}$ and $r_{0}$ are the central position, strength
and size of the PW, respectively. Similar to a magnetic field, the attracting
potential providing an electrostatic force ${\bf F}=-\nabla U(r)$ alters the
moving trajectory of Weyl fermions and forces them into closed orbits.
However, different from the electrons in normal 2D electron gas or the Dirac
fermions on the surface of 3D topological insulators Xu2 , Weyl fermions can
not be confined on the WSM surfaces due to the unidirectional propagation of
the topological surface states (see the yellow and green segments in Fig.
1(b)).

For the Weyl fermion residing on the top surface states with a group velocity
pointing in -$x$ direction, once it is trapped inside the quantum-dot like PW
(see Fig. 1(a)), the attractive electrostatic force tries to change its real
space trajectory (the top red curve) and meanwhile the -$z$ component pushes
it into the bulk states in the momentum space (the left black half circle in
Fig. 1(b)). The velocity is then tilted in the -$z$ direction and the real
space trajectory is deflected to the left. The wave-function which firstly
localizes on the top surface starts to permeate into the body region, and the
central position of the Weyl fermion moves down as shown with the dotted black
line in Fig. 1(a) (see Fig. S1 for the wave-function evolutions in the Weyl
orbits) SI .

If there exists only one PW on the top surface, the Weyl fermion is merely
pumped into the body region and stops moving downward to affect the bottom
surface. While if the bottom PW is placed as well, the electrostatic force
would continue driving the Weyl fermion already in the body and pumps them
into the bottom surface. The wave-function is then pushed from the bulk states
into the bottom surface states with a group velocity in $x$ direction (the
bottom red curve in Fig. 1(a)). For the Weyl fermion pumped from the bottom
surface into the top one, a similar statement can be made. Thus the Weyl
fermion accomplishes half a circular motion in both top/bottom surfaces under
the drive of PWs. This results in a nonlocal correlation which is formed
between the top and bottom PWs and mediated by the bulk states inside the Weyl
orbit.

Figure 2: (Color online) Nonlocal DOS response in WSM nanoribbons. (a) to (d)
Four infinitely long WSM nanoribbons (along $x$ direction) with no PWs (a),
top PW (b), bottom PW (c) and top/bottom PWs (d) on the surface. (e) The DOS
$\rho_{\rm no}$ and $\rho_{\rm b}$ of the center section (the bright region in
the subfigure) for the nanoribbon in (a) and (b). Here the site i counts from
the left to the right, and then from the bottom to the top. (f) The DOS
$\rho_{\rm t}$ and $\rho_{\rm tb}$ of the center section for the nanoribbon in
(c) and (d). (g) The DOS response $\rho_{\rm b}-\rho_{\rm no}$ and $\rho_{\rm
tb}-\rho_{\rm t}$. Throughout the paper we set:
$N_{y}=20,N_{z}=61,\gamma=0,t_{x}=t_{y}=t_{z}=t,r_{0}=10a,\omega_{x}=0.004t/a^{2},\omega_{z}=0.01t/a^{2}$.
The central position of the PWs is set on the center of the WSM slab surface.
For (e-g) we set the energy $E=0.2t$.

To verify the nonlocal correlation discussed above, we conceive four
infinitely long WSM nanoribbon systems as shown in Fig. 2(a-d). In Fig. 2(a) a
pure nanoribbon with no PWs is shown, but in Fig. 2(b, c, d) the PW is put on
the bottom, top and top/bottom surfaces, respectively. To eliminate that the
surface states bypass the side walls and reach the other surface due to their
helicity, we set the nanoribbons infinitely long in $x$ direction and restrict
them in $y$ and $z$ directions with $N_{y}$ and $N_{z}$ layers. Figure 2(e)
gives the DOS $\rho_{i}$ of the central section with and without the bottom PW
(see Supplementary Material). Here the energy $E$ is chosen to be $0.2t$ above
the Weyl node and the Fermi surface is a single Weyl orbit with no second bulk
states mixing in. As can be seen from the no-potential case (the blue curve
$\rho_{\rm no}$), the DOS dominates in the top/bottom surfaces, and decays
quickly into a stable distribution away about three layers. The nonzero
$\rho_{i}$ in the body arises from the few bulk states in the Weyl orbit.
After the bottom PW is placed, the DOS $\rho_{\rm b}$ (the red dashed curve in
Fig. 2(e)) shows a sharp peak in the central position of the PW and the peak
height decays exponentially with the increased layer. Away from the bottom PW,
the DOS $\rho_{\rm b}$ almost overlaps with $\rho_{\rm no}$, indicating that
the influence of the single PW is just local Hosur .

Figure 2(f) shows the DOS $\rho_{\rm t}$ and $\rho_{\rm tb}$ of the central
section with top and top/bottom PWs, respectively. The DOS $\rho_{\rm t}$ is
symmetric with $\rho_{\rm b}$ about the central site $i_{c}=(N_{z}N_{y})/2+1$
because the systems in Fig. 2(b, c) can be transformed into each other by a
combined operation of the time-reversal operator $\hat{T}$ and the two-fold
rotation $C_{2}^{x}$ along the $x$ axis. When the bottom PW is further placed,
in addition to a peak in the bottom PW, the DOS inside the top PW shows a
dramatic decrease. Meanwhile the DOS of the bottom PW in turn is also
influenced as compared with the peak height of $\rho_{\rm b}$ in Fig. 2(e).

In order to clearly show the nonlocal correlation, we define the DOS
differences $\rho_{\rm b}-\rho_{\rm no}$ and $\rho_{\rm tb}-\rho_{\rm t}$ as
the responses of the applied bottom PW and plot them in Fig. 2(g) (also see
Fig. S2). While without the top PW, $\rho_{\rm b}-\rho_{\rm no}$ is almost
zero when deviated from the center of the bottom PW by $3a$. However, when in
the existence of the top PW, the DOS in top surface has a giant response to
the bottom PW and $\rho_{\rm tb}-\rho_{\rm t}$ shows deep valley at the center
of the top PW. The bottom PW itself is mutually influenced as a feedback to
the top nonlocal response. The giant nonlocal response in the double PWs case
indicates that there does exist a remote nonlocal correlation between the top
and bottom PWs. Since the Weyl fermions can not take a side circle to reach
the opposite surface, and the coupling between the top and bottom surface
states is set to be weak and unable to induce such a giant nonlocal response
(see Supplementary Material), the only reason that can cause this remote
correlation is the Weyl orbit regime we discussed above.

Next we show that except for the DOS response, an enhanced nonlocal transport
in a six-terminal system could also be induced by the nonlocal correlation. On
the basis of Fig. 2(a, c, d), we additionally connect four electrodes on the
top and bottom surfaces as probing terminals Hou2 ; Datta . The electrode is
composed of two-dimensional normal metals which can be described by a tight-
binding Hamiltonian $H_{\rm
NM}=\sum_{i,\sigma}E_{0}a_{i\sigma}^{\dagger}a_{i\sigma}+t_{\rm
NM}\sum_{<ij>,\sigma}a_{i\sigma}^{\dagger}a_{j\sigma}$. Here
$a_{i\sigma}^{\dagger}(a_{i\sigma})$ is the creation (annihilation) operator
at site $i$ with spin $\sigma$, $E_{0}$ is the on-site energy, and $<ij>$
denotes the nearest sites with a hopping amplitude $t_{\rm NM}$. Figure 3(a)
shows the right side view of the six-terminal system refitted from Fig. 2(d).
The top and bottom red arrows denote the topological surface states, and the
yellow circles represent the PWs.

We first give the transmission coefficients $T_{m1}$ from terminal 1 to
terminal $m$ ($m$=2, 3, 4, front, and back) without the PWs in Fig. 3(b). For
the incident energy $E<0.21t$, the transmission coefficients $T_{41}$ and
$T_{b1}$ have large values while the other transmissions are quite small with
their probability being nearly zero. This indicates that almost all the Weyl
fermions propagate along the surfaces and the inter-surface transport can
hardly happen. Further increasing $E$, $T_{21}$ and $T_{f1}$ have small step-
like increases because the second bulk states above the Weyl orbit begin to
participate.

Figure 3: (Color online) Nonlocal transport in the six-terminal system. (a)
The six-terminal system composed of an infinitely long WSM nanoribbon and four
normal metal electrodes which are connected on the top and bottom surfaces.
The red arrows denote the propagation of the surface states, and the yellow
circles mark the positions of the surface PWs. The dashed line between the PWs
shows the nonlocal correlation. (b) The transmission coefficients $T_{i1}$ as
functions of energy $E$ in the six-terminal system without PWs. (c) The
nonlocal transmission coefficient $T_{21}$ for the six-terminal system with no
PWs, top PW, and top/bottom PWs, respectively. (d) The DOS profile of the
central section in an infinitely long WSM nanoribbon with top/bottom PWs. The
energy is set to $E=0.186t$.

In Fig. 3(c) we restrict the transport into the Weyl orbit regime by setting
the energy range to $[0,0.2t]$ and investigate the influence of the surface
PWs on the nonlocal transmission $T_{21}$. In the no-PWs case, $T_{21}$ has
small values in most of the energy regions (see the blue curve) except for
showing sharp peaks in the Van Hove singularities (VHS) of the quasi-1D sub-
bands DasSarma  which arises from the finite size effect of the nanoribbon.
Putting a PW on the top surface, $T_{21}$ increases slightly because few
electrons in the bulk states are scattered into terminal 2. Further adding a
PW on the bottom surface, the transmission coefficient $T_{21}$ shows a
dramatic enhancement in the whole energy range $[0,0.2t]$ (see the dotted
yellow curve) with its amplitude serval times larger than the single PW case.
This enhancement can not originate from the normal scattering of the bulk
states by two PWs and is a prominent sign for the nonlocal correlation between
the top/bottom surfaces. Once the correlation is established by virtue of the
Weyl orbit, Weyl fermions propagating in the bottom surface can be pumped into
the bulk and then run into the the top surface, which constructs a
transmission channel enhancing the nonlocal inter-surface transport. In Fig.
S3, we design an experimental device for detecting this nonlocal transport
signal.

At a critical energy $E_{c}=0.186t$, we note a remarkable conductance
(transmission) peak with its value about 7 times of the single PW case and 13
times of the no PWs case (see Fig. S3). Although the peak position is close to
the VHS of the nearest sub-band, its appearance can not be explained by the
large DOS from the VHS because the height and broadening are both larger than
the usual peaks of $E<E_{c}$. The conductance peak observed in a mesoscopic
system always implies the existence of a novel quantum state. To validate
this, we plot the DOS profile of the central section in the infinitely long
nanoribbon at energy $E_{c}$ in Fig. 3(d). The bright fringes in the top and
bottom layers of Fig. 3(d) indicate the surface states, and the bright spots
mark the positions of the two PWs. Between the two PWs, there exists a novel
bright DOS stripe with a bridge-like shape linking the top and bottom
surfaces. _This suggests a new resonant state formed by two spatially
separated confining potentials with the wave-function spreading to the non-
confined region._ We here name it the “Weyl bridge state”. Its appearance
arises from the constructive interference of Weyl fermions taking a round trip
between the top and bottom layers driven by the PWs. In Supplementary Material
we make a semi-classical analysis on the real space trajectory of the Weyl
bridge state (see Fig.S4).

In the Weyl bridge state the correlation between the PWs gets to its maximum
which contributes to the large nonlocal transmission. This bridge state is
quite similar to the worm-hole effect in a 3D topological insulator Rosenberg
and the 3D quantum Hall state in the topological semimetals Zhang3 . The most
important distinction between them is that the Weyl bridge state does not rely
on any globally applied magnetic or electric field, and is an intrinsic
generation of the exotic topological band structure of WSMs. So the
nonlocality discovered here should be unique in topological semimetals with a
Weyl orbit structure. Finally, we want to make more discussions on the Weyl
bridge state by comparing it with a natural phenomenon of lightning (Movie
S1): if there accumulate enough free charges in the cloud, putting a conductor
in the ground would induce an electric breakdown through the atmosphere. In
our WSM systems, the top/bottom PWs play the role of the cloud and the
conductor respectively, and the Weyl bridge state is just a breakdown of Weyl
fermions through the WSM body.

_Conclusion._ We have proven a nonlocal correlation effect in WSMs as a result
of the Weyl orbit. As a unique property of topological semimetals, this
nonlocal correlation is totally mediated by the quantum states residing on the
Weyl orbit and can be simply generated by two surface PWs. The movable PWs
allow a flexible and precise engineer on this nonlocal correlation and provide
a new way to investigate the surface-bulk relationship in topological
semimetals. The present result opens a new vision on the nonlocality in
quantum physics, and has guiding significance in designing new electric device
with fancy functions, such as the impurity detection and remote control in
mesoscopic and micro-scales.

_Acknowledgement_ Z. H. thanks helpful discussions with Hua Jiang. This work
was financially supported by National Key R and D Program of China (Grant No.
2017YFA0303301), NSF-China (Grant No. 11921005), the Strategic Priority
Research Program of Chinese Academy of Sciences (Grant No. DB28000000), and
Beijing Municipal Science & Technology Commission (Grant No.
Z181100004218001).

## References

  * (1) A. Einstein, B. Podolsk, and N. Rosen, Can quantum-mechanical description of physical reality be considered complete? Phys. Rev. 47, 777 (1935).
  * (2) R. Horodecki, P. Horodecki, M. Horodecki, and K. Horodecki, Quantum entanglement. Rev. Mod. Phys. 81, 865 (2009).
  * (3) J. S. Bell, On the Einstein Podolsky Rosen paradox. Physics 1, 195 (1964).
  * (4) N. Gisin, Bell s inequality holds for all non-product states. Phys. Lett. A 154, 201 (1991).
  * (5) S. Popescu, Nonlocality beyond quantum mechanics. Nat. Phys. 10, 264 (2014).
  * (6) C. H. Bennett and S. J. Wiesner, Communication via one- and two-particle operators on Einstein-Podolsky-Rosen States. Phys. Rev. Lett. 69, 2881 (1992).
  * (7) C. H. Bennett, G. Brassard, C. Crépeau, R. Jozsa, A. Peres, and W. K. Wootters, Teleporting an unknown quantum state via dual classical and Einstein-Podolsky-Rosen channels. Phys. Rev. Lett. 70, 1895 (1993).
  * (8) A.Y. Kitaev, Unpaired Majorana fermions in quantum wires. Phys. Usp. 44, 131 (2001).
  * (9) L. Fu and C.L. Kane, Superconducting proximity effect and Majorana fermions at the surface of a topological insulator. Phys. Rev. Lett. 100, 096407 (2008).
  * (10) R. M. Lutchyn, J. D. Sau, and S. Das Sarma, Majorana fermions and a topological phase transition in semiconductor-superconductor heterostructures. Phys. Rev. Lett. 105, 077001 (2010).
  * (11) J. Klinovajia and D. Loss, Fractional charge and spin states in topological insulator constrictions. Phys. Rev. B 92, 121410(R) (2015).
  * (12) J. Alicea, Y. Oreg, G. Refael, F. von Oppen, and M. P. A. Fisher, Non-Abelian statistics and topological quantum information processing in 1D wire networks. Nat. Phys. 7, 412 (2011).
  * (13) Y. Wu, H. Liu, J. Liu, H. Jiang, and X. C. Xie, Double-frequency Aharonov-Bohm effect and non-Abelian braiding properties of Jackiw-Rebbi zero-mode. arXiv 1901.06138 [cond-mat.mes.-hall] (18 January 2019).
  * (14) A. L. Fetter and J. D. Walecka, Screening in an electron gas in _Quantum Theory of Many-Particle Systems_ (McGraw-Hill Book Company, New York, 1971), pp. 175-180.
  * (15) P. J. W. Moll et al., Transport evidence for Fermi-arc-mediated chirality transfer in the Dirac semimetal Cd3As2. Nature 535, 266 (2016).
  * (16) C. Zhang et al., Evolution of Weyl orbit and quantum Hall effect in Dirac semimetal Cd3As. Nat. Commun. 8, 1272 (2017).
  * (17) C. Zhang et al., Quantum Hall effect based on Weyl orbits in $\rm Cd_{3}As_{2}$. Nature 565, 331 (2019).
  * (18) X. Wan et al., Topological semimetal and Fermi-arc surface states in the electronic structure of pyrochlore iridates. Phys. Rev. B 83, 205101 (2011).
  * (19) A. A. Burkov and L. Balents, Weyl semimetal in a topological insulator multilayer. Phys. Rev. Lett. 107, 127205 (2011).
  * (20) S.-Y. Xu et al., Discovery of a Weyl fermion semimetal and topological Fermi arcs. Science 349, 613 (2015).
  * (21) B. Q. Lv et al., Experimental discovery of Weyl semimetal TaAs. Phys. Rev. X 5, 031013 (2015).
  * (22) S.-M. Huang et al., A Weyl fermion semimetal with surface Fermi arcs in the transition metal monopnictide TaAs class. Nat. Commun. 6, 7373 (2015).
  * (23) Z. Hou, Y.-F. Zhou, N.-X. Yang, and Q.-F. Sun, Chirality-dependent electron transport in Weyl semimetal p-n-p junctions. Commun. Phys. 2, 86 (2019).
  * (24) N.-X. Yang, Y.-F. Zhou, Z. Hou, and Q.-F. Sun, Anomalous spin Nernst effect in Weyl semimetals. J. Phys.: Condens. Matter 31, 435301 (2019).
  * (25) C. M. Wang, H.-P. Sun, H.-Z. Lu, and X. C. Xie, 3D Quantum Hall Effect of Fermi Arcs in Topological Semimetals. Phys. Rev. Lett. 119, 136806 (2017).
  * (26) Y. Xu et al., Observation of topological surface state quantum Hall effect in an intrinsic three-dimensional topological insulator. Nat. Phys. 10, 956 (2014).
  * (27) See supplementary material for materials and methods.
  * (28) P. Hosur, Friedel oscillations due to Fermi arcs in Weyl semimetals. Phys. Rev. B 86, 195102 (2012).
  * (29) Z. Hou, Y. Xing, A.-M Guo, and Q.-F. Sun, Crossed Andreev effects in two-dimensional quantum Hall systems. Phys. Rev. B 94, 064516 (2016).
  * (30) S. Datta, _Electronic Transport in Mesoscopic Systems_ (Cambridge University Press, Cambridge, 1995).
  * (31) S. DasSarma and X. C. Xie, Calculated transport properties of ultrasubmicrometer quasi-one-dimensional inversion lines. Phys. Rev. B 35, 9875(R) (1987).
  * (32) G. Rosenberg, H.-M. Guo, and M. Franz, Wormhole effect in a strong topological insulator. Phys. Rev. B 82, 041104(R) (2010).
  * (33) Y. Zhang, D. Bulmash, P. Hosur, A. C. Potter, and A. Vishwanath. Quantum oscillations from generic surface Fermi arcs and bulk chiral modes in Weyl semimetals. Sci. Rep. 6, 23741 (2016).


# A Novel Approach to Queue-Reactive Models: The Importance of Order Sizes

Hamza Bodor<EMAIL_ADDRESS>Laurent Carlier
<EMAIL_ADDRESS>
###### Abstract

In this article, we delve into the applications and extensions of the queue-
reactive model for the simulation of limit order books. Our approach
emphasizes the importance of order sizes, in conjunction with their type and
arrival rate, by integrating the current state of the order book to determine,
not only the intensity of order arrivals and their type, but also their sizes.
These extensions generate simulated markets that are in line with numerous
stylized facts of the market. Our empirical calibration, using futures on
German bonds, reveals that the extended queue-reactive model significantly
improves the description of order flow properties and the shape of queue
distributions. Moreover, our findings demonstrate that the extended model
produces simulated markets with a volatility comparable to historical real
data, utilizing only endogenous information from the limit order book. This
research underscores the potential of the queue-reactive model and its
extensions in accurately simulating market dynamics and providing valuable
insights into the complex nature of limit order book modeling.

###### keywords:

Limit order book , Market simulators , Stylized facts , Market microstructure

☑✓ ✗✗

[inst1]organization=Université Paris 1 Panthéon-Sorbonne, Centre d’Economie de
la Sorbonne,addressline=106 Boulevard de l’Hôpital, city=Paris Cedex 13,
postcode=75642, country=France

[inst2]organization=BNP Paribas Corporate and Institutional Banking, Global
Markets Data & Artificial Intelligence Lab,addressline=20 boulevard des
Italiens, city=Paris, postcode=75009, country=France

##  1 Introduction

Limit order book (LOB) simulation is a critical area of study with far-
reaching implications in the field of Finance. It serves as the backbone for
understanding and predicting market dynamics, thereby aiding in the
development of effective trading strategies (Cartea et al., 2015). The primary
utility of LOB simulation is manifested in its application to the development
and assessment of algorithmic trading strategies. By simulating various market
scenarios and order book dynamics, LOB models provide crucial insights that
enable traders to formulate strategies that are responsive to different market
conditions. This not only aids in understanding emerging market trends but
also plays a significant role in enhancing the decision-making processes in
algorithmic trading (Cohen and Szpruch, 2012; Balch et al., 2019). Moreover,
it plays a pivotal role in risk management, where understanding the LOB can
help mitigate potential losses due to market volatility (Byrd et al., 2019).
LOB simulation aims to accurately replicate the complex dynamics of financial
markets at the microscopic scale. This involves capturing the intricate
interactions between different market participants, the impact of external
events on market behavior, and the inherent uncertainty and randomness in
market movements. However, achieving a high degree of accuracy in LOB
simulation is a complex task. The intrinsic complexity of financial markets,
coupled with the vast amount of data involved and the need for high
computational efficiency, makes this a challenging problem.

In the realm of limit order book simulation, a variety of methods have been
developed, each offering distinct advantages and challenges. Traditional
approaches often rely on theoretical models like Markovian Poisson processes
(Avellaneda and Stoikov, 2008; Cont et al., 2010; Huang et al., 2015) and
Hawkes processes (Toke, 2011; Lu and Abergel, 2018; Fosset et al., 2022; Wu et
al., 2022). These methods are primarily used to model the pattern and impact
of order arrivals in LOBs. Rooted in the concept of regular patterns in order
flows, these models leverage statistical assumptions to yield solutions that
are concise and analytically tractable.

Alternatively, agent-based methods model the limit order book as a system of
interacting agents, each following a set of rules or strategies. These methods
can capture the heterogeneity of market participants and the strategic
interactions between them. Abides (Byrd et al., 2019), an interactive agent-
based discrete event simulation environment, has encouraged multiple works on
agent-based limit order book generation (Vyetrenko et al., 2020; Amrouni et
al., 2021). More recently, Lussange et al. (2021) propose an agent composed of
two parts: one part predicts the asset price based on the state of the order
book, and the second component uses this prediction to decide on the action to
take. They demonstrate that a configuration with agents of this type, along
with Zero Intelligence (ZI) agents, allows for the production of a simulated
market faithful to multiple stylized facts of the equity market.

Recently, machine learning techniques, particularly Generative Adversarial
Networks (GANs), have also been applied to limit order book simulation
(Sirignano and Cont, 2019; Zhang et al., 2019; Coletta et al., 2023). These
methods can learn complex patterns from data and generate realistic order book
dynamics. They have shown promising results and good capability of taking into
account different aspects of the limit order book.

##  2 Methodology

###  2.1 Queue-Reactive Model

The Queue-Reactive (QR) model of Huang et al. (2015) has gained significant
attention from numerous researchers due to its simple yet effective simulation
framework for the limit order book. This model verifies multiple stylized
facts and addresses various issues that other basic simulators may encounter
(Mariotti et al., 2023).

In Huang et al. (2015), the limit order book is viewed as a $2K$-dimensional
vector, where $K$ represents the depth of the LOB. The reference price,
often111The reference price typically matches the mid-price when the spread is
an even multiple of the tick size. Otherwise, it is adjusted to either half a
tick above or below the mid-price, depending on the specific market
conditions. For more details on this mechanism and its implications on market
dynamics, see subsequent paragraphs. the mid-price, divides the LOB into two
halves: the bid side and the ask side. Each side is populated with an
aggregation of limit orders, represented by $Q_{\pm i}$, where $i$ denotes the
distance from the reference price in terms of ticks, the smallest price
increment change an asset can make on the market.

The quantities at each limit are subject to changes due to three types of
orders: _Limit_ (L) orders, which increase the size of the queue; _Cancel_ (C)
orders, which remove an existing order; and _Market_ (M) orders, which
immediately match with an existing limit quantity, resulting in trades (hence,
this type may also be referred to as _Trade_ orders).

In the queue-reactive model, the arrival of order flow at a given price level
is modeled as a heterogeneous Poisson process. The intensity of this process
depends solely on the current state of the order book, particularly the queue
sizes. The intensities of limit, cancel, and market order arrivals are denoted
by $\lambda^{L}$, $\lambda^{C}$, and $\lambda^{M}$ respectively, and these
intensities are functions of the queue size at each limit.

Figure 1: Illustration of queue size evolution. This diagram represents the
sequential updates in the size of a given queue over time. Initially, the
queue size is $q_{k-1}$. It is then updated: following a time interval of
$\Delta t_{k-1}$, the event $e_{k-1}$, a cancellation, reduces the queue size.
Subsequently, after a time interval $\Delta t_{k}$, a new event, $e_{k}$
(characterized as a limit order), occurs, injecting additional liquidity and
increasing the queue size by $q_{k+1}-q_{k}$. These events are identified by
their corresponding order types: $\eta_{k-1}=C$ for cancellation and
$\eta_{k}=L$ for a limit order, affecting the queue’s liquidity state.

Figure 1 shows how the queue size of a given level is updated. where event
$k-1$ decreases the size of the queue from $q_{k-1}$ to $q_{k}$ after $\Delta
t_{k-1}$. Similarly, the size of the queue increases to $q_{k+1}$ after
$\Delta t_{k}$ when event $k$ occurs.

While this explains how the queue sizes evolve under the queue-reactive model,
it is also important to understand how the price changes under this model. In
fact, Figure 2 displays a scenario of price change under the QR model where
this occurs when the best bid or the best ask price is totally consumed
(either after a cancel or a market order). When this happens and the spread
value is an even multiplier of the tick size, the reference price $p_{ref}$
will now have two possible values: $p_{mid}+0.5~{}\textrm{tick}$ or
$p_{mid}-0.5~{}\textrm{tick}$. The model postulates then that with probability
$\theta$, the reference price changes value, allowing the opposite side to
fill the new price. Otherwise, the reference price stays at its value, and the
levels do not change. The parameter $\theta$ is calibrated in order to
reproduce the mean reversion ratio, as discussed in Huang et al. (2015).

Figure 2: Transition of the reference price $p_{ref}$ in response to spread
changes. When the spread is an even multiple of the tick size, $p_{ref}$ can
adjust by $p_{mid}\pm 0.5$ tick. The figure illustrates this behavior, showing
that with probability $\theta$, the reference price shifts, enabling order
fulfillment at the new price level. $\mathcal{B}(\theta)$ indicates a
Bernoulli trial with success probability $\theta$, determining the reference
price change scenario.

For a set $\mathcal{E}$ of events at a given queue $i$,
$\mathcal{E}=(e_{k}=(\eta_{k},q_{k},\Delta t_{k}))_{k}$, the likelihood of the
QR model can be expressed as:

$\mathcal{L}(\mathbf{\lambda}|\mathcal{E})=\prod_{k}e^{-\lambda^{\text{global}}(q_{k})\Delta
t_{k}}\lambda^{\eta_{k}}(q_{k})$ (1)

where:

  * 1.

$\eta_{k}$ is the type of the $k^{th}$ event,
$\eta\in\mathcal{T},\mathcal{T}=\\{L,C,M\\}$ is the set of possible event
types.

  * 2.

$q_{k}$ is the size of the queue just before the $k^{th}$ event.

  * 3.

$\Delta t_{k}$ is the time difference between events $k-1$ and $k$.

  * 4.

$\lambda^{\text{global}}(q_{k})=\sum_{\eta\in\mathcal{T}}\lambda^{\eta}(q_{k})$,
the intensity of arrival of an order when the size of the queue is $q_{k}$.

For a given size $n$ of the queue, measured in Average Event Size222This is
the average of all events sizes occurring regardless of their type. Huang et
al. (2015) propose also the usage of Average Trade Size (ATS) as a normalizing
factor, which accounts only for the size of trade events. (AES), maximum
likelihood estimators are known in closed form:

$\displaystyle\hat{\Lambda}(n)$ $\displaystyle=\left(\frac{1}{\\#\\{k\mid
q_{k}=n\\}}\sum_{\\{k\mid q_{k}=n\\}}\Delta t_{k}\right)^{-1}$
$\displaystyle\hat{\lambda}^{L}(n)$
$\displaystyle=\hat{\Lambda}(n)\frac{\\#\\{e_{k}\in\mathcal{E}~{}|~{}\eta_{k}=L,q_{k}=n\\}}{\\#\\{e_{k}\in\mathcal{E}~{}|~{}q_{k}=n\\}}$
$\displaystyle\hat{\lambda}^{C}(n)$
$\displaystyle=\hat{\Lambda}(n)\frac{\\#\\{e_{k}\in\mathcal{E}~{}|~{}\eta_{k}=C,q_{k}=n\\}}{\\#\\{e_{k}\in\mathcal{E}~{}|~{}q_{k}=n\\}}$
$\displaystyle\hat{\lambda}^{M}(n)$
$\displaystyle=\hat{\Lambda}(n)\frac{\\#\\{e_{k}\in\mathcal{E}~{}|~{}\eta_{k}=M,q_{k}=n\\}}{\\#\\{e_{k}\in\mathcal{E}~{}|~{}q_{k}=n\\}}$

In our study, we base our experimental results on Bund futures data. Bund
futures are well-known for their high liquidity and the abundance of available
data, making them ideal for detailed analysis. The Bund is a large tick asset,
for which the spread is usually of one tick size, which helps reduce the
impact of large spreads on model behavior. This is quite different from the
assets used in Huang et al. (2015), allowing us to test the model in a market
setting with different characteristics.

Figure 3: Results of queue-reactive model on Bund futures - Intensities as
functions of queue size (best levels).

Figure 3 shows the calibrated intensities on our data, one can conclude that:

  * 1.

Intensities vary significantly with the size of the queue, which justifies the
choice of this model over a simple Poisson model.

  * 2.

The limit intensity generally decreases, except when queue sizes are _tiny_.
This phenomenon can be understood by considering the behavior of Market
Participants (MPs). MPs may perceive it as risky to enter a level with an
extremely low volume, until a certain threshold (typically between 7 AES to 10
AES) is reached. At this threshold, there is a noticeable spike in limit
orders, which can be attributed to the MPs’ pursuit of priority. However, as
the size of the queue grows, this pursuit diminishes.

  * 3.

The same behavior is observed for cancel orders. However, it has the opposite
interpretation; this may be because MPs who do not want to be executed
withdraw their orders when the size of the queue is around 10 AES. After that,
market participants generally prefer to be executed rather than cancel and
move to a deeper level. This contradiction in behavior with the limit
intensity may indicate that there are two major types of market participants:
one type always seeks to be executed as quickly as possible without paying the
spread (without being aggressive), and others who prefer to seek priority
without being executed, possibly waiting for better opportunities. This
hypothesis of a fragmented market has also been adopted by other research
works such as MacDonald (1994) or Lux and Marchesi (1998).

  * 4.

The trading rate decays with the available volume at the best level. This
phenomenon is easily explained by market participants _rushing for liquidity_
when it is scarce, and _waiting for a better price_ when liquidity is
abundant, similar results are noticed on the equity markets in the Huang et
al. (2015) study.

Figure 4: Intensities of limit, cancel, and market orders after the dataset is
divided into two equal time segments.

Moreover, the dataset, encompassing an entire year’s worth of data, was evenly
split into two halves, ensuring the chronological sequence of trading days was
preserved. This deliberate division was undertaken to examine the stability
and reliability of the model calibrations over two distinct time periods.
After calibrating the models for each interval separately, we analyzed the
resulting order flow intensities for limit, cancel, and market orders relative
to queue sizes. The outcomes, illustrated in Figure 4, reveal a noteworthy
consistency across both periods, thereby validating the temporal stability of
the calibrated intensities within the span of the examined period.

Now that we have the intensity of arrival of orders, and their types, the last
and only missing information about the order is its size. We explore two
methods of orders sizing:

  * 1.

Each order is assigned a uniform size, determined by rounding up the average
event size ($\lceil AES_{i}\rceil$) for each level $i$. This approach,
referred to as QRU (QR model with Unitary order sizes), aligns with the model
used in (Huang et al., 2015).

  * 2.

The size of each order is sampled from the stationary distribution of order
sizes for its respective level, denoted as QR.

###  2.2 Queue-Reactive Model Extensions

In addition to these two variants of the QR model, we propose two other
extensions heavily inspired by it. The goal is to examine the impact of
incorporating a more sophisticated model for order sizes on the simulator
performances. The first extension involves a modified QR model that includes
five possible event types, by adding two additional types, cancel_all and
market_all, that model the complete consumption of the queue. This aims to
simulate instances where the size of a consumption order precisely matches the
queue size. The other extension involves incorporating the order size as a
variable in determining the nature of order arrival, similar to its influence
by event type.

####  2.2.1 Five-Type Queue-Reactive model (FTQR)

As mentioned earlier, and as will be discussed in section 3, one of the main
limitations of the queue-reactive model is its inability to reproduce
historical market volatility. This is primarily due to a slower progression of
the LOB when the size of one of the best queues is small, compared to
historical progression. In fact, as the sizes of the orders are drawn from the
stationary distribution of sizes, which is very skewed towards small sizes, it
leads to a very slow decrease of the queue, hence, fewer price changes and low
volatility.

In fact, Figure 5 presents the stationary distribution of order sizes
truncated at 30 lots 333In the context of BUND Futures, a queue comprises
multiple contracts (lots), each valued at 100,000 Euros.($\approx 5AES$)
(stationary), and the conditional distribution of order sizes when the size of
the queue is 30 lots (conditional). It is an example of the state of order
sizes when the queue size is relatively low. This illustration highlights the
disparity between the two distributions, with the conditional distribution
leaning towards larger values and, correspondingly, greater queue consumption.

In particular, there is a pronounced surge of orders that deplete 100% of the
queue, predominantly for market orders. The frequency of such orders is more
than 40 times higher than that in the stationary queue size consumption ratio.
This significant discrepancy provides an explanation as to why queue
consumption is markedly slower when stationary sizes are utilized. Especially
after showing the link between queue consumption and price dynamics in Figure
2. This, in turn, elucidates why the simulated price exhibits considerably
less volatility compared to actual prices.

(a) Cancel orders.

(b) Market orders.

Figure 5: Comparison of stationary distribution of order sizes up to 30 (in
gray) with conditional distribution when queue size equals 30 (in red),
presented on a logarithmic scale.

These insights motivate a first extension of the QR model that take into
account the phenomenon of complete queue consumption when the size of the
queue is small. In such situations, two specific events tend to occur:

  * 1.

Market participants seize the available quantity and aggressively execute it
(using market orders).

  * 2.

Market participants who possess these scarce quantities feel uncomfortable
being alone in the queue and cancel all their orders. Another way of seeing it
is that when market participants are alone in the queue, it is _free_ to
cancel their order without worrying about losing priority in the queue.

In this context, we propose adding two additional events that correspond to
the complete consumption of available quantity, either through trades or
cancellations. The likelihood of this new model is similar to that presented
in Equation 1 by including _market_all_ and _cancel_all_ events in
$\mathcal{T}$, the set of possible event types.

Figure 6: Comparison of QR and FTQR intensities on Bund future data. Figure 7:
Calibrated intensity functions for additional types of the FTQR model applied
to Bund future data, truncated at values less than 25 lots. Intensities for
larger queue sizes are negligible.

Figures 6 and 7 above show the obtained intensities of the five-types queue-
reactive model and compare them with the vanilla QR model intensities,
revealing the following observations:

  * 1.

Limit orders: The limit order intensities are identical in both models, as the
adjustment does not affect limit orders.

  * 2.

Cancel orders: The cancel order intensities are also the same, except for
smaller queue sizes where the intensities of the adjusted model are lower.
This difference originates from the additional cancel type that consumes the
entire level. Nevertheless, the impact of such cancel orders on overall
liquidity is negligible compared to other cancel orders.

  * 3.

Market orders: This order type is the most influenced by the adjustment. We
notice that the new intensities are lower than the original ones, particularly
for smaller queue sizes, as these situations provide more opportunities for
rarity chasing.

  * 4.

Market_all orders: The new intensities compensate for the earlier ones,
particularly for small queue values, reflecting a significant value. This
finding indicates that many orders consume liquidity at a high rate when this
liquidity is scarce. The intensity of market_all orders decays with the queue
size.

  * 5.

Cancel_all orders: This order type exhibits a behavior similar to the previous
one. It is mainly interpreted as the intensity of arrival of events where a
market participant, who is alone in the queue, completely withdraws from the
market.

####  2.2.2 Size-Aware Queue-Reactive model (SAQR)

In this section, we explore the impact of order sizes on the dynamics of limit
order book simulations. While research has extensively analyzed the rates and
types of event arrivals, the size of orders has received less attention.

To build upon this dimension, we propose an adaptation of the QR model that
takes into account the size of orders. We achieve this by changing the
possible values of $\mathcal{T}$ in Equation (1) to a 2-dimensional order
type, $\mathcal{T}=\\{L,C,T\\}\times\mathbb{N}^{*}$. This means that the
intensity nature not only specifies the type of the order but also its size.
This adjustment does not change the formulation of the likelihood. And the
solution, for $(\eta,s)$ and a queue size $n$, is:

$\hat{\lambda}^{(\eta,s)}(n)=\hat{\Lambda}(n)\frac{\\#\\{e_{k}\in\mathcal{E}~{}|~{}\eta_{k}=\eta,s_{k}=s,q_{k}=n\\}}{\\#\\{e_{k}\in\mathcal{E}~{}|~{}q_{k}=n\\}}$

where $s_{k}$ is an additional model variable and represents the size of the
$k^{th}$ order.

(a) Limit orders

(b) Cancel orders

(c) Market orders

Figure 8: Matrix heatmap of normalized calibrated intensities On Bund future
data (logarithmic scale). Normalization is done such that intensities sum to
one for each row.

Figure 8 presents the intensity matrix for the three types of orders: limit,
cancel and market orders (averaged for each queue size). It is evident from
the figure that while the distribution of limit order sizes remains relatively
consistent across the different queue size, the same cannot be said for cancel
and market orders.

In particular, the diagonal values in the intensity matrix, which represent
orders that consume the entire available liquidity, are noticeably dominant.
This dominance is especially pronounced for market orders, confirming the
tendency of these orders to consume the entire quantity available at the best
levels.

Interestingly, this phenomenon is more prevalent when the available quantity
is relatively low. This observation highlights the fact that the distribution
of order sizes is not uniform, but rather varies significantly depending on
the size of the queue. This nuanced understanding of order size distribution
is crucial for accurately simulating and analyzing the dynamics of the limit
order book.

Furthermore, one can notice in the conditional market intensities some
particularly _salient_ intensities at 9 AES, 17 AES and 26 AES. These are
buckets including orders of sizes 50 lots, 100 lots and 150 lots,
respectively. These round numbers are generally more occurrent in historical
data than other order sizes, hence the high intensity of arrival of orders of
these sizes.

###  2.3 Hawkes Models

Hawkes processes, named after Alan Hawkes (Hawkes, 1971), are a class of self-
exciting and mutually-exciting point processes. Unlike Poisson processes, they
dispense with the independence-of-increments property, allowing past events to
influence the likelihood of future events. This characteristic makes them
particularly adapted for modeling clustering and contagion phenomena, which
are prevalent in social media networks, seismology, and financial markets.

Mathematically, a Hawkes process is defined by its intensity function,
$\lambda(t)$, which depends on the entire history of the process. The
intensity function is often formulated as:

$\lambda(t)=\mu+\int_{-\infty}^{t}h(t-s)dN(s)$

where $\mu$ is the baseline intensity, $h$ is a memory kernel describing the
influence of past events, and $N(s)$ is the count of events up to time $s$.

The complexity of financial markets, with multiple interacting variables,
necessitates a more comprehensive approach to capture these dynamics
accurately. This leads us to the multivariate Hawkes process, which extends
the basic concept to encompass multiple types of interdependent events.

Mathematically, a multivariate Hawkes process is defined by a set of intensity
functions, $\lambda_{i}$, each depending on the entire history of the process.
The intensity function for the $i^{th}$ component of the process at time $t$
can be expressed as:

$\lambda_{i}(t)=\mu_{i}+\sum_{j=1}^{n}\int_{-\infty}^{t}h_{ij}(t-s)dN_{j}(s)$

where $\mu_{i}$ is the baseline intensity for the $i$-th component, $h_{ij}$
is a memory kernel describing the influence of past events of type $j$ on the
current rate of type $i$, and $N_{j}(s)$ is the count of type $j$ events up to
time $s$. This model captures not only the self-excitations ($j=i$) of events,
but also the cross-excitation ($j\neq i$) between events.

In the context of limit order book (LOB) modeling, Hawkes processes have been
instrumental in describing the flow of orders. These models view different
order types as distinct events in a multivariate point process, capturing the
dynamics of market reactions to these orders and their interrelations.
Significant developments in LOB modeling include contributions from Zhao
(2010) and Toke (2011), who have expanded upon the Cont et al. (2010)
framework. Zhao’s empirical analysis of crude oil futures led to the rejection
of the assumption that market event inter-arrival times were independent and
identically distributed, challenging the conventional use of independent
Poisson processes. Instead, Zhao (2010) introduced a Hawkes process, as
suggested by Bauwens and Hautsch (2009), to model the arrival rate of market
events as a function of recent order arrival rates and quantities. This
approach better represented periods of clustered high and low arrival rates,
aligning more closely with empirical observations (Ellul, 2001; Hall and
Hautsch, 2006). Similarly, Toke (2011) replaced the Poisson processes in the
Cont et al. (2010) model with multivariate mutually-exciting Hawkes processes,
tailoring them for each type of market event. More recently, Wu et al. (2022)
integrated the queue-reactive model with Hawkes processes, developing a model
that captures both the excitation dynamics and adjusts intensities in
accordance with the queue size.

In our study, we follow a methodology akin to that of Zhao (2010) and Toke
(2011), aiming to compare it to the different versions of the queue-reactive
model. To this end, we calibrate a Hawkes process to develop our simulation
framework. Our parametrization strategy is as follows:

  * 1.

We treat all events occurring at the best limits as a single, six-dimensional
multivariate process. This approach accounts for three events at each best
price, allowing us to capture not only the excitation between events within
the same queue but also the interdependencies between the bid and ask sides of
the market.

  * 2.

Exponential kernels are employed to model the time-dependent influence of past
events on future occurrences. Which means, the memory kernel functions are of
the form $h_{ij}(t)=\alpha_{ij}e^{-\beta_{ij}t}$.

  * 3.

We adhere to the same price change paradigm as established in the queue-
reactive model, ensuring consistency and comparability in our analyses.

For the calibration of this model, we utilized the tick package of Bacry et
al. (2017). Once calibrated, the same package facilitated the simulation of
event sequences and interarrival times. These simulations allowed to construct
the order flow, which was subsequently incorporated into our matching engine
to generate the final simulated order book.

In parallel with the queue-reactive model, we examine two variants of the
Hawkes-based model in relation to order sizes. The first variant, referred to
as Hawkes_U, operates on the assumption of unitary order sizes. The second
variant, dubbed Hawkes_S, utilizes the stationary distribution of order sizes.
This bifurcation allows for a nuanced examination of the models’ performance
under different assumptions about market order size models.

##  3 Results

###  3.1 Data & Calibration Paradigm

The data used for calibration is from Euro-Bund Futures (FGBL)444Euro-Bund
Futures are derivative contracts based on 10-year German government bonds
(Bund) and traded on the Eurex Exchange. Each futures contract, with a
notional coupon rate of 6%, has a remaining term to maturity of 8.5 to 10.5
years on the delivery day. covering active trading days on the year 2021. Data
is filtered into periods from 9 AM to 6 PM. The maximum depth considered is
five levels on each side.

Level |  | #L
---
$(\times 10^{6})$
| #C
---
$(\times 10^{6})$
| #M
---
$(\times 10^{4})$
AES |  | AIT
---
(ms)
|
1 | 119.6 | 115.6 | 756 | 5.92 | 56 |  |
2 | 24.09 | 23.84 | 2.88 | 5.15 | 272 |  |
3 | 15.76 | 15.29 | 0 | 4.87 | 370 |  |
4 | 10.72 | 10.45 | 0 | 5.98 | 468 |  |
5 | 11.12 | 10.76 | 0 | 4.11 | 523 |  |
Table 1: Descriptive statistics about events for each level.

Table 1 provides statistics about the data utilized for model calibration. It
displays the number of observations for each event type, the Average Event
Size (AES), and the Average Inter-event Time (AIT) for each queue. Notably,
the table underscores the high frequency of events impacting queues close to
the mid-price, especially the best prices.

Additionally, our dataset includes updates to the limit order book and the
daily trade history. However, the models necessitate training with order flow,
representing the sequence of events in the limit order book. To obtain this
order flow, we carefully constructed it by discerning whether each update
indicates liquidity provision (via a limit order) or liquidity consumption
(through cancel and market orders). Further details on order flow construction
are available in Bodor and Carlier (2024). Additionally, a matching engine was
developed to pair orders produced by each model, subsequently reflecting them
in a simulated order book.

This data is further pre-processed following the methodology in Huang et al.
(2015). Specifically:

  * 1.

Queue sizes are quantized based on the average event size. If the size of
queue $i$ is denoted by $q_{i}$ and its average event size is $AES_{i}$, then
the updated queue size is given by
$q_{i}\leftarrow\lceil\frac{q_{i}}{AES_{i}}\rceil$.

  * 2.

The data is divided into periods with a constant reference price. The
interarrival time counter for events is reset to 0 for all queues when there
is a change in the reference price.

  * 3.

The parameter that controls the reference price movements is $\theta=0.7$.
This probability has been tweaked to reproduce the observed ratio of price
continuations to price alternations.

For model intensity calibration and subsequent LOB simulations, we:

  * 1.

Utilize the calibrated intensities to determine both the time of arrival and
the type of new orders.

  * 2.

For modeling price changes, we adhere to the aforementioned paradigm. However,
we refrain from redrawing the LOB from its stationary distribution to
exclusively replicate the endogenous market behavior.

  * 3.

Upon a change in the reference price, the size of the queue for newly
appearing prices is sampled from their respective stationary distributions.

###  3.2 Stylized Facts

In our quest to assess the efficiency of different market simulators, we
advocate for their comparison against a set of established stylized facts from
the market. These facts span micro to macro-level empirical findings,
particularly tailored for the Bund futures market. As detailed in Bodor and
Carlier (2024), we delve into an intricate analysis of stylized facts
pertinent to German bond futures, spotlighting Schatz, Bobl, Bund, and Buxl,
grounded on exhaustive limit order book data. This exploration discerns both
commonalities and unique traits across these futures, shedding light on market
dynamics and pointing towards the refinement of market simulation tools. Our
primary focus is to benchmark our models against the stylized facts observed
in the Bund futures market (Bodor and Carlier, 2024), as detailed in the
following sections.

####  3.2.1 Price dynamics and volatility

A fundamental feature of an effective simulator is its ability to simulate
realistic price dynamics and replicate the volatility inherent to the target
market. Figure 9 illustrates the price dynamics of both historical data and
simulations generated using the different models.

It becomes apparent that the queue-reactive model, particularly the unitary
sized variant, exhibits slower dynamics than the actual market, coupled with
reduced volatility. Conversely, the FTQR model seems to engender a markedly
noisier market environment. Meanwhile, the SAQR model’s performance closely
aligns with real price patterns. Regarding the Hawkes-based simulation models,
they display characteristics similar to those of the standard queue-reactive
model. Notably, these simulators also tend to mirror the slower market
dynamics observed in the queue-reactive model, contrasting with the quicker
pace of real-world markets.

(a)

(b)

Figure 9: Mid-price dynamics and annualized volatility plots (1 second
sampling frequency).

More generally, Table 2 presents a comprehensive comparison of price
volatility across different days and distinct periods of the day. This
comparison employs two metrics:

  * 1.

Relative difference:
$\text{mean}\left(100\%\times\frac{\sigma_{s}-\sigma_{r}}{\sigma_{r}}\right)$

  * 2.

Quadratic error: $\text{mean}\left((\sigma_{s}-\sigma_{r})^{2}\right)$

Here, $\sigma_{s}$ and $\sigma_{r}$ denote the annualized volatility of the
simulated and real markets, respectively.

The subsequent table delineates the comparative volatility of simulated
markets vis-à-vis the real market. Notably, the SAQR model persistently
outperforms other models across all periods of the day, reflecting the
smallest average relative difference and quadratic error. This sharply
contrasts the QRU model, which exhibits considerably diminished volatility
relative to the real market.

The performance trajectories of the QR and FTQR models warrant special
attention, as they fluctuate based on the time of day. The QR model exhibits
closer adherence to real market dynamics during the 10 AM to 2 PM bracket,
generally perceived as the “calm” periods of the day. Conversely, the FTQR
model manifests superior performance from 3 PM to 6 PM, periods known for high
market activity driven by events like closing auctions and options expiry.

| 9 AM - 6 PM | 10 AM - 2 PM | 3 PM - 6 PM
---|---|---|---
|  | Relative
---
Difference
(%)
| Quadratic
---
Error
| Relative
---
Difference
(%)
| Quadratic
---
Error
| Relative
---
Difference
(%)
| Quadratic
---
Error
QRU | -91 | 18.37 | -90 | 14.97 | -91 | 21.50
QR | -52 | 7.34 | -48 | 5.17 | -54 | 9.56
FTQR | 69 | 8.88 | 82 | 10.71 | 49 | 8.12
SAQR | 1 | 1.35 | 8 | 0.82 | -3 | 2.33
Hawkes_U | -72 | 12.46 | -70 | 9.64 | -73 | 15.16
Hawkes_S | -69 | 11.70 | -67 | 8.96 | -71 | 14.33
Table 2: Statistics about volatility of simulated markets compared with real
market on different periods of the day.

####  3.2.2 Transactions magnitude and market activity

An essential aspect of limit order book simulators is their ability to
accurately capture the magnitude of transactions and orders. This entails the
capability of the simulator to replicate a similar order of magnitude of
traded volumes within a given time window. As illustrated in Figure 10, which
shows that the SAQR model exhibits the highest proficiency in replicating
traded volumes during typical trading periods of the day. Both the QR and QRU
models demonstrate a relatively low transactions volume. However, the FTQR
model exhibits a significantly higher rate. The Hawkes-based models are also
quite accurate at reproducing a relatively similar market activity, showcasing
the importance of excitation between events for the verification of this
stylized fact.

Figure 10: Traded volumes on windows of size 10 minutes. | 9 AM - 6 PM | 10 AM - 2 PM | 3 PM - 6 PM
---|---|---|---
|  | Relative
---
Difference
(%)
| Quadratic
---
Error
($\times 10^{9}$)
| Relative
---
Difference
(%)
| Quadratic
---
Error
($\times 10^{9}$)
| Relative
---
Difference
(%)
| Quadratic
---
Error
($\times 10^{9}$)
QRU | -80 | 2.02 | -76 | 1.32 | -82 | 2.81
QR | -70 | 1.55 | -65 | 0.95 | -74 | 2.26
FTQR | 49 | 0.85 | 76 | 1.38 | 30 | 0.63
SAQR | -26 | 0.36 | -13 | 0.14 | -36 | 0.62
Hawkes_U | -40 | 0.60 | -28 | 0.33 | -48 | 0.97
Hawkes_S | -32 | 0.41 | -19 | 0.21 | -41 | 0.80
Table 3: Statistics about traded volumes of simulated during different periods
of the day.

Table 3 provides a comprehensive comparison between the traded volumes of the
simulated market and the real market. Consistently, the SAQR model emerges as
the closest to real data, exhibiting the smallest quadratic error across all
day periods, followed by the Hawkes-based models, that show also quite good
results, with a notable improvement of the Hawkes_S model over Hawkes_U.
Interestingly, the QR model outperforms the FTQR model during the periods from
10 AM to 2 PM, while the FTQR model appears to be better suited for
replicating the high market activity observed from 3 PM to 6 PM.

####  3.2.3 Order book queue size distribution

An additional noteworthy empirical observation pertains to the distribution of
available liquidity in the market. Numerous studies have not only acknowledged
the consistency of this distribution but also reported the successful fit of
the Gamma distribution to it (Abergel et al., 2016; Vyetrenko et al., 2020;
Bodor and Carlier, 2024).

In Figure 11, the empirical distribution of queue sizes at the most favorable
levels is effectively captured. All models fit accurately the Gamma
distribution, highlighting their capability in modeling queue dynamics.
Notably, the SAQR and FTQR models demonstrate a particularly strong alignment
with actual market data, underscoring their accuracy in reflecting real market
conditions. In contrast, the Hawkes-based models do not mimic the queue size
distribution as accurately. They tend to simulate many limit orders,
regardless of the queue’s actual size. This results in an unrealistic pile-up
of orders, which does not match what we usually see in real markets. These
observations highlight the strengths of the queue-reactive model. It adjusts
the flow of orders – both incoming and outgoing – based on the size of the
queue. This empirical validation of the queue-reactive model’s ergodicity,
which was theoretically established for unitary sizes in Huang et al. (2015)
and extended to a general order size model in Lehalle et al. (2021),
underscores its efficacy in replicating actual market dynamics, especially in
terms of queue size distribution.

Figure 11: Distribution of order book volumes and Gamma law fit.

In order to compare these distributions across multiple days of data, we adopt
the Kolmogorov-Smirnov statistic as a distance between distributions. The
Kolmogorov-Smirnov (KS) test is a nonparametric test used to determine the
similarity between two probability distributions. The statistic of the KS test
is defined as:

$D_{n}=\sup_{x}|F_{n}(x)-F(x)|$

where $F_{n}(x)$ is the empirical cumulative distribution function (CDF) of
the observed data and $F(x)$ is the CDF of the reference distribution or
theoretical model. The value of $D_{n}$ lies between 0 and 1. The smaller the
value of $D_{n}$, the more similar the observed data distribution is to the
reference distribution. Therefore, the KS statistic can be used to measure the
similarity between two distributions.

Table 4 presents the descriptive statistics of the Kolmogorov-Smirnov
statistic, comparing the real and simulated queue size distributions from
different days of calibration data. The table indicates that both the SAQR and
QR models perform well on average, despite the fact that the volume
distributions generated by these two models are quite distinct. This
discrepancy can be attributed to the variability in market behavior, which
aligns more closely with the QR distribution on some days and the SAQR
distribution on others. Regardless, both models can be considered to validate
the stylized fact. In contrast, Hawkes-based models exhibit very wide
distribution distances from the real distribution, further confirming that
these models suffer from the queue size divergence problem, where,
intrinsically, there is no guarantee of a return to equilibrium, unlike the
queue reactive model.

| QRU | QR | FTQR | SAQR | Hawkes_U | Hawkes_S
---|---|---|---|---|---|---
mean | 0.32 | 0.22 | 0.28 | 0.20 | 0.63 | 0.69
std | 0.09 | 0.11 | 0.08 | 0.06 | 0.05 | 0.04
min | 0.23 | 0.14 | 0.05 | 0.08 | 0.57 | 0.64
25% | 0.27 | 0.15 | 0.26 | 0.15 | 0.60 | 0.66
50% | 0.30 | 0.18 | 0.30 | 0.21 | 0.62 | 0.68
75% | 0.31 | 0.22 | 0.33 | 0.26 | 0.65 | 0.71
max | 0.68 | 0.61 | 0.37 | 0.28 | 0.78 | 0.81
Table 4: Kolmogorov-Smirnov statistic of the distribution of queue sizes in
simulated markets compared to the real market, description over multiple days
of calibration data.

####  3.2.4 Signature plot

The signature plot of a price is defined as the variance of price increments
at a given lag value $h$, divided by $h$. In other words, it is defined as
$\sigma^{2}_{h}=\frac{\mathbb{V}\left(P(t+h)-P(t)\right)}{h}$, $P(t)$ being
the mid-price of the asset at time $t$. It is interesting to observe the
evolution of this metric with the value of $h$, as this highlights the
transition of price movements at different time scales. Figure 12 shows this
transition for an example of one real day and simulated markets. While all the
markets seem to exhibit the same shape of decrease in volatility as the time
lag value increases, the SAQR model is the best at capturing the real
signature plot at different time scales.

Figure 12: Signature plot.

More generally, the following table provides more statistics regarding the
comparison of the signature plot of real data from different days to the
signature plot of simulated data, the distance between curves is defined as
the sum of difference of values for each lag $h$. The table shows that the
QRU, QR, Hawkes_U and Hawkes_S models generally underestimate the volatility
while the FTQR heavily overestimates it. The SAQR is the model with the least
quadratic error and least bias, confirming its ability to capture the real
volatility at different scales. However, it is important to note that while
the SAQR model excels in terms of quadratic error and bias, it still struggles
to perfectly capture the exact decaying profile of the real data, especially
at smaller time scales.

|  | Relative
---
Difference
(%)
| $\mathbf{L^{2}}$
---
Norm
$\times 10^{-5}$
QRU | -9.80 | 2.89
QR | -7.55 | 2.28
FTQR | 41.32 | 9.99
SAQR | 4.43 | 1.05
Hawkes_U | -8.57 | 2.09
Hawkes_S | -7.91 | 1.95
Table 5: Statistics over multiple days of calibration data - signature plot
values.

####  3.2.5 Distribution of returns

(a)

(b)

Figure 13: Distribution of logarithmic returns $\log(\frac{P(t+\tau)}{P(t)})$
of real and simulated data.

Figure 13 displays a comparison of daily return distributions across various
simulated markets and the real market. Part (a) of the figure compares the
distributions of all models against the real market, while part (b) focuses on
the best-performing models for enhanced clarity. This comparison reveals that
the SAQR model’s distribution of returns closely resembles that of the real
market. In part (a), it is noticeable that the QR, QRU, Hawkes_U, and Hawkes_S
models generate returns that are too concentrated around zero, reflecting a
slower volatility than what is seen in the actual market. On the other hand,
part (b) illustrates how the SAQR model’s return distribution closely mirrors
the real market. This finding is corroborated by Table 6, which indicates that
the SAQR and FTQR models have return distributions most similar to those of
the real market, with the SAQR model being the closest.

| QRU | QR | FTQR | SAQR | Hawkes_U | Hawkes_S
---|---|---|---|---|---|---
mean | 0.52 | 0.35 | 0.21 | 0.16 | 0.43 | 0.42
std | 0.04 | 0.05 | 0.05 | 0.03 | 0.03 | 0.03
min | 0.43 | 0.25 | 0.12 | 0.11 | 0.39 | 0.37
25% | 0.49 | 0.31 | 0.17 | 0.14 | 0.41 | 0.40
50% | 0.51 | 0.35 | 0.20 | 0.15 | 0.43 | 0.42
75% | 0.55 | 0.37 | 0.23 | 0.17 | 0.46 | 0.44
max | 0.60 | 0.44 | 0.35 | 0.25 | 0.50 | 0.48
Table 6: Kolmogorov-Smirnov statistic of the distribution of returns in
simulated markets compared to the real market. Description over multiple days
of calibration data

####  3.2.6 Order sizes distribution

Figure 14: Distribution of order sizes (the y-axis is in logarithmic scale).

Figure 14 provides insights into the distribution of order sizes in both real
and simulated markets. Several observations can be made:

  * 1.

As anticipated, the QRU and Hawkes_U models perform poorly in reproducing this
stylized fact. All orders in these models are unitary and clustered around the
average event size (AES), resulting in a limited range of order sizes, and
completely different distribution.

  * 2.

The QR and Hawkes_S models exhibit a good fit to the historical distribution
of order sizes. This outcome is expected since these models sample, by
construction, order sizes from the stationary distribution of real market
order sizes.

  * 3.

Conversely, the SAQR model demonstrates a strong capability in reproducing the
distribution of order sizes. It is important to emphasize that, in contrast to
the QR model, this capability is not inherent to the SAQR model’s design. The
order sizes in the SAQR model are drawn based on a distribution conditional on
the queue size. This feature allows the SAQR model to capture the empirical
order size distribution more accurately, accounting for the realism of the
queue sizes dynamics of this model.

  * 4.

The FTQR model, while similar to the QR model in terms of order size picking
paradigm, presents a distribution slightly shifted towards larger sizes. This
shift is primarily due to the introduction of the two event types _cancel_all_
and _market_all_ which encourage the placement of larger orders. This could
also justify the increased volatility and activity observed in this model.

####  3.2.7 Order book shape

Figure 15: Order book mean volumes (averaged).

Figure 15 illustrates the average volume available at each queue up to the
fifth level, normalized by the average of all queue sizes. This figure
indicates that the limit order book is generally more populated at deeper
levels, a phenomenon also observed on other markets (Abergel et al., 2016;
Chakraborti et al., 2011), especially on the Bund (see Bodor and Carlier
(2024)).

This stylized fact is particularly evident in the FTQR and SAQR models.
However, it appears less applicable to the QRU and QR models. This difference
underscores the superior ability of the FTQR and SAQR models to accurately
simulate queue size distributions at various price depths. The Hawkes-based
models are absent in the figure as they simulate only best prices order flow.

####  3.2.8 Summary

Additional stylized facts are detailed in A. Tables 7 and 3.2.8 present a
comprehensive list of stylized facts examined in our study, along with a
performance summary for each model concerning these facts. We have used a
color-coded system for clarity: green signifies that a model successfully
captures a stylized fact, orange denotes partial accuracy, red indicates a
miss, and black dash means irrelevant. The SAQR model emerges as a clear
frontrunner, corroborating the most stylized facts and thereby underscoring
its proficiency in mimicking real market dynamics. Conversely, the QRU and
Hawkes_U models lag behind, emphasizing the significance of incorporating
order sizes in market simulations. To sum up, accurately modeling order sizes
is instrumental in enhancing the authenticity of simulated markets and
bolstering the efficacy of the queue-reactive model.

[HTML]00C057 Index |  [HTML]00C057Stylized fact name
---|---
1 | Distribution of order sizes
2 | Power-law of order sizes distribution
3 | Signature plot
4 | Distribution of available volumes in the queue
5 | Price dynamics and volatility
6 | Long range dependency
7 | Distribution of returns
8 | Order book shape
9 | Absence of autocorrelation
10 | Traded volumes in a fixed window
11 | Weibull fit of interarrival time of trades
12 | Excitation between events
Table 7: Table of stylized facts.

[HTML]00C057 [HTML]00C057Model | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12
---|---|---|---|---|---|---|---|---|---|---|---|---
[HTML]00C057QRU | ✗ | ✗ | ✗ | ☑ | ✗ | ✗ | ✗ | ✗ | ☑ | ✗ | ☑ | ✗
[HTML]00C057QR | ☑ | ☑ | ✗ | ☑ | ✗ | ☑ | ✗ | ✗ | ☑ | ✗ | ☑ | ✗
[HTML]00C057FTQR | ✗ | ☑ | ✗ | ☑ | ✗ | ☑ | ✗ | ☑ | ☑ | ✗ | ☑ | ✗
[HTML]00C057SAQR | ☑ | ☑ | ☑ | ☑ | ☑ | ☑ | ☑ | ☑ | ☑ | ☑ | ☑ | ✗
[HTML]00C057Hawkes_U | ✗ | ✗ | ✗ | ✗ | ✗ | ☑ | ✗ | - | ☑ | ☑ | ☑ | ☑
[HTML]00C057Hawkes_S | ☑ | ☑ | ✗ | ✗ | ✗ | ☑ | ✗ | - | ☑ | ☑ | ☑ | ☑

Table 8: Summary of stylized facts verified by each model.

##  4 Conclusion and Perspectives

In this paper, we explore adaptations of the queue-reactive model to
accurately simulate the dynamics of limit order books, using empirical data,
with a focus on German bond futures. Our research uniquely incorporates and
demonstrates the critical importance of order sizes, alongside order types and
arrival rates, in realistically modeling limit order book flows. This
methodological enhancement is validated by aligning closely with numerous
market stylized facts, emphasizing the model’s ability to accurately replicate
real market behaviors. Moving forward, we aim to further refine the model by
integrating additional dynamics such as event excitation and intra-day
seasonality, among others, ensuring a balance between simplicity and
comprehensive market representation.

## Acknowledgements

We extend our sincere thanks to the first author’s PhD supervisor, Pr. Olivier
Guéant, for his invaluable guidance and insights throughout this project. Our
gratitude is also due to the Vietnam Symposium in Banking and Finance 2023 for
the opportunity to present and discuss our work with esteemed researchers.
Finally, we would like to express our appreciation to Pr. Charles-Albert
Lehalle and Pr. Mathieu Rosenbaum for the constructive discussions about their
paper and model, and for the precious advice on future research directions.

## Disclosure statement

No potential conflict of interest was reported by the authors.

## Declaration of Funding

This work was conducted as part of a CIFRE thesis in collaboration with the
Sorbonne Economics Centre, with financial support from BNP Paribas Arbitrage.
The contribution of BNP Paribas Arbitrage has been pivotal to the research,
ensuring access to essential data, providing necessary funding, and supplying
resources and materials to conduct the research work. Despite these
contributions, BNP Paribas Arbitrage had no influence on the study’s design,
analysis, interpretation of the results, or the decision-making process
regarding this publication.

## References

  * Cartea et al. (2015) Á. Cartea, S. Jaimungal, J. Penalva, Algorithmic and High-Frequency Trading, Cambridge University Press, 2015.
  * Cohen and Szpruch (2012) S. N. Cohen, L. Szpruch,  A limit order book model for latency arbitrage,  Mathematics and Financial Economics 6 (2012) 211–227.
  * Balch et al. (2019) T. H. Balch, M. Mahfouz, J. Lockhart, M. Hybinette, D. Byrd,  How to evaluate trading strategies: Single agent market replay or multiple agent interactive simulation?,  arXiv preprint arXiv:1906.12010 (2019).
  * Byrd et al. (2019) D. Byrd, M. Hybinette, T. H. Balch,  Abides: Towards high-fidelity market simulation for ai research,  arXiv preprint arXiv:1904.12066 (2019).
  * Avellaneda and Stoikov (2008) M. Avellaneda, S. Stoikov,  High-frequency trading in a limit order book,  Quantitative Finance 8 (2008) 217–224.
  * Cont et al. (2010) R. Cont, A. Kukanov, S. Stoikov,  Stochastic models for order book dynamics,  Operations research 58 (2010) 549–563.
  * Huang et al. (2015) W. Huang, C.-A. Lehalle, M. Rosenbaum,  Simulating and analyzing order book data: The queue-reactive model,  Journal of the American Statistical Association 110 (2015) 107–122.
  * Toke (2011) I. M. Toke,  “market making” in an order book model and its impact on the spread,  Econophysics of Order-driven Markets: Proceedings of Econophys-Kolkata V (2011) 49–64.
  * Lu and Abergel (2018) X. Lu, F. Abergel,  High-dimensional hawkes processes for limit order books: modelling, empirical analysis and numerical calibration,  Quantitative Finance 18 (2018) 249–264.
  * Fosset et al. (2022) A. Fosset, J.-P. Bouchaud, M. Benzaquen,  Non-parametric estimation of quadratic hawkes processes for order book events,  The European Journal of Finance 28 (2022) 663–678.
  * Wu et al. (2022) P. Wu, M. Rambaldi, J.-F. Muzy, E. Bacry,  Queue-reactive hawkes models for the order flow,  Market microstructure and liquidity (2022).
  * Vyetrenko et al. (2020) S. Vyetrenko, D. Byrd, N. Petosa, M. Mahfouz, D. Dervovic, M. Veloso, T. Balch,  Get real: Realism metrics for robust limit order book market simulations,  Proceedings of the First ACM International Conference on AI in Finance (2020) 1–8.
  * Amrouni et al. (2021) S. Amrouni, A. Moulin, J. Vann, S. Vyetrenko, T. Balch, M. Veloso,  Abides-gym: gym environments for multi-agent discrete event simulation and application to financial markets,  Proceedings of the Second ACM International Conference on AI in Finance (2021) 1–9.
  * Lussange et al. (2021) J. Lussange, I. Lazarevich, S. Bourgeois-Gironde, S. Palminteri, B. Gutkin,  Modelling stock markets by multi-agent reinforcement learning,  Computational Economics 57 (2021) 113–147.
  * Sirignano and Cont (2019) J. Sirignano, R. Cont,  Universal features of price formation in financial markets: perspectives from deep learning,  Quantitative Finance 19 (2019) 1449–1459.
  * Zhang et al. (2019) P. Zhang, H. Chen, Z. Liu, Z. Liu, Y. Yan,  Integrating the data augmentation scheme with various classifiers for acoustic scene modeling,  arXiv preprint arXiv:1907.06639 (2019).
  * Coletta et al. (2023) A. Coletta, J. Jerome, R. Savani, S. Vyetrenko,  Conditional generators for limit order book environments: Explainability, challenges, and robustness,  arXiv preprint arXiv:2306.12806 (2023).
  * Mariotti et al. (2023) T. Mariotti, F. Lillo, G. Toscano,  From zero-intelligence to queue-reactive: limit-order-book modeling for high-frequency volatility estimation and optimal execution,  Quantitative Finance 23 (2023) 367–388.
  * MacDonald (1994) R. MacDonald,  Exchange Rate Theory: Chaotic of Foreign Exchange Markets,  The Economic Journal 104 (1994) 966–967. URL: https://doi.org/10.2307/2235001. doi:10.2307/2235001. arXiv:https://academic.oup.com/ej/article-pdf/104/425/966/27040459/ej0966.pdf.
  * Lux and Marchesi (1998) T. Lux, M. Marchesi,  Scaling and criticality in a stochastic multi-agent model of a financial market,  Nature 397 (1998) 498–500.
  * Hawkes (1971) A. G. Hawkes,  Spectra of some self-exciting and mutually exciting point processes,  Biometrika 58 (1971) 83–90.
  * Zhao (2010) L. Zhao, A model of limit order book dynamics and a consistent estimation procedure, Ph.D. thesis, Carnegie Mellon University, 2010.
  * Bauwens and Hautsch (2009) L. Bauwens, N. Hautsch,  Modelling financial high frequency data using point processes,  Handbook of Financial Time Series (2009) 953–979.
  * Ellul (2001) A. Ellul,  Determinants of order choice on the new york stock exchange,  Academy of Management Journal 44 (2001) 81.
  * Hall and Hautsch (2006) A. Hall, N. Hautsch,  Order aggressiveness and order book dynamics,  Empirical Economics 30 (2006) 973–1005.
  * Bacry et al. (2017) E. Bacry, M. Bompaire, S. Gaïffas, S. Poulsen,  tick: a Python library for statistical learning, with a particular emphasis on time-dependent modeling,  arXiv preprint arXiv:1707.03003 (2017).
  * Bodor and Carlier (2024) H. Bodor, L. Carlier,  Stylized facts and market microstructure: An in-depth exploration of german bond futures market,  arXiv preprint arXiv:2401.10722 (2024).
  * Abergel et al. (2016) F. Abergel, M. Anane, A. Chakraborti, A. Jedidi, I. M. Toke, Limit order books, Cambridge University Press, 2016.
  * Lehalle et al. (2021) C.-A. Lehalle, O. Mounjid, M. Rosenbaum,  Optimal liquidity-based trading tactics,  Stochastic Systems 11 (2021) 368–390.
  * Chakraborti et al. (2011) A. Chakraborti, I. M. Toke, M. Patriarca, F. Abergel,  Econophysics review: I. empirical facts,  Quantitative Finance 11 (2011) 991–1012.
  * Cont (2001) R. Cont,  Empirical properties of asset returns: stylized facts and statistical issues,  Quantitative Finance 1 (2001) 223–236.
  * Bouchaud et al. (2018) J.-P. Bouchaud, J. Bonart, J. Donier, M. Gould, Trades, quotes and prices: financial markets under the microscope, Cambridge University Press, 2018.

##  Appendix A Complementary Stylized Facts Results

###  A.1 Autocorrelation of price returns

Figure 16 displays the autocorrelation of returns for different sampling
frequencies and for the six models, it shows that there is generally no
autocorrelation between returns at different time scales.

Figure 16: Autocorrelation of returns for different sampling frequencies.

###  A.2 Long Range Dependence

While a realistic market should show no correlation between returns, it is
commonly observed in different markets that the correlation of absolute
returns shows a slow decay (Cont, 2001). Moreover, this decay follows a power-
law distribution. Figure 17 shows how simulated markets exhibit this slow
decay of absolute returns correlation, with a pretty good fit of power-law
distribution, except for the QRU and Hawkes_U models.

Figure 17: Long range dependence.

###  A.3 Weibull Fit of Interarrival Time of Trades

Another known stylized fact is the Weibull fit of interarrival times of market
orders. Results in Figure 18 show that all the models verify this stylized
fact, as the Weibull distribution is the closest fit. Moreover, this result
aligns with previous empirical observations of the Bund (Bodor and Carlier,
2024), and is similar to findings for other assets such as equity (Abergel et
al., 2016). It is worth noting that this stylized fact is an important
consideration in the analysis of financial markets.

Figure 18: Distribution of interarrival times of market orders.

###  A.4 Excitation between events

Figure 19 presents the matrix of conditional probabilities of an event
$E_{t^{+}}$ (columns) given a past event $E_{t^{-}}$ (rows), denoted as
$\text{I\kern-1.49994ptP}\left(E_{t^{+}}|E_{t^{-}}\right)$. The matrices
reveal that the queue-reactive models inadequately replicate the empirical
matrices, as evidenced by nearly uniform rows in the matrix. This indicates
that the occurrence of a past event does not significantly influence the
likelihood of future events in these models. In contrast, Hawkes-based models
demonstrate a notable accuracy in reproducing the excitation between events,
attributable to their design and the inclusion of a memory component.

Figure 19: Transition probabilities between events.

###  A.5 Power-law of Order Sizes

The power-law fit of order sizes is also one of the stylized facts that is
verified in multiple markets (Chakraborti et al., 2011; Abergel et al., 2016;
Bouchaud et al., 2018). In particular, Bund futures market exhibits this
behavior (Bodor and Carlier, 2024). Figure 20 shows that all the models
satisfy this stylized fact, except for the QRU and Hawkes_U models, mainly
because of the unitary nature of the order sizes.

Figure 20: Distribution of trade order sizes and power-law fit.