11institutetext:  Astrophysics Research Centre, University of KwaZulu-Natal,
Durban, 3696, South Africa 22institutetext: School of Mathematics, Statistics
$\&$ Computer Science, University of KwaZulu-Natal, Westville Campus, Durban
4041, South Africa 33institutetext: Wits Centre for Astrophysics, School of
Physics, University of the Witwatersrand, Private Bag 3, 2050, Johannesburg,
South Africa 44institutetext: INAF-Istituto di Radioastronomia, via Gobetti
101, 40129 Bologna, Italy 55institutetext: Centre for Radio Astronomy
Techniques and Technologies, Department of Physics and Electronics, Rhodes
University, P.O. Box 94, Makhanda 6140, South Africa 66institutetext: South
African Radio Astronomy Observatory, 2 Fir Street, Black River Park,
Observatory, Cape Town 7925, South Africa 77institutetext: European Southern
Observatory, Karl-Schwarzschild-Str. 2, 85748 Garching, Germany
88institutetext: Department of Astronomy, Haverford College, Haverford, PA
19041, USA 99institutetext: Instituto de Física, Pontificia Universidad
Católica de Valparaíso, Casilla 4059, Valparaíso, Chile 1010institutetext:
Instituto de Astrofísica and Centro de AstroIngeniería, Facultad de Física,
Pontificia Universidad Católica de Chile, Av. Vicuña Mackenna 4860, 7820436
Macul, Santiago, Chile 1111institutetext: NASA Goddard Space Flight Center,
8800 Greenbelt Rd, Greenbelt, MD 20771, USA

# The MeerKAT Massive Distant Clusters Survey: A Radio Halo in a Massive
Galaxy Cluster at z = 1.23

S. P. Sikhosana<EMAIL_ADDRESS>M. Hilton 3322 G. Bernardi
445566 K. Kesebonye 33 D. Y. Klutse 1122 K. Knowles 115566 K. Moodley 1122 T.
Mroczkowski 77 B. Partridge 88 C. Sifón 99 C. Vargas 1010 E. Wollack 1111

(Accepted XXX. Received YYY; in original form ZZZ)

In the current paradigm, high redshift radio halos are expected to be scarce
due to inverse Compton energy losses and redshift dimming, which cause them to
be intrinsically faint. This low occurrence fraction is predicted by cosmic
ray electron turbulent re-acceleration models. To date, only a handful of
radio halos have been detected at redshift z ¿ 0.8. We report the MeerKAT
detection of a radio halo hosted by a galaxy cluster ACT-CL J0329.2-2330 at
$z=1.23$, making it the highest redshift halo detected thus far. Using
$L$-band and $UHF$-band observations, we derive a radio halo spectral index of
$\alpha^{1.3GHz}_{0.8GHz}$ = 1.3 $\pm$ 0.4 and a radio power of
P${}_{\mathrm{1.4\leavevmode\nobreak\ GHz}}=(4.4\pm 1.5)\times 10^{24}$ W
Hz-1. This result further confirms that there is rapid magnetic field
amplification in galaxy clusters at high redshift.

###### Key Words.:

galaxies: clusters: intracluster medium –- radio continuum: galaxies –-
X-rays: galaxies: clusters

##  1 Introduction

Galaxy clusters, which are at the top of the hierarchical structure formation,
grow by accreting smaller systems during cluster mergers (Press & Schechter,
1974). These merger events result in the release of enormous amounts of energy
($\sim$1064 ergs), which is partially dispersed through shocks and turbulence
in the intracluster medium (ICM). The turbulent activity induces magnetic
field amplification and the (re)acceleration of cosmic ray electrons (CRes) to
relativistic energies ($\gamma\sim 10^{3}$ – $10^{4}$). Such environments are
conducive to the formation of cluster-wide synchrotron radiation that is not
attributed to individual galaxies (like in the case of emission from active
galactic nuclei). This Mpc-scale non-thermal diffuse emission is generally
categorised as radio halos and relics. The origin of halos and relics has long
been debated by proponents of Leptonic and Hadronic models, with observations
favouring Leptonic models (see Brunetti & Jones, 2014; van Weeren et al.,
2019, for reviews).

Radio relics are highly polarised (up to $\sim$70$\%$) (Stuardi et al., 2022;
Di Gennaro et al., 2023), elongated radio sources located at the periphery of
galaxy clusters. Recent observations indicate that relics also exhibit
filamentary morphologies (Rajpurohit et al., 2020; de Gasperin et al., 2022).
The origin of relics is linked to the shock waves that propagate in the ICM
during cluster mergers. The electrons scattering upstream and downstream of
the shock region are accelerated to relativistic energies via diffusive shock
acceleration (DSA; Fermi, 1949; Drury, 1983; Hoeft & Brüggen, 2007). X-ray
studies reveal that the Mach numbers of cluster shocks are too low
($\mathcal{M}$ $<$ 3) to reproduce the observed relic structures and
luminosities purely from accelerating electrons in the thermal pool
(Markevitch et al., 2005; Botteon et al., 2016; Di Gennaro et al., 2019).
Thus, the current scenario is that relics are formed by the DSA of a pre-
accelerated CRe population (Pinzke et al., 2013; Bonafede et al., 2014;
Botteon et al., 2020a). Indeed, a few recent studies have connected relics to
head-tail radio galaxies, which are candidate sources of the seed electrons
(e.g. van Weeren et al., 2017; Rajpurohit et al., 2021).

Radio halos are located in the central region of galaxy clusters and have a
regular morphology, which tends to trace the X-ray emitting ICM (Govoni et
al., 2001). Unlike relics, halos have low polarisation percentages. From
current observations, the origin of halos strongly aligns with the Leptonic
model, where pre-existing seed relativistic electrons are re-accelerated by
merger-induced turbulence via stochastic Fermi II-type processes (Brunetti et
al., 2001; Cassano & Brunetti, 2005; Feretti et al., 2012; Donnert & Brunetti,
2014). Observations have linked radio halos to host clusters with high masses
(M500 ¿ 4$\times 10^{14}$M⊙) and X-ray and/or optical merger signatures
(Markevitch et al., 2005; Cassano et al., 2013; Shimwell et al., 2014; Lindner
et al., 2014; Kale et al., 2015; Knowles et al., 2019; Raja et al., 2021). The
presence of radio halos with ultra-steep radio spectra111 Sν $\propto$
$\nu^{-\alpha}$, where Sν is the source flux density at frequency $\nu$ and
$\alpha$ is the spectral index used in this work. ($\alpha$ ¿ 1.5) further
supports the Leptonic model. It has been theoretically shown that the re-
acceleration of secondary electrons arising from proton-proton collisions may
also contribute to synchrotron radiation (Brunetti & Lazarian, 2011; Pinzke et
al., 2017), even with the $\gamma$-ray limits of the Fermi-LAT observations
(Adam et al., 2021). However, this will not be a significant contribution,
unlike the case for Hadronic models (Dennison, 1980; Dolag & Enßlin, 2000).
The sensitivity of low-frequency instruments, such as the Low-Frequency Array
(LOFAR; van Haarlem et al., 2013), has resulted in the discovery of new
classes of diffuse emission. The recent discovery of mega halos (Cuciti et
al., 2022) has further challenged our understanding of CRe transportation in
galaxy clusters.

At high redshifts, the surface brightness of radio halos and relics is
affected by the inverse Compton (IC) energy losses of the CRe (CRe IC
scattering lifetimes scale as (1 + z)4) and redshift dimming effects (Brunetti
& Jones, 2014; Di Gennaro et al., 2021). It is therefore expected that the
radio halo occurrence rate decreases at higher redshifts. Lindner et al.
(2014) were the first to detect a radio halo at a relatively high redshift of
$z=0.870$. The host cluster is the massive ACT-CL J0102-4915 (‘El Gordo’;
Menanteau et al. 2012, $M_{\rm 500c}=(1.17\pm 0.17)\times 10^{15}$ M⊙). For
some time, this had been the only radio halo detected at $z\gtrsim 0.8$. Di
Mascolo et al. (2021) reported a tentative detection of diffuse emission in
SPT J2016–5844, which is at redshift $z$= 1.13. However, they could not claim
a definite detection due to the limited resolution of the Australian Square
Kilometre Array Pathfinder and the lack of sensitivity of the Australia
Telescope Compact Array.

Due to the scarcity of high redshift halos and relics, statistical studies
have been constrained to nearby galaxy cluster samples ($z\sim 0.1-0.4$;
Cassano et al., 2013; Cuciti et al., 2021). Recently, Di Gennaro et al. (2021)
observed a Planck$-$selected sample of high redshift (0.6 ¡ $z$ ¡ 0.9)
clusters (Planck Collaboration et al., 2016). In their study, they reported a
halo occurrence fraction of 50$\%$ and that the halos exhibited radio
luminosities similar to those hosted by nearby clusters ($z\sim 0.2$). They
hence concluded that observing such luminosities at higher redshifts would
mean that magnetic fields are amplified faster during the first phases of
cluster formation. Therefore, this implies that the magnetic field strengths
in high redshift clusters are similar to those in galaxy clusters at lower
redshift.

In this letter, we report the discovery of the highest-redshift radio halo now
known, detected in deep MeerKAT $L$ and $UHF$-band observations, hosted by
ACT-CL J0329.2-2330 (herein ACT-CL J0329). ACT-CL J0329 is reported in the SZ
cluster catalogue compiled as part of the Atacama Cosmology Telescope (ACT)
Data Release 5 (Hilton et al., 2021), and also appears in the SPTpol Extended
Cluster Survey catalogue (Bleem et al., 2020). It is at a redshift $z$ = 1.23
and has a mass of M${}_{\rm 500c}=9.7^{+1.7}_{-1.6}\times 10^{14}$ M⊙ (Hilton
et al., 2021). The ACT SZ cluster properties are summarised in Table
LABEL:tab:j0329. This radio halo detection gives us insight into non-thermal
astrophysical phenomena that occurred when the universe was less than half its
current age ($\sim$8.5 billion years ago). The letter is organised as follows.
We outline the MeerKAT observations and data reduction in Section 2. In
Section 3, we present the MeerKAT results, and in Section 4, we present our
conclusions. We adopt a $\Lambda$CDM flat cosmology with $H_{0}=70$ km s-1
Mpc-1, $\Omega_{m}=0.3$, and $\Omega_{\Lambda}=0.7$. At the redshift of ACT-CL
J0329 ($z$ = 1.23), the luminosity distance is 8535.6 Mpc, and 1$\arcsec$
corresponds to 8.32 kpc.

Table 1: Properties of ACT-CL J0329.2-2330. R.A.J2000 (hh:mm:ss.s) | 03:29:16.54
---|---
Dec.J2000 (dd:mm:ss.s) |  $-$23:30:07.56
redshift | 1.23
$M_{\rm 500cCal,SZ}\;(10^{14}\;M_{\odot})$ | 9.7${}^{+1.7}_{-1.6}$
$M_{\rm 500cUncorr,SZ}\;(10^{14}\;M_{\odot})$ | 8.0${}^{+1.3}_{-1.1}$
Signal-to-Noise ratio | 22.1
Notes. SZ values are from the ACT DR5 catalogue (Hilton et al., 2021). $M_{\rm
500cCal,SZ}$ is the ACT richness-based weak-lensing calibrated cluster mass
(thought to be the best estimate of mass). $M_{\rm 500cUncorr,SZ}$ is the ACT
cluster mass estimate assuming the Arnaud et al. (2010) scaling relation,
without Eddington bias correction, which is used for comparison with masses in
the Planck PSZ2 cluster catalogue.

##  2 MeerKAT Observations and Data Reduction

MeerKAT is a 64-antenna array operated by the South African Radio Astronomy
Observatory (SARAO; Jonas, 2009; Camilo et al., 2018). MeerKAT’s longest
baseline is 8 km, resulting in 6$\arcsec$ resolution at $L$-band (900 $-$ 1670
MHz). Approximately three-quarters of the antennas are located within a 1 km
radius, with a minimum baseline of 29 m resulting in a $UHF$-band maximum
scale of 52$\arcmin$, making MeerKAT an ideal instrument for observing faint
extended emission.

ACT-CL J0329 was initially observed as part of the MeerKAT Massive Distant
Cluster Survey (MMDCS) project proposal. This survey ultimately aims to obtain
$L/UHF$-band observations of the 30 most massive clusters at $z$ ¿ 1 in the
ACT DR5 catalogue (Hilton et al., 2021) with the purpose of studying diffuse
radio emission, star formation, and AGN activity. The selected clusters also
fall within the Dark Energy Camera Legacy Survey (DECaLS; Dey et al., 2019)
optical/infrared survey footprint. The MeerKAT observations were carried out
at $L$-band with a total on-target time of 3.5 hours, using a dump rate of 8
seconds and 4096 channels. J0408$-$6545 and J0409$-$1757 were observed for the
amplitude and phase calibration. The detection of a radio halo in ACT-CL
J0329, which was part of the sample, prompted the submission of a follow-up
proposal at $UHF$-band (580 $-$ 1015 MHz) (PI: S.P. Sikhosana). For the $UHF$
observations, the target was observed for 5 hours using the 4k channel mode,
and the dump rate was 4 seconds. A summary of the observations is in Table
LABEL:tab:observations. We detail the methods used to calibrate and image the
data below.

For both $L$ and $UHF$-band observations, we used the Containerized Automated
Radio Astronomy Calibration (CARACal222https://github.com/caracal-
pipeline/caracal) pipeline (Józsa et al., 2020) to calibrate and image the raw
visibilities. CARACal uses the Stimela333https://github.com/ratt-ru/Stimela-
classic (Makhathini, 2018) framework to execute tasks from multiple astronomy
software tools. The pipeline follows the standard radio interferometric data
reduction steps, which include flagging, primary calibration, self-
calibration, and imaging. We started by flagging all the frequency bands known
to be affected by radio frequency interference (RFI) using the Common
Astronomy Software Application CASA Team et al. (CASA;
2022)444https://casa.nrao.edu/ flagdata task, these flags account for the high
flagged percentages reported in Table LABEL:tab:observations. We used
Tricolour (Hugo et al., 2022) to excise the remaining low-level RFI. We used
the casa tasks for producing delay, bandpass, and gain calibration solutions
for the primary calibration step. We used setjy to correct for the absolute
flux density scale of J0408$-$6545\. Self-calibration was carried out using
CubiCal (Kenyon et al., 2018). The imaging was done using WSClean (Offringa et
al., 2014) with multiscale and wideband deconvolution algorithms enabled and
the Briggs robust (Briggs, 1995) set to 0. The details of the full-resolution
(FR) images are presented in Table LABEL:tab:observations, and the primary
beam-corrected images are shown in Figures 8 and 9.

Table 2: MeerKAT $L$ and $UHF$-band observations and full resolution image properties. | $L$ | $UHF$ |
---|---|---|---
Central frequency (MHz) | 1283 | 815.9 |
$\Delta\nu$ before flagging (MHz) | 856 | 435 |
No. of antennas | 64 | 60 |
Observing date (Y-M-D) | 2023-01-27 | 2023-10-12 |
Amplitude calibrator | J0408$-$6545 | J0408$-$6545 |
Phase calibrator | J0409$-$1757 | J0409$-$1757 |
On-source time (hrs) | 3.5 | 5 |
$t_{\rm int}$ (sec) | 8 | 4 |
Flagged ($\%$) | 55 | 36 |
$\theta_{\rm synth}$ (″$\times$ ″) | 8.0 $\times$ 6.7 | 10.5 $\times$ 8.6 |
p.a. (°) | 162.4 | 142.7 |
$\sigma_{\rm rms}$ ($\upmu$Jy/beam) | 6.6 | 13.2 |

In Figure 1, faint diffuse emission is visible in the central region of the
cluster in FR images; however, there are bright compact sources embedded. We
then implemented further imaging steps to model the compact sources, remove
them in the visibility plane, and image at low resolution. These steps were
carried out for both $L$ and $UHF$-band datasets. We began by modelling the
compact sources in the field by producing a high-resolution (HR) image using
WSclean’s inner-tukey taper set to 10 k$\lambda$, which corresponds to 173 kpc
at the cluster redshift. The chosen taper sufficiently captures the compact
sources while resolving the diffuse emission component. We used the model
column generated by the HR images to subtract the compact sources in the
visibility plane; this was done by the casa task uvsub. We then produced a FR
compact-source-subtracted image and used it to quantify the contamination from
the removal of compact sources (see Figures 6 and 7). We used the method
described in van Weeren et al. (2021) to quantify the contamination. Finally,
we produced a low-resolution (LR) image by applying an outer taper of 12
k$\lambda$. The flux density measurements of the diffuse radio sources were
extracted from the LR images using polygon regions for each source, which are
guided by the 3$\sigma$ contour levels. The uncertainty associated with the
flux density of the diffuse radio source, measured at frequency $\nu$ is given
by

$\Delta S_{\nu}=\sqrt{(\delta S_{\nu}\times
S_{\nu})^{2}+N_{beams}\sigma_{rms}^{2}+\sigma_{sub}^{2}}\;,$ (1)

where $\delta S_{\nu}$ is the flux calibration uncertainty, which we assume to
be 5$\%$ for $L$-band (Knowles et al., 2022) and 14$\%$ for $UHF$-band
(Obtained from cross-matching with NVSS sources within the primary beam, using
$\alpha$=0.7 to scale the flux densities.), $N_{beams}$ is the number of beams
within the region that the flux density was measured, and $\sigma_{rms}$ is
the local rms noise of the image. $\sigma_{sub}$ is the uncertainty due to
compact-source-subtraction (van Weeren et al., 2021), it is given by

$\sigma_{sub}=\sum_{i}N_{beam,s_{i}}\sigma_{rms}^{2}\;,$ (2)

where the sum is taken over the number of compact sources ($i$) subtracted
within the polygon region of the diffuse radio source and $N_{beam,s_{i}}$ is
the number of beams covering the $i^{th}$ compact source.

We used UV- and resolution-matched images, which are centred at 1.28 GHz and
815.9 MHz, to create a spectral index map. To ensure that both images cover
the same UV range at the shortest baselines, we applied a minuv cut of 100
$\lambda$. We matched the resolution of the images by re-imaging the $L$-band
data with beam parameters matching the $UHF$ LR image. The spectral index map
was produced using Broadband Radio Astronomy Tools (BRATS, Harwood et al.,
2013, 2015). We set the cutoff to 3$\sigma$ and the calibration uncertainty to
5$\%$; the rest of the parameters were left to default settings.

##  3 Results and Discussion

The MeerKAT $L$ and $UHF$-band LR images reveal extended emission at the
centre of ACT-CL J0329. This emission has a largest linear size of 1.1 Mpc at
1.28 GHz and has regular morphology. The morphology traces the X-ray emission
as shown in Figure 2, where the X-ray image is obtained from archival Chandra
ACIS-I observations (observation IDs: 18282 and 18882) with an exposure time
of 47 ksec, using the 0.5-4 keV energy band. Thus, we categorise the source as
a radio halo. There is a nearby elongated source to the east that connects to
the halo emission; in the FR image, this source appears to be a head-tail (HT)
radio galaxy. After subtracting the compact emission of the HT galaxy, there
is residual emission in the LR images, which we measure to have flux densities
of 0.29 $\pm$ 0.05 mJy and 0.78 $\pm$ 0.3 mJy at $L$ and $UHF$-band,
respectively. However, the emission of the elongated source is too faint in
the LR images to produce a spectral index map; hence, we cannot unambiguously
conclude whether it is a radio relic. The presence of the HT galaxy suggests
that the seed electrons may be deposited into the ICM as proposed by the re-
acceleration model (Lusetti et al., 2024).

Figure 1: Left:$L$-band full-resolution image of the radio halo. The local rms
of the image is 1$\sigma$ = 6.6 $\upmu$Jy/beam; the beam size is indicated by
the red ellipse. The overlaid contours are from the low-resolution image,
which has a local rms of 1$\sigma$ = 7.1 $\upmu$Jy/beam. The contour levels
are $\sigma\times$[3,6,10] Right: $UHF$-band full-resolution image of the
radio halo. The local rms of the image is 1$\sigma$ = 13.2 $\upmu$Jy/beam; the
beam size is indicated by the red ellipse. The overlaid contours are from the
low-resolution image, which has a local rms of 1$\sigma$ = 17.6
$\upmu$Jy/beam. The contour levels are $\sigma\times$[3,6,10]

Using the 3$\sigma$ contours of the LR images (see Figure 1) as boundaries, we
find that the radio halo has a linear scale of 1.1$\times$0.9 Mpc at 1.28 GHz
and 900$\times$950 kpc at 815.9 MHz. The radio halo has a smooth, regular
morphology that traces the thermal bremsstrahlung emission of the ICM, as seen
in Figure 2. We extract the flux densities in the same region; the radio halo
has a flux density of 3.44 $\pm$ 0.2 mJy and 6.11 $\pm$ 0.9 mJy at $L$ and
$UHF$-band, respectively. Using the multi-frequency flux density measurements,
we derive the integrated radio halo spectral index to be
$\alpha^{1.3GHz}_{0.8GHz}$ = 1.3 $\pm$ 0.4. The best-fit value of the
integrated spectral index places the halo as a moderately steep spectrum radio
halo. However, the large uncertainty associated with the spectral index limits
our discussion. Bootstrapping the $L$-band flux density, we infer a halo radio
power of P${}_{\rm 1.4GHz}=(4.4\pm 1.5)\times 10^{24}$ W Hz-1. The
uncertainties associated with the integrated spectral index and the radio
power are obtained from with 1000 Monte Carlo realisations. We compare this
luminosity to the samples observed in Di Gennaro et al. (2021), (Knowles et
al., 2021), and the radio halo in the El Gordo cluster (Lindner et al., 2014).
Plotting the P${}_{\mathrm{1.4\leavevmode\nobreak\ GHz}}-$M500 scaling
relation, we find that the halo radio power falls within the scatter of both
the high redshift and nearby radio halos. This is indicated in Figure 3, where
the best-fit line is from scaling relation studies of nearby radio halos
($z\sim$0.2) carried out by Cassano et al. (2013). The radio luminosity of the
halo falls within the correlation, indicating that the halo is as luminous as
nearby ($z\sim$ 0.2) halos. This result further confirms that there is rapid
magnetic field amplification in galaxy clusters at high redshifts, as
discussed in Di Gennaro et al. (2021).

Figure 2: Archival Chandra image of ACT-CL J0329, observation IDs: 18282 and
18882. The cluster was observed for 47 ksec, using ACIS-I in VFAINT mode, and
the image covers the 0.5-4 keV energy band. The image was smoothed to a FWHM
of 12$\arcsec$. The overlaid contours are from the MeerKAT $UHF$-band LR
images; the contour levels are $\sigma\times$[3,4,8], where 1$\sigma$ = 17.6
$\upmu$Jy/beam. The beam size is indicated by the red ellipse.  Figure 3:
P${}_{\mathrm{1.4\leavevmode\nobreak\ GHz}}-$M500 correlation plot of radio
halos. To account for the discrepancies between the ACT and Planck mass
derivations, we use the uncorrected M500c values for the ACT clusters (Hilton
et al., 2021). The solid line shows the best-fit relation, with parameters
extracted from the low-redshift sample ($z\sim$0.2) presented in Cassano et
al. (2013). The (Di Gennaro et al., 2021) sample consists of radio halos at z
¿ 0.6, and the (Knowles et al., 2021) sample covers 0.22 ¡ z ¡ 0.65. We also
include El Gordo, which is at z = 0.87 (Lindner et al., 2014).

Spectral index studies of high redshift radio halos are vital for
understanding the CRe re-acceleration mechanisms in the early universe. Figure
4 shows that the spectral index distribution varies from 0.5 to 3.3, with the
eastern region of the halo having steeper spectra while the western region has
flatter spectra. From the map, we extract a median spectral index of $\alpha$
= 1.7, which is steeper in comparison to the integrated spectral index value.
The variations in spectral index maps can sometimes be due to measurement
errors. According to Vacca et al. (2014) and Pearce et al. (2017), if the
variations are due to measurement errors, the median error value should be
comparable to the standard deviation. The spectral index distributed for the
maps is indicated in Figure 5. The median error of the distribution is 0.05,
and the standard deviation is 0.98, thus implying that the fluctuations in the
spectral index map are intrinsic. A few studies have reported varying spectral
index distribution in radio halos (Botteon et al., 2020b; Hoang et al., 2021;
Rajpurohit et al., 2023). According to Brunetti et al. (2001) and Brunetti &
Lazarian (2007), fluctuations in the spectral index over the radio halo
regions arise due to the non-uniformity of the magnetic field distribution and
different acceleration efficiencies. The fluctuations observed in the ACT-CL
J0329 spectral index map are peculiar in that, the steeper spectral index
values are concentrated in the eastern region, which does not overlap with a
compact source. Thus, the physical interpretation could be that the turbulent
energy is not homogeneously dissipated in the halo volume.

Figure 4: Spectral index map of the radio halo and the corresponding error
map. The overlaid black and magenta contours are from the MeerKAT $UHF$-band
and $L$-band’s high-resolution images, respectively. For both images, the
contour levels are $\sigma\times$[3,9,18,27]. 1$\sigma$ is 12.4 $\upmu$Jy/beam
for the $UHF$-band image and 4.9 $\upmu$Jy/beam for the $L$-band image. Figure
5: Histogram of the spectral index distribution for the halo in ACT-CL J0329.
The black solid line indicates the median spectral index value $\alpha$=1.7.
The dashed black lines show the standard deviation around the median,
$\sigma$=0.98. The dot-dashed red line is the integrated spectral index.

##  4 Conclusion

In this letter, we have presented MeerKAT $L$ and $UHF$-band observations of
ACT-CL J0329.2-2330, a galaxy cluster at $z$=1.23. The low-resolution images
reveal a radio halo in the cluster; this is the highest redshift halo reported
to date. The halo has an integrated spectral index of $\alpha$ = 1.3 $\pm$ 0.4
and a radio power of 4.4 $\pm$ 1.5$\times$1024 W Hz-1. This indicates that the
halo is as luminous as the halos found in nearby massive galaxy clusters. The
spectral index map reveals varying spectral indices across the halo region.

## Acknowledgements

The MeerKAT telescope is operated by the South African Radio Astronomy
Observatory, which is a facility of the National Research Foundation, an
agency of the Department of Science and Innovation. We acknowledge the
assistance of the South African Radio Astronomy Observatory (SARAO) science
commissioning and operations team (led by Sharmila Goedhart). We acknowledge
the financial assistance of the SARAO towards this research. KM and MH
acknowledge support from the National Research Foundation of South Africa. CS
acknowledges support from the Agencia Nacional de Investigación y Desarrollo
(ANID) through Basal project FB210003.

## Data Availability

The data derived in this study are available upon reasonable request from the
corresponding author.

## References

  * Adam et al. (2021) Adam, R., Goksu, H., Brown, S., Rudnick, L., & Ferrari, C. 2021, A&A, 648, A60
  * Arnaud et al. (2010) Arnaud, M., Pratt, G. W., Piffaretti, R., et al. 2010, A&A, 517, A92
  * Bleem et al. (2020) Bleem, L. E., Bocquet, S., Stalder, B., et al. 2020, ApJS, 247, 25
  * Bonafede et al. (2014) Bonafede, A., Intema, H. T., Brüggen, M., et al. 2014, ApJ, 785, 1
  * Botteon et al. (2020a) Botteon, A., Brunetti, G., Ryu, D., & Roh, S. 2020a, A&A, 634, A64
  * Botteon et al. (2020b) Botteon, A., Brunetti, G., van Weeren, R. J., et al. 2020b, ApJ, 897, 93
  * Botteon et al. (2016) Botteon, A., Gastaldello, F., Brunetti, G., & Dallacasa, D. 2016, MNRAS, 460, L84
  * Briggs (1995) Briggs, D. S. 1995, PhD thesis, New Mexico Institute of Mining and Technology
  * Brunetti & Jones (2014) Brunetti, G. & Jones, T. W. 2014, International Journal of Modern Physics D, 23, 1430007
  * Brunetti & Lazarian (2007) Brunetti, G. & Lazarian, A. 2007, MNRAS, 378, 245
  * Brunetti & Lazarian (2011) Brunetti, G. & Lazarian, A. 2011, MNRAS, 410, 127
  * Brunetti et al. (2001) Brunetti, G., Setti, G., Feretti, L., & Giovannini, G. 2001, MNRAS, 320, 365
  * Camilo et al. (2018) Camilo, F., Scholz, P., Serylak, M., et al. 2018, ApJ, 856, 180
  * CASA Team et al. (2022) CASA Team, Bean, B., Bhatnagar, S., et al. 2022, PASP, 134, 114501
  * Cassano & Brunetti (2005) Cassano, R. & Brunetti, G. 2005, MNRAS, 357, 1313
  * Cassano et al. (2013) Cassano, R., Ettori, S., Brunetti, G., et al. 2013, ApJ, 777, 141
  * Cuciti et al. (2021) Cuciti, V., Cassano, R., Brunetti, G., et al. 2021, A&A, 647, A51
  * Cuciti et al. (2022) Cuciti, V., de Gasperin, F., Brüggen, M., et al. 2022, Nature, 609, 911
  * de Gasperin et al. (2022) de Gasperin, F., Rudnick, L., Finoguenov, A., et al. 2022, A&A, 659, A146
  * Dennison (1980) Dennison, B. 1980, in BAAS, Vol. 12, 471
  * Dey et al. (2019) Dey, A., Schlegel, D. J., Lang, D., et al. 2019, AJ, 157, 168
  * Di Gennaro et al. (2023) Di Gennaro, G., Brüggen, M., van Weeren, R. J., et al. 2023, A&A, 675, A51
  * Di Gennaro et al. (2019) Di Gennaro, G., van Weeren, R. J., Andrade-Santos, F., et al. 2019, ApJ, 873, 64
  * Di Gennaro et al. (2021) Di Gennaro, G., van Weeren, R. J., Brunetti, G., et al. 2021, Nature Astronomy, 5, 268
  * Di Mascolo et al. (2021) Di Mascolo, L., Mroczkowski, T., Perrott, Y., et al. 2021, A&A, 650, A153
  * Dolag & Enßlin (2000) Dolag, K. & Enßlin, T. A. 2000, A&A, 362, 151
  * Donnert & Brunetti (2014) Donnert, J. & Brunetti, G. 2014, MNRAS, 443, 3564
  * Drury (1983) Drury, L. O. 1983, Reports on Progress in Physics, 46, 973
  * Feretti et al. (2012) Feretti, L., Giovannini, G., Govoni, F., & Murgia, M. 2012, A&A Rev., 20, 54
  * Fermi (1949) Fermi, E. 1949, Physical Review, 75, 1169
  * Gaia Collaboration et al. (2023) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2023, A&A, 674, A1
  * Govoni et al. (2001) Govoni, F., Enßlin, T. A., Feretti, L., & Giovannini, G. 2001, A&A, 369, 441
  * Harwood et al. (2015) Harwood, J. J., Hardcastle, M. J., & Croston, J. H. 2015, MNRAS, 454, 3403
  * Harwood et al. (2013) Harwood, J. J., Hardcastle, M. J., Croston, J. H., & Goodger, J. L. 2013, MNRAS, 435, 3353
  * Hilton et al. (2021) Hilton, M., Sifón, C., Naess, S., et al. 2021, ApJS, 253, 3
  * Hoang et al. (2021) Hoang, D. N., Zhang, X., Stuardi, C., et al. 2021, A&A, 656, A154
  * Hoeft & Brüggen (2007) Hoeft, M. & Brüggen, M. 2007, MNRAS, 375, 77
  * Hugo et al. (2022) Hugo, B. V., Perkins, S., Merry, B., Mauch, T., & Smirnov, O. M. 2022, in Astronomical Society of the Pacific Conference Series, Vol. 532, Astronomical Society of the Pacific Conference Series, ed. J. E. Ruiz, F. Pierfedereci, & P. Teuben, 541
  * Jonas (2009) Jonas, J. L. 2009, IEEE Proceedings, 97, 1522
  * Józsa et al. (2020) Józsa, G. I. G., White, S. V., Thorat, K., et al. 2020, CARACal: Containerized Automated Radio Astronomy Calibration pipeline, Astrophysics Source Code Library, record ascl:2006.014
  * Kale et al. (2015) Kale, R., Venturi, T., Giacintucci, S., et al. 2015, A&A, 579, A92
  * Kenyon et al. (2018) Kenyon, J. S., Smirnov, O. M., Grobler, T. L., & Perkins, S. J. 2018, MNRAS, 478, 2399
  * Knowles et al. (2019) Knowles, K., Baker, A. J., Bond, J. R., et al. 2019, MNRAS, 486, 1332
  * Knowles et al. (2022) Knowles, K., Cotton, W. D., Rudnick, L., et al. 2022, A&A, 657, A56
  * Knowles et al. (2021) Knowles, K., Pillay, D. S., Amodeo, S., et al. 2021, MNRAS, 504, 1749
  * Lindner et al. (2014) Lindner, R. R., Baker, A. J., Hughes, J. P., et al. 2014, ApJ, 786, 49
  * Lusetti et al. (2024) Lusetti, G., de Gasperin, F., Cuciti, V., et al. 2024, MNRAS, 528, 141
  * Makhathini (2018) Makhathini, S. 2018, PhD thesis, Rhodes University, Drosty Rd, Grahamstown, 6139, Eastern Cape, South Africa, available via http://hdl.handle.net/10962/57348
  * Markevitch et al. (2005) Markevitch, M., Govoni, F., Brunetti, G., & Jerius, D. 2005, ApJ, 627, 733
  * Menanteau et al. (2012) Menanteau, F., Hughes, J. P., Sifón, C., et al. 2012, ApJ, 748, 7
  * Offringa et al. (2014) Offringa, A. R., McKinley, B., Hurley-Walker, N., et al. 2014, MNRAS, 444, 606
  * Pearce et al. (2017) Pearce, C. J. J., van Weeren, R. J., Andrade-Santos, F., et al. 2017, ApJ, 845, 81
  * Pinzke et al. (2013) Pinzke, A., Oh, S. P., & Pfrommer, C. 2013, MNRAS, 435, 1061
  * Pinzke et al. (2017) Pinzke, A., Oh, S. P., & Pfrommer, C. 2017, MNRAS, 465, 4800
  * Planck Collaboration et al. (2016) Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2016, A&A, 594, A27
  * Press & Schechter (1974) Press, W. H. & Schechter, P. 1974, ApJ, 187, 425
  * Raja et al. (2021) Raja, R., Rahaman, M., Datta, A., et al. 2021, MNRAS, 500, 2236
  * Rajpurohit et al. (2020) Rajpurohit, K., Hoeft, M., Vazza, F., et al. 2020, A&A, 636, A30
  * Rajpurohit et al. (2023) Rajpurohit, K., Osinga, E., Brienza, M., et al. 2023, A&A, 669, A1
  * Rajpurohit et al. (2021) Rajpurohit, K., Wittor, D., van Weeren, R. J., et al. 2021, A&A, 646, A56
  * Shimwell et al. (2014) Shimwell, T. W., Brown, S., Feain, I. J., et al. 2014, MNRAS, 440, 2901
  * Stuardi et al. (2022) Stuardi, C., Bonafede, A., Rajpurohit, K., et al. 2022, A&A, 666, A8
  * Vacca et al. (2014) Vacca, V., Feretti, L., Giovannini, G., et al. 2014, A&A, 561, A52
  * van Haarlem et al. (2013) van Haarlem, M. P., Wise, M. W., Gunst, A. W., et al. 2013, A&A, 556, A2
  * van Weeren et al. (2019) van Weeren, R. J., de Gasperin, F., Akamatsu, H., et al. 2019, Space Sci. Rev., 215, 16
  * van Weeren et al. (2017) van Weeren, R. J., Ogrean, G. A., Jones, C., et al. 2017, ApJ, 835, 197
  * van Weeren et al. (2021) van Weeren, R. J., Shimwell, T. W., Botteon, A., et al. 2021, A&A, 651, A115

##  Appendix A Ancillary images of ACT-CL J0329.2-2330

Figure 6: The L-band full-resolution compact-source-subtracted image with
contour overlays of the high-resolution image. The contour levels are
$\sigma\times$[3,12,24,30], where 1$\sigma$ is 4.9 $\upmu$Jy/beam. The beam of
the high-resolution image is 6.5$\arcsec\times$5.5$\arcsec$, 166$\degree$ pa.
Figure 7: The UHF-band full-resolution compact-source-subtracted image with
contour overlays of the high-resolution image. The contour levels are
$\sigma\times$[3,12,24,30], where 1$\sigma$ is 12.4 $\upmu$Jy/beam. The beam
of the high-resolution image is 10.7$\arcsec\times$8.1$\arcsec$, 146 $\degree$
pa. Figure 8:  The MeerKAT L-band full field of view image of ACT-CL J0329.
The rms noise of the image is 6.6 $\upmu$Jy/beam. Figure 9:  The MeerKAT UHF-
band full field of view image of ACT-CL J0329. The rms noise of the image is
13.2 $\upmu$Jy/beam. Figure 10: Optical grz image of ACT-CL J0329 from DECaLS
(Dey et al. 2019) The overlaid contours are from the MeerKAT $UHF$-band LR
images, the contour levels are $\sigma\times$[3,6,10]; where 1$\sigma$ = 17.6
$\upmu$Jy/beam. The white source at the centre is a foreground star located at
a distance of $\approx 1500$ pc as recorded in Gaia DR3 (Gaia Collaboration et
al. 2023).


# Discovering Useful Compact Sets of Sequential Rules in a Long Sequence

Erwan Bourrand1,2, Luis Galárraga3, Esther Galbrun4, Elisa Fromont1, Alexandre
Termier1  1 Univ Rennes, IRISA UMR 6074, Rennes, France
2 Advisor_SLA, Inc.
3 Inria RBA Rennes, France
4 University of Eastern Finland
{luis.galarraga, elisa.fromont<EMAIL_ADDRESS>
###### Abstract

We are interested in understanding the underlying generation process for long
sequences of symbolic events. To do so, we propose Cossu, an algorithm to mine
small and meaningful sets of sequential rules. The rules are selected using an
MDL-inspired criterion that favors compactness and relies on a novel rule-
based encoding scheme for sequences. Our evaluation shows that Cossu can
successfully retrieve relevant sets of closed sequential rules from a long
sequence. Such rules constitute an interpretable model that exhibits
competitive accuracy for the tasks of next-element prediction and
classification.

##  I Introduction

Long sequences of symbols are ubiquitous. Examples include DNA sequences,
server logs, traces of network packets, and even long texts. Discovering
regularities in such sequences allows for a better understanding of the
sequence generation process, and can be useful e.g., for diagnostic and
prediction.

Some models, such as LSTMs [1], excel at detecting those regularities and use
them for accurate prediction. On the downside, such models are hardly
understandable by human users. Models based on pattern mining, in contrast,
provide high-level human-readable descriptions of the underlying structure of
the data. However, this interpretability typically comes at the cost of a
restricted predictive power. We focus on the latter category of models, since
we are interested in prediction but also in knowledge discovery.

In the realm of pattern mining, sequential rules are a well-established model
for understanding and predicting sequential data. Sequential rules take the
form $A\Rightarrow C$, where both the antecedent $A$ and the consequent $C$
are sequences of items, e.g., events, nucleotides, or words. A _confidence_
score is often attached to the rule to measure how likely it is to observe $C$
after $A$ in the sequence. For instance, the rule
$A_{1}A_{2}A_{3}\Rightarrow\textit{PowerFail}$ with confidence $80\%$ tell us
that there is a four in five chance of undergoing a power failure after having
observed the sequence of alarms $A_{1}$, $A_{2}$ and $A_{3}$. This can be used
both for prediction, i.e., to issue an early warning, and for diagnostic,
i.e., to identify the link between the alarms and the power failure.

While sequential rule mining has been extensively studied in the literature
[2, 3, 4, 5, 6], two major issues remain. First, the vast majority of
approaches consider a database of short sequences rather than a long sequence.
While a long sequence can always be split into a database, this incurs a loss
of information at the boundaries. Second, to the best of our knowledge, no
approach for sequential rule mining addresses the so-called pattern explosion,
the fact that millions of sequential rules may be produced due to the
combinatorial search space.

We propose the first sequential rule mining approach that takes as input a
long sequence and outputs a compact set of rules. The selection of the rules
is based on the principle of Minimum Description Length (MDL), a paradigm that
has proved effective for model selection in pattern mining [7, 8]. MDL
approaches rely on an encoding scheme. Ours is presented in Section IV, along
with our algorithm to mine rules, in Section V. Experiments, presented in
Section VI, demonstrate the relevance of the rules found, as well as their
predictive power for the tasks of next-element prediction, and classification
for long sequences.

##  II Background

Our work stands at the crossroads of sequential and MDL-based pattern mining.
Below we give brief overviews of these two fields, in turn.

Sequential Pattern Mining. Mining sequential patterns on data is a well-
studied problem that has given rise to a plethora of methods accounting for
different notions of patterns. The vast majority of approaches assume a
database of sequences as input, where sequences are ordered collections of
itemsets. We refer the reader to the survey in [2] for a detailed account.

Among the most recent works, Fumarola et al. [5] mine frequent closed
sequences, i.e. sequences of maximal length such that any extension will
inevitably have lower support. Closed sequences can alleviate pattern
explosion to some extent as they are a subset of frequent sequences. Other
approaches [9, 6, 4] mine sequential rules of the form $A\Rightarrow C$ where
$A$ and $C$ may be itemsets [9, 6] or sub-sequences [4] with the constraint
that $C$ occurs after $A$. Itemsets are generally easier to mine than
sequences, however they are less expressive because they do not account for
order in the appearance of items. Still, none of these methods is directly
portable to long sequences as they all assume the data has been partitioned
into small sequences (and hence often rely on a different definition of
support). A few methods [3, 10] can natively mine rules on long sequences.
Unfortunately, such methods are not resilient to pattern explosion and are
limited to association rules between itemsets, thus they do not capture the
order between items within the antecedent and the consequent.

Mining patterns with an MDL criterion. The Minimum Description Length (MDL)
principle [11] is a criterion rooted in information theory that stipulates
that the best model for a dataset is the model that compresses it best. For a
dataset $D$ and a family of models $\mathcal{H}$, the best model for $D$,
according to the (two-parts) MDL, is the one that minimizes the description
length

$H^{*}=\arg\min_{H\in\mathcal{H}}\;\operatorname{L}(H)+\operatorname{L}(D|H)\,,$
(1)

where $\operatorname{L}$ denotes code length in bits. Equation (1) strikes a
balance between model complexity, as measured by $\operatorname{L}(H)$, and
fitness to the data, as measured by $\operatorname{L}(D|H)$. When applying MDL
to pattern mining, a model is a collection of patterns, e.g., itemsets [7],
sequences [8] or, in our case, sequential rules. One of the main ingredients
of any MDL-inspired pattern mining approach is the _encoding scheme_ , namely,
the protocol to encode the data with the patterns and to encode the patterns
themselves. Once an encoding mechanism is in place, we can generate candidate
sets of patterns, calculate the corresponding code lengths, and select the set
resulting in the shortest code length.

Krimp [7] was one of the pioneering efforts to apply MDL to pattern mining,
more specifically to itemset mining. Given a (large) collection of itemsets
mined from the data, Krimp selects a small representative subset. The selected
itemsets were also empirically shown to be effective for classification.
Recently, Fischer and Vreeken [12] proposed an approach to mine compact sets
of rules from data without the sequential dimension. The approaches proposed
in [13, 14] extract sequences with gaps (but not rules) from univariate and
multivariate long sequences. Unlike Krimp [7], they combine the generation and
selection of candidate patterns, which boosts efficiency and pattern quality.
Shokoohi-Yekta et al. [15] draw inspiration from the MDL principle to mine
sequential rules from time-series discretized into an integer domain. A mining
approach tailored for the streaming setting is proposed in [16], but its
compliance with the MDL principle is debatable.

##  III Definitions and Notation

A _sequence_ $S$ of length $n$ over an alphabet $\Sigma$ is an ordered
collection of $n$ occurrences of symbols from $\Sigma$. We denote by $S[i]$
the $i$-th element in $S$ ($1\leq i\leq n$), and let $S[i,j]$ denote the
contiguous sequence $\langle S[i]\dots S[j]\rangle$ when $1\leq i\leq j\leq n$
or the empty sequence $\emptyset$ otherwise. For simplicity, we write the
sequence $\langle\sigma_{1},\sigma_{2}\dots\sigma_{k}\rangle$ as
$\sigma_{1}\sigma_{2}\dots\sigma_{k}$. We denote the concatenation of two
sequences by $\oplus$.

We say that sequence $S^{\prime}$ is a _subsequence_ of $S$ (respectively $S$
is a _supersequence_ of $S^{\prime}$), denoted as $S^{\prime}\sqsubseteq S$,
iff there exist $i$ and $j$ in $[1,n]$ such that $S^{\prime}=S[i,j]$. Then, we
call $S[i,j]$ a _match_ of $S^{\prime}$ in $S$. More specifically, we denote
by $S^{\prime}\sqsubseteq_{i}S$ the fact that $S^{\prime}$ is a subsequence of
$S$ with a match starting at position $i$, i.e. $S^{\prime}=S[i,i+\left\lvert
S^{\prime}\right\rvert-1]$, and by $S^{\prime}\sqsubseteq^{j}S$ the fact that
$S^{\prime}$ is a subsequence of $S$ with a match ending at position $j$, i.e.
$S^{\prime}=S[j-\left\lvert S^{\prime}\right\rvert+1,j]$.

The _support_ of a subsequence $S^{\prime}$ in $S$, denoted by
$\operatorname{supp}_{S}(S^{\prime})$, is the number of distinct matches of
$S^{\prime}$ in $S$, that is,
$\operatorname{supp}_{S}(S^{\prime})=\left\lvert\\{i\in[1,n]:\;S^{\prime}\sqsubseteq_{i}S\\}\right\rvert=\left\lvert\\{j\in[1,n]:\;S^{\prime}\sqsubseteq^{j}S\\}\right\rvert$.
Finally, a sequence $S^{\prime}$ is _closed_ in $S$ if there is no
supersequence $S^{\prime\prime}\sqsupset S^{\prime}$ in $S$ such that
$\operatorname{supp}_{S}(S^{\prime\prime})=\operatorname{supp}_{S}(S^{\prime})$.

###### Example 1

Sequences.  Given the sequences $S=abceabcadeab$ and $S^{\prime}=abc$,
$S^{\prime}\sqsubseteq_{1}S$, $S^{\prime}\sqsubseteq_{5}S$,
$S^{\prime}\sqsubseteq^{3}S$ and $S^{\prime}\sqsubseteq^{7}S$ are all true,
hence $\operatorname{supp}_{S}(S^{\prime})=2$.

A _rule_ is a pair of sequences $(A,C)$ over alphabet $\Sigma$, such that
$C\neq\emptyset$. It is denoted as $A\Rightarrow C$. $A$ and $C$ are called
the _antecedent_ and the _consequent_ of the rule, respectively. A rule such
that $A=\emptyset$ and $\left\lvert C\right\rvert=1$ is called a _singleton
rule_. The set of singleton rules over alphabet $\Sigma$ is
$\Psi_{\Sigma}=\\{\emptyset\Rightarrow\langle\sigma\rangle:\;\sigma\in\Sigma\\}$.

A rule $R:A\Rightarrow C$ _triggers_ on sequence $S$ at position $i$, denoted
as $R\vdash_{i}S$, iff the antecedent has a match in $S$ ending at position
$i$, i.e., iff $A\sqsubseteq^{i}S$. The rule $R$ _applies_ on sequence $S$ at
position $i$, denoted by $R\Vdash_{i}S$, iff the antecedent has a match in $S$
ending at position $i$ and the consequent has a match in $S$ starting at
position $i+1$, i.e. if $A\sqsubseteq^{i}S$ and $C\sqsubseteq_{i+1}S$.
Clearly, $R\Vdash_{i}S$ implies $R\vdash_{i}S$.

The _support_ of a rule is the number of times it applies in the sequence,
whereas its _confidence_ is the ratio of the number of times the rule applies
to the number of times the rule triggers:

$\operatorname{supp}_{S}(R)=\left\lvert\\{i\in[1,n]:\;R\Vdash_{i}S\\}\right\rvert$

and

$\operatorname{conf}_{S}(R)=\frac{\left\lvert\\{i\in[1,n]:\;R\Vdash_{i}S\\}\right\rvert}{\left\lvert\\{i\in[1,n]:\;R\vdash_{i}S\\}\right\rvert}\,.$

###### Example 2

Rules.  In sequence $S$ from Example 1, rule $R:ab\Rightarrow c$ has support
$\operatorname{supp}_{S}(R)=2$ and confidence
$\operatorname{conf}_{S}(R)=2/3$, because $R$ triggers at positions $1$, $5$
and $11$ but only applies at positions $1$ and $5$.

Given a sequence $S$ up to position $m$ and a rule $R:A\Rightarrow C$, we say
that $R$ is _active_ at stage $j$ ($0\leq j<\left\lvert C\right\rvert$) and
_predicts_ $C[j+1]$, if $A\oplus C[1,j]=S[m-\left\lvert A\right\rvert-j,m]$.
This corresponds to cases where the antecedent of $R$ is followed by the first
$j$ elements of its consequent (empty in case $j=0$) in $S$. Formally, we
define the predicates $\operatorname{active}(R,S,j)$, to indicate that rule
$R$ is active on sequence $S$ at stage $j$, and
$\operatorname{predict}(R,S,j)$, to return $C[j+1]$ if
$\operatorname{active}(R,S,j)$ is true.

##  IV Encoding Scheme

To illustrate our encoding scheme, we resort to a sender-receiver metaphor,
where the sender first transmits a set of rules and their corresponding
weights, followed by a sequence of code words, one for each element of the
sequence, in such a way that the receiver can reconstruct the original
sequence.

Encoding a Sequence via Sequential Rules. Assume we want to predict the next
element after $cabba$, based on two rules $R_{1}:abba\Rightarrow b$ and
$R_{2}:ba\Rightarrow cd$. If the confidence of $R_{1}$ is higher than the
confidence of $R_{2}$, we would guess that $b$ is more likely than $c$ to
appear next. We rely on this simple principle to encode a sequence element by
element using a set of sequential rules. Assuming the receiver knows the set
of rules and their confidence, as represented by weights, as well as the
history of the sequence decoded so far, they can compute a probability
distribution over the possible next elements. Hence, the sender can represent
the next element as a code word chosen according to this probability
distribution. If the set of rules effectively produces a probability
distribution concentrated on the element that indeed occurs next, the code
word for that element will be short, allowing for a concise representation of
the input sequence.

Given a sequence $S$ of $m$ elements seen so far, and a set of weighted rules
$(A_{R}\Rightarrow C_{R},w_{R})$, such that $A_{R}$, $C_{R}$, and $w_{R}$ are
respectively the antecedent, the consequent, and the weight of rule $R$, our
goal is to compute the probability of the next element given this set of rules
– denoted by $\mathcal{R}$. This amounts to computing a probability
distribution over $\Sigma$. For each possible next element $\sigma\in\Sigma$,
we sum the weights of the active rules that predict $\sigma$ and divide by the
sum of weights of all active rules:

$P_{\mathcal{R},S}(\sigma)=\frac{\sum_{(R,j)\in\mathcal{A}_{S,\sigma}}w_{R}}{\sum_{(R,j)\in\mathcal{A}_{S}}w_{R}}$

where $\mathcal{A}_{S}=\\{(R,j),(R,w_{R})\in\mathcal{R}\text{
s.t.}\,\operatorname{active}(R,S,j)\\}$, and
$\mathcal{A}_{S,\sigma}=\\{(R,j)\in\mathcal{A}_{S}\text{
s.t.}\,\operatorname{predict}(R,S,j)=\sigma\\}$. A full sequence $S$ of length
$n$ is transmitted element by element. At each step $m$, we encode the next
symbol $S[m]$ with a code word chosen according to the probability assigned to
it based on $\mathcal{R}$ and the portion of the sequence seen so far,
$S[1,m-1]$. Hence, the overall code length for the full sequence is

$\operatorname{L}(S|\mathcal{R})=\sum_{m\in[1,n]}-\log_{2}\big{(}P_{\mathcal{R},S[1,m-1]}(S[m])\big{)}.$
(2)

The receiver, having the same collection of rules $\mathcal{R}$ and previous
elements of the sequence, can perfom the same computation to dynamically
reconstruct the code at his end, and decode the transmitted element, using an
agreed canonical order on the alphabet to break ties if necessary.

To ensure comparability, the models must achieve _lossless_ compression. For
this reason, by construction our models always contain the singleton rule set
$\Psi_{\Sigma}$. This way, every symbol of the alphabet receives a non-zero
probability of occurrence at any step, and can hence be transmitted. In fact,
the singleton rule set constitutes the basis of the simplest possible model,
what we call the _empty model_ and denote by $\mathcal{R}^{\emptyset}$, since
it does not contain any proper rule. To each rule $R_{\sigma}\in\Psi_{\Sigma}$
we associate a weight equal to the _background probability_ of the predicted
element, estimated as its frequency in the sequence,
$f_{\sigma}=\left\lvert\\{i\in[1,n]:S[i]=\sigma\\}\right\rvert/n$.

Encoding the Rules. For the receiver to have access to the rules
$\mathcal{R}$, the sender must transmit the antecedents, consequents, and
weights at the start of the exchange. The antecedent and consequent are
subsequences, hence they can be encoded simply by stating the length of the
sequence then listing their elements using codes of length
$-\log_{2}(f_{\sigma})$, in order. The weight, on the other hand, is a real
number of bounded precision in $(0,1)$. Each weight is written as a finite
list of decimals that we encode by applying universal coding to the integer
value obtained by listing these decimals in reverse order. Put differently, we
encode weight $w=0.d_{1}d_{2}d_{3}...d_{k}$ with a code word of length
$\operatorname{L}_{D}(w)=\operatorname{L}_{\mathbb{N}}(d_{k}...d_{3}d_{2}d_{1})$.
$\operatorname{L}_{\mathbb{N}}(z)$ is the length of the code word assigned by
the universal code to integer $z$, which is such that
$\operatorname{L}_{\mathbb{N}}(z)=\log_{2}^{*}(z)+\log_{2}(c_{0})$, with
$c_{0}$ a constant adjusted to ensure the Kraft inequality is satisfied. This
penalizes weights by the number of significant digits they contain rather than
by their value. Overall, the code length for a rule $R$ is the sum of the
lengths of the code words for the antecedent, the consequent and the weight,
respectively:

$\operatorname{L}(R)=\left[\operatorname{L}_{\mathbb{N}}(\left\lvert
A_{R}\right\rvert+1)+\sum_{\sigma\in A_{R}}-\log_{2}(f_{\sigma})\right]$
$\qquad+\left[\operatorname{L}_{\mathbb{N}}(\left\lvert
C_{R}\right\rvert)+\sum_{\sigma\in
C_{R}}-\log_{2}(f_{\sigma})\right]+\operatorname{L}_{D}(w_{R})\,.$

We encode a rule table by stating the number of rules and then listing them

$\operatorname{L}(\mathcal{R})=\operatorname{L}_{\mathbb{N}}(\left\lvert\mathcal{R}\right\rvert)+\sum_{R\in\mathcal{R}}\operatorname{L}(R)\,.$
(3)

##  V The Cossu algorithm

We now introduce an algorithm to mine a set of sequential rules that
compresses the input sequence well under the encoding scheme presented in
Section IV. Finding such a set of rules is intractable in practice, thus we
resort to heuristics.

Our Cossu (COmpact Sets of Sequential rUles) algorithm is outlined in
Algorithm 1. Given a sequence $S$ over alphabet $\Sigma$ as input, Cossu
returns a set of rules $\mathcal{R}$ in two phases. First, _rule construction_
(lines 3–4) generates a collection of candidate rules $\mathcal{C}$, which are
then evaluated during _rule selection_ (lines 5–12).

Rule Construction. Cossu starts by extracting closed frequent subsequences by
applying an off-the-shelf sequence miner [17] to the input sequence $S$ with a
minimum support threshold of $2$. Sequential rules are then generated from the
sequences by considering all possible partitionings into an antecedent and a
non-empty consequent. We use the _compression gain_ as a rough estimate of the
individual ability of a rule to compress the sequence:

$\operatorname{gain}(R)=\operatorname{conf}_{S}(R)\cdot\operatorname{supp}_{S}(R)\cdot\operatorname{L}(C_{R})-(\operatorname{L}(A_{R})+\operatorname{L}(C_{R}))\;.$
(4)

This score puts in balance the potential benefit and cost of adding the rule
to the rule set. The benefit depends on how often the consequent is correctly
predicted (first term), whereas the cost is the code length of the rule’s
antecedent and consequent (second term). The latter term constitutes a lower
bound on the cost, with a weight yet undetermined. At the end of this first
phase, the set of candidates consists of those rules that have a strictly
positive compression gain.

Rule Selection. Cossu resorts to a greedy strategy to build the rule set
$\mathcal{R}$. Initially, $\mathcal{R}$ consists of the singleton rules (line
5). The algorithm then processes the candidate rules in $\mathcal{C}$ by order
of decreasing compression gain, i.e. from most promising to least promising
(line 6). The candidate rule is tentatively added to the rule set and the
weights are re-adjusted (line 7). If the resulting rule set yields better
compression than the current one, the current set is replaced (line 9) and the
algorithm goes into a pruning loop. Otherwise, the rule is discarded and the
algorithm moves on to the next candidate. The pruning loop (lines 10–12)
checks whether the newly added rule makes any of the previously incorporated
non-singleton rules obsolete. To do so, it tentatively removes each of the
rules in turn and checks whether the compression improves as a result, in
which case the rule is permanently removed from $\mathcal{R}$.

Algorithm 1 Cossu: finding a compact set of sequential rules

1:A long sequence $S$ of elements over an alphabet $\Sigma$

2:A compact set of sequential rules $\mathcal{R}$

3:$\mathcal{S}\leftarrow\text{MineClosedSequences}(S)$

4:$\mathcal{C}\leftarrow\\{R\in\text{PrepareRules}(\mathcal{S}),R\not\in\Psi_{\Sigma}\text{
and }\operatorname{gain}(R)>0\\}$

5:$\mathcal{R}\leftarrow\mathcal{R}^{\emptyset}$, with adjusted weights
$\hat{\vec{w}}=\arg\min_{\vec{w}}\operatorname{L}(S|\mathcal{R}^{\emptyset})$

6:for $R\in\mathcal{C}\text{, ordered by decreasing }\operatorname{gain}(R)$
do

7: $\mathcal{R}^{\prime}\leftarrow\mathcal{R}\cup\\{R\\}$,
with re-adjusted weights
$\hat{\vec{w}}=\arg\min_{\vec{w}}\operatorname{L}(S|\mathcal{R}^{\prime})$

8: if
$\operatorname{L}(\mathcal{R}^{\prime},S)<\operatorname{L}(\mathcal{R},S)$
then

9: $\mathcal{R}\leftarrow\mathcal{R}^{\prime}$

10: for $R^{\prime}\in\mathcal{R}\setminus\mathcal{R}^{\emptyset}$ do

11: if
$\operatorname{L}(\mathcal{R}\setminus\\{R^{\prime}\\},S)\leq\operatorname{L}(\mathcal{R},S)$
then

12: $\mathcal{R}\leftarrow\mathcal{R}\setminus\\{R^{\prime}\\}$

13:return $\mathcal{R}$

Adjusting Rule Weights. Given a set of selected rules, Cossu must determine a
suitable vector of weights $\vec{w}$ associated to the rules, so as to
optimize the code length of both the model $\operatorname{L}(\mathcal{R})$
(complexity) and the sequence $\operatorname{L}(S|\mathcal{R})$ (fit). Because
optimizing both aspects concurrently is very challenging, we assume that the
code length of the weights is fixed, by fixing the floating point precision of
the weights, and focus on the problem of minimizing
$\operatorname{L}(S|\mathcal{R})$, i.e. solving

$\arg\min_{\vec{w}}\sum_{m\in[1,n]}-\log_{2}\big{(}P_{\mathcal{R},S[1,m-1]}(S[m])\big{)}\;.$

We do so greedily, by adjusting the weight of each rule in turn, while keeping
the value of the other weights fixed. The minimization resorts to the golden
section search algorithm [18], initializing the weights to $1$. In the end, we
scale all weights to ensure that their values lie in the interval $(0,1)$.
Designing an algorithm that is able to optimize both
$\operatorname{L}(\mathcal{R})$ and $\operatorname{L}(S|\mathcal{R})$ is a
major direction for future work.

##  VI Evaluation

We evaluate our proposed algorithm on both synthetic and real-world datasets.
In particular, we evaluate the quality of the rules learned by Cossu on real-
world data, by using them as features in two classical machine learning tasks
on long sequences: next event prediction and classification. While Cossu is
not designed specifically for these tasks, we show it still yields acceptable
performance and reports rules that capture useful regularities in the data.

###  VI-A Evaluation on synthetic data

To study the behaviour of our algorithm in controlled settings, we generate
synthetic random sequences of symbols and inject patterns that match a set of
hand-crafted sequential rules. Then, we verify whether Cossu is able to
recover the planted rules.

Data generation and experimental setting. The data generation process
constructs an initial sequence $S$ by drawing $n$ symbols from a distribution
$\mathscr{D}$ on an alphabet $\Sigma$. We then insert a set $\mathscr{T}$ of
non-trivial rules (target) in the data. For each rule $R\in\mathscr{T}$, and
with a probability equal to $\mathit{ip}$, we insert the consequent of the
rule after every match of the antecedent in $S$. This makes the sequence $S$
longer, thus we truncate it to keep only the $n$ first elements.

By default $n=5000$, $\Sigma=\\{A,B,C,D,E\\}$, $\mathscr{R}=\\{A\rightarrow
B\\}$, with insertion probability $\mathit{ip}=50\%$, and $\mathscr{D}$ is an
equiprobable distribution on the elements of $\Sigma$. In each experiment we
vary one of these parameters and leave the others at their default value.

We highlight that because of the data and rule generation strategies, the
insertion probability of a rule is not equivalent to its actual confidence $c$
in the final sequence. Their relationship is given by the expression
$c=\mathit{ip}+(1-\mathit{ip})\cdot P(B|A)$, where $P(B|A)$ is the confidence
of $A\rightarrow B$ in the original sequence.

Our evaluation metric is the _hit rate_ , defined as the percentage of times
Cossu retrieves the exact set of rules that were introduced, which we call the
gold standard. We report the hit rate across $100$ different sequences.

TABLE I: Hit rate and average number of candidate rules $\left\lvert\mathcal{C}\right\rvert$ for different insertion probabilities ($\mathit{ip}$) of the single rule $A\rightarrow B$. $\mathit{ip}$ | 5% | 10% | 20% | 30% | 40% | 50% | 60% | 70% | 80% | 90% | 100%
---|---|---|---|---|---|---|---|---|---|---|---
Hit rate (%) | 3 | 75 | 100 | 100 | 100 | 100 | 100 | 100 | 100 | 94 | 100
$\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$) | 13.7 | 13.7 | 13.8 | 14.0 | 14.2 | 14.5 | 14.8 | 15.1 | 15.5 | 15.8 | 16.2

Different rule insertion probability. In Table I we show the hit rate and
average number of rules found by Cossu for different insertion probabilities
when inserting a single rule. Except when $\mathit{ip}<20\%$, Cossu always
achieves a hit rate of between $90\%$ and $100\%$. For low insertion rates,
the confidence of the target rule is closer to the background probabilities of
the 5 symbols. This makes Cossu unable to distinguish the rules from the
background and explains the low hit rates ($\leq 75\%$).

TABLE II: Hit rate and average number of candidate rules $\left\lvert\mathcal{C}\right\rvert$ for different alphabet sizes ($\left\lvert\Sigma\right\rvert$), no rule/one rule. $\left\lvert\Sigma\right\rvert$ | 2 | 3 | 5 | 10 | 50 | 100
---|---|---|---|---|---|---
Hit rate (%) | 100 / | 42 | 100 / | 100 | 100 / | 100 | 100 / | 100 | 100 / | 100 | 100 / | 100
$\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$) | 54.2 / | 58.3 | 26.1 / | 28.0 | 13.6 / | 14.5 | 7.0 / | 7.3 | 3.3 / | 3.3 | 1.9 / | 1.9

Different alphabet size. In Table II, we change the size of the alphabet
$\Sigma$ when generating the original sequence. We can see that when no rule
is inserted, Cossu does not return any rule either even if the number of
possible candidate rules $\left\lvert\mathcal{C}\right\rvert$ is often huge
(and would have been output without the MDL constraints). The hit rate is only
$42\%$ when there are only 2 symbols in the alphabet and we insert one single
rule ($A\rightarrow B$). However, when $\left\lvert\Sigma\right\rvert>2$,
Cossu is able to find the inserted rule (and only this rule).

TABLE III: Hit rate and average number of candidate rules $\left\lvert\mathcal{C}\right\rvert$ for different sequence lengths, no rules/one rule. $\left\lvert S\right\rvert$ | 100 | 500 | 1000 | 2000 | 5000 | 10000
---|---|---|---|---|---|---
Hit rate (%) | 100 / | 16 | 100 / | 99 | 100 / | 100 | 100 / | 100 | 100 / | 100 | 100 / | 100
$\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$) | 0.1 / | 0.1 | 1.0 / | 1.0 | 2.2 / | 2.3 | 5.8 / | 5.1 | 13.7 / | 14.5 | 29.8 / | 31.5

Different sequence length. In Table III we show the hit rate of Cossu when the
sequence length varies. We can observe that if the sequence is too short
($\left\lvert S\right\rvert\leq 500$) the rule $A\rightarrow B$ inserted is
not occurring enough to be retrievable by our MDL-based algorithm: On average,
for $\left\lvert S\right\rvert=100$ the antecedent $A$ appears $20$ times so
the subsequence $AB$ is present around 10 times in the sequence and hence does
not constitute a useful regularity. When the sequence becomes long enough
($\left\lvert S\right\rvert>500$), Cossu correctly finds the inserted rule
only.

TABLE IV: Hit rate and average number of candidate rules $\left\lvert\mathcal{C}\right\rvert$ for different inserted rule size (each time, only one rule is inserted in the generated sequence, with $\mathit{ip}=0.5$). Injected rule | $A\rightarrow B$ | $A\rightarrow BC$ | $A\rightarrow BCD$ | $AB\rightarrow C$
---|---|---|---|---
Hit rate (%) | 100 | 100 | 100 | 100
$\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$) | 14.475 | 15.510 | 16.500 | 13.798
Injected rule | $AB\rightarrow CD$ | $AB\rightarrow CDE$ | $ABC\rightarrow D$ | $ABC\rightarrow DE$
Hit rate (%) | 100 | 100 | 43 | 87
$\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$) | 13.956 | 14.155 | 13.655 | 13.703

Different rule size. In Table IV, we study the behaviour of Cossu when there
are more than one symbol in the antecedent and/or consequent parts of the
inserted rule. Note that the longer the antecedent is, the smaller the chances
for the rule to occur in the sequence are, and the less likely Cossu will need
the rule to compress the sequence. On the contrary, the longer the consequent
is (which is inserted each time the antecedent is present in the original
data), the more likely Cossu will select the rule to compress the data (esp.
for short antecedents). The hit rate is $100\%$ for all rules where the size
of the antecedent is lower than $3$ symbols. However, it is lower ($43\%$)
when the antecedent size is greater or equal to $3$ symbols but again better
($87\%$) when the consequent of the rule is longer.

Different target rule set size. In Table V we study the behaviour of Cossu
when more than one rule is inserted in the synthetic sequence. Note that
$A\rightarrow B\rightarrow C$ means that rules $A\rightarrow B$ and
$B\rightarrow C$ were inserted. Recall also that the lower the number of
symbols $\left\lvert\Sigma\right\rvert$, the more difficult it is for Cossu to
find the rule in the data. We can see that Cossu is stable, as it is able to
retrieve any number of inserted rules in any configuration. In the difficult
case where $\left\lvert\Sigma\right\rvert=3$ (with 2 rules), Cossu adds an
additional rule $\emptyset\rightarrow BC$ in 4 of the runs as a way to assign
a more precise weight to the rules.

TABLE V: Hit rate and average number of candidate rules $\left\lvert\mathcal{C}\right\rvert$ for different alphabet size and different inserted rules (with $\mathit{ip}=0.5$ for all the inserted rules). Injected rules | Hit rate (%) | $\left\lvert\Sigma\right\rvert$ |  $\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$)
---|---|---|---
$A\rightarrow B\,;\,C\rightarrow D$ | 100 | 4 | 19.865
$A\rightarrow B\,;\,C\rightarrow D\,;\,E\rightarrow F$ | 100 | 6 | 12.957
$A\rightarrow B\rightarrow C$ | 96 | 3 | 30.405
$A\rightarrow B\rightarrow C\rightarrow D$ | 100 | 4 | 21.880
$A\rightarrow B\rightarrow C\rightarrow D\rightarrow E$ | 100 | 5 | 17.647
$A$ $\rightarrow$ | $100$ | $7$ | $12.977$
$B\rightarrow E$ $\rightarrow$ |  |  |
$D\rightarrow F\rightarrow G$ |  |  |
$C$ $\rightarrow$ |  |  |

Different alphabet distribution. Finally, in Table VI we evaluate the
behaviour of Cossu when the discrete probability distribution of the symbols
in the generated data (before the rule insertion) is not uniform. We consider
distributions where a symbol (the antecedent $A$ or the consequent $B$ of the
inserted rule $A\rightarrow B$) is significantly more (denoted “strong”) or
less (denoted “weak”) represented compared to the other symbols. In the
“increasing” case, the probabilities of the two symbols A and B concerned by
the inserted rule are lower than the probabilities of the other symbols: this
makes it easier for Cossu to retrieve the rule (hence the $99\%$ hit rate). In
the “decreasing” case, the base probability of B is already high and the
additional occurrences produced by the rule insertions in the data are often
negligible compared to the original number of occurrences of B (hence the
$72\%$ hit rate). We can observe that, for the most unbalanced distribution
named “strong antecedent” (resp. “strong consequent”) in the table, the hit
rate drops to $0\%$ (resp. $42\%$). In the first case the rule systematically
retrieved is $\emptyset\rightarrow AB$ instead of $A\rightarrow B$. In the
second case the base probability of the symbol $B$ is already high and, as
discussed before, the rule cannot be easily spotted. Overall, this shows that
Cossu is sensitive to distributions with a particularly high entropy. This is
further explored on real-world data in the following sections.

TABLE VI: Hit rate and average number of candidate rules $\left\lvert\mathcal{C}\right\rvert$ for different alphabet distributions when inserting rule $A\rightarrow B$.  Distribution | P(A) | P(B) | P(C) | P(D) | P(E) | Hit rate (%) |  $\left\lvert\mathcal{C}\right\rvert$ ($10^{3}$)
---|---|---|---|---|---|---|---
Uniform | $1/\phantom{1}5$ | $1/\phantom{1}5$ | $1/\phantom{1}5$ | $1/\phantom{1}5$ | $1/\phantom{1}5$ | $100$ | 14.5
Increasing | $1/15$ | $2/15$ | $3/15$ | $4/15$ | $5/15$ | $99$ | 15.8
Decreasing | $5/15$ | $4/15$ | $3/15$ | $2/15$ | $1/15$ | $72$ | 17.3
Weak ante. | $1/41$ | $10/41$ | $10/41$ | $10/41$ | $10/41$ | $98$ | 16.2
Weak cons. | $10/41$ | $1/41$ | $10/41$ | $10/41$ | $10/41$ | $84$ | 16.0
Strong ante. | $10/14$ | $1/14$ | $1/14$ | $1/14$ | $1/14$ | $0$ | 26.2
Strong cons. | $1/14$ | $10/14$ | $1/14$ | $1/14$ | $1/14$ | $42$ | 27.4

In summary, we observed in this evaluation that Cossu is able to recover the
planted rules under good conditions and that performance degrades as noise
increases.

###  VI-B Evaluation on a prediction task

In this experiment, our dataset consists of a life log posted on the site
quantifiedawesome.com111http://quantifiedawesome.com/records that records the
daily activities of its owner, an enthusiast of the “Quantified Self”
movement, since 2011. We downloaded the publicly available data for the period
11/2012–08/2019 (46445 events). This dataset strikes an interesting balance
between regularity in the activities (so that patterns can be learnt) and
variations in the daily patterns (so that learning is far from trivial). Also,
the interpretation of the rules does not require special domain expertise.

Cossu: $0.7134$ |  HMM 15: $0.7040$
---|---
Bigram: $0.6979$ |  HMM 15: $0.7166$

Figure 1: ROC curves and corresponding Area Under Curve (AUC) comparing the prediction performance of Cossu on the Quantified Awesome dataset, against two HMM with 5 and 15 hidden states, and a 2-gram predictor with different minimum confidence thresholds (left). Most frequent rules found by Cossu on Quantified Awesome dataset (right). TABLE VII: Results for next event prediction on the Quantified Awesome dataset. Note that the baseline prediction does not depend on a confidence threshold ($\tau$). Method |  $\tau$ |  00 % |  020 % |  040 % |  060 % |  080 %
---|---|---|---|---|---|---
Cossu |  Precision (%) |  $\phantom{0}26.67$ |  $\phantom{0}38.57$ |  $\phantom{0}56.93$ |  $\phantom{0}73.69$ |  $\phantom{0}85.05$
|  Recall (%) |  $100.00$ |  $\phantom{0}69.91$ |  $\phantom{0}41.03$ |  $\phantom{0}18.16$ |  $\phantom{00}9.17$
|  F1 score (%) |  $\phantom{0}42.11$ |  $\phantom{0}49.71$ |  $\phantom{0}47.69$ |  $\phantom{0}29.14$ |  $\phantom{0}16.56$
HMM 5 |  Precision (%) |  $\phantom{0}24.30$ |  $\phantom{0}26.28$ |  $\phantom{0}39.54$ |  $\phantom{0}58.49$ |  $\phantom{0}85.87$
|  Recall (%) |  $100.00$ |  $\phantom{0}96.51$ |  $\phantom{0}56.73$ |  $\phantom{0}21.13$ |  $\phantom{00}4.83$
|  F1 score (%) |  $\phantom{0}39.10$ |  $\phantom{0}41.31$ |  $\phantom{0}46.60$ |  $\phantom{0}31.04$ |  $\phantom{00}9.15$
HMM 15 |  Precision (%) |  $\phantom{0}27.08$ |  $\phantom{0}28.34$ |  $\phantom{0}39.99$ |  $\phantom{0}57.40$ |  $\phantom{0}78.23$
|  Recall (%) |  $100.00$ |  $\phantom{0}98.45$ |  $\phantom{0}71.18$ |  $\phantom{0}27.60$ |  $\phantom{0}11.98$
|  F1 score (%) |  $\phantom{0}42.62$ |  $\phantom{0}44.01$ |  $\phantom{0}51.21$ |  $\phantom{0}37.28$ |  $\phantom{0}20.78$
Bigram |  Precision (%) |  $\phantom{0}25.99$ |  $\phantom{0}32.26$ |  $\phantom{0}52.24$ |  $\phantom{0}67.72$ |  $\phantom{0}82.23$
|  Recall (%) |  $100.00$ |  $\phantom{0}81.68$ |  $\phantom{0}35.95$ |  $\phantom{0}15.69$ |  $\phantom{0}11.74$
|  F1 score (%) |  $\phantom{0}41.26$ |  $\phantom{0}46.25$ |  $\phantom{0}42.59$ |  $\phantom{0}25.48$ |  $\phantom{0}20.55$
Baseline |  Precision (%) $11.53$ Recall (%) $100.00$ F1 score (%) $20.68$
TABLE VIII: Basic statistics of our experimental datasets for classification (left) and accuracy of the different methods (right). The best performing method is in bold. _Training sequence_ is the length of the sequence used by Cossu to learn the rules for the target classes $c_{1}$ and $c_{2}$. $|\Sigma|$ is the number of different words after pre-processing the data. A dash (-) indicates that the classifier is not applicable to the data. Dataset | Nb. of | Training sequence |  $|\Sigma|$ | Cossu | SVM | HMM | BERT
---|---|---|---|---|---|---|---
| instances | $c_{1}$ | $c_{2}$ |  |  |  |  |
Presidential Debate | 555 | 17658 | 20457 | 2414 | 0.798 | 0.779 | 0.803 | 0.731
Newsgroups | 180 | 13373 | 10895 | 5235 | 0.96 | 0.96 | 0.95 | 0.98
Critics | 1800 | 8722 | 8891 | 5558 | 0.707 | 0.697 | 0.707 | 0.705
Quant. Awesome | 1095 | 14562 | 5698 | 109 | 0.72 | 0.59 | 0.69 | -
Stylometry | 2000 | 14273 | 17331 | 3966 | 0.72 | 0.7 | 0.75 | 0.75
Stylometry (PoS) | 1200 | 11342 | 10629 | 31 | 0.64 | 0.59 | 0.605 | -

Experimental setup. We divide the data into multiple train/test setups with
one year (365 consecutive days) of training data and the immediately following
week as testing data. These train/test couples are built for each of the 350
such test weeks (starting on Mondays) for which the data is available. Each
test week is a sequence of events and a prediction is made for every event of
the sequence, considering as input the sequence of events up to the event to
predict (excluded).

To predict the next event in the sequence we used the rules mined by Cossu to
compute the probability of apparition of every symbol of the alphabet (as the
inverse of the log of the symbol’s code length). The symbol with the highest
probability is predicted if its probability is higher than a given minimum
confidence threshold $\tau$.

Competitors. We compare Cossu to methods that can both predict next events and
provide some understanding about the data. These methods are:

Baseline

predicts the same event as exactly one year before at this time;

Bigram

predicts as follows: all pairs of events $aX$ (bigrams) are retrieved from the
training data, and the frequencies of the different symbols $X$ following an
$a$ are computed. The highest frequency symbol is predicted if its frequency
is higher than the minimum confidence threshold $\tau$.

Hidden Markov Models (HMM)

are trained on a one-year sequence and used to compute, for each symbol of the
alphabet, the probability to be the next one. The highest probability symbol
is predicted if its probability is above $\tau$. We test HMMs with 5 and 15
hidden states.

Results. We present the results at different confidence thresholds in Table
VII.

The F1 performance of the baseline is poor (20.68%), and Cossu is always
better, except for $\tau=80\%$ where very few rules from Cossu are allowed to
make a prediction. This makes the recall plummet but gives one of the best
precision values of the experiments. The table also shows that Cossu
outperforms bigrams, and is comparable to HMM with 15 hidden states. This is
confirmed by the representation of the same results as a ROC curve shown in
Figure 1 (which also provides the AUC of the curves). Overall, HMMs have
higher recall, while Cossu has higher precision. This illustrates well the
difference of the two approaches: an HMM learns a global model of the data
that can cover a wider range of situations and is more robust to noise. On the
downside, this global model is not interpretable. In contrast, Cossu’s rules
can be seen as a set of interpretable local models that are combined for
prediction: they might do very well on some portions of the sequence but be
less suitable on noisy portions.

As an illustration, Figure 1 shows the top-10 rules most frequently selected
by Cossu during the training phases. For example, the rule
$\textit{Cook}\rightarrow\textit{Dinner}$ has appeared in Cossu’s output in
86.86% of the training phases. A point is displayed for each rule on the test
weeks where this rule was discovered during training (one year immediately
before). We note that among rules identified by Cossu some are applicable
across the entire sequence and represent very mundane activities, such as
$\textit{Cook}\rightarrow\textit{Dinner}$ or
$\textit{Dinner}\rightarrow\textit{Clean\;the\;kitchen}$, whereas others are
more specific to particular time periods. For example, in the year preceding
week 200 in the graphic, the log’s author gave birth to her first child. This
drastically changed her daily habits, as can be seen for instance by the
disappearance of the rule $\textit{Read\;fiction}\rightarrow\textit{Sleep,
Routine}$, and the emergence of the rule $\textit{Subway,
Childcare}\rightarrow\textit{Subway}$.

###  VI-C Evaluation on a classification task

In this experiment, we use the following datasets:

Quantified Awesome

is the same dataset as in Section VI-B, for which we predict whether a
sequence of events happened on a weekday or during the weekend;

Presidential debates

     222https://www.kaggle.com/mrisdal/2016-us-presidential-debates

consists of 555 sentences uttered by Donald Trump and Hillary Clinton in the
context of the 2016 US election; the goal is to predict the speaker given a
sentence;

Newsgroups

     333https://kdd.ics.uci.edu/databases/20newsgroups/20newsgroups.html

is a dataset of 180 posts about the topics _electronic_ and _religion_ ;

Film critics

     444https://github.com/clairett/pytorch-sentiment-classification/tree/master/data/SST2

is a dataset of 1800 movie reviews used for sentiment analysis;

Stylometry

     555https://bit.ly/3wBqITA

is a database of 2000 writings from H.P. Lovecraft or E.A. Poe, where the goal
is to predict the author based on its writing style;

Stylometry (PoS)

is the sequence of part-of-speech tags from the previous dataset.

For all datasets (except Quantified Awesome), our alphabet $\Sigma$ consists
of all the words present in the corpus after removing special characters and
stop words, and stemming of the remaining words. All datasets were split into
training, validation, and test subsets. Table VIII shows basic statistics of
the datasets.

Experimental setup. In line with the MDL-based classifier proposed in [7], we
can use Cossu to classify sequences as follows: For each class $c$, we build a
_training sequence_ by concatenating all the sequences labeled with $c$. We
then extract sets of rules with Cossu on the resulting long sequence. This set
of rules defines a model $H_{c}$ that characterizes the generation process for
sequences in class $c$. To predict the class of an input sequence $S$, we
compute the description length $L(S|H_{c})$ for each class $c$ and return the
class that compresses $S$ best, i.e., the class that minimizes $L(S|H_{c})$.

We compare the Cossu classifier to (i) a SVM classifier trained on a bag-of-
words representation of the text instances, (ii) an HMM-based classifier, and
(iii) logistic regression trained on top of the BERT language model [19].

Results. Classification accuracies on the test sets are reported in Table
VIII. For all methods, the reported performance was obtained with
hyperparameters tuned for best performance in cross-validation. Despite not
being a native classifier, Cossu achieves the best performance on three of the
experimental datasets, and exhibits comparable performance to HMM and state-
of-the-art text classifiers. Furthermore, and unlike its competitors, Cossu is
inherently interpretable because its rules serve as an explicit explanation
for the outcome of classification. As anecdotal examples, we show a few rules
characterizing the classes of the Stylometry dataset:

E.A. Poe: | rue $\rightarrow$ morge
---|---
| ourang $\rightarrow$ outang
| little old $\rightarrow$ gentleman
H.P. Lovecraft: | pnakotic $\rightarrow$ manuscript
| catch eight $\rightarrow$ clock coach arkham
| necronomicon mad $\rightarrow$ arab abdul alhazred

##  VII Conclusion

We have presented Cossu, an algorithm to mine sequential rules from a long
sequence of symbolic events. Our approach is inspired by the MDL principle and
retrieves a compact and relevant set of rules. Experiments on real-world data
show that, in addition to providing an interpretable model for sequential
data, the retrieved rules achieve competitive accuracy on the tasks of next
event prediction and a classification, as compared to other symbolic methods.

As future work, we would like to design a procedure to optimize the weights of
the rules while dynamically adjusting their precision. We would also like to
allow our method to handle gaps, perhaps by taking inspiration from [8]. An
implementation of the method is available at
https://gitlab.inria.fr/lgalarra/cossu.

Acknowledgements. The authors would like to thank Nikolaj Tatti for providing
a fast implementation of closed gapless sequential patterns miner.

## References

  * [1] S. Hochreiter and J. Schmidhuber, “Long Short-Term Memory,” _Neural Computing_ , vol. 9, no. 8, pp. 1735–1780, 1997.
  * [2] P. Fournier-Viger, J. C.-W. Lin, R. U. Kiran, Y. S. Koh, and R. Thomas, “A Survey of Sequential Pattern Mining,” _Data Science and Pattern Recognition_ , vol. 1, no. 1, pp. 54–77, 2017.
  * [3] B. Cule and B. Goethals, “Mining Association Rules in Long Sequences,” in _Advances in Knowledge Discovery and Data Mining_ , 2010, pp. 300–309.
  * [4] P. Fournier-Viger, T. Gueniche, and V. Tseng, “Using Partially-Ordered Sequential Rules to Generate More Accurate Sequence Prediction,” in _Advanced Data Mining and Applications_ , 2012, pp. 431–442.
  * [5] F. Fumarola, P. F. Lanotte, M. Ceci, and D. Malerba, “CloFAST: Closed Sequential Pattern Mining Using Sparse and Vertical Id-lists,” _Knowledge and Information Systems_ , vol. 48, no. 2, pp. 429–463, 2016.
  * [6] P. Fournier-Viger, T. Gueniche, S. Zida, and V. S. Tseng, “ERMiner: Sequential Rule Mining Using Equivalence Classes,” in _Proc. of the Int. Symp. on Intelligent Data Analysis_ , 2014, pp. 108–119.
  * [7] J. Vreeken, M. van Leeuwen, and A. Siebes, “Krimp: Mining Itemsets that Compress,” _Data Mining and Knowledge Discovery_ , vol. 23, no. 1, pp. 169–214, 2011.
  * [8] N. Tatti and J. Vreeken, “The Long and the Short of It: Summarising Event Sequences with Serial Episodes,” in _Proc. of ACM KDD_ , 2012, pp. 462–470.
  * [9] P. Fournier-Viger, R. Nkambou, and V. S.-M. Tseng, “RuleGrowth: Mining Sequential Rules Common to Several Sequences by Pattern-growth,” in _Proc. of ACM SAC_ , 2011, pp. 956–961.
  * [10] N. Gupta, N. Mangal, K. Tiwari, and P. Mitra, _Mining Quantitative Association Rules in Protein Sequences_ , 2006, pp. 273–281.
  * [11] P. Grünwald, _The Minimum Description Length Principle_. The MIT Press, 2007.
  * [12] J. Fischer and J. Vreeken, “Sets of Robust Rules, and How to Find Them,” in _Proc. of ECML-PKDD_ , 2019.
  * [13] R. Bertens, J. Vreeken, and A. Siebes, “Keeping it Short and Simple: Summarising Complex Event Sequences with Multivariate Patterns,” in _Proc. of ACM KDD_ , 2016, pp. 735–744.
  * [14] A. Bhattacharyya and J. Vreeken, “Efficiently Summarising Event Sequences with Rich Interleaving Patterns,” in _Proc. of SIAM ICDM_ , 2017, pp. 795–803.
  * [15] M. Shokoohi-Yekta, Y. Chen, B. Campana, B. Hu, J. Zakaria, and E. Keogh, “Discovery of meaningful rules in time series,” in _Proc. of ACM KDD_ , 2015, pp. 1085–1094.
  * [16] Y. Yan, L. Cao, S. Madden, and E. A. Rundensteiner, “SWIFT: Mining Representative Patterns from Large Event Streams,” _Proc. of the VLDB Endow._ , vol. 12, no. 3, pp. 265–277, 2018.
  * [17] N. Tatti and B. Cule, “Mining closed strict episodes,” _Data Min. Knowl. Discov._ , vol. 25, no. 1, pp. 34–66, 2012.
  * [18] J. Kiefer, “Sequential Minimax Search for a Maximum,” _Proc. of the American Mathematical Society_ , vol. 4, no. 3, pp. 502–506, 1953.
  * [19] J. Devlin, M.-W. Chang, K. Lee, and K. Toutanova, “BERT: Pre-training of Deep Bidirectional Transformers for Language Understanding,” in _Proc. of the NAACL Conference_ , 2019, pp. 4171–4186.