# Single-pixel transmission matrix recovery via 2-photon fluorescence

Shupeng Zhao,1,2 Bernhard Rauer,1 Lorenzo Valzania,1 Jonathan Dong,3
Ruifeng Liu,2 Fuli Li,2 Sylvain Gigan,1 Hilton B. de Aguiar,∗1

1Laboratoire Kastler Brossel, ENS-Université PSL, CNRS, Sorbonne Université,
Collège de France. 24 rue Lhomond, 75005 Paris, France.
2Shaanxi Province Key Laboratory for Quantum Information and Quantum
Optoelectronic Devices,
and Department of Applied Physics, School of Science, Xi’an Jiaotong
University,
Xi’an 710049, China
3Biomedical Imaging Group, Ecole polytechnique fédérale de Lausanne (EPFL).
∗E-mail<EMAIL_ADDRESS>
> Imaging at depth in opaque materials has long been a challenge. Recently,
> wavefront shaping has enabled significant advance for deep imaging.
> Nevertheless, most non-invasive wavefront shaping methods require cameras,
> lack the sensitivity for deep imaging under weak optical signals, or can
> only focus on a single "guidestar". Here, we retrieve the transmission
> matrix (TM) non-invasively using two-photon fluorescence exploiting a
> general single-pixel detection framework, allowing to achieve single-target
> focus on multiple guidestars spread beyond the memory effect range. In
> addition, if we assume memory effect correlations exist in the transmission
> matrix, we are able to significantly reduce the number of measurements
> needed.

## Introduction

Fluorescence-based microscopy is the workhorse imaging modality in biology and
medicine. However, imaging opaque biological specimens at-depth is a current
challenge in optics (?, ?) as light passing through biological tissue
experiences an environment of heterogeneous refractive index, leading to
aberrations and scattering. To overcome these effects, feedback-based
wavefront shaping techniques have emerged as a new paradigm to inverse the
effect of scattering and focus light to a diffraction-limited spot (?). This
is achieved by iteratively modulating the input wavefront with a spatial light
modulator (SLM) aided by a feedback signal from the targeted point. However,
this correction is only valid for one focus position and the process has to be
repeated for another position. Retrieving the transmission matrix (TM), which
connects input and output planes of the scattering medium, brings more freedom
to focus light at any position of the target plane (?). Despite significant
advances in imaging by characterizing the transmission matrix of a system,
there are many applications scenarios where this cannot be neatly done by
simply imaging the output target plane. Non-invasive settings required in
biological imaging, typically using epi-geometry, are more challenging as they
require an undisturbed feedback from within the medium for optimization- or
TM-based schemes. In such cases, wavefront optimization requires the
assistance of beacons or so-called "guidestars" in the target plane (?).
Photoacoustic feedback (?, ?), and acousto-optic tagging (?, ?, ?, ?, ?, ?)
are compelling methods because of acoustic waves’s resilience to scattering
(?), which enables the reconstruction of a TM but requires complex super-
resolution methods to enhance the resolution that is limited by the sound
wavelength (?, ?, ?, ?). Another powerful method relies on the measurement of
a time-gated matrix in reflection (?, ?), but is based on retro-reflected
ballistic excitation photons, hence limited in depth because of their
exponential attenuation.

Fluorescence-based feedback signals can be easily separated from the
excitation wavelength, however reach focusing at the optical diffraction limit
has been a major challenge (?). Because the fluorescence signal is strongly
scattered by the complex medium, the main difficulty in wavefront shaping
experiments is how to assess whether the excitation energy was delivered to
only one guidestar. Optical non-linear fluorescence (?), or other
computational non-linear methods (?, ?, ?) have been proposed to solve this
problem. Once achieved focusing on a single guidestar, memory effect (ME) of
the scattering medium (?, ?) can be used for raster-scanning imaging (?).
Nevertheless, the imaging range is limited by the ME range. A potential
solution to overcome this limitation is to recover the TM using camera-based
methods, which have been used to demix the emission signal from the different
guidestars (?, ?). Summarizing current status (?), previous fluorescence based
strategies can either focus at one single position using single-pixel
detectors (?, ?, ?, ?, ?, ?, ?) (or cameras (?, ?)), or multiple positions by
reconstructing the TM using cameras (?, ?, ?).Thus, a method based on single-
pixel detection that can reconstruct the TM is urgently needed, as camera-
based methods lack sensitivity at-depth (due to weaker signal generation).

Here, we propose a general model of fluorescence microscopy to perform non-
invasive focusing and imaging at multiple target positions through the
retrieval of the coherent TM, using only a single-pixel detector. We start by
establishing a general framework and show that the forward model is only
solvable with non-linear fluorescence signals. In particular, we use a
multiplexed phase retrieval scheme using gradient-descent-based algorithm to
extract the TM. We demonstrate focusing on multiple fluorescing targets beyond
the memory effect. To speed up the measurement process, we consider the memory
effect which significantly decreases the amount of measurements needed.

## Principle

Figure 1: Experimental and reconstruction principles. (a) Simplified schematic
view of a general fluorescence microscope with a single-pixel detector. Random
wavefront generated by a SLM (each wavefront represents one row of matrix
$\mathbf{A}$) impinges on a scattering medium with a transmission matrix
$\mathbf{T}$ and excites 2-photon fluorescence (2PF) of the extended object
$\mathbf{o}$. The fluorescence signal $\mathbf{y}$ is collected by a single-
pixel detector (photo-multiplier tube (PMT)). Retrieval of the TM is achieved
in two different extreme cases 1 and 2, where tailored algorithms are used for
solving different models. (b) In case 1, $\mathbf{T}$ is assumed to have
completely uncorrelated elements, i.e. a random matrix. We solve for
$\mathbf{T}$ using as input the measurements $\mathbf{y}$ and the wavefront
modulation embedded in the user-input matrix $\mathbf{A}$. The reconstruction
of $\mathbf{T}$ allows to focus on each guidestar in the object. (c) In case
2, $\mathbf{T}$ has correlations, i.e. the infinite ME range is considered,
and is given by the pixel-wise product of a rank-1 matrix
$\mathbf{T}^{{}^{\prime}}$ and a "tilt matrix" $\mathbf{F}^{T}$, where
$\mathbf{F}$ is a matrix corresponding to 2D Fourier transform. An iterative
method is used where, in each iteration, we solve for
$\mathbf{T}^{{}^{\prime}}$ at first, then the wavefronts estimation are used
for raster-scanning the sample to update $\mathbf{o}$ until convergence.

We are interested in the situation described in Fig 1a: an extended incoherent
object hidden behind a scattering medium is illuminated by a series of random
speckle patterns, created by modulating an incident laser via spatial light
modulator. The total fluorescence is integrated by a single-pixel detector.
That is, we do not resolve the spatial fluctuations of the fluorescence
pattern. We can model this in a general matrix-form phase-retrieval framework,
where the detected signal $\mathbf{y}\in\mathbb{R}^{N_{\rm{m}}\times 1}$ is
described by

$\mathbf{y}=\lvert\mathbf{A}\mathbf{T}\rvert^{2p}\mathbf{o},$ (1)

where $\mathbf{A}\in\mathbb{C}^{N_{\rm{m}}\times N_{\rm{slm}}}$ is a known
measurement matrix, with each of the $N_{\rm{m}}$ row representing a
realization of modulated wavefront defined by $N_{\rm{slm}}$ pixels,
$\mathbf{T}\in\mathbb{C}^{N_{\rm{slm}}\times N_{\rm{obj}}}$ is a field
transmission matrix from the SLM to the fluorescent object plane,
$\mathbf{o}\in\mathbb{R}^{N_{\rm{obj}}\times 1}$, and $p$ the order of the
$p$-photon excited fluorescence signal. In what follows, we consider two
extreme cases: thick scattering medium with negligible ME range (the columns
in $\mathbf{T}$ are fully decorrelated, case 1 in Fig.1b) and thin scattering
medium with infinite ME range (the columns in $\mathbf{T}$ are fully
correlated, case 2 in Fig.1c).

The object being sparse, only the positions where the fluorescent response of
the object is non-zero contribute to the signal. Therefore we can consider
only some columns of $\mathbf{T}$ where the fluorescent response is non-zero.
To do that, we merge the vector $\mathbf{o}$ into the matrix $\mathbf{T}$ and
Eq. (1) becomes:

$\mathbf{y}=\sum_{i=1}^{N^{{}^{\prime}}_{\rm{obj}}}\lvert\mathbf{A}\mathbf{t}_{i}\rvert^{2p},$
(2)

where $N^{{}^{\prime}}_{\rm{obj}}<N_{\rm{obj}}$ represents the number of
columns left after merged with $\mathbf{o}$. $\mathbf{t}_{i}$ represent the
columns of the transmission matrix $\mathbf{T}$ weighted by the corresponding
elements in $\mathbf{o}^{\frac{1}{2p}}$. Note that $p=1$ corresponds to linear
fluorescence, and in this case the solution of Eq.(2) is not unique (?),
because there are infinite matrices $\mathbf{T}$ that satisfy the same
measurement $\mathbf{y}$ (see the Supplementary Information for detail).
Conversely, for $p=2$ the ambiguity is lifted and Eq.(2) has unique solution
provided that a certain number of measurements is fulfilled (see the
Supplementary Information for a mathematical proof). Therefore, as long as we
detect nonlinear fluorescence signals and collect enough random wavefront
measurements, the transmission matrix $\mathbf{T}$ can be uniquely determined.
Note that $\mathbf{A}$ is, in our case, obtained by displaying a series of
random patterns on a SLM, but could be achieved by any other deterministic
modulation scheme (galvanometric mirrors, digital micromirror devices etc).

In order to reconstruct $\mathbf{T}$ with 2-photon fluorescence signal, we use
a gradient-descent-based algorithm to solve Eq.(2) when $p=2$ (see methods for
further details on the Algorithm 1). Since we are solving a multiplexed phase
retrieval problem, the number of wavefronts measured $N_{\rm{m}}$ should scale
with the number of controlled modes $N_{\rm{slm}}$ (the dimension of unknown
signal) and the number of excitable speckles grains at the object plane where
the fluorescence response is nonzero: $N^{{}^{\prime}}_{obj}$ (the number of
unknown signals). It means the required measurement number depends on the
complexity of the object.

While this method allows to find a solution, including the presence of memory
effect into our forward model allows us to reduce the measurement number. As
showing in Fig.1c, under the infinite memory effect condition, there is only
phase ramp difference between columns of $\mathbf{T}$. These phase ramp can be
modeled by a "tilt matrix" $\mathbf{F}^{T}$, where $\mathbf{F}$ is a matrix
corresponding to 2D Fourier transform (Multiplying $\mathbf{F}$ with a vector
is equivalent to first reshaping the vector into a matrix, then doing 2D
Fourier transform, and finally reshaping it back into a vector.). We can then
rewrite Eq.(1) to

$\mathbf{y}=\lvert\mathbf{A}(\mathbf{F}^{T}\odot\mathbf{T}^{{}^{\prime}})\rvert^{4}\mathbf{o},$
(3)

where $\mathbf{T}^{{}^{\prime}}$ is a matrix that every column is the same
vector $\mathbf{t}$ that corresponds to a wavefront ($\odot$ stands for
element-wise product). In order to retrieve both $\mathbf{T}^{{}^{\prime}}$
and $\mathbf{o}$ in Eq.3, we propose a synergistic algorithm to ensure
convergence. As schematically shown in Fig. 1c: (1) collect measurement vector
$\mathbf{y}$ with the measurement matrix $\mathbf{A}$, and initialize
$\mathbf{o}$ by raster-scanning the sample, (2) update $\mathbf{t}$ by solving
Eq.(3) numerically, (3) update $\mathbf{o}$ by raster-scanning the sample with
the compensate wavefront $\rm{arg}(\mathbf{t}^{*})$, (4) repeat (2) and (3)
until a strong 2-photon signal is detected. (see methods for further details
on the Algorithm 2). Since the matrix $\mathbf{T}^{{}^{\prime}}$ has only one
vector $\mathbf{t}$ as unknown, the multiplexed phase retrieval problem
simplifies to a normal phase retrieval problem with only one column to be
reconstructed. Hence, the measurement number decreases considerably going from
$N_{\rm{m}}\propto N_{\rm{slm}}\times N^{{}^{\prime}}_{\rm{obj}}$ to typically
$N_{\rm{m}}\propto N_{\rm{slm}}$. Note here that instead of only using an SLM
to generate a measurement matrix $\mathbf{A}$, we can use a scan element and
SLM together to generate $\mathbf{A}$: this allows to speed up measurements
significantly (e.g. using galvanometric mirrors).

## Results

Figure 2: $\mathbf{T}$ retrieval beyond the ME. (a) 2PF signal detected
($\mathbf{y}$) vs. number of random wavefronts. Representative wavefronts used
are shown in the insets. The two red dashed line represents the average level
of the first and last 100 signal samples, respectively. (b) Inspected focii
(control CCD) for each column after retrieval of $\mathbf{T}$. The insets show
line profile. (c) Sum of all images for each focus demonstrating unique single
focus for each column retrieved. (d) Brightfield image taken with control
camera confirms the number of sources. (e) Imaging after analysis on
reconstructed $\mathbf{T}$ using correlations between the different columns
(The red circle in (c) depicts the memory effect range).

To demonstrate the feasibility of the methods, we have performed initial
experiments on sparse fluorescence samples. We modulate the wavefront sequence
using a SLM to define measurement matrix $\mathbf{A}$ (see method for details
on the setup). The SLM is located conjugate to the back focal plane of the
illumination objective. The fluorescence samples are drop-cast on a glass
coverslip whose bottom surface, facing the illumination objective, was sand-
blasted to introduce a single strongly scattering layer (?). After the
scattering medium and the sample, the 2PF signal $\mathbf{y}$ is collected by
a photomultiplier tube (PMT). An additional charge-coupled device (CCD) camera
is used only to inspect the quality of the focusing performance. A typical
realization of the detected 2PF signal is shown in Fig. 2a where the
fluctuations are due to the random wavefront realizations: on short-term
averages these fluctuations are constant not showing any bleaching, and are
furthermore well contrasted. We then solve Eq.(2) when $p=2$, from which we
retrieve $\mathbf{T}$ to find the wavefronts $\rm{arg}(\mathbf{t}^{*}_{i})$ to
focus on individual multiple targets (Fig.2b), or focii summed (Fig. 2c) as
inspected by the CCD camera. Bright field images of the sources match well
with the results (Fig. 2d). Residual correlations in $\mathbf{T}$ allows for
imaging without raster-scanning the focus (Fig. 1d) and beyond the memory
effect range (the red dashed circle in Fig.2c represents the memory effect
range). For this imaging example, we calculate the correlation of each two
columns of $\mathbf{T}$. Computing successively these pairwise correlations
between close targets allows retrieving the full object, well beyond a single
ME range (see Supplementary Information for further details). Importantly, we
did not use any assumption on the ME to retrieve the transmission matrix, but
only to retrieve a posteriori the relative position between sources.

Figure 3: $\mathbf{T}$ retrieval assuming ME correlations. (a) 2PF imaging by
raster scanning the sample using the speckles depicted in (b), themselves
captured with the CCD before each iteration. The images are low-resolution in
order to speed up acquisition. High resolution images before (c1) and after
two iterations (c2) show the cluster, and its morphology confirmed by the
brightfield camera (c3).

To speed up measurement and reconstruction time, we then studied a system that
has ME in $\mathbf{T}$, and explicitly take it into account in the algorithm.
Fig.3 shows the experimental result when 8 random patterns are used for
raster-scanning one by one to generate the measurement matrix $\mathbf{A}$.
The initialization of $\mathbf{o}$ is (Fig. 3a1) obtained by raster-scanning
the object without correcting for the scattering layer. In each iteration, we
solve Eq.(3) first, from which we retrieve $\mathbf{T}^{{}^{\prime}}$, itself
containing the wavefront $\mathbf{t}$. Then, $\rm{arg}(\mathbf{t}^{*})$ is
displayed on the SLM to correct the scattering layer for raster-scanning to
obtain an update of $\mathbf{o}$ (Fig. 3a2). This process repeats until a
desired stop criteria is achieved (in our case, a strong signal). Typically 3
iterations are needed to converge to a single focus as shown by the inspection
camera at each iteration (Fig.3b1-3). To ensure fast convergence, the
dimension of $\mathbf{o}$ is set equal to $N_{\rm{slm}}$. Comparison of a
high-resolution imaging taken after the optimization procedure (Fig. 3c2)
shows excellent agreement with the inspection camera (Fig. 3c3).

Figure 4:  Retrieval of $\mathbf{T}$ for extended objects. The samples in (a)
and (b) are fibres that coming from lens tissue marked by a fluorescence pen.
The samples in (c) is a bunch of 1$\mu$m beads. (a1), (b1) and (c1) are the
imaging result by scanning in $64\times 64$ steps without wavefront
correction. (a2), (b2) and (c2) are the corresponding results after the
wavefront correction. Here, we use 16 random patterns scanning to generate the
measurement matrix $\mathbf{A}$, and put $\mathbf{A}$ and $\mathbf{y}$ as
input for Alg.2. The samples in (a), (b) take 3 iterations to converge, and
(c) takes 5 iterations. (a3), (b3) and (c3) are the corresponding bright field
images.

Inspection at each iteration shows that the speckle converges to a single foci
despite the size of the object. Fig.4 shows the experimental results for
extended objects. Fig.4a-b are fluorescence-dyed fibres into orthogonal
directions spanning the whole speckle envelope. Before optimization, clear
replica of the object are seen because of the sparse 2PF speckle excitation.
After optimization, convergence to a single focus is ensured. Even in extreme
cases where the sample extend quasi homogeneously the speckle envelope
(Fig.4c), there is clearly a convergence to a single focus as shown by
comparison with brightfield images (Fig.4 c3, control camera). In these
examples, because we reconstruct an extended complex $\mathbf{o}$, to speed up
convergence, we use 16 random patterns raster-scanning to generate
$\mathbf{A}$ which weakly scales with the object complexity. The three results
showing here take 3, 3 and 5 iterations to converge, respectively.

## Discussion

We presented a method to non-invasively extract the transmission matrix of a
buried fluorescent object, using a single-pixel detector. To validate the
principle, we conducted an experiment by placing a single-pixel detector on
one side of the object. However, in practice, it can also be placed on the
other side to enable epi-detection mode. In both scenarios, the forward model
of the system remains the same. Based on a phase-retrieval formalism, we
exploit two extreme cases on the level of correlations embedded in
$\mathbf{T}$. If no correlations exists (as it would be the case for a very
optically-thick medium), we show that retrieval of $\mathbf{T}$ is unique for
nonlinear optical processes, such as two-photon fluorescence signals. The
downside of this approach is that it requires a number of measurements that
scales with the object complexity, which leads to lengthy acquisitions. We can
largely circumvent this limitation if $\mathbf{T}$ displays correlations, in
the form of memory effect, by introducing a tilt matrix $\mathbf{F}^{T}$.
Rigorous solution is achieved by updating $\mathbf{o}$ at each new estimation
of $\mathbf{T}$. In this case, the number of wavefront measurements does not
scale with object sparsity, so measurements can be significantly faster.

Two kind of methods exist for two-photon imaging through scattering medium.
The first is known as wavefront shaping based method (?, ?), which need fast
SLMs in order to ultimately outpace the speckle decorrelation time imposed by
living tissue. The second is focus scanning holographic aberration probing
(F-SHARP) (?) and can be viewed as an interferometric method. Both, and works
that follows (?, ?), are optimization based methods and can only image within
the memory effect range. We propose a method for imaging beyond the memory
effect range, although it needs more measurements than optimization based
methods, the data processing is not needed between the measurements which
leaves a great possibility to improve the measurement speed. In addition, due
to the random illumination strategy, the proposed method can prevent bleach
compared to optimization based methods. A particular advantage here is that,
when the imaging object is smaller than the ME range, the proposed method can
realize faster focusing speed compared to wavefront shaping methods as in
principle galvanometric and digital micromirror devices can be used, and have
better stability compared to interferometric methods.

## Methods

### Experimental setup and sample preparation

The experimental setup is sketched in Fig.8. The source used is an OPO
(Coherent Mira OPO-X) pumped by a pulsed Ti-sapphire laser (Coherent Chameleon
Ultra II) with 140 fs pulse duration and repetition rate of 80 MHz. For 2P
excitation, the OPO is tuned to 1050 nm. To modulate the wavefront, we use a
liquid crystal phase-only SLM (Meadowlarks HSP512L) which is imaged on the
back focal plane of a water immersion objective (Zeiss W Plan-Apochromat
40x/1.0 DIC M27). After illuminating the sample, a second microscope objective
(Zeiss EC Plan-NEOFLUAR 40x/1.3) images the excitation light pattern onto a
CCD camera (Basler acA1300-30$\rm{\mu}$m) and collects the fluorescence light
which is steered towards a PMT (Hamamatsu H7422P-40) by a dichrocic mirror
(Semrock Di03-R785-t1-25x36). Additional filters in front of the PMT (Thorlabs
FESH0600, FESH0700) ensure that only nonlinear fluorescent light is detected.

For the scattering medium, a rough surface is introduced on a #1.5 cover slip
by sandblasting its bottom surface with 220 grit sand. The beads (Fluospheres
carboxylate, 1.0$\rm{\mu}$m, red (580nm/605nm)) are drop-casted on the top
surface of the coverslip, about 170$\rm{\mu}$m above the scattering layer.
They are covered by the NOA65 glue to create a near-index matched environment.
Then, a second #1.5 cover slip is placed on top of the beads. Finally, the
glue is cured under UV light, the curing process is kept short to avoid
bleaching the beads. The fibres are provided by Thorlab lens tissue marked by
a orange fluorescence pen, they are directly covered by a second #1.5 cover
slip on top.

### Reconstruction details

input : number of gradient descent, $nn$; step size, $\Delta$; measurement
matrix, $\mathbf{A}$; measurement vector, $\mathbf{y}$; column number,
$N^{{}^{\prime}}_{\rm{obj}}$.

output : columns of transmission matrix, $\mathbf{t}_{i}$ ($i=1,2\cdots
N^{{}^{\prime}}_{\rm{obj}}$).

1 for  _count $\leftarrow 1$ to $nn$_ do

2

3 $\mathbf{t}_{i}=\mathbf{t}_{i}-\Delta\frac{\frac{\partial
f}{\partial\mathbf{t}^{*}_{i}}}{\lVert\frac{\partial
f}{\partial\mathbf{t}^{*}_{i}}\rVert}_{2}$

4

5 end for

Algorithm 1 Multiplexed phase retrieval

In order to reconstruct $\mathbf{t}_{i}$ in Eq.(2) when $p=2$. The Alg.1 is
proposed to solve the following optimized problem:

$\rm{min}_{\mathbf{t}_{i}}f,\quad\rm{with}\quad
f=\lVert\mathbf{y}-\sum_{i=1}^{N^{{}^{\prime}}_{obj}}\lvert\mathbf{A}\mathbf{t}_{i}\rvert^{4}\rVert^{2}_{2},$
(4)

where $N^{{}^{\prime}}_{\rm{obj}}$ is chosen to be equal or larger than the
actual number of sources behind the scattering medium. Some of the columns of
the reconstructed $\mathbf{T}$ will be repeated if
$N^{{}^{\prime}}_{\rm{obj}}$ is chosen larger (see Supplementary Information
for more details). For the initialization, we randomly initialize each column
$\mathbf{t}_{i}$ of the transmission matrix. According to (?, ?, ?), the
derivative of the complex quadratic cost function with respect to
$\mathbf{t}^{*}_{i}$, i.e. $\frac{\partial f}{\partial\mathbf{t}^{*}_{i}}$, is
necessary for updating $\mathbf{t}_{i}$ in each iteration, and $\frac{\partial
f}{\partial\mathbf{t}^{*}_{i}}$ can be easily caculated according to the
Wirtinger derivatives as

$\frac{\partial
f}{\partial\mathbf{t}^{*}_{i}}=4\mathbf{A}^{{\dagger}}[(\sum_{i=1}^{N^{{}^{\prime}}_{\rm{obj}}}\lvert\mathbf{A}\mathbf{t}_{i}\rvert^{4}-\mathbf{y})\odot{\lvert\mathbf{A}\mathbf{t}_{i}\rvert}^{2}\odot(\mathbf{A}\mathbf{t}_{i})].$
(5)

For the result in Fig.2, $16\times N^{{}^{\prime}}_{\rm{obj}}\times
N_{\rm{slm}}$ random patterns are displayed on the SLM to obtain $\mathbf{y}$,
where $N_{\rm{slm}}=16\times 16$ and $N^{{}^{\prime}}_{\rm{obj}}=7$. The
reconstruction of Alg.1 takes less than 10 seconds on GPU.

input : number of gradient descent, $nn$; step size, $\Delta$; measurement
matrix, $\mathbf{A}$; measurement vector, $\mathbf{y}$.

output : aberration, $\mathbf{t}$; object, $\mathbf{o}$.

1 while  _2-photon signal is low_ do

2 for  _count $\leftarrow 1$ to $nn$_ do

3

4 $\mathbf{t}=\mathbf{t}-\Delta\frac{\frac{\partial
f}{\partial\mathbf{t}^{*}}}{\lVert\frac{\partial
f}{\partial\mathbf{t}^{*}}\rVert_{2}}$

5

6 end for

7

8 Update $\mathbf{o}$ by raster-scanning the sample with the compensate
wavefront $\rm{arg}(\mathbf{t}^{*})$

9 end while

Algorithm 2 Phase retrieval + Scanning

In order to reconstruct $\mathbf{T}^{{}^{\prime}}$ concatenated by vector
$\mathbf{t}$ in Eq.(3), the Alg.2 is proposed to solve the following optimized
problem:

$\rm{min}_{\mathbf{t},\mathbf{o}}f,\quad\rm{with}\quad
f=\sum^{N_{m}}_{q=1}(y_{q}-\mathbf{o}^{T}\lvert\mathbf{F}(\mathbf{a}^{*}_{q}\odot\mathbf{t})\rvert^{4})^{2},$
(6)

where $\mathbf{a}^{{\dagger}}_{q}$ represents the $q$th row of $\mathbf{A}$,
$y_{q}$ represents the $q$th element of $\mathbf{y}$. For the initialization,
$\mathbf{t}$ is randomly initialized, and $\mathbf{o}$ is obtained by raster-
scanning the object with a scanning device (applying phase ramps on the SLM in
our case). The gradient $\frac{\partial f}{\partial\mathbf{t}^{*}}$ is derived
as:

$\begin{split}\frac{\partial
f}{\partial\mathbf{t}^{*}}&=4\sum^{N_{m}}_{q=1}(\mathbf{o}^{T}\lvert\mathbf{F}(\mathbf{a}^{*}_{q}\odot\mathbf{t})\rvert^{4}-y_{q})\mathbf{a}_{q}\\\
&\odot\\{\mathbf{F}^{{\dagger}}[\mathbf{o}\odot\lvert\mathbf{F}(\mathbf{a}^{*}_{q}\odot
t)\rvert^{2}\odot(\mathbf{F}(\mathbf{a}^{*}_{q}\odot\mathbf{t}))]\\}.\end{split}$
(7)

However, the auto-gradient derived by Pytorch are used in our code for faster
computation speed. In each iteration, $nn$ times of gradient descent are
performed to reconstruct the $\mathbf{t}$, then the $\mathbf{o}$ is updated by
raster-scanning the sample with a scanning device and displaying compensate
pattern $\rm{arg}(\mathbf{t}^{*})$ on the SLM simultaneously (applying a
complete set of phase ramps on top of $\rm{arg}(\mathbf{t}^{*})$ on the SLM in
our case). For the results in Fig.3 and Fig.4, $8\times N_{\rm{slm}}$ and
$16\times N_{\rm{slm}}$ measurements are taken to get $\mathbf{y}$
respectively. They are obtained by using 8 and 16 random patterns on the SLM
and raster-scanning the sample with them one by one, $N_{\rm{slm}}=16\times
16$ measurements for each raster-scanning. The sampling strategy here can
significantly reduce the sampling time if the galvo-mirror are used as
scanning device and placed on the conjugate plane of the SLM. The higher
resolution images in Fig.3 and Fig.4 are obtained by raster-scanning the
sample with SLM in $64\times 64$ steps. In each iteration, the reconstruction
of $\mathbf{t}$ take less than 10 seconds on GPU.

The computational resource we used: GPU, 2 Nvidia GeForce RTX 2080 Ti; RAM, 64
GB.

## Code available

Simulation codes are available at
https://github.com/comediaLKB/Single_pixel_TM_recovery.

## Acknowledgements

L.V. acknowledges funding from the Swiss National Science Foundation (project
P400P2_199329). S.Z. acknowledges funding from the China Scholarship Council
(CSC) (202006280163). This project was funded by the European Union’s Horizon
2020 research and innovation program under the FET-Open grant No. 863203
(Dynamic). The authors acknowledge the assistance of Eric Bezzam in optimizing
the code.

## Competing interests

The authors declare no competing interests.

## References

  * 1. J. Bertolotti, O. Katz, Nature Physics 18, 1008 (2022).
  * 2. S. Gigan, et al., Journal of Physics: Photonics 4, 042501 (2022).
  * 3. I. M. Vellekoop, A. Mosk, Optics letters 32, 2309 (2007).
  * 4. S. M. Popoff, et al., Physical review letters 104, 100601 (2010).
  * 5. R. Horstmeyer, H. Ruan, C. Yang, Nature photonics 9, 563 (2015).
  * 6. F. Kong, et al., Optics letters 36, 2053 (2011).
  * 7. T. Chaigne, et al., Nature Photonics 8, 58 (2014).
  * 8. J. W. Tay, P. Lai, Y. Suzuki, L. V. Wang, Scientific reports 4, 1 (2014).
  * 9. K. Si, R. Fiolka, M. Cui, Nature photonics 6, 657 (2012).
  * 10. H. Ruan, M. Jang, B. Judkewitz, C. Yang, Scientific reports 4, 1 (2014).
  * 11. Y. M. Wang, B. Judkewitz, C. A. DiMarzio, C. Yang, Nature communications 3, 928 (2012).
  * 12. J. Tang, R. N. Germain, M. Cui, Proceedings of the National Academy of Sciences 109, 8434 (2012).
  * 13. X. Xu, H. Liu, L. V. Wang, Nature photonics 5, 154 (2011).
  * 14. O. Katz, F. Ramaz, S. Gigan, M. Fink, Nature communications 10, 717 (2019).
  * 15. D. B. Conkey, et al., Nature communications 6, 7902 (2015).
  * 16. B. Judkewitz, Y. M. Wang, R. Horstmeyer, A. Mathy, C. Yang, Nature photonics 7, 300 (2013).
  * 17. K. Si, R. Fiolka, M. Cui, Scientific reports 2, 1 (2012).
  * 18. A. Badon, et al., Science advances 2, e1600370 (2016).
  * 19. S. Jeong, et al., Nature Photonics 12, 277 (2018).
  * 20. O. Katz, E. Small, Y. Guan, Y. Silberberg, Optica 1, 170 (2014).
  * 21. A. Boniface, B. Blochet, J. Dong, S. Gigan, Optica 6, 1381 (2019).
  * 22. D. Aizik, I. Gkioulekas, A. Levin, arXiv preprint arXiv:2205.06215 (2022).
  * 23. D. Aizik, A. Levin, arXiv preprint arXiv:2301.11421 (2023).
  * 24. S. Feng, C. Kane, P. A. Lee, A. D. Stone, Physical review letters 61, 834 (1988).
  * 25. I. Freund, M. Rosenbluh, S. Feng, Physical review letters 61, 2328 (1988).
  * 26. A. Boniface, J. Dong, S. Gigan, Nature communications 11, 1 (2020).
  * 27. A. d’Arco, F. Xia, A. Boniface, J. Dong, S. Gigan, Optics Express 30, 30845 (2022).
  * 28. Y. Liu, et al., Nature communications 6, 1 (2015).
  * 29. I. N. Papadopoulos, J.-S. Jouhanneau, J. F. Poulet, B. Judkewitz, Nature Photonics 11, 116 (2017).
  * 30. I. N. Papadopoulos, et al., Light: Science & Applications 9, 1 (2020).
  * 31. M. A. May, et al., Nature Communications 12, 1 (2021).
  * 32. A. Daniel, D. Oron, Y. Silberberg, Optics Express 27, 21778 (2019).
  * 33. Y. Baek, H. B. de Aguiar, S. Gigan, arXiv preprint arXiv:2301.13140 (2023).
  * 34. J. Dong, F. Krzakala, S. Gigan, ICASSP 2019-2019 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP) (IEEE, 2019), pp. 4963–4967.
  * 35. B. Rauer, H. B. de Aguiar, L. Bourdieu, S. Gigan, Optics Letters 47, 6233 (2022).
  * 36. W. Bridges, et al., Applied Optics 13, 291 (1974).
  * 37. R. Remmert, Theory of complex functions, vol. 122 (Springer Science & Business Media, 1991).
  * 38. R. F. Fischer, Precoding and signal shaping for digital transmission (John Wiley & Sons, 2005).
  * 39. L. Bian, et al., Optics express 23, 4856 (2015).
  * 40. L. Zhu, et al., Nature communications 13, 1 (2022).
  * 41. F. Soldevila, et al., arXiv preprint arXiv:2302.06519 (2023).

## Supplementary Information:

### The uniqueness of the solution depending on $p$

To aid in understanding, we consider a scenario that only 2 beads under the
scattering medium. We aim to demonstrate that, given a sufficient number of
random measurements, the solution to Eq.(2) in the main text is non-unique
when $p=1$, while it is unique when $p=2$. It‘s not hard to see that the
conclusion can be generalized to the case of more than 2 beads.

The corresponding 2 columns of $\mathbf{T}\in\mathbb{C}^{N_{\rm{slm}}\times
2}$ are $\mathbf{t}_{1}$ and $\mathbf{t}_{2}$, respectively.
$\mathbf{a}^{{\dagger}}\in\mathbb{C}^{1\times N_{\rm{slm}}}$ represents one
row of the measurement matrix $\mathbf{A}$. and the corresponding fluorescence
signal $y$ is:

$y=\lvert\mathbf{a}^{{\dagger}}\mathbf{t}_{1}\rvert^{2p}+\lvert\mathbf{a}^{{\dagger}}\mathbf{t}_{2}\rvert^{2p},$
(8)

Assuming $\mathbf{t}^{{\dagger}}_{1}\mathbf{t}_{2}=0$, which means the columns
of transmission matrix are orthogonal to each other (since the transmission
matrix are complex random matrix, this a weak assumption). We can then
represent inner product in a canonical basis starting with

$\mathbf{l}_{1}=\begin{pmatrix}1\\\ 0\\\ 0\\\
\vdots\end{pmatrix}\quad\rm{and}\quad\mathbf{l}_{2}=\begin{pmatrix}0\\\ 1\\\
0\\\ \vdots\end{pmatrix}.$ (9)

Eq.(8) in this basis becomes:

$\begin{split}y&=\lvert(\mathbf{R}\mathbf{a})^{{\dagger}}\mathbf{R}\mathbf{t}_{1}\rvert^{2p}+\lvert(\mathbf{R}\mathbf{a})^{{\dagger}}\mathbf{R}\mathbf{t}_{2}\rvert^{2p}\\\
&=\lVert\mathbf{t}_{1}\rVert^{2p}\lvert\mathbf{\hat{a}}^{{\dagger}}\mathbf{l}_{1}\rvert^{2p}+\lVert\mathbf{t}_{2}\rVert^{2p}\lvert\mathbf{\hat{a}}^{{\dagger}}\mathbf{l}_{2}\rvert^{2p}\\\
&=\lVert\mathbf{t}_{1}\rVert^{2p}\lvert\hat{a}_{1}\rvert^{2p}+\lVert\mathbf{t}_{2}\rVert^{2p}\lvert\hat{a}_{2}\rvert^{2p},\end{split}$
(10)

where $\mathbf{R}$ represents a unitary basis rotation matrix, and
$\hat{a}_{i}$ represents the $i$th element of the vector $\mathbf{\hat{a}}$.

The question now is: Are there any
$(\mathbf{u}_{1},\mathbf{u}_{2})\neq(\mathbf{l}_{1},\mathbf{l}_{2})$ that
satisfy the following formula?

$\forall\mathbf{\hat{a}}^{{\dagger}}\quad\lVert\mathbf{t}_{1}\rVert^{2p}\lvert\hat{a}_{1}\rvert^{2p}+\lVert\mathbf{t}_{2}\rVert^{2p}\lvert\hat{a}_{2}\rvert^{2p}=\lvert\mathbf{\hat{a}}^{{\dagger}}\mathbf{u}_{1}\rvert^{2p}+\lvert\mathbf{\hat{a}}^{{\dagger}}\mathbf{u}_{2}\rvert^{2p}.$
(11)

In the $p=1$ case, for any $2\times 2$ unitary matrix $\mathbf{S}$ we can find
a set of solution:

$\begin{pmatrix}\mathbf{u}_{1}\\\ \mathbf{u}_{2}\\\
\end{pmatrix}=\mathbf{S}\begin{pmatrix}\lVert\mathbf{t}_{1}\rVert\mathbf{l}_{1}\\\
\lVert\mathbf{t}_{2}\rVert\mathbf{l}_{2}\\\ \end{pmatrix}$ (12)

satisfy Eq.(11). For example,
$(\mathbf{u}_{1},\mathbf{u}_{2})=(\frac{1}{\sqrt{2}}(\lVert\mathbf{t}_{1}\rVert\mathbf{l}_{1}+\lVert\mathbf{t}_{2}\rVert\mathbf{l}_{2}),\frac{1}{\sqrt{2}}(\lVert\mathbf{t}_{1}\rVert\mathbf{l}_{1}-\lVert\mathbf{t}_{2}\rVert\mathbf{l}_{2}))$.
Thus, the Eq.(2) in the main text has infinite solution when $p=1$.

In the $p=2$ case, Eq.(11) becomes

$\forall\mathbf{\hat{a}}^{{\dagger}}\quad\lVert\mathbf{t}_{1}\rVert^{4}\lvert\hat{a}_{1}\rvert^{4}+\lVert\mathbf{t}_{2}\rVert^{4}\lvert\hat{a}_{2}\rvert^{4}=\lvert\mathbf{\hat{a}}^{{\dagger}}\mathbf{u}_{1}\rvert^{4}+\lvert\mathbf{\hat{a}}^{{\dagger}}\mathbf{u}_{2}\rvert^{4}.$
(13)

Let $\hat{a}_{1}=\hat{a}_{2}=0$ and $\hat{a}_{i}=\forall,i\neq 1\mkern
8.0mu\rm{or}\mkern 8.0mu2$, we have

$\mathbf{u}_{1}=\begin{pmatrix}\alpha\\\ \beta\\\ 0\\\ 0\\\ \vdots\\\
\end{pmatrix}\quad\mathbf{u}_{2}=\begin{pmatrix}\gamma\\\ \delta\\\ 0\\\ 0\\\
\vdots\\\ \end{pmatrix}.$ (14)

Let $\hat{a}_{1}=\forall$, $\hat{a}_{2}=0$, we have

$\lvert\alpha\rvert^{4}+\lvert\gamma\rvert^{4}=\lVert\mathbf{t}_{1}\rVert^{4}.$
(15)

Let $\hat{a}_{1}=0$, $\hat{a}_{2}=\forall$, we have

$\lvert\beta\rvert^{4}+\lvert\delta\rvert^{4}=\lVert\mathbf{t}_{2}\rVert^{4}.$
(16)

Let $\hat{a}_{1}=\hat{a}_{2}=\forall$, we have

$\lvert\alpha+\beta\rvert^{4}+\lvert\gamma+\delta\rvert^{4}=\lVert\mathbf{t}_{1}\rVert^{4}+\lVert\mathbf{t}_{2}\rVert^{4}.$
(17)

Let $\hat{a}_{1}=-\hat{a}_{2}=\forall$, we have

$\lvert\alpha-\beta\rvert^{4}+\lvert\gamma-\delta\rvert^{4}=\lVert\mathbf{t}_{1}\rVert^{4}+\lVert\mathbf{t}_{2}\rVert^{4}.$
(18)

Eq.(17) plus Eq.(18) and using Eq.(15) and Eq.(16) we have

$(\alpha\beta^{*})^{2}+(\alpha^{*}\beta)^{2}+4\lvert\alpha^{*}\beta\rvert^{2}+(\gamma\delta^{*})^{2}+(\gamma^{*}\delta)^{2}+4\lvert\gamma^{*}\delta\rvert^{2}=0.$
(19)

Defining

$\alpha\beta^{*}=m+in,\quad\alpha^{*}\beta=m-in,\quad\gamma\delta^{*}=p+iq,\quad\gamma^{*}\delta=p-iq,$
(20)

into Eq.(19), where $i$ represents the imaginary unit and $m,n,p,q$ are real,
we have

$6m^{2}+2n^{2}+6p^{2}+2q^{2}=0.$ (21)

Thus, $m=n=p=q=0$. It means $\alpha\beta^{*}=\gamma\delta^{*}=0$. So
$\alpha\mkern 8.0mu\rm{or}\mkern 8.0mu\beta=0$ and $\gamma\mkern
8.0mu\rm{or}\mkern 8.0mu\delta=0$. Therefore, we conclude that

$(\mathbf{u}_{1},\mathbf{u}_{2})=(\mathbf{l}_{1},\mathbf{l}_{2})\mkern
8.0mu\rm{or}\mkern 8.0mu(\mathbf{l}_{2},\mathbf{l}_{1}).$ (22)

Thus, Eq.(2) in the main text has unique solution when $p=2$.

###  $\mathbf{T}$ reconstruction quality for using different
$N^{{}^{\prime}}_{\rm{obj}}$

In order to estimate the $\mathbf{T}$ reconstruction quality using different
$N^{{}^{\prime}}_{\rm{obj}}$, three kind of random matrix $T$ are generated
for ground truth, the number of columns they contain is 3, 7 and 11
respectively. Different $N^{{}^{\prime}}_{\rm{obj}}$ ranging from 1 to 14 are
using in Alg.1 to reconstruct $\mathbf{T}$. The reconstruction fidelity is
shown in Fig.5, the red, blue and green are corresponding to the ground truth
$\mathbf{T}$ with 3, 7 and 11 columns respectively.

The fidelity of each column $\mathbf{t}_{i}$ is calculated between the
reconstructed $\mathbf{t}_{i}$ and every column of the ground-truth
$\mathbf{T}$ then pick the maximum. The fidelity of the transmission matrix
$\mathbf{T}$ is the average fidelity of all the column
$\mathbf{t}_{i},i=1,2,\cdots,N^{{}^{\prime}}_{\rm{obj}}$. As showing in the
result, when the $N^{{}^{\prime}}_{\rm{obj}}$ is chosen smaller than the
sources number, the reconstructed column $\mathbf{t}_{i}$ is not guaranteed to
match one of the column in the ground-truth matrix $\mathbf{T}$. On the
contrary, when the $N^{{}^{\prime}}_{\rm{obj}}$ is chosen bigger than the
sources number, each reconstructed $t_{i}$ match one of the column in the
ground-truth $\mathbf{T}$. Actually, in this case some of the columns of the
reconstructed $\mathbf{T}$ will be repeated. In the simulation, $N_{m}=65536$
and $N_{\rm{slm}}=16\times 16$. The each solid line represents the average of
10 realization, and the shadow represents the fidelity range distribution.

### Transmission matrix based image reconstruction

Previous works that realized imaging beyond the memory effect (ME) by
analysing the retrieved distance map between all the emitters (?, ?, ?, ?),
but they all used a camera to retrieve all the emission fluorescent patterns
first, which is not sensitive enough for non-linear microscopy. Here, we
retrieved the relative position between every two emitters by using the
transmission matrix $\mathbf{T}$. Fig.6 shows the image reconstruction work
flow by analysing the $T$ we retrieved in result Fig.2e in the manuscript. If
the two emitters are really close, the corresponding two columns in
$\mathbf{T}$ should be highly correlated and provide relative position
information. The relative position map of the two emitters is calculated by:

$o_{i,k}=\lvert\mathcal{F}_{2d}(\frac{\mathbf{t}^{*}_{i}}{\lvert\mathbf{t}^{*}_{i}\rvert}\odot\frac{\mathbf{t}_{k}}{\lvert\mathbf{t}_{k}\rvert})\rvert^{2},$
(23)

where $\mathcal{F}_{2d}$ represent the 2d Fourier transform, $\odot$ means
element wise production. To infer whether the columns $\mathbf{t}_{i}$ and
$\mathbf{t}_{k}$ are within or beyond the ME range, we study
$\alpha=\frac{max\\{o_{i,k}\\}}{max\\{o_{k,k}\\}}$ as a function of relative
distance, where $max\\{o_{i,k}\\}$ stands for maximum value of $o_{i,k}$. A
given threshold, $\alpha_{tres}$ of $\alpha$ is introduced to evaluate it. For
example, if $\alpha$ is greater than $\alpha_{tres}$, the two columns belong
to the same ME range. Otherwise, they belong to different ME ranges, the
relative position map $o_{i,k}$ is reassigned to a zero matrix. The pairwise
relative position map $o_{i,k}$ after thresholding is showing in Fig.6.

Then the sub images of the object in the vicinity of the emitter $k$ are
obtained by $O_{k}=\sum_{i=1}^{N^{{}^{\prime}}_{\rm{obj}}}o_{i,k}$. In order
to merge the sub images to a global image, we also need the relative position
vector $\mathbf{r}_{i,k}$ of the two emitters $i$ and $k$. If emitters $i$ and
$k$ are in the same ME range, $\mathbf{r}_{i,k}$ is a vector in the map
$o_{i,k}$ that origin at the center and end at the narrow delta-like peak. If
emitters $i$ and $k$ are beyond the ME range but emitter $j$ is between them,
we can still calculate the $\mathbf{r}_{i,k}$ between them as
$\mathbf{r}_{i,k}=\mathbf{r}_{i,j}+\mathbf{r}_{j,k}$. Therefore, even if the
ensemble of emitters expands well beyond the ME range, the full spatial
distribution can be recovered if the different isoplanatic patches are
“connected" by emitters.

Finally, we use the method proposed in (?) to merge all the sub images, and
realize imaging beyond ME range. In short, the global reconstruction of
$O_{global}$ can be obtained by composing all the sub images $O_{k}$ into one
images, taking into account their relative positions with respect to the first
emitter, ${\mathbf{r}}_{k,1}$:

$O_{global}=\sum_{k=1}^{N^{{}^{\prime}}_{\rm{obj}}}O_{k}(\mathbf{r}-\mathbf{r}_{k,1}).$
(24)

Therefore, a image beyond the ME range is obtained. The red dashed circle
showing in the reference focusing image in Fig.6 represents the ME range.
Another result with 8 fluorescence beads is shown in Fig.7.

Figure 5: The fidelity of retrieval $\mathbf{T}$ with different
$N^{{}^{\prime}}_{\rm{obj}}$: Three different kind of random matrix
$\mathbf{T}$ are generated for ground truth, they have 3 (red data), 7 (blue
data) and 11 (green data) columns respectively. $\mathbf{T}$ reconstruction
fidelity are evaluated when using different $N^{{}^{\prime}}_{\rm{obj}}$. The
each solid line represents the average of 10 realization, and the shadow
represents the fidelity range distribution. Figure 6: Imaging procedure of
Fig.2e in the manuscript: At first, for each column $\mathbf{t}_{i}$, the
relative position map $o_{i,k}$ with column $\mathbf{t}_{k}$ is calculated and
thresholded. Then, the sub images $O_{i}$ centered by the emitter
corresponding to $\mathbf{t}_{i}$ are obtained by summing all the
corresponding relative position maps. Finally, the global image is merged by
all the sub images. Figure 7: Image merge strategy based on relative position
vector $\mathbf{r}_{i,k}$. (a) and (b) are relative position maps $o_{1,1}$
and $o_{1,6}$. (c) 3 sub images $O_{1},O_{6},O_{2}$ are stitched by 2 relative
position vector $\mathbf{r}_{1,6}$ and $\mathbf{r}_{6,2}$. (d) The
$O_{global}$ after image merge. (e) Brightfield image taken with control
camera confirms the number of sources. (f) Sum of all images for each focus
demonstrating unique single focus for each column retrieved. Figure 8:
Experimental setup: SLM, spatial light modulator; Sample, sample; SM,
scattering medium; DBS, dichroic beam splitters; CCD, charge-coupled device;
PMT, photo-multiplier tube


# Submission and Formatting Instructions for
International Conference on Machine Learning (ICML 2022)

Firstname1 Lastname1  Firstname2 Lastname2  Firstname3 Lastname3  Firstname4
Lastname4  Firstname5 Lastname5  Firstname6 Lastname6  Firstname7 Lastname7
Firstname8 Lastname8  Firstname8 Lastname8

###### Abstract

This document provides a basic paper template and submission guidelines.
Abstracts must be a single paragraph, ideally between 4–6 sentences long.
Gross violations will trigger corrections at the camera-ready phase.

Machine Learning, ICML

##  1 Electronic Submission

Submission to ICML 2022 will be entirely electronic, via a web site (not
email). Information about the submission process and LaTeX templates are
available on the conference web site at:

http://icml.cc/

The guidelines below will be enforced for initial submissions and camera-ready
copies. Here is a brief summary:

  * •

Submissions must be in PDF.

  * •

New to this year: If your paper has appendices, submit the appendix together
with the main body and the references as a single file. Reviewers will not
look for appendices as a separate PDF file. So if you submit such an extra
file, reviewers will very likely miss it.

  * •

Page limit: The main body of the paper has to be fitted to 8 pages, excluding
references and appendices; the space for the latter two is not limited. For
the final version of the paper, authors can add one extra page to the main
body.

  * •

Do not include author information or acknowledgements in your initial
submission.

  * •

Your paper should be in 10 point Times font.

  * •

Make sure your PDF file only uses Type-1 fonts.

  * •

Place figure captions _under_ the figure (and omit titles from inside the
graphic file itself). Place table captions _over_ the table.

  * •

References must include page numbers whenever possible and be as complete as
possible. Place multiple citations in chronological order.

  * •

Do not alter the style template; in particular, do not compress the paper
format by reducing the vertical spaces.

  * •

Keep your abstract brief and self-contained, one paragraph and roughly 4–6
sentences. Gross violations will require correction at the camera-ready phase.
The title should have content words capitalized.

###  1.1 Submitting Papers

Paper Deadline: The deadline for paper submission that is advertised on the
conference website is strict. If your full, anonymized, submission does not
reach us on time, it will not be considered for publication.

Anonymous Submission: ICML uses double-blind review: no identifying author
information may appear on the title page or in the paper itself. Section 2.3
gives further details.

Simultaneous Submission: ICML will not accept any paper which, at the time of
submission, is under review for another conference or has already been
published. This policy also applies to papers that overlap substantially in
technical content with conference papers under review or previously published.
ICML submissions must not be submitted to other conferences and journals
during ICML’s review period. Informal publications, such as technical reports
or papers in workshop proceedings which do not appear in print, do not fall
under these restrictions.

Authors must provide their manuscripts in PDF format. Furthermore, please make
sure that files contain only embedded Type-1 fonts (e.g., using the program
pdffonts in linux or using File/DocumentProperties/Fonts in Acrobat). Other
fonts (like Type-3) might come from graphics files imported into the document.

Authors using Word must convert their document to PDF. Most of the latest
versions of Word have the facility to do this automatically. Submissions will
not be accepted in Word format or any format other than PDF. Really. We’re not
joking. Don’t send Word.

Those who use LaTeX should avoid including Type-3 fonts. Those using latex and
dvips may need the following two commands:

    dvips -Ppdf -tletter -G0 -o paper.ps paper.dvi
    ps2pdf paper.ps

It is a zero following the “-G”, which tells dvips to use the config.pdf file.
Newer TeX distributions don’t always need this option.

Using pdflatex rather than latex, often gives better results. This program
avoids the Type-3 font problem, and supports more advanced features in the
microtype package.

Graphics files should be a reasonable size, and included from an appropriate
format. Use vector formats (.eps/.pdf) for plots, lossless bitmap formats
(.png) for raster graphics with sharp lines, and jpeg for photo-like images.

The style file uses the hyperref package to make clickable links in documents.
If this causes problems for you, add nohyperref as one of the options to the
icml2022 usepackage statement.

###  1.2 Submitting Final Camera-Ready Copy

The final versions of papers accepted for publication should follow the same
format and naming convention as initial submissions, except that author
information (names and affiliations) should be given. See Section 2.3.2 for
formatting instructions.

The footnote, “Preliminary work. Under review by the International Conference
on Machine Learning (ICML). Do not distribute.” must be modified to
“Proceedings of the $\mathit{39}^{th}$ International Conference on Machine
Learning, Baltimore, Maryland, USA, PMLR 162, 2022. Copyright 2022 by the
author(s).”

For those using the LaTeX style file, this change (and others) is handled
automatically by simply changing $\mathtt{\backslash
usepackage\\{icml2022\\}}$ to

$\mathtt{\backslash usepackage[accepted]\\{icml2022\\}}$

Authors using Word must edit the footnote on the first page of the document
themselves.

Camera-ready copies should have the title of the paper as running head on each
page except the first one. The running title consists of a single line
centered above a horizontal rule which is $1$ point thick. The running head
should be centered, bold and in $9$ point type. The rule should be $10$ points
above the main text. For those using the LaTeX style file, the original title
is automatically set as running head using the fancyhdr package which is
included in the ICML 2022 style file package. In case that the original title
exceeds the size restrictions, a shorter form can be supplied by using

`\icmltitlerunning{...}`

just before $\mathtt{\backslash begin\\{document\\}}$. Authors using Word must
edit the header of the document themselves.

##  2 Format of the Paper

All submissions must follow the specified format.

###  2.1 Dimensions

The text of the paper should be formatted in two columns, with an overall
width of 6.75 inches, height of 9.0 inches, and 0.25 inches between the
columns. The left margin should be 0.75 inches and the top margin 1.0 inch
(2.54 cm). The right and bottom margins will depend on whether you print on US
letter or A4 paper, but all final versions must be produced for US letter
size. Do not write anything on the margins.

The paper body should be set in 10 point type with a vertical spacing of 11
points. Please use Times typeface throughout the text.

###  2.2 Title

The paper title should be set in 14 point bold type and centered between two
horizontal rules that are 1 point thick, with 1.0 inch between the top rule
and the top edge of the page. Capitalize the first letter of content words and
put the rest of the title in lower case.

###  2.3 Author Information for Submission

ICML uses double-blind review, so author information must not appear. If you
are using LaTeX and the icml2022.sty file, use `\icmlauthor{...}` to specify
authors and `\icmlaffiliation{...}` to specify affiliations. (Read the TeX
code used to produce this document for an example usage.) The author
information will not be printed unless accepted is passed as an argument to
the style file. Submissions that include the author information will not be
reviewed.

####  2.3.1 Self-Citations

If you are citing published papers for which you are an author, refer to
yourself in the third person. In particular, do not use phrases that reveal
your identity (e.g., “in previous work (langley00), we have shown …”).

Do not anonymize citations in the reference section. The only exception are
manuscripts that are not yet published (e.g., under submission). If you choose
to refer to such unpublished manuscripts (anonymous), anonymized copies have
to be submitted as Supplementary Material via CMT. However, keep in mind that
an ICML paper should be self contained and should contain sufficient detail
for the reviewers to evaluate the work. In particular, reviewers are not
required to look at the Supplementary Material when writing their review (they
are not required to look at more than the first $8$ pages of the submitted
document).

####  2.3.2 Camera-Ready Author Information

If a paper is accepted, a final camera-ready copy must be prepared. For
camera-ready papers, author information should start 0.3 inches below the
bottom rule surrounding the title. The authors’ names should appear in 10
point bold type, in a row, separated by white space, and centered. Author
names should not be broken across lines. Unbolded superscripted numbers,
starting 1, should be used to refer to affiliations.

Affiliations should be numbered in the order of appearance. A single footnote
block of text should be used to list all the affiliations. (Academic
affiliations should list Department, University, City, State/Region, Country.
Similarly for industrial affiliations.)

Each distinct affiliations should be listed once. If an author has multiple
affiliations, multiple superscripts should be placed after the name, separated
by thin spaces. If the authors would like to highlight equal contribution by
multiple first authors, those authors should have an asterisk placed after
their name in superscript, and the term “*Equal contribution” should be placed
in the footnote block ahead of the list of affiliations. A list of
corresponding authors and their emails (in the format Full Name
<EMAIL_ADDRESS>can follow the list of affiliations. Ideally only one or
two names should be listed.

A sample file with author names is included in the ICML2022 style file
package. Turn on the [accepted] option to the stylefile to see the names
rendered. All of the guidelines above are implemented by the LaTeX style file.

###  2.4 Abstract

The paper abstract should begin in the left column, 0.4 inches below the final
address. The heading ‘Abstract’ should be centered, bold, and in 11 point
type. The abstract body should use 10 point type, with a vertical spacing of
11 points, and should be indented 0.25 inches more than normal on left-hand
and right-hand margins. Insert 0.4 inches of blank space after the body. Keep
your abstract brief and self-contained, limiting it to one paragraph and
roughly 4–6 sentences. Gross violations will require correction at the camera-
ready phase.

###  2.5 Partitioning the Text

You should organize your paper into sections and paragraphs to help readers
place a structure on the material and understand its contributions.

####  2.5.1 Sections and Subsections

Section headings should be numbered, flush left, and set in 11 pt bold type
with the content words capitalized. Leave 0.25 inches of space before the
heading and 0.15 inches after the heading.

Similarly, subsection headings should be numbered, flush left, and set in 10
pt bold type with the content words capitalized. Leave 0.2 inches of space
before the heading and 0.13 inches afterward.

Finally, subsubsection headings should be numbered, flush left, and set in 10
pt small caps with the content words capitalized. Leave 0.18 inches of space
before the heading and 0.1 inches after the heading.

Please use no more than three levels of headings.

####  2.5.2 Paragraphs and Footnotes

Within each section or subsection, you should further partition the paper into
paragraphs. Do not indent the first line of a given paragraph, but insert a
blank line between succeeding ones.

You can use footnotes111Footnotes should be complete sentences. to provide
readers with additional information about a topic without interrupting the
flow of the paper. Indicate footnotes with a number in the text where the
point is most relevant. Place the footnote in 9 point type at the bottom of
the column in which it appears. Precede the first footnote in a column with a
horizontal rule of 0.8 inches.222Multiple footnotes can appear in each column,
in the same order as they appear in the text, but spread them across columns
and pages if possible.

Figure 1: Historical locations and number of accepted papers for International
Machine Learning Conferences (ICML 1993 – ICML 2008) and International
Workshops on Machine Learning (ML 1988 – ML 1992). At the time this figure was
produced, the number of accepted papers for ICML 2008 was unknown and instead
estimated.

###  2.6 Figures

You may want to include figures in the paper to illustrate your approach and
results. Such artwork should be centered, legible, and separated from the
text. Lines should be dark and at least 0.5 points thick for purposes of
reproduction, and text should not appear on a gray background.

Label all distinct components of each figure. If the figure takes the form of
a graph, then give a name for each axis and include a legend that briefly
describes each curve. Do not include a title inside the figure; instead, the
caption should serve this function.

Number figures sequentially, placing the figure number and caption _after_ the
graphics, with at least 0.1 inches of space before the caption and 0.1 inches
after it, as in Figure 1. The figure caption should be set in 9 point type and
centered unless it runs two or more lines, in which case it should be flush
left. You may float figures to the top or bottom of a column, and you may set
wide figures across both columns (use the environment figure* in LaTeX).
Always place two-column figures at the top or bottom of the page.

###  2.7 Algorithms

If you are using LaTeX, please use the “algorithm” and “algorithmic”
environments to format pseudocode. These require the corresponding stylefiles,
algorithm.sty and algorithmic.sty, which are supplied with this package.
Algorithm 1 shows an example.

Algorithm 1 Bubble Sort

Input: data $x_{i}$, size $m$

repeat

Initialize $noChange=true$.

for $i=1$ to $m-1$ do

if $x_{i}>x_{i+1}$ then

Swap $x_{i}$ and $x_{i+1}$

$noChange=false$

end if

end for

until $noChange$ is $true$

###  2.8 Tables

You may also want to include tables that summarize material. Like figures,
these should be centered, legible, and numbered consecutively. However, place
the title _above_ the table with at least 0.1 inches of space before the title
and the same after it, as in Table 1. The table title should be set in 9 point
type and centered unless it runs two or more lines, in which case it should be
flush left.

Table 1: Classification accuracies for naive Bayes and flexible Bayes on
various data sets.
Data set | Naive | Flexible | Better?
---|---|---|---
Breast |  95.9$\pm$ 0.2 |  96.7$\pm$ 0.2 | $\surd$
Cleveland |  83.3$\pm$ 0.6 |  80.0$\pm$ 0.6 | $\times$
Glass2 |  61.9$\pm$ 1.4 |  83.8$\pm$ 0.7 | $\surd$
Credit |  74.8$\pm$ 0.5 |  78.3$\pm$ 0.6 |
Horse |  73.3$\pm$ 0.9 |  69.7$\pm$ 1.0 | $\times$
Meta |  67.1$\pm$ 0.6 |  76.5$\pm$ 0.5 | $\surd$
Pima |  75.1$\pm$ 0.6 |  73.9$\pm$ 0.5 |
Vehicle |  44.9$\pm$ 0.6 |  61.5$\pm$ 0.4 | $\surd$

Tables contain textual material, whereas figures contain graphical material.
Specify the contents of each row and column in the table’s topmost row. Again,
you may float tables to a column’s top or bottom, and set wide tables across
both columns. Place two-column tables at the top or bottom of the page.

###  2.9 Theorems and such

The preferred way is to number definitions, propositions, lemmas, etc.
consecutively, within sections, as shown below.

######  Definition 2.1.

A function $f:X\to Y$ is injective if for any $x,y\in X$ different, $f(x)\neq
f(y)$.

Using Definition 2.1 we immediate get the following result:

######  Proposition 2.2.

If $f$ is injective mapping a set $X$ to another set $Y$, the cardinality of
$Y$ is at least as large as that of $X$

###### Proof.

Left as an exercise to the reader. ∎

Lemma 2.3 stated next will prove to be useful.

######  Lemma 2.3.

For any $f:X\to Y$ and $g:Y\to Z$ injective functions, $f\circ g$ is
injective.

######  Theorem 2.4.

If $f:X\to Y$ is bijective, the cardinality of $X$ and $Y$ are the same.

An easy corollary of Theorem 2.4 is the following:

######  Corollary 2.5.

If $f:X\to Y$ is bijective, the cardinality of $X$ is at least as large as
that of $Y$.

######  Assumption 2.6.

The set $X$ is finite.

######  Remark 2.7.

According to some, it is only the finite case (cf. 2.6) that is interesting.

###  2.10 Citations and References

Please use APA reference format regardless of your formatter or word
processor. If you rely on the LaTeX bibliographic facility, use natbib.sty and
icml2022.bst included in the style-file package to obtain this format.

Citations within the text should include the authors’ last names and year. If
the authors’ names are included in the sentence, place only the year in
parentheses, for example when referencing Arthur Samuel’s pioneering work
(Samuel59). Otherwise place the entire reference in parentheses with the
authors and year separated by a comma (Samuel59). List multiple references
separated by semicolons (kearns89; Samuel59; mitchell80). Use the ‘et al.’
construct only for citations with three or more authors or after listing all
authors to a publication in an earlier reference (MachineLearningI).

Authors should cite their own work in the third person in the initial version
of their paper submitted for blind review. Please refer to Section 2.3 for
detailed instructions on how to cite your own papers.

Use an unnumbered first-level section heading for the references, and use a
hanging indent style, with the first line of the reference flush against the
left margin and subsequent lines indented by 10 points. The references at the
end of this document give examples for journal articles (Samuel59), conference
publications (langley00), book chapters (Newell81), books (DudaHart2nd),
edited volumes (MachineLearningI), technical reports (mitchell80), and
dissertations (kearns89).

Alphabetize references by the surnames of the first authors, with single
author entries preceding multiple author entries. Order references for the
same authors by year of publication, with the earliest first. Make sure that
each reference includes all relevant information (e.g., page numbers).

Please put some effort into making references complete, presentable, and
consistent, e.g. use the actual current name of authors. If using bibtex,
please protect capital letters of names and abbreviations in titles, for
example, use {B}ayesian or {L}ipschitz in your .bib file.

## Accessibility

Authors are kindly asked to make their submissions as accessible as possible
for everyone including people with disabilities and sensory or neurological
differences. Tips of how to achieve this and what to pay attention to will be
provided on the conference website http://icml.cc/.

## Software and Data

If a paper is accepted, we strongly encourage the publication of software and
data with the camera-ready version of the paper whenever appropriate. This can
be done by including a URL in the camera-ready copy. However, do not include
URLs that reveal your institution or identity in your submission for review.
Instead, provide an anonymous URL or upload the material as “Supplementary
Material” into the CMT reviewing system. Note that reviewers are not required
to look at this material when writing their review.

## Acknowledgements

Do not include acknowledgements in the initial version of the paper submitted
for blind review.

If a paper is accepted, the final camera-ready version can (and probably
should) include acknowledgements. In this case, please place such
acknowledgements in an unnumbered section at the end of the paper. Typically,
this will include thanks to reviewers who gave useful comments, to colleagues
who contributed to the ideas, and to funding agencies and corporate sponsors
that provided financial support.

##  Appendix A You _can_ have an appendix here.

You can have as much text here as you want. The main body must be at most $8$
pages long. For the final version, one more page can be added. If you want,
you can use an appendix like this one, even using the one-column format.


# Dynamic Deep Convolutional Candlestick Learner

Jun-Hao Chen
Department of Computer Science and
Information Engineering
National Taiwan University
&Yun-Cheng Tsai
Department of Data Science
Soochow University
<EMAIL_ADDRESS>

###### Abstract

Candlestick pattern is one of the most fundamental and valuable graphical
tools in financial trading that supports traders observing the current market
conditions to make the proper decision. This task has a long history and, most
of the time, human experts. Recently, efforts have been made to automatically
classify these patterns with the deep learning models. The GAF-CNN model is a
well-suited way to imitate how human traders capture the candlestick pattern
by integrating spatial features visually. However, with the great potential of
the GAF encoding, this classification task can be extended to a more
complicated object detection level. This work presents an innovative
integration of modern object detection techniques and GAF time-series encoding
on candlestick pattern tasks. We make crucial modifications to the
representative yet straightforward YOLO version 1 model based on our time-
series encoding method and the property of such data type. Powered by the deep
neural networks and the unique architectural design, the proposed model
performs pretty well in candlestick classification and location recognition.
The results show tremendous potential in applying modern object detection
techniques on time-series tasks in a real-time manner.

_K_ eywords Financial Vision $\cdot$ Object Detection $\cdot$ You Only Look
Once (YOLO) $\cdot$ Time-Series Encoding $\cdot$ Gramian Angular Field (GAF)
$\cdot$ Convolutional Neural Network (CNN) $\cdot$ Candlestick $\cdot$ Trading

##  1 Introduction

Object detection plays one of the central roles in computer vision. It has a
wide variety of applications [1, 2, 3]. Object detection tasks need to
classify each object and precisely recognize their locations within the image.
The successful implementation of such functions may benefit different fields
in computer vision. There are several candidates, and artificial neural
networks (ANN) is one of them [1]. Artificial neural networks are universal
function approximators. A single hidden layer neural network should be able to
approximate any computable function [4]. However, the width of this hidden
layer may exceed the limited computational resources and make it impractical.
The development of convolutional neural networks (CNN) [5, 6, 7, 8, 9, 10]
revolutionized computer vision and the whole machine learning community. With
CNN, it is possible to design deep neural networks capable of capturing the
spatial relationships in the data and are with human-level or superhuman
performances in classification accuracies. Moreover, several variants of CNN,
such as RCNN [11] and YOLO [3], begin to show up since $2014$ and bring us new
possibilities in the application of CNN. The latest refined deep learning-
based models still keep outperforming the state-of-the-art even now, not to
mention the flourishing versatile applications in the real world. Impressed by
the performance of object detection models, in this work, we aim to harness
the strength of location recognition to time-series-related tasks in the
financial field. Time-series data is a sequence of data points indexed in the
temporal order. It is a widely used data type in domains like the financial
and medical industries. Meanwhile, time-series data are continuously
accumulating in these sectors, creating an urgent need for well-designed
machine learning models to analyze them. It is highly desirable in the
financial industry if there is a model that can predict the future based on
previous data. For example, in the stock market, past patterns are commonly
used to predict potential future outcomes. Hence, the classification task
plays a crucial role in time-series analytics. Our previous GAF-CNN model with
a simple structure could quickly solve the multi-variate time-series
classification problem in financial candlestick data [12]. Nevertheless, the
classification model is limited to predicting the classes based on fixed
windows size inputs.

The candlestick data is a specific genre of multi-variate time-series data
that simultaneously contains open, high, low, and close time series. In
financial trading, the candlestick is a frequently used tool supporting
traders to observe the financial market pattern and make decisions. To
automatically capture the candlestick pattern under different time scales, the
model must classify the classes and recognize the locations according to
traders’ practical needs. Both demands, especially the latter, can potentially
be solved by the strength of object detection models mentioned earlier. To
harness the power of the modern object detection deep learning models, we need
to process the raw time-series data to fit CNN’s input format first. A simple
way to achieve this is to use time-series encoding methods converting one-
dimensional time-series data into a two-dimensional matrix. Gramian Angular
Field (GAF) is one of the encoding techniques which represents time-series
data with polar coordinates and converts them into the symmetric matrix [13].
Due to its simple operations and rich information, we use GAF as our time-
series preprocessing method.

The next thing to consider with the input data encoded as a matrix is the
object detection model. The object detection model based on deep learning
mainly contains two genres: (1) one-stage and (2) two-stage. The former has
the advantage of speed, while the latter is better in performance. To be
convenient and straightforward, a one-stage YOLO-$v1$ model is in this study
[3]. The YOLO model is one of the most remarkable and admired object detection
models consistently evolving since $2014$. Despite the recent development in
more powerful models, we found that the simple original version $1$ is
sufficiently competent in our purpose. Moreover, based on the GAF matrix and
time-series use cases, several adjustments will be applied to make the model
more concise and practical in time-series applications.

We create a custom dataset for research regarding the data since there is no
existing candlestick dataset for object detection tasks. we use foreign
exchange EUR/USD $1$-minute data, from January $1$, $2000$, to January $1$,
$2020$, as our raw data. As mentioned before, we label these data as eight
candlestick pattern classes referring to the Major Candlestick Signals [14].
We choose the top $10$ for each class as the targets. Next, we use these top
10 targets and apply the Dynamic Time Warping (DTW) method [15] to collect the
entire dataset with different window sizes. The processed final dataset
contains eight candlestick pattern classes with window sizes from $5$ to $16$.
More detail of our rule-based and data collection method is in Section 2.

Finally, we adopt both accuracy and confusion matrix to evaluate the
performance of our candlestick detecting model. The innovative time-series
object detection approach we proposed can reach around $88.35$% classification
accuracy and $75.55$% accuracy on window sizes determination in the testing
set. These results have sufficiently verified the object detection model. It
is competent at time-series-related tasks and would pave the way to financial
pattern detection in real-time scenarios.

##  2 Background

###  2.1 Candlestick Data

Entry and exit points detection play an essential role in financial trading.
Crucial information for strategic decision-making often hides in the depths of
unpredictable market prices. To capture this information more effectively,
Munehisa Homma proposed an intuitive visualized tool, called _Candlestick
Chart_ , to better observe the rice market’s trading circumstances in the
$18$th century [16]. A Candlestick chart is a visual representation that
contains four-dimensional time-series data, the open, high, low, and close
prices (OHLC) for a specific period. The period can be arbitrarily customized,
usually depending on the transaction strategy. The OHLC data is defined as
follows:

  1. 1.

Open: the first price of a specific period.

  2. 2.

High: the highest price of a specific period.

  3. 3.

Low: the lowest price of a specific period.

  4. 4.

Close: the last price of a specific period.

Apart from the original OHLC data, _Real Body_ and _Shadows_ , operated by
OHLC, are also practical features helping traders to make trading judgments
well. The definitions are listed below. The real body will be colorized with
black or red if the open price is higher than the close price. Otherwise, it
will be colorized with white or green. The terminologies above are commonly
used in professional traders’ discussions. Figure 1 illustrates the components
of the candlestick chart in detail.

  1. 1.

Real Body: the length between the open and close price of a specific period.

  2. 2.

Upper Shadow: the length between the high price and $max(open,close)$ of a
specific period.

  3. 3.

Lower Shadow: the length between the low price and $min(open,close)$ of a
specific period.

Figure 1: Candlestick chart provides useful information for trading, such as
open, high, low, and close prices.

To make profits from financial trading, one of the critical skills is reading
out the current market status from the candlestick data. It is challenging yet
essential to extract meaningful information from candlestick data. There are
certain kinds of _candlestick pattern_ which provide valuable information
about the dynamics of the market. The candlestick pattern can help capture the
critical circumstances among the supply and demand behind the financial market
[17]. For a detailed description of the basic candlestick patterns, we refer
to the book _The Major Candlesticks Signals_ [14]. Each pattern introduced in
the book is composed of two parts:

  1. 1.

Trend pattern: A specific pattern within the prices of a period. The uptrend
and downtrend are used in the reference to depict a pattern.

  2. 2.

Three bars pattern: The significant characteristics of the pattern appear in
the last three bars/candles, which are crucial for human traders to make
decisions.

Besides, referring to [14], eight of these patterns are adopted in this study,
including the Morning Star, Evening Star, Bullish Engulfing, Bearish
Engulfing, Shooting Star, Inverted Hammer, Bullish Harami, and Bearish Harami
patterns. Each of them will be well-illustrated in the Appendix.

###  2.2 GAF Encoding

Gramian Angular Field (GAF) is a time-series encoding method proposed by Wang
and Oates [13]. The process encodes one-dimensional time-series data into a
two-dimensional symmetric matrix with the matrix elements defined on the time-
series represented in polar coordinates. The process contains three steps.
First, normalize the raw time-series data $X$ to $[0,1]$ with a simple Min-max
normalization as equation 1, where notation $\widetilde{x}_{i}$ denotes each
normalized time-series in the entire normalized dataset $\widetilde{X}$.

$\displaystyle\widetilde{x}_{i}$
$\displaystyle=\frac{x_{i}-\min(X)}{\max(X)-\min(X)}$ (1)

Second, represent $\widetilde{x}_{i}$ into a polar coordinate system and
calculate the cosine angle and radius with equation 2. After this step, all of
the normalized time-series data will be transformed into the angle
representation with the arccosine function.

$\displaystyle\phi$
$\displaystyle=\arccos(\widetilde{x}_{i}),0\leq\widetilde{x}_{i}\leq
1,\widetilde{x}_{i}\in\widetilde{X}$ (2) $\displaystyle r$
$\displaystyle=\frac{t_{i}}{N},t_{i}\in\mathbb{N}$

Finally, use the cosine function to construct the symmetric matrix with the
combination of the angles $\phi_{i}+\phi_{j}$. The elements in the GAF matrix
contain all permutations of any two tips of time-series data. In other words,
the entire GAF matrix includes rich information of temporal correlations of
each two data points. Overall, the processed final dataset contains eight
candlestick pattern classes with window sizes from 5 to 16.

GAF $\displaystyle=\cos(\phi_{i}+\phi_{j})$ (3)
$\displaystyle=\widetilde{X}^{T}\cdot\widetilde{X}-\sqrt{I-\widetilde{X}^{2}}^{T}\cdot\sqrt{I-\widetilde{X}^{2}}$

After the encoding process, it can be found out that the diagonal elements of
the GAF matrix represent the encoded information of original time-series data.
The order of the timeline starts from the top left to the bottom right in the
GAF matrix. Moreover, the non-diagonal elements represent the correlation
between the two time-series data encoded information. In addition, the GAF
matrix also contains an essential property. Mapping the normalized time-series
data to the GAF matrix is bijective only if $\phi\in[0,\pi]$. This allows the
GAF matrix to have one and only one result of time-series data with the unique
inverse transformation through its diagonal elements.

###  2.3 GAF-CNN in Candlestick Classification

The GAF-CNN is literally a two-step process to train the CNN model with time-
series data. However, before harnessing the strength of CNN models, time-
series data need to be encoded into matrix format. Our previous study shows
two valuable features to encode candlestick patterns: (1) open, high, low, and
close (OHLC) and (2) close, upper shadow, lower shadow, and real-body (CULR).
Firstly, the OHLC features are just the raw candlestick data representation,
as Section 2.1 mentions. And then, the CULR features are inspired by the
observation habit of human traders. Human traders are used to capturing the
subtle features through visual length rather than numerical value. Both parts
perform well under a simple Lenet-5 structure. The following Figure 2 is a GAF
encoding example of both features of the evening star pattern.

Figure 2: A GAF encoding example of OHLC and CULR features of the evening star
pattern. The left-hand side shows the OHLC features, and the right-hand side
is the CULR features.

However, the length of window sizes is fixed in the simple classification
task. It cannot fit more general financial applications—for example, different
time scale detection. Besides, GAF encoding has excellent potential in
integrating deep learning object detection techniques to develop a dynamic
window length detection model. Therefore, the rest of this study will present
a method to realize object detection with GAF time-series encoding.

##  3 Data Preparation

Since there is no existing object detection dataset for the candlestick data,
we need to create a custom dataset for the research. A three-step rule-based
process is employed to prepare the dataset used in this study. First, we label
the eight candlestick patterns according to the book _The Major Candlesticks
Signals_ [14] and target each class’s top $10$ patterns. According to the
abstract depiction of the eight candlestick patterns in the book, we handmade
both rules of (1) trend pattern and (2) three bars pattern rules with our
handmade conditions. The trend pattern is defined as the top $20$ percentiles
of positive and negative slopes. The value below the $20$ percentile will
regard as no trend existing. Furthermore, we choose the patterns with the top
$10$ steep slopes as the targets, including ten positive and ten negative
slopes. Additionally, the three bars pattern definitions are more complex than
the trend pattern; therefore, we present an Evening Star case as an example
below. More details of the rules will show in our released code.

In Figure 3, several essential required rules are depicted to form an Evening
Star pattern, and most of the conditions can rule literally. For instance, the
first condition requires the bar/candle’s body to be white so that the rule
will be the (close - open) value of the first bar should be greater than zero.
However, not all depiction can rule quickly; some of the descriptions, like
"enough," is very abstract, which is also very subjective to different
traders. This study uses the value greater than $75$ and lower than $25$
percentile to represent the "enough" case. Still, the actual application is
supposed to label with traders themselves.

Figure 3: The three bars pattern depiction of the Evening Star pattern. The
textual description illustrates the three bars pattern. The trend pattern is
only present by the arrow.

Second, we then apply the Dynamic Time Warping [18] method to compare the raw
foreign exchange data. Each class’s top $10$ patterns only collect the data if
it is lower than the DTW distance threshold. Here, the DTW distance means the
similarity of the compared data to the top $10$ patterns. The $20$ percentile
of the distance can collect enough data to verify our idea. After this step,
we will get plenty of similar patterns to the top $10$. However, the data
still lacks the clear and definite window sizes to be detected by the model.

The rule-based process’s last step is to define these window sizes precisely.
Reasonably, we assume that traders can recognize their patterns’ start and the
end. In the language of object detection, this means that the bounding box’s
position should be well-defined and not an ambiguous label. Hence, we adopt
the price reversion of trend representing the patterns’ start so that the
window size will be the distance between the start and the end of each
pattern. The price reversion is defined as the candlestick’s color change in
contrast to the trend after the $5th$ bar. For example, Figure 4 is the case
of the Evening Star pattern containing a positive trend. Therefore, the price
reversion here means the color changes from white to black. Based on the
difficulty of collecting real data, our rule-based process simulates traders’
labeled patterns. All rule-based data are required to satisfy the strictness
and consistency among the same class and window size.

Figure 4: The illustration of the window size. According to the real-time
usage characteristic of time-series data, we only need to define the start of
the pattern.

##  4 Methods

The typical design of the YOLO-$v1$ model is detecting objects in the visual
images, such as Microsoft COCO or CIFAR-$10$ dataset. However, it is a
significantly different data form in the case of the time-series GAF encoding
matrix. Moreover, according to the real-time usage property of time-series
data, the bounding box can also be simplified. Hence, we make the model’s
predictions more concise and compelling in the case of the GAF matrix. We make
three crucial modifications from the original model based on the properties of
the GAF matrix and the real-time usage situation as below.

  1. 1.

Symmetric matrix:

The original YOLO-$v1$ model adopts the bounding box with the shape as
$(x,y,h,w)$. The notation $(x,y)$ represents the position of the $x$ and
$y$-axis, and $(w,h)$ means the height and width predictions of the bounding
box. Based on the property of the symmetric matrix, we only need the
information of one side of the off-diagonal elements. In other words, we can
remove the redundant side of the bounding box prediction; for example, we only
need the information of $x$ and $h$. Therefore, the predictions $(x,y,h,w)$ of
the bounding box in the original YOLO-$v1$ model can simplify to $(x,h)$.

  2. 2.

Application scenario of time-series model:

The position $x$ in the YOLO-$v1$ model can be anywhere corresponding to the
image, but it is not the same in the GAF matrix case. The moving-window
framework needs to be adopted in a real-time usage situation. It needs to keep
inferring the latest results and constantly to update time-series data,
especially in financial trading. With the moving-window framework, we do not
need to localize the object’s position. The proper location should always fix
at the latest time-series data point, which is the bottom right corner of the
GAF matrix. For instance, the position $x$ should always be $16$ under a
$16\times 16$ size GAF matrix. Figure 5 illustrates the specific position. In
Figure 5, we use a bounding box with a window size equal to 16 as an instance.
The candlestick data form is on the left-hand side with the timeline from left
to right. The corresponding GAF matrix with the timeline from top-left to
bottom-right is on the right-hand side. As the timeline shows, $x$ equals $15$
is the latest position of this time-series data. With the rightest side as the
starting point, the object can be framed by the width based on the above-
mentioned symmetric property.

  3. 3.

Remaining modification:

The original grid system in the YOLO-$v1$ model divides the input image into
an $S\times S$ grid, where $S=7$ in YOLO-$v1$, but the grid system needs to be
modified to correspond to the GAF matrix. In this research, we are building a
model to output candlestick pattern labels and the size of the time window.
What we consider is only the length counting from the bottom-right corner. We
do not consider other starting points. Hence, we can simplify the grid system
from $S\times S$ to $1\times 1$. In other words, the grid system is not needed
in our time-series moving-window framework. Figure 5 also helps to figure out
this point.

Figure 5: As mentioned in Section 2.2, the bottom right of the GAF matrix
represents the latest time-series data. Instead of the original grid system,
the start point of GAF matrix only needs to define the location $(x,y)$ at the
bottom right corner. And the model will define the window size $w$.

Figure 6 shows the entire architecture of our modified YOLO-$v1$ model. Each
$conv$-layer contains a convolutional layer, a batch normalization layer, and
a leaky-relu activation. Since input with $16\times 16$ sizes, we design the
max-pooling layer four times to make the feature map one-dimensional in the
$conv6$-layer. We also reduce each layer’s parameters to avoid using too many
parameters in a more straightforward time-series task. With all the
modifications, the last layer’s predictions contain $12$ values: two pairs of
$(w,c)$ and eight classes.

Figure 6: The structure of our YOLO-$v1$ model. The modified YOLO-$v1$ model
includes six $conv$-layers, followed by a dense-layer. Each $conv$-layer
consists of a convolutional layer, a batch normalization layer, and a leaky-
relu activation.

##  5 Results

We train the modified time-series object detection model with the dataset
labeled with eight classes and window size range from $5$ to $16$. The entire
dataset is pre-split into training, validation, and testing sets with the
$80/20/20$% ratio. Besides, the process of splitting data is arranged to avoid
using future information in time-series data. Detailed information about the
experiment parameters setting in this work is shown in Table 1.

Parameters | Values
---|---
optimizer | adadelta
learning rate | 0.001
rho | 0.95
decay | 0.0005
epochs | 4000
Table 1: The parameters used in our training process.

The testing results in Table 2 show that the model can achieve at least
$80.1$% accuracy in each class. The overall classification accuracy reaches
$88.35$% on average. The average performance is very close to $90.7$%, which
is the original GAF-CNN model’s performance solving only the simple
classification problem [19]. Moreover, Figure 7 shows the confusion matrix of
testing on window size. The prediction results of the model are mostly on the
diagonal, implying that the model performance is pretty well. Concretely, the
window size prediction can achieve almost $75.55$% correct on average in the
testing set.

Label | Number of Samples | Accuracy (%)
---|---|---
1 | 789 | 98.9
2 | 803 | 98.1
3 | 855 | 90.3
4 | 918 | 85.3
5 | 641 | 82.4
6 | 582 | 80.1
7 | 892 | 82.8
8 | 838 | 88.9
Avg |  | 88.35
Table 2: The results of testing accuracies. The accuracy of each label means
the one versus all accuracy. The accuracy on average is calculated by the
arithmetic mean. Figure 7: The confusion matrix of testing on window sizes.
The row and column represent the true label and predicted label, respectively.
The background color represents the accuracy of each window label. The results
show stable performance of the model on determining the window size.

Under the limited dataset, both the classification and location recognition
results are optimistic. Based on the invertible transformation property of the
GAF matrix, we can convert the detection results from GAF matrices back to the
more intuitive candlestick chart. Figure 8 presents a concrete detection
example of the Evening Star pattern. The right-hand side is the bounding box
of the stacked OHLC inputs, and the left-hand side is the candlesticks case.
The predicted window size is labeled in the red box. Figure 9 shows more
examples of our detection results.

Figure 8: An example of detection results for both candlestick chart and GAF
matrices. Figure 9: Further examples of our detection results. The patterns’
names and scores are labelled in each subfigure.

Our results show that the object detection model with simple architecture
could solve both the classification and location recognition well in time-
series data. It also helps detect patterns with dynamic length and different
periods according to traders’ practical needs. Furthermore, with the moving-
window framework, real-time time-series object detection can easily be
realized, benefiting professional traders in financial industries.

##  6 Discussions

###  6.1 Advanced applications with modern object detection models

As demonstrated above, the GAF encoding method integrates the time-series data
with object detection models. It has superior potential for solving dynamic
length detection of time-series data. In this work, we present an instance of
the YOLO-$v1$ model. However, there are still numerous variants of modern
object detection models that can be further applied. For example, the two-
stage models like RCNN [11] may perform better inaccuracy, not to mention the
state-of-the-art approaches. Besides, it is interesting to implement other
computer vision techniques like semantic segmentation [20] and instance
segmentation [21] to the time-series domain. We leave these potential
directions for future research.

###  6.2 Robustness of dynamic models

For a machine learning model to be confidently used in critical application
scenarios, such models must be robust to adversarial noises. Deep learning
models are known to be susceptible to such attacks [22]. In the field of
objection detection, there are known attacks, for example [23]. In the
previous work [24], the authors studied the potential attacks on the deep
candlestick learner, which operates on static inputs. We build the model on
top of objection detection models in this work, which can process real-time
data. It is interesting to study the effects of such attacks on our new model
and the potential defense mechanisms.

##  7 Conclusions

This work presents an innovative computer vision application combining object
detection and time-series data processing in the financial trading field. We
consider the actual time-series application scenario and adjust to simplify
the original YOLO-$v1$ model. With these crucial modifications, the modified
time-series object detection model becomes more efficient and can apply in the
real-time detecting framework. Even with limited labeled data and computing
power, we can still achieve $88.35$% classification accuracy and almost
$75.55$% window size prediction accuracy on average in the testing set. These
results are very positive in both classification and location recognition.
From the perspective of trading, patterns are generally too subjective to be
captured by a general algorithm. However, there is tremendous potential in
applying the object detection model to solve these emotional time-series
problems. The help of expert labeling and the universal approximation
capability of deep learning. We expect the proposed framework would benefit
the financial sectors and other fields dealing with sophisticated time-series
data.

## References

  * [1] Li Liu, Wanli Ouyang, Xiaogang Wang, Paul Fieguth, Jie Chen, Xinwang Liu, and Matti Pietikäinen.  Deep learning for generic object detection: A survey.  International journal of computer vision, 128(2):261–318, 2020\.
  * [2] Tsung-Yi Lin, Piotr Dollár, Ross Girshick, Kaiming He, Bharath Hariharan, and Serge Belongie.  Feature pyramid networks for object detection.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 2117–2125, 2017.
  * [3] Joseph Redmon, Santosh Divvala, Ross Girshick, and Ali Farhadi.  You only look once: Unified, real-time object detection.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 779–788, 2016.
  * [4] Kurt Hornik, Maxwell Stinchcombe, Halbert White, et al.  Multilayer feedforward networks are universal approximators.  Neural networks, 2(5):359–366, 1989.
  * [5] Alex Krizhevsky, Ilya Sutskever, and Geoffrey E Hinton.  Imagenet classification with deep convolutional neural networks.  In Advances in neural information processing systems, pages 1097–1105, 2012.
  * [6] Christian Szegedy, Wei Liu, Yangqing Jia, Pierre Sermanet, Scott Reed, Dragomir Anguelov, Dumitru Erhan, Vincent Vanhoucke, and Andrew Rabinovich.  Going deeper with convolutions.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 1–9, 2015.
  * [7] Karen Simonyan and Andrew Zisserman.  Very deep convolutional networks for large-scale image recognition.  arXiv preprint arXiv:1409.1556, 2014.
  * [8] Yann LeCun, Yoshua Bengio, and Geoffrey Hinton.  Deep learning.  nature, 521(7553):436–444, 2015.
  * [9] Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner.  Gradient-based learning applied to document recognition.  Proceedings of the IEEE, 86(11):2278–2324, 1998.
  * [10] Ian Goodfellow, Yoshua Bengio, Aaron Courville, and Yoshua Bengio.  Deep learning, volume 1.  MIT press Cambridge, 2016.
  * [11] Ross Girshick, Jeff Donahue, Trevor Darrell, and Jitendra Malik.  Rich feature hierarchies for accurate object detection and semantic segmentation.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 580–587, 2014.
  * [12] Yun-Cheng Tsai, Jun-Hao Chen, and Chun-Chieh Wang.  Encoding candlesticks as images for patterns classification using convolutional neural networks.  arXiv preprint arXiv:1901.05237, 2019.
  * [13] Zhiguang Wang and Tim Oates.  Imaging time-series to improve classification and imputation.  arXiv preprint arXiv:1506.00327, 2015.
  * [14] Stephen W. Bigalow.  The Major Candlesticks Signals.  The Candlestick Forum LLC, 2014.
  * [15] Donald J Berndt and James Clifford.  Using dynamic time warping to find patterns in time series.  In KDD workshop, number 16 in 10, pages 359–370. Seattle, WA, USA:, 1994.
  * [16] Steve Nison.  Japanese candlestick charting techniques: a contemporary guide to the ancient investment techniques of the Far East.  Penguin, 2001.
  * [17] Chiung-Hon Leon Lee, Alan Liu, and Wen-Sung Chen.  Pattern discovery of fuzzy time series for financial prediction.  IEEE Transactions on Knowledge and data Engineering, 18(5):613–625, 2006.
  * [18] Meinard Müller.  Dynamic time warping.  Information retrieval for music and motion, pages 69–84, 2007.
  * [19] Jun-Hao Chen and Yun-Cheng Tsai.  Encoding candlesticks as images for pattern classification using convolutional neural networks.  Financial Innovation, 6:1–19, 2020.
  * [20] Liang-Chieh Chen, George Papandreou, Iasonas Kokkinos, Kevin Murphy, and Alan L Yuille.  Deeplab: Semantic image segmentation with deep convolutional nets, atrous convolution, and fully connected crfs.  IEEE transactions on pattern analysis and machine intelligence, 40(4):834–848, 2017.
  * [21] Kaiming He, Georgia Gkioxari, Piotr Dollár, and Ross Girshick.  Mask r-cnn.  In Proceedings of the IEEE international conference on computer vision, pages 2961–2969, 2017.
  * [22] Alexey Kurakin, Ian Goodfellow, and Samy Bengio.  Adversarial machine learning at scale.  arXiv preprint arXiv:1611.01236, 2016.
  * [23] Kevin Eykholt, Ivan Evtimov, Earlence Fernandes, Bo Li, Amir Rahmati, Chaowei Xiao, Atul Prakash, Tadayoshi Kohno, and Dawn Song.  Robust physical-world attacks on deep learning visual classification.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 1625–1634, 2018.
  * [24] Jun-Hao Chen, Samuel Yen-Chi Chen, Yun-Cheng Tsai, and Chih-Shiang Shur.  Adversarial robustness of deep convolutional candlestick learner.  arXiv preprint arXiv:2006.03686, 2020.