# BST1047+1156: A (Failing) Ultradiffuse Tidal Dwarf in the Leo I Group

J. Christopher Mihos  Department of Astronomy, Case Western Reserve
University, Cleveland OH 44106, USA  Patrick R. Durrell  Department of
Physics, Astronomy, Geology, and Environmental Sciences, Youngstown State
University, Youngstown, OH 44555 USA  Aaron E. Watkins  Centre for
Astrophysics Research, University of Hertfordshire, College Lane, Hatfield
AL10 9AB, UK  Stacy S. McGaugh  Department of Astronomy, Case Western Reserve
University, Cleveland OH 44106, USA  John Feldmeier  Department of Physics,
Astronomy, Geology, and Environmental Sciences, Youngstown State University,
Youngstown, OH 44555 USA

###### Abstract

We use deep Hubble Space Telescope imaging to study the resolved stellar
populations in BST1047+1156, a gas-rich, ultradiffuse dwarf galaxy found in
the intragroup environment of the Leo I galaxy group. While our imaging
reaches approximately two magnitudes below the tip of the red giant branch at
the Leo I distance of 11 Mpc, we find no evidence for an old red giant
sequence that would signal an extended star formation history for the object.
Instead, we clearly detect the red and blue helium burning sequences of its
stellar populations, as well as the fainter blue main sequence, all indicative
of a recent burst of star formation having taken place over the past 50–250
Myr. Comparing to isochrones for young metal-poor stellar populations, we
infer this post-starburst population to be moderately metal poor, with
metallicity [M/H] in the range $-1$ to $-1.5$. The combination of a young,
moderately metal-poor post starburst population and no old stars motivates a
scenario in which BST1047 was recently formed during a weak burst of star
formation in gas that was tidally stripped from the outskirts of the
neighboring massive spiral M96. BST1047’s extremely diffuse nature, lack of
ongoing star formation, and disturbed HI morphology all argue that it is a
transitory object, a “failing tidal dwarf” in the process of being disrupted
by interactions within the Leo I group. Finally, in the environment
surrounding BST1047, our imaging also reveals the old, metal-poor
([M/H]$=-1.3\pm 0.2$) stellar halo of M96 at a projected radius of 50 kpc.

Dwarf Galaxies — Galaxy evolution — Galaxy Interactions — Low Surface
Brightness Galaxies

††software:  astropy (The Astropy Collaboration, 2013, 2018, 2022), DOLPHOT
(Dolphin, 2000), numpy (Harris et al., 2020), matplotlib (Hunter et al.,
2007), scipy (Virtanen et al., 2020),

##  1 Introduction

The properties of extreme low surface brightness (LSB) galaxies continue to
challenge models of galaxy formation and evolution. While much attention has
been focused recently on the “ultra-diffuse galaxies” found in dense galaxy
clusters, gas-rich LSBs found in the field and group environments (e.g.,
McGaugh & Bothun, 1994; Cannon et al., 2015; Leisman et al., 2017) may have a
less complicated evolutionary path, and better probe mechanisms driving galaxy
formation at the lowest densities. For example, the high gas fractions and low
metallicities of LSB galaxies (McGaugh & Bothun, 1994; Ellison et al., 2008;
Pilyugin et al., 2014) argue that they have converted little of their baryonic
mass into stars. This is likely due to their extremely low gas densities (van
der Hulst et al., 1993; van Zee et al., 1997; Wyder et al., 2009), which
result in a sputtering and inefficient star formation history (Schombert et
al., 2001; Schombert & McGaugh, 2014, 2015). Thus, these galaxies raise
questions both macro and micro: how galaxy formation is linked to the global
environment, and how stars form on smaller scales within galaxies.

The recent discovery of the extreme LSB galaxy BST1047+1156 (Mihos et al.,
2018, hereafter BST1047) is particularly notable in this context. With an HI
velocity that places it unambiguously within the Leo I galaxy group (D=11 Mpc;
Graham et al., 1997; Lee & Jang, 2016), BST1047 has the lowest surface
brightness of any known star forming galaxy ($\mu_{B,peak}$=28.8 mag
arcsec-2), an isophotal radius of $R_{30}\approx$ 2 kpc, and a total gas mass
of $4.5\times 10^{7}$ $M_{\odot}$(Mihos et al., 2018). The object’s peak HI
column density ($1.4\times 10^{20}\ {\rm cm}^{-2}$) is well below that in
which stars typically form (Bigiel et al., 2008, 2010; Krumholz et al., 2009;
Clark & Glover, 2014), yet its extremely blue optical colors ($B-V=0.14\pm
0.09$) and GALEX far-UV emission both argue for the presence of young stars
(Mihos et al., 2018). BST1047’s combination of extraordinarily high gas
fraction ($f_{g}\approx 0.99$), extremely blue optical colors, and vanishingly
low surface brightness makes it the most extreme gas-rich LSB object known to
date.

Exactly how BST1047 formed, and what has triggered its recent star formation,
remains unclear. The Leo I group is awash in extended HI, including the large
“Leo HI Ring” surrounding NGC 3379 to the north (Schneider, 1985; Schneider et
al., 1986), likely a remnant of past tidal interactions (Michel-Dansac et al.,
2010; Corbelli et al., 2021). BST1047 itself is embedded in a low density HI
stream connecting the Ring to the spiral galaxy M96. This, plus the fact that
BST1047 sports a pair of HI tidal tails of its own, suggests the object may be
an extremely diffuse LSB galaxy recovering from a weak burst of tidally-
triggered star formation. Alternatively, BST1047 may be a “tidal dwarf
galaxy,” (e.g., Duc et al., 2000; Lelli et al., 2015), spawned directly from
tidally compressed gas, with the young stars marking its formation age. Since
tidal dwarfs should be free of dark matter (Barnes & Hernquist, 1992;
Elmegreen et al., 1993) and perhaps only tenuously bound, under this scenario
BST1047 may be a short-lived object — a “failing” tidal dwarf caught in the
throes of tidal disruption in the group environment.

Figure 1: Optical and HI imaging of the Leo I Group and BST1047+1156, taken
from Mihos et al. (2018). The upper left panel shows the deep wide-field
B-band imaging of Watkins et al. (2014), while the upper right panel shows the
HI map from Oosterloo et al. (2010). In these panels, the white box shows the
location of BST1047 in the Leo I group, and is blown up in the lower panels.
The lower left panel shows the B-band image, while the lower center panel
shows the B-image after being masked of compact sources and rebinned in
9$\times$9 pixel boxes to show low surface brightness emission. The lower
right panel shows the HI map on the same scale. In the lower panels, the
yellow circle is 70″ in radius, twice the size of the $R_{30}$ isophote.

Either of these scenarios has important ramifications for issues surrounding
theories of formation and evolution of low mass galaxies. If BST1047 is a
diffuse but long-lived LSB galaxy, with an established, old stellar
population, it challenges star formation models which posit that stars should
not form at such low gas densities. Under such models, where has the older
population come from? How can galaxies this diffuse sustain such prolonged
star formation histories? Conversely, if BST1047 is a disrupting tidal dwarf,
it would provide insight into the evolutionary link between tidal
interactions, formation and disruption of dwarf galaxies, and the deposition
of young stars into the intragroup medium. Key to resolving the question of
BST1047’s origin is an understanding of its stellar populations — in
particular, does it have a well established old red giant branch sequence,
indicative of an long-lived star forming history, or are the stellar
populations exclusively young, as might be expected BST1047 was recently
formed during a tidal encounter?

To answer these questions, we use deep Hubble Space Telescope ACS imaging to
study the stellar populations of BST1047. Using the F606W and F814W filters,
our imaging extends roughly two magnitudes below the expected tip of the red
giant branch at the Leo I distance, allowing us to detect and characterize
stellar populations across a range of ages, including any red giant branch
stars, red and blue helium burning stars, and potentially even upper main
sequence stars. These various populations give constraints on the ages and the
metallicities of both young and old stellar populations, providing strong
constraints on the extended star formation history in BST1047.

##  2 Observational Data

###  2.1 Imaging and Reduction

We imaged BST1047 using the Wide Field Channel (WFC) of the Advanced Camera
for Surveys (ACS) on the Hubble Space Telescope (HST) under program GO-16762.
The imaging field, shown in the left panel of Figure 2, places BST1047 in the
eastern side of the ACS field of view, avoiding nearby bright stars and
leaving the western side blank for background estimation.

The field was imaged over 8 orbits in F606W and 7 orbits in F814W; each orbit
consisted of two 1185s exposures, yielding total exposure times of 16458s and
16590s in F606W and F814W, respectively (one F606W exposure was cut short due
to guide star loss). Each visit made use of a small ($\sim$ 3.5–4.5 pixel)
custom four point box dither pattern to aid in sub-pixel sampling of the ACS
images and to also avoid placing any objects on bad or hot pixels. The
different visits were further shifted in slightly larger (20 pixel) offsets to
avoid other artifacts and facilitate effective cosmic ray removal in our long
(1/2 orbit) exposures. As BST1047 is small enough to fit into a single WFC
chip, the galaxy was centered on the WFC1 chip, and no attempt was made to
cover the ACS chip gap.

Point-source photometry is carried out on the individual, CTE-corrected flc
images using DOLPHOT (described below), which requires a sufficiently deep
drizzled image to use as an astrometric reference. To create this image, the
individual images from different visits needed to be precisely aligned. We
found that images from the three visits (including eight F606W images and two
F814W images) that were astrometrically calibrated with the GSC v2.4.2 catalog
were slightly offset ($\sim 0.5$pixel in F606W; $\sim 0.2$pixels in F814W)
from the remaining 20 images calibrated to the newer Gaia eDR3 catalog. To
improve the relative image alignments, we used the drizzlepac/tweakreg package
to adjust the image world coordinate systems based on point source positions
on the individual flc images measured using Source-Extractor (Bertin &
Arnouts, 1996). After these corrections, we used drizzlepac/astrodrizzle to
create stacked deep F606W and F814W images of the ACS field; the F814W image
is shown in the left panel of Figure 2.

Figure 2: Left panel: Stacked F814W ACS image of BST1047+1156, with a total
exposure time of 14$\times$1215s. North is up, east is to the left, and the
field of view is 202″$\times$202″. Right panel: Spatial distribution of point
sources detected in the ACS imaging. All point sources within the solid blue
ellipse are used for the analysis of BST1047, while point sources found
outside the dashed red ellipse are considered part of the background or
surrounding environment.

###  2.2 Point Source Photometry and Artificial Star Tests

With such an extremely low surface brightness
($\langle\mu_{B}\rangle_{e}=28.9$ mag arcsec-2; Mihos et al., 2018), the
integrated light from BST1047 is too faint to show up unambiguously in our ACS
imaging; instead, we only detect it through its resolved stellar populations.
We use the software package DOLPHOT (an updated version of HSTPhot; Dolphin,
2000) to perform point-source photometry of objects on the individual CTE-
corrected flc images using pre-computed Tiny Tim PSFs (Krist, 1995). We
performed object detection and photometry on all 30 individual images (16 in
F606W; 14 in F814W) images simultaneously, using the deep F814W drizzled image
stack created above for the reference image. We used the Nov 2019 version of
DOLPHOT 2.0 (available at http://americano.dolphinsim.com/dolphot/) to pre-
process the raw flc images, applying bad-pixel masks and pixel-area masks
(acsmask), splitting the images into the individual WFC1/2 chip images
(splitgroups), and constructing an initial background sky map for each
chip/image from each image (calcsky).

Photometry with DOLPHOT is very dependent on the choice of input parameters
(see Williams et al., 2014), so we experimented with a number of the
parameters, finally settling on values similar to those used in previous deep
photometric studies with ACS (e.g., Williams et al., 2014; Mihos et al., 2018;
Shen et al., 2021) and/or suggested by the DOLPHOT/ACS User’s Guide. As our
ACS field is relatively uncrowded, we adopted a photometric aperture RAper=4.0
pix, a PSF fitting region of RPSF=10 pix, and the FITSKY=1 option for
derivation of the sky background The only changes we made to the usual DOLPHOT
workflow were the derivation of the aperture corrections on each chip/image.
With so few bright stellar objects in our frames, some of the individual
DOLPHOT-computed aperture corrections (and thus the individual F606W/F814W
magnitudes) could be affected, even for brighter stars. To improve this, we
input our own visually-selected list of 53 isolated stellar objects over the
entire field which DOLPHOT could use to compute aperture corrections. The
final aperture corrections for each chip/image/filter were based on anywhere
from 6 to 29 measured stars. Finally, the instrumental magnitudes were
converted to the VEGAMAG HST photometric system. We used updated zeropoints
(at the time of observations, using the ACS zeropoint calculator
(https://acszeropoints.stsci.edu/) of 26.398 for F606W and 25.502 for F814W.
We present all photometry in the VEGAMAG system unless explicitly stated
otherwise.

To ensure the most accurate point source photometry, we apply the following
selection parameters to the photometric catalog. We start by selecting only
those objects with DOLPHOT object TYPE=1 (“good star”) and signal-to-noise
S/N$>$3.5 in both the F606W and F814W filters. We also only select sources
that are uncrowded (CROWD$<$0.25) and have a goodness of fit value of
CHI$<$2.4 in both filters; these values are based both on visual inspection of
bright stars and galaxies in our images, and the results from the artificial
stars detailed below. At fainter magnitudes (F814W$>$26), contamination from
unresolved background galaxies becomes problematic. To reduce this
contamination, we also make a magnitude-dependent cut on the DOLPHOT SHARP
parameter, using $|SHARP|<0.04+0.3e^{(m-m_{\rm crit})}$, with $m_{\rm crit}=$
29.5 and 28.7 in F606W and F814W, respectively. This function is similar to
that used in the our previous HST studies of stellar populations in M101 and
the Virgo Cluster (Mihos et al., 2018, 2022), and the function parameters are
chosen based both on the observed photometric catalog and on our artificial
star analysis. We have also checked that the sources rejected under our SHARP
criteria do not show stellar population-like patterns in the color-magnitude
diagram that might suggest we are over-aggressively rejecting actual stars in
the Leo I group environment. The spatial distribution of point sources
selected in this fashion are shown in the right panel of Figure 2.

To assess the photometric completeness and bias of the photometry in our ACS
imaging, we use DOLPHOT to insert and measure 100,000 artificial stars over
the magnitude range $22<{\rm F606W}<30$ and color range $-0.5<{\rm
F606W-F814W}<2.0$. We process the artificial stars using the same photometric
selection criteria used for the actual data, and plot in Figure 3 the
completeness fraction and shift in magnitude and color (defined as input minus
measured) as a function of F814W magnitude and F606W$-$F814W color. Because of
our joint selection in F606W and F814W, completeness is a function of both
magnitude and color, with 50% completeness at F814W=28.2 in the blue (at
F606W$-$F814W=0.0), and rising to F814W=27.8 in the red (at
F606W$-$F814W=1.0). At magnitudes brighter than F814W=27.0 we see little
systematic shift in either magnitude or color, but at fainter magnitudes
shifts in both are evident at the $\approx$ 0.1 mag level, consistent with our
previous analysis of ACS data in Mihos et al. (2018). In our analysis that
follows, we always plot magnitudes and colors as measured, correcting only for
foreground extinction ($A_{\rm F606W}=0.062,A_{\rm F814W}=0.038$, Schlafly &
Finkbeiner, 2011), and use the results of the artificial star tests to adjust
the theoretical stellar isochrones to account for these systematic effects
when interpreting our photometric results.

Figure 3: Results of artificial star tests. Left panel: Completeness, Middle
panel: F814W magnitude shift, Right panel: F606W$-$F814W color shift. Shifts
are measured as input minus measured values, such that a positive magnitude
shift corresponds to a star being measured systematically too bright, and a
positive color shift corresponds to a star being measured systematically too
blue. In each panel, the white line shows the 50% completeness limit.

##  3 Analysis

The right panel of Figure 2 shows the spatial distribution of point sources in
our ACS field; an excess of sources corresponding to the stellar population of
BST1047 can clearly be seen on the eastern half of the FOV. The distribution
of point sources appears slightly elongated roughly along the north-south
axis, and shows small scale clumpiness as well. We construct a color-magnitude
diagram (CMD) for BST1047 by extracting all point sources within an ellipse
(determined by eye and shown in Figure 2) centered at
$(\alpha,\delta)_{J2000}$=(10:47:43.59, 11:55:47.0), and having an ellipticity
of 0.85, semimajor axis of 50″, and position angle of 17°. The center of this
ellipse is approximately 14.3″ south of BST1047’s center coordinate originally
reported in Mihos et al. (2018). For comparison, we also construct a
background CMD by extracting sources that lie outside a 350 pixel (17.5″)
buffer around the BST1047 ellipse, shown as the dotted red ellipse in Figure
2. The extracted CMDs for each region (BST1047 and background) are shown in
the top panels of Figure 4.

While the background region is meant as a control for the BST1047 field, it
has a much larger area (by a factor of 2.87), and thus over-represents the
potential contamination to BST1047’s CMD. The lower right panel corrects for
this difference in area by randomly subsampling sources in the background
region by a factor of 2.87 to match the area of the BST1047 field, thus acting
as a more representative control sample for BST1047.

Figure 4: Color-magnitude diagrams for point sources in the ACS imaging. The
top panels shows the extracted CMDs for BST1047 (left) and the background
(right). The lower left panel shows the background CMD randomly subsampled to
match the area of BST1047, while the lower right panel shows the point source
CMD for the Abell 2744 Flanking Field data (taken from the analysis of Mihos
et al. 2018), randomly subsampled to match the area of the background field.
Each CMD in the lower panels thus acts as a control field for the CMD
immediately above it. In each subpanel, the right axis shows the apparent
magnitude, the left axis shows the absolute magnitude at the adopted 11.0 Mpc
distance of the Leo I Group, and the red dotted line shows the expected RGB
tip magnitude. Typical errorbars as a function of magnitude are shown in
black, and the 50% completeness limit for the BST1047 imaging is shown as the
dashed black line.

Of course, the background region itself is not a pure background. As BST1047
resides within the Leo I group, and also sits projected only 15′ (48 kpc)
northeast of the luminous spiral galaxy M96, our ACS pointing samples not only
background sources, but also stars in M96’s extended stellar halo, as well as
any potential Leo I intragroup stars (Watkins et al., 2014; Ragusa et al.,
2022). To estimate a cleaner background CMD, we turn to the deep HST imaging
of the Abell 2744 Flanking Field (Lotz et al., 2017). That imaging used the
same filters used here, and in Mihos et al. (2018) we extracted point source
photometry for that imaging using the same techniques as described above.
Thus, it acts as a reasonable control field for our background region here. In
the lower right panel of Figure 4 we show the Abell 2744 Flanking Field
photometry, using the same selection criteria as used in this study, and
subsampled down by a factor of 1.62 to match the area of our background
region.

We start with a discussion of the CMD in the full background region (upper
right panel of Figure 4), comparing it to its control field, the subsampled
Abell 2744 FF field directly below it. The most striking feature of the
background region is the clear signature of a metal-poor red giant branch
population, terminating at the expected RGB tip at F814W=26.2. Brighter than
this, there are a number of red stars in the field, possibly AGB stars or true
background contaminants. At these magnitudes, and over the small ACS field of
view, foreground contamination from Milky Way stars should be small; comparing
to the TRILEGAL models of Girardi et al. (2005); Girardi (2016) we would
expect only a handful of objects brighter than the observed RGB tip. At
fainter magnitudes (F814W$>$26) we also see a swarm of bluer sources with
colors $0.0<{\rm F606W-F814W}<0.5$, but these sources appear comparable in
number to those seen in the Abell 2744 Flanking Field, and are likely
unresolved background sources. Finally, we also see a handful of brighter
sources in this bluer color range, but not obviously in excess of the
background expectation.

Turning to the CMD for BST1047 itself, we again see a clear red sequence of
stars, but one that is distinctly bluer than that in the background region.
Whereas the sequence in the background region reaches a color of
F606W$-$F814W$\approx$1.03 when it reaches the RGB tip, the sequence in
BST1047 has a color of F606W$-$F814W$\approx$0.78 at a comparable brightness,
and continues on to brighter magnitudes above the expectation for the RGB tip.

We demonstrate this color difference in Figure 5, which shows the color
distribution of point sources of all colors, in the magnitude range $m_{\rm
tip}\leq{\rm F814W}\leq m_{\rm tip}+0.75$ (i.e., within 0.75 magnitudes of the
expected RGB tip). The left panel shows the relative color distribution in
each region, where the color difference between the red sequences in each
region is clear. The right panel of Figure 5 shows the surface density of
sources as a function of color – in other words, the color distribution
normalized by area. Here too, the difference in the red sequences is dramatic:
not only are they different in color, the density of red stars is much higher
in BST1047 than in the background. The two sequences are clearly tracing
different populations of stars.

Figure 5: Color distribution of point sources in the magnitude range $m_{\rm
tip,F814W}<{\rm F814W}<m_{\rm tip,F814W}+0.75$. The left panel shows the
relative fraction of sources in each field, while the right panel shows the
actual number density of sources.

At magnitudes brighter than the RGB tip, the CMDs in Figure 4 show an excess
of stars in BST1047 both at red colors (F606W$-$F814W $\approx$ 0.8–1.0) and
in the blue (F606W$-$F814W $\approx$ 0.0–0.4) compared to the background
region. The morphology and color of these bright red and blue sequences
suggest they are helium burning sequences from evolving massive stars,
signatures of recent star formation in BST1047. At fainter magnitudes (F814W
$\sim$ 26–28) we also see an excess population of stars with very blue colors
of F606W$-$F814W $<$ 0.0 compared to the background. These sources can also be
seen in the color distributions shown in Figure 5, and may represent massive
stars still on the upper main sequence.

In Figure 6 we overlay isochrones for old stellar populations of varying
metallicities onto the CMDs for both BST1047 and the background regions. We
use the PARSEC 1.2S isochrones (Bressan et al., 2012; Marigo et al., 2017),
with a fixed age of 10 Gyr, and with a range of metallicities spanning
[M/H]=$-2$ to $-0.7$. We adjust the tracks to reflect the systematic
photometric shifts discussed in §2.2, but note that in this portion of the CMD
the shifts are negligible at the RGB tip, and always $<$0.02 mag even down at
the 50% completeness limit. Looking at the background region, which likely
includes populations in M96’s outer stellar halo, the red sequence seen there
is well-matched by old RGB tracks with metallicity [M/H] in the range $-1$ to
$-1.5$. Presuming this is M96’s halo we are seeing, and adopting a typical
halo alpha abundance of $[\alpha/{\rm Fe}]=+0.3$, the metallicity corresponds
to [Fe/H] $\approx-1.2$ to $-1.7$ (Salaris et al., 1993; Streich et al.,
2014). These metallicities are similar to those found in the outer halos of
nearby spirals in the GHOSTS project (Monachesi et al., 2016), again arguing
these stars belong to the old halo population of M96. In contrast, the old
isochrones provide a poor match for the red sequence in BST1047; the most
metal-poor isochrone ([M/H]$=-2.0$) only reaches a color of F606W$-$F814W =
0.9, significantly redder than the mean color of the red sequence in BST1047
(F606W$-$F814W = 0.76). This red sequence in BST1047 likely then consists of
young red helium burning stars or very metal poor intermediate age RGB stars.

Figure 6: CMDs for BST1047 (left) and background (right) regions, overlaid
with PARSEC 1.2S isochrones for 10 Gyr old populations of varying
metallicities. These isochrones have been adjusted to reflect the small
systematic photometric shifts in magnitude and color in our data as determined
from the artificial star tests (see Section 2.2), but in the red these shifts
are negligible at magnitudes brighter than F814W $\approx$ 27.5.

We compare the CMD of BST1047 to younger isochrones in Figure 8, which
overplots the PARSEC 1.2S isochrones for a range of young and intermediate
ages, using metallicities of [M/H] = $-2.0$, $-1.3$, and $-0.7$. We again
adjust the isochrones for the systematic photometric shifts. At bright
magnitudes (F814W$<$27) the shifts remain negligible, but in the blue at
fainter magnitudes (F606W$-$F814W$<$0, F814W$>$27, the systematic blueward
shift becomes more noticeable, shifting the isochrones bluer by $\approx
0.05-0.1$ mag near the 50% completeness limit and leading to the slight
“bluish bulge” of the tracks in this region. With these effects in mind, these
tracks show that the most luminous stars, 1–2 magnitudes brighter than the RGB
tip, are consistent with blue and red helium burning sequences arising from
massive stars younger than a few hundred million years old. At fainter
magnitudes, the population of objects with very blue F606W$-$F814W colors
$<0.0$ may be tracing massive main sequence stars as young as 50 Myr. Looking
at the intermediate age RGB tracks, even at younger ages, RGB sequences are
still generally too red to match the red sequence we see in BST1047, except
perhaps at the very lowest metallicity ([M/H]$=-2$) and with relatively young
($<$2.5 Gyr) RGB populations. However, given the clear detection of BHeB stars
in the same region, the most natural explanation for the red sequence in BST
1047 is that it is the associated RHeB sequence, with little evidence for a
significant population of old RGB stars in the field.

In terms of metallicity, the sparseness of the stellar populations and the
tight spacing of tracks within the helium burning sequences make it hard to
give tight metallicity constraints for the population. Nonetheless, the stars
are clearly metal poor, with [M/H] $\approx-1.0$ or somewhat lower. More
metal-rich than this, the red helium burning sequences turn much redder than
observed in BST1047 where the sequence remains bluer than F606W$-$F814W = 1.0.
At the most metal-poor extreme, [M/H]=$-2.0$, both the red and blue helium
burning tracks start to shift bluer than seen in the observed CMDs, making it
unlikely that the populations are this metal-poor.

Figure 7: The CMD location (top) and spatial position (bottom) of the bright
variable star detected in our imaging. In both plots, the red symbol shows the
location of the variable. The CMD shows all point sources within the BST1047
region, which is shown by the blue oval in the lower panel. Lines and symbols
in the CMD are as described in Figure 4, and the lower figure shows the
202″$\times$202″ACS field of view, with North up and East to the left.

In our photometry, we also find one variable star with properties potentially
consistent with being a luminous Cepheid variable. The source is located 20″
north of the center of BST1047, at $(\alpha,\delta)_{J2000}$ = (10:47:43.32,
+11:56:08.4). It has a mean magnitude of roughly F814W $\approx$ 24.4 and
shows variability at the level of 0.35 mag in the individual ACS images,
significantly larger than the single-image relative magnitude uncertainty of
0.05 mag at that magnitude. Because of the sparse cadence of our observations,
secure photometry is difficult, but calculated from our two most concurrent
F606W and F814W images (separated by 4.5 days), the object has a color of
roughly F606W$-$F814=0.7, which would put it near the Cepheid instability
strip in the F606W/F814W color-magnitude diagram (see, e.g., McCommas et al.,
2009). While our data lack the proper cadence for accurate phasing, if the
object is a Cepheid in BST1047, then with an absolute magnitude of $M_{\rm
F814}=-5.8$ and using the F814W period-luminosity of Riess et al. (2019), the
object would have a period of 26 days, roughly twice the time span of our
imaging data, and consistent with the time variability we see in the source.
Without proper imaging cadence, it is difficult to place strong constraints on
the properties of the object, but if it is a Cepheid, that would also be
consistent with the other signatures of young massive stars that we observe in
BST1047.

Figure 8: CMDs for BST1047 overlaid with PARSEC isochrones of stellar
populations with varying ages and metallicities of [M/H]=$-2.0$ (left), $-1.3$
(middle), and $-0.7$ (right). These isochrones have been adjusted to reflect
the small systematic photometric shifts in magnitude and color in our data as
determined from the artificial star tests (see Section 2.2). These shifts in
the isochrones are most noticeable as a slight blue “bulge” in the tracks at
F814$>27.5$ and F606W$-$F814W $<0$.

So far in our analysis we have selected sources spatially, by region (sources
within BST1047 vs those in the surrounding background region) and considered
the CMDs of the regions separately. A complementary approach is to create a
CMD for the entire ACS field, sub-select sources by their location within this
full-field CMD, and ask where these subsamples are located spatially in the
field. We show such an analysis in Figure 9. The panel at the top of the
figure shows the CMD for the whole ACS image, where we have defined regions in
the CMD that highlight the different putative stellar populations discussed so
far. In particular, we highlight the blue and red helium burning sequences
(‘BHeB’ and ‘RHeB’, respectively), the blue main sequence region (‘MS’), and
the location of the older red giant branch population (‘RGB’). Selecting stars
that fall in these regions in the CMD, we then plot in the bottom panels the
spatial location of these CMD-selected sources.

Figure 9 clearly demonstrates that sources selected from CMD regions
corresponding to young populations — the MS, BHeB, and RHeB regions — are
preferentially found within BST1047. In contrast, objects drawn from the older
RGB region are spread much more evenly across the ACS field, with no
preferential clustering in or near BST1047. These spatial population patterns
are consistent with a scenario in which BST1047 is dominated by recent star
formation, with little or no evidence for an older stellar population. The
smoothly distributed old RGB stars in the field are much more likely to come
from M96’s stellar halo, with perhaps some additional contribution from Leo I
intragroup stars (although any such contribution is likely to be small;
Watkins et al., 2014; Ragusa et al., 2022). Additionally, there is a hint of a
weak gradient in the spatial distribution of RGB stars across the ACS field.
Comparing RGB counts on the western and eastern halves of the image, we find
54 RGB stars on the west side and 36 on the eastern side, roughly a $2\sigma$
difference. With M96 located 15′ southwest of the field, this gradient could
be tracing the radial dropoff in M96’s halo population, or just a signature of
patchiness in M96’s halo and/or the intragroup starlight in the region. For
now we leave additional analysis and a more detailed discussion of the
properties of the M96 halo population to a future paper.

Figure 9: Spatial distribution of point sources in our ACS field, selected by
their position on the color-magnitude diagram. The top panel shows the CMD for
the full ACS field, with regions color coded by their presumed evolutionary
stages. The bottom panels show the ACS field of view, mapping the point
sources corresponding to each selected CMD region. Field orientation and scale
is the same as in Figure 2, and the dotted oval shows the region containing
BST1047.

##  4 The Origin of BST1047+1156

The detection of massive young stars in the blue and red helium burning
sequences confirms a recent burst of star formation in BST1047, as originally
inferred from the very blue broad-band colors of the galaxy’s integrated light
(Mihos et al., 2018). The most luminous stars in these sequences have absolute
magnitudes of $M_{\rm F814W}\approx-5$, consistent with massive stars
($M_{*}\approx 6-7$ $M_{\odot}$) with lifetimes $<$ 100 Myr, but the sequences
continue down to fainter magnitudes ($M_{\rm F814W}\approx-3$) arguing for an
extended phase of star formation extending to at least 300 Myr ago. However,
this burst must have been relatively short lived; the lack of detected
H$\alpha$ emission in BST1047 (Donahue et al., 1995) sets an upper limit on
the present-day star formation rate of $\lesssim 5\times 10^{-5}$
$M_{\odot}$yr-1 (Mihos et al., 2018). While we defer a full population
modeling analysis to a future paper, we note here that not only are the red
and blue population sequences consistent with a recent short burst of
moderately metal-poor star formation, such a scenario is also quantitatively
consistent with both the star counts seen in the HST imaging and the
integrated light properties reported in Mihos et al. (2018). For example,
using the PARSEC stellar population modeling tools (Bressan et al., 2012;
Marigo et al., 2017, available at http://stev.oapd.inaf.it/cmd), a short (50
Myr) Gaussian burst of star formation of age $\sim 150$ Myr, mass
$M_{*}=2\times 10^{10}$ $M_{\odot}$, [M/H]=-1.5 and a Kroupa (2001) initial
mass function yields $\approx$ 40 bright stars in the range $-5\leq
M_{F814W}\leq-4$ and $\approx$ 180 stars in the fainter range $-4\leq
M_{F814W}\leq-3$. Performing a similar census of stars in the observed CMD and
correcting for background as shown in 4 yields 40 $\pm$ 7 stars and 160 $\pm
15$ stars in the brighter and fainter magnitude range, respectively. This
model also yields a total integrated B magnitude of $M_{B}=-10.1$ and color of
$B-V=0.08$, compared to the measured values of $M_{B}=-10.2\pm 0.14$ and
$B-V=0.14\pm 0.09$ (Mihos et al., 2018). Thus, all extant data are consistent
with a recent, fading post starburst population in BST1047 that formed within
the past few hundred million years.

Aside from the presence of high mass stars in BST1047, the other notable
feature of its CMD is the lack of any prominent red giant branch population.
While we cannot rule out a modest number of intermediate age ($<$ 5 Gyr) RGB
stars, they would need to be very metal poor, with metallicities
[M/H]$\sim-2$, to be hidden within the younger RHeB sequence. This lack of an
old stellar population in BST1047 is in marked contrast to stellar populations
of other types of diffuse star forming galaxies. The population of dwarf
irregulars in the Local Group, while showing a wide variation of star
formation histories, typically show old populations indicative of extended
star formation histories (e.g., Grebel, 1997; Weisz et al., 2014), including
even the extremely faint and diffuse dwarfs such as Leo T (Weisz et al., 2012)
or Leo P (McQuinn et al., 2015). Another natural comparison would be to the
population of blue, star-forming field low surface brightness galaxies.
Resolved stellar population work in field LSBs has identified the helium
burning sequences from evolving young stars (e.g., Schombert & McGaugh, 2014,
2015) as well as red giant branch stars that trace the older stellar
populations (Schombert & McGaugh, 2021). The lack of RGB stars in BST1047 thus
stands in contrast to the resolved populations in field LSBs, and more
generally the integrated colors of field LSB galaxies are typically much
redder than those of BST1047. For example, the $B-V$ color of BST1047 is 0.14
$\pm$ 0.09 (Mihos et al., 2018), compared to colors of $B-V\approx$ 0.3 – 0.6
for field LSBs (McGaugh & Bothun, 1994). Recent studies show that reddening
from dust in LSBs is typically quite low (Junais et al., 2023), arguing that
the redder colors of these field LSB galaxies indicate a substantial
contribution of light from old stars relative to what we observe with BST1047.
The lack of old stars in BST1047 likely then rules out scenarios where the
galaxy is simply a extremely low surface brightness outlier in the population
of star-forming field LSB galaxies.

The metallicities of the young stellar populations in BST1047 also argue
against a model in which the object is a pre-existing low mass dwarf galaxy.
Given the low inferred stellar mass for BST1047 (2–4 $\times
10^{5}$$M_{\odot}$, this work and Mihos et al. 2018), placing it on the mass-
metallicity relationship for dwarf galaxies (e.g., Kirby et al., 2013) would
predict a metallicity of [Fe/H]$=-2$, appreciably lower than the metallicity
inferred from the analysis shown in Figure 8.

A more likely scenario for BST1047 is that it formed during a tidal
interaction between galaxies within the Leo I Group. Tidal interactions can
strip gas from the gas-rich outer disks of spiral galaxies, expelling that gas
into the surrounding environment. Concurrently with the stripping the gas in
the tidal debris can be collisionally compressed, leading to a burst of star
formation and, potentially, to the formation of a tidal dwarf galaxy (e.g.,
Duc et al., 2000; Bournaud & Duc, 2006; Lelli et al., 2015). In this aspect,
BST1047 may be most similar to (albeit fainter and much more diffuse than)
tidal dwarf candidates found in the M81 group (Durrell et al., 2004; Mouhcine
& Ibata, 2009; Chiboucas et al., 2013), which also appear to lack old stellar
populations. Because tidal dwarfs form from pre-enriched material stripped
from a larger host galaxy, these objects should also be elevated in
metallicity compared to regular dwarf galaxies of the same mass (Duc &
Mirabel, 1998; Weilbacher et al., 2003), just as we find for BST1047. Indeed,
the metallicity of the young stars in BST1047 is comparable to that found in
the outskirts of large spirals (e.g., Zaritsky et al., 1994; van Zee et al.,
1998; Berg et al., 2020), and is also distinct from the higher, solar-like
metallicities found in the Leo Ring to the north (Corbelli et al., 2021).
Thus, the most likely origin for BST1047 is from gas that was stripped from
the outer disk of the spiral galaxy M96, as also suggested by the distorted
tidal HI morphology of M96’s outer disk (Figure 1 and Oosterloo et al., 2010).

Such a scenario would also explain the presence of massive young stars forming
in an object with such low gas density, well below that more typically found
in star forming environments (e.g., Bigiel et al., 2008, 2010; Wyder et al.,
2009). If BST1047 formed during a tidal encounter, the initial compression of
gas in the tidal caustics would have led to much higher gas densities capable
of driving a weak starburst now traced by the young populations in BST1047.
Subsequent tidal or ram-pressure stripping of the object, coupled perhaps with
energy input from stellar winds and supernovae from the evolving starburst
population, could have then dissociated the molecular gas and left BST1047
with a very diffuse ISM. The peak column density in BST1047 today is very low
($\sim$ 1 $M_{\odot}$ pc-2 Mihos et al., 2018), and likely incapable of
fueling any additional star formation. While the relatively large beam size of
the HI data leaves open the possibility of pockets of high density gas on
small scales, recent CO observations of BST1047 have failed to detect
molecular gas in the system (Corbelli et al., 2023), although those studies
did not survey BST1047’s full HI extent. Nonetheless, there is no evidence for
any current, ongoing star formation in BST1047 today.

The ultimate fate of BST1047 remains unclear. Its HI morphology (Figure 1)
shows streamers of HI extending to the southeast towards M96, likely a
signature of tidal stripping of BST1047, or perhaps ram-pressure stripping
from hot gas in M96’s halo or the group environment. BST1047’s very low
density makes it susceptible to stripping processes, particularly if it is a
tidal dwarf with no cocooning halo of dark matter to keep it bound. Mihos et
al. (2018) used the observed HI kinematics to show that BST1047’s dynamical
mass and baryonic mass were comparable ($\sim 5-6\times 10^{7}$ $M_{\odot}$),
providing support for the idea that the object lacks dark matter, as expected
for tidal dwarfs. If BST1047 is in the process of being disrupted by the
environment of the Leo I Group, we may be catching this object in a very
transitory phase, with its fading post-starburst population soon to be
stripped and expelled into the group environment. As such, the galaxy may be a
prime example of a “failing” tidal dwarf, born in the tidal debris of a recent
encounter, but lacking sufficient mass to overcome the destructive dynamical
processes found within the group environment.

##  5 Summary

We have used deep Hubble ACS imaging in F606W/F814W to study the resolved
stellar populations in the gas-rich ultradiffuse object BST1047+1156 in the
Leo I Group. At zero color, our photometry reaches limiting magnitudes of
F606Wlim=28.7 and F814Wlim=28.2, extending two magnitudes down the red giant
branch at the 11.0 Mpc distance of the Leo I Group. We clearly detect the
stellar population associated with BST1047, identifying the red and blue
helium burning sequences expected from an evolving population of massive
stars. We also find an excess of fainter blue stars likely to be slightly less
massive stars still on the main sequence. The distribution of color and
luminosity of stars in BST1047 are consistent with a modestly metal poor
stellar population ([M/H] $\sim-1.0$ to $-1.5$) with ages of a few hundred
million years, consistent with the integrated colors and surface brightness
measured in ground based imaging (Mihos et al., 2018).

However, we find no trace of a red giant branch sequence in the stellar
populations of BST1047 despite going sufficiently deep to detect such stars.
This lack of an old stellar population argues strongly against scenarios in
which BST1047 is a long lived LSB galaxy that has merely had a weak burst of
star formation due to interactions within the group environment. Instead, the
combination of its exclusively young and moderately metal-poor stellar
populations, its diffuse nature, and its disturbed HI morphology argues that
we are seeing a transient object, likely formed from gas recently stripped
from the outer disk of M96 due to tidal forces at work within the group
environment. These tidal forces continue to strip gas and stars away from
BST1047 today, feeding the the intragroup stellar population of the Leo I
Group. BST1047 is thus likely to be a failing tidal dwarf, formed from the
tidal debris of M96 but with such low density that it is destined to
ultimately disperse into the intragroup population of the Leo I group.

Finally, in the environment surrounding BST1047, we also clearly detect red
giant stars in the stellar halo of M96. These stars are distributed fairly
uniformly across the ACS field of view, showing no spatial correlation with
the location of BST1047. From the location of the red giant sequence on the
color-magnitude diagram we infer a moderately low stellar metallicity of [M/H]
$\approx-1.3\pm 0.2$. These data probe the stellar populations in the galaxy’s
outer halo at the extremely large projected radial distance of 50 kpc, and we
plan a future paper incorporating the data from the adjacent WFC3 parallel
field to study the properties of M96’s outer stellar halo in more detail.

The authors would like to thank Christian Soto and Norman Grogin for their
help with planning and refining the HST observations, as well as the anonymous
referee for helpful suggestions that improved the presentation of our results.
This research is based on observations made with the NASA/ESA Hubble Space
Telescope for program #GO-15258 and obtained at the Space Telescope Science
Institute (STScI). STScI is operated by the Association of Universities for
Research in Astronomy, Inc., under NASA contract NAS5-26555. Support for this
program was provided by NASA through grants to J.C.M. and P.R.D. from STScI.
A.E.W. acknowledges support from the STFC through grants ST/S00615X/1 and
ST/X001318/1. HST (ACS). The Hubble Space Telescope imaging data used in this
study can be accessed at the Mikulski Archive for Space Telescopes (MAST) at
the Space Telescope Science Institute via https://doi.org/10.17909/300a-ns91
(catalog DOI: 10.17909/300a-ns91).

## References

  * The Astropy Collaboration (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33. doi:10.1051/0004-6361/201322068
  * The Astropy Collaboration (2018) The Astropy Collaboration 2018, AJ, 156, 123. doi:10.3847/1538-3881/aabc4f
  * The Astropy Collaboration (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167. doi:10.3847/1538-4357/ac7c74
  * Barnes & Hernquist (1992) Barnes, J. E. & Hernquist, L. 1992, Nature, 360, 715. doi:10.1038/360715a0
  * Bertin & Arnouts (1996) Bertin, E. & Arnouts, S. 1996, A&AS, 117, 393.
  * Berg et al. (2020) Berg, D. A., Pogge, R. W., Skillman, E. D., et al. 2020, ApJ, 893, 96. doi:10.3847/1538-4357/ab7eab
  * Bigiel et al. (2010) Bigiel, F., Leroy, A., Walter, F., et al. 2010, AJ, 140, 1194. doi:10.1088/0004-6256/140/5/1194
  * Bigiel et al. (2008) Bigiel, F., Leroy, A., Walter, F., et al. 2008, AJ, 136, 2846. doi:10.1088/0004-6256/136/6/2846
  * Bournaud & Duc (2006) Bournaud, F. & Duc, P.-A. 2006, A&A, 456, 481. doi:10.1051/0004-6361:20065248
  * Bressan et al. (2012) Bressan, A., Marigo, P., Girardi, L., et al. 2012, MNRAS, 427, 127. doi:10.1111/j.1365-2966.2012.21948.x
  * Cannon et al. (2015) Cannon, J. M., Martinkus, C. P., Leisman, L., et al. 2015, AJ, 149, 72. doi:10.1088/0004-6256/149/2/72
  * Chiboucas et al. (2013) Chiboucas, K., Jacobs, B. A., Tully, R. B., et al. 2013, AJ, 146, 126. doi:10.1088/0004-6256/146/5/126
  * Clark & Glover (2014) Clark, P. C. & Glover, S. C. O. 2014, MNRAS, 444, 2396. doi:10.1093/mnras/stu1589
  * Corbelli et al. (2021) Corbelli, E., Cresci, G., Mannucci, F., et al. 2021, ApJ, 908, L39. doi:10.3847/2041-8213/abdf64
  * Corbelli et al. (2023) Corbelli, E., Thilker, D., Mannucci, F., et al. 2023, A&A, 671, A104. doi:10.1051/0004-6361/202244941
  * Dolphin (2000) Dolphin, A. E. 2000, PASP, 112, 1383. doi:10.1086/316630
  * Donahue et al. (1995) Donahue, M., Aldering, G., & Stocke, J. T. 1995, ApJ, 450, L45. doi:10.1086/316771
  * Duc et al. (2000) Duc, P.-A., Brinks, E., Springel, V., et al. 2000, AJ, 120, 1238. doi:10.1086/301516
  * Duc & Mirabel (1998) Duc, P.-A. & Mirabel, I. F. 1998, A&A, 333, 813
  * Durrell et al. (2004) Durrell, P. R., Decesar, M. E., Ciardullo, R., et al. 2004, Recycling Intergalactic and Interstellar Matter, 217, 90. doi:10.48550/arXiv.astro-ph/0311130
  * Ellison et al. (2008) Ellison, S. L., Patton, D. R., Simard, L., et al. 2008, ApJ, 672, L107. doi:10.1086/527296
  * Elmegreen et al. (1993) Elmegreen, B. G., Kaufman, M., & Thomasson, M. 1993, ApJ, 412, 90. doi:10.1086/172903
  * Girardi et al. (2005) Girardi, L., Groenewegen, M. A. T., Hatziminaoglou, E., et al. 2005, A&A, 436, 895. doi:10.1051/0004-6361:20042352
  * Girardi (2016) Girardi, L. 2016, Astronomische Nachrichten, 337, 871. doi:10.1002/asna.201612388
  * Graham et al. (1997) Graham, J. A., Phelps, R. L., Freedman, W. L., et al. 1997, ApJ, 477, 535. doi:10.1086/303740
  * Grebel (1997) Grebel, E. K. 1997, Reviews in Modern Astronomy, 10, 29
  * Harmsen et al. (2017) Harmsen, B., Monachesi, A., Bell, E. F., et al. 2017, MNRAS, 466, 1491. doi:10.1093/mnras/stw2992
  * Harris et al. (2020) Harris, C.R., Millman, K.J., van der Walt, S.J. et al. 2020, Nature, 585, 357
  * Hunter et al. (2007) Hunter, J.D., 2007, Computing in Science and Engineering, 9:3, 90.
  * Jang & Lee (2017) Jang, I. S. & Lee, M. G. 2017, ApJ, 836, 74. doi:10.3847/1538-4357/836/1/74
  * Junais et al. (2023) Junais, Małek, K., Boissier, S., et al. 2023, A&A, 676, A41. doi:10.1051/0004-6361/202346528
  * Kirby et al. (2013) Kirby, E. N., Cohen, J. G., Guhathakurta, P., et al. 2013, ApJ, 779, 102. doi:10.1088/0004-637X/779/2/102
  * Krist (1995) Krist, J. 1995, Astronomical Data Analysis Software and Systems IV, 77, 349
  * Kroupa (2001) Kroupa, P. 2001, MNRAS, 322, 231. doi:10.1046/j.1365-8711.2001.04022.x
  * Krumholz et al. (2009) Krumholz, M. R., McKee, C. F., & Tumlinson, J. 2009, ApJ, 699, 850. doi:10.1088/0004-637X/699/1/850
  * Lee et al. (2006) Lee, H., Skillman, E. D., Cannon, J. M., et al. 2006, ApJ, 647, 970. doi:10.1086/505573
  * Lee & Jang (2016) Lee, M. G. & Jang, I. S. 2016, ApJ, 822, 70. doi:10.3847/0004-637X/822/2/70
  * Lee et al. (2017) Lee, M. G., Kang, J., Lee, J. H., et al. 2017, ApJ, 844, 157. doi:10.3847/1538-4357/aa78fb
  * Leisman et al. (2017) Leisman, L., Haynes, M. P., Janowiecki, S., et al. 2017, ApJ, 842, 133. doi:10.3847/1538-4357/aa7575
  * Lelli et al. (2015) Lelli, F., Duc, P.-A., Brinks, E., et al. 2015, A&A, 584, A113. doi:10.1051/0004-6361/201526613
  * Lotz et al. (2017) Lotz, J. M., Koekemoer, A., Coe, D., et al. 2017, ApJ, 837, 97. doi:10.3847/1538-4357/837/1/97
  * Marigo et al. (2017) Marigo, P., Girardi, L., Bressan, A., et al. 2017, ApJ, 835, 77. doi:10.3847/1538-4357/835/1/77
  * McCommas et al. (2009) McCommas, L. P., Yoachim, P., Williams, B. F., et al. 2009, AJ, 137, 4707. doi:10.1088/0004-6256/137/6/4707
  * McGaugh & Bothun (1994) McGaugh, S. S. & Bothun, G. D. 1994, AJ, 107, 530. doi:10.1086/116874
  * McQuinn et al. (2015) McQuinn, K. B. W., Skillman, E. D., Dolphin, A., et al. 2015, ApJ, 812, 158. doi:10.1088/0004-637X/812/2/158
  * Michel-Dansac et al. (2010) Michel-Dansac, L., Duc, P.-A., Bournaud, F., et al. 2010, ApJ, 717, L143. doi:10.1088/2041-8205/717/2/L143
  * Mihos et al. (2018) Mihos, J. C., Carr, C. T., Watkins, A. E., et al. 2018, ApJ, 863, L7. doi:10.3847/2041-8213/aad62e
  * Mihos et al. (2018) Mihos, J. C., Durrell, P. R., Feldmeier, J. J., et al. 2018, ApJ, 862, 99. doi:10.3847/1538-4357/aacd14
  * Mihos et al. (2022) Mihos, J. C., Durrell, P. R., Toloba, E., et al. 2022, ApJ, 924, 87. doi:10.3847/1538-4357/ac35d9
  * Monachesi et al. (2016) Monachesi, A., Bell, E. F., Radburn-Smith, D. J., et al. 2016, MNRAS, 457, 1419. doi:10.1093/mnras/stv2987
  * Mouhcine & Ibata (2009) Mouhcine, M. & Ibata, R. 2009, MNRAS, 399, 737. doi:10.1111/j.1365-2966.2009.15135.x
  * Oosterloo et al. (2010) Oosterloo, T., Morganti, R., Crocker, A., et al. 2010, MNRAS, 409, 500. doi:10.1111/j.1365-2966.2010.17351.x
  * Pilyugin et al. (2014) Pilyugin, L. S., Grebel, E. K., Zinchenko, I. A., et al. 2014, AJ, 148, 134. doi:10.1088/0004-6256/148/6/134
  * Ragusa et al. (2022) Ragusa, R., Mirabile, M., Spavone, M., et al. 2022, Frontiers in Astronomy and Space Sciences, 9, 852810. doi:10.3389/fspas.2022.852810
  * Riess et al. (2019) Riess, A. G., Casertano, S., Yuan, W., et al. 2019, ApJ, 876, 85. doi:10.3847/1538-4357/ab1422
  * Salaris et al. (1993) Salaris, M., Chieffi, A., & Straniero, O. 1993, ApJ, 414, 580. doi:10.1086/173105
  * Schlafly & Finkbeiner (2011) Schlafly, E. F. & Finkbeiner, D. P. 2011, ApJ, 737, 103. doi:10.1088/0004-637X/737/2/103
  * Schneider (1985) Schneider, S. 1985, ApJ, 288, L33. doi:10.1086/184416
  * Schneider et al. (1986) Schneider, S. E., Salpeter, E. E., & Terzian, Y. 1986, AJ, 91, 13. doi:10.1086/113975
  * Schombert & McGaugh (2014) Schombert, J. & McGaugh, S. 2014, PASA, 31, e036. doi:10.1017/pasa.2014.32
  * Schombert & McGaugh (2015) Schombert, J. & McGaugh, S. 2015, AJ, 150, 72. doi:10.1088/0004-6256/150/3/72
  * Schombert & McGaugh (2021) Schombert, J. & McGaugh, S. 2021, AJ, 161, 91. doi:10.3847/1538-3881/abd54d
  * Schombert et al. (2001) Schombert, J. M., McGaugh, S. S., & Eder, J. A. 2001, AJ, 121, 2420. doi:10.1086/320398
  * Shen et al. (2021) Shen, Z., van Dokkum, P., & Danieli, S. 2021, ApJ, 909, 179. doi:10.3847/1538-4357/abdd29
  * Streich et al. (2014) Streich, D., de Jong, R. S., Bailin, J., et al. 2014, A&A, 563, A5. doi:10.1051/0004-6361/201220956
  * van der Hulst et al. (1993) van der Hulst, J. M., Skillman, E. D., Smith, T. R., et al. 1993, AJ, 106, 548. doi:10.1086/116660
  * van Zee et al. (1997) van Zee, L., Haynes, M. P., Salzer, J. J., et al. 1997, AJ, 113, 1618. doi:10.1086/118379
  * van Zee et al. (1998) van Zee, L., Salzer, J. J., Haynes, M. P., et al. 1998, AJ, 116, 2805. doi:10.1086/300647
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E. et al. 2020, Nature Methods, 17, 261
  * Watkins et al. (2014) Watkins, A. E., Mihos, J. C., Harding, P., et al. 2014, ApJ, 791, 38. doi:10.1088/0004-637X/791/1/38
  * Weilbacher et al. (2003) Weilbacher, P. M., Duc, P.-A., & Fritze-v. Alvensleben, U. 2003, A&A, 397, 545. doi:10.1051/0004-6361:20021522
  * Weisz et al. (2012) Weisz, D. R., Zucker, D. B., Dolphin, A. E., et al. 2012, ApJ, 748, 88. doi:10.1088/0004-637X/748/2/88
  * Weisz et al. (2014) Weisz, D. R., Dolphin, A. E., Skillman, E. D., et al. 2014, ApJ, 789, 147. doi:10.1088/0004-637X/789/2/147
  * Williams et al. (2014) Williams, B. F., Lang, D., Dalcanton, J. J., et al. 2014, ApJS, 215, 9. doi:10.1088/0067-0049/215/1/9
  * Wyder et al. (2009) Wyder, T. K., Martin, D. C., Barlow, T. A., et al. 2009, ApJ, 696, 1834. doi:10.1088/0004-637X/696/2/1834
  * Zaritsky et al. (1994) Zaritsky, D., Kennicutt, R. C., & Huchra, J. P. 1994, ApJ, 420, 87. doi:10.1086/173544


# Learning Functions to Study the Benefit of Multitask Learning

Gabriele Bettgenhäuser
Spoken Language Systems (LSV)
Saarland University
Saarbrücken, 66123 Germany
<EMAIL_ADDRESS>
&Michael A. Hedderich
Spoken Language Systems (LSV)
Saarland University
Saarbrücken, 66123 Germany
<EMAIL_ADDRESS>
&Dietrich Klakow
Spoken Language Systems (LSV)
Saarland University
Saarbrücken, 66123 Germany
<EMAIL_ADDRESS>

###### Abstract

We study and quantify the generalization patterns of multitask learning (MTL)
models for sequence labeling tasks. MTL models are trained to optimize a set
of related tasks jointly. Although multitask learning has achieved improved
performance in some problems, there are also tasks that lose performance when
trained together. These mixed results motivate us to study the factors that
impact the performance of MTL models. We note that theoretical bounds and
convergence rates for MTL models exist, but they rely on strong assumptions
such as task relatedness and the use of balanced datasets. To remedy these
limitations, we propose the creation of a task simulator and the use of
Symbolic Regression to learn expressions relating model performance to
possible factors of influence. For MTL, we study the model performance against
the number of tasks (T), the number of samples per task (n) and the task
relatedness measured by the adjusted mutual information (AMI). In our
experiments, we could empirically find formulas relating model performance
with factors of $\sqrt{n}$, $\sqrt{T}$, which are equivalent to sound
mathematical proofs in [1], and we went beyond by discovering that performance
relates to a factor of $\sqrt{AMI}$.

##  1 Introduction

Deep learning methods have achieved state-of-the-art performance in many
applications when learning to solve a single problem. In this domain, current
efforts are partly focused on engineering new neural architectures or
exploring novel methods to improve the accuracy of these single-task models.
There is another strain of research that explores learning models that share a
certain degree of knowledge among different related tasks. A common approach
is to use multitask learning (MTL).

In this setting, two or more related tasks are jointly trained such that it is
expected that the model achieves better generalization for all the tasks. The
main assumption is that related tasks might contain complementary information,
which should act as an inductive bias to guide the model towards a better
optimum as compared to training one model for each task in isolation [2].

Multitask learning has been applied to many problems over time, particularly
in NLP. Collobert and Weston [3] proposed a unified architecture to mutually
learn different sequence labeling tasks, a language modeling task and a
semantically related words problem. McCann et al. [4] used an MTL-based model
to solve 10 different tasks in natural language (including semantic parsing,
machine translation, summarization and sentiment analysis) by casting them as
question answering. Crichton et al. [5] proposed an MTL-based Named Entity
Recognition (NER) in biomedical text mining. In Biomedicine, there are a
variety of named entities datasets for the different sub-domains (e.g. genes,
proteins, chemicals and species names). Some of these datasets are very small,
which harms the performance of a model trained to classify named entities.
However, obtaining labeled datasets is usually done manually by experts and
are expensive to develop. Therefore, Crichton et al. [5] cast domain-specific
NER datasets as the different tasks and proposed a model to mutually learn
them to improve generalization. Their results disclose yet another benefit of
MTL models, by learning multiple tasks the model capitalizes from more data,
such that tasks with small training sets can generalize better.

Even though great average improvements were achieved, the results are often
mixed when evaluating each task in isolation and comparing them with the
state-of-the-art performance of single-task models [3, 5, 6, 7, 8]. The
challenge of selecting a set of tasks, which guarantees better generalization
for all tasks in multitask settings, remains partially unsolved. This is, to a
certain extent, due to the few theoretical investigations that clearly expose
the conditions under which multitask learning leads to better performance.

The few works that do try to theoretically understand generalization in MTL,
however, rely mostly on strong assumptions about how tasks are related or are
only applied to simple scenarios, such as settings where all tasks contribute
with datasets of the same size. Baxter [9] proposed generalization bounds for
MTL using VC-dimension and covering numbers, however, he assumes that task
relatedness is given a priori. Maurer [10] used the Rademacher complexity to
analyze linear multitask methods. In a later work, Maurer et al. [1] have
expanded their evaluation to nonlinear methods and assumed that similar tasks
shared a common feature representation. They proposed a bound based on the
Gaussian average while evaluating their theoretical results for noiseless
binary classification tasks.

The mixed results obtained in MTL motivate to study models that are more
theoretically grounded while avoiding being trapped into prior assumptions and
constrained problems. In this work, we skip strictly deriving bounds for a
well-defined domain and rather follow the scientific methodology that inspired
Johannes Kepler to be the first to correctly explain the elliptical orbit of
planets amongst other riddles involving planetary motion. In his work, Kepler
applied the scientific methodology of using observational data to discover
functions relating variables of interest. The laws empirically obtained by
Kepler were later mathematically proven by the theoretical works of Isaac
Newton.

Following an approach similar to Kepler’s, we will try to learn formulas from
data by using a method called Symbolic Regression. The empirically obtained
expressions are able to explain the generalization performance of multitask
models with respect to given parameter(s) (e.g. the number of tasks and task
relatedness). With that, we also expect to find a more elementary way to group
tasks in multitask models and elucidate the benefits of MTL.

#### Our Contributions

  * •

We developed a task simulator to create a multitask environment of sequence
labeling tasks. This simulator enables the user to control the number of
tasks, the level of relatedness between the tasks, and the number of samples
per task.

  * •

We used a novel approach based on Symbolic Regression to study the convergence
patterns of machine learning algorithms according to parameters of interest.
More specifically, we apply this methodology to the non-trivial case of
multitask learning models for sequence labeling tasks.

  * •

We discovered simple formulas that quantify the performance of MTL models as a
function of the controlled parameters enabled by the task simulator.

##  2 Background and Related Work

In the following, we supply the background information and related literature
to formally introduce the main concepts essential to this work.

#### Multitask Learning

Multitask learning is a method that aims to achieve better generalization by
learning $T$ related tasks jointly. The improvement in MTL happens because if
tasks are related, features that are learned for one task might benefit the
learning of some other related task. Multitask learning also has an implicit
data augmentation effect, assuming that different tasks might contribute with
their own data set. MTL will only induce better generalization if the
objectives of each task are not conflicting, which is usually the case when
tasks are related. Therefore, MTL is not robust when facing outlier tasks,
where performance can degrade quickly [11].

Caruana [2] describes multitask learning as a form of inductive bias. He
argues that learning one task at a time can be counterproductive because it
misuses the additional knowledge that training signals from different tasks in
the same domain usually have. Training signals from one task can improve a
model by introducing a bias when learning some other task which will cause the
whole model to prefer better hypotheses. The inductive bias behaves similarly
to regularizers. It prevents the model from overfitting since it is not trying
to optimize a single goal, but multiple objective functions for the different
tasks.

Maurer et al. [1] used the Empirical Risk Minimization (ERM) framework to
derive their theoretical results. A general formulation of ERM for multitask
learning settings adds the notion of optimizing for multiple tasks jointly.
Given $T$ tasks $\\{t_{1},t_{2},...,t_{T}\\}$, where each task has its own
dataset with samples of the form $z_{t}=(x_{t,i},y_{t,i})$. Suppose we are
given $n$ independent and identically distributed samples for each task and
these samples are generated from unknown distributions $p_{t}(z)$. Let
$Q(z_{t},w_{t})$ be a set of loss functions of the form
$L(y_{t},f(x_{t},w_{t}))$ and $w_{t}\in\Omega$. Then, the goal is:

$\min_{w_{1},...w_{T}}\frac{1}{T}\sum_{t=1}^{T}\frac{1}{n}\sum_{i=1}^{n}L(y_{t,i},f(x_{t,i},w_{t}))$
(1)

#### Bounds and Convergence Rates for Multitask Learning

Baxter [9] derived theoretical bounds for MTL based on the VC-dimension and
covering numbers theory under the assumption that tasks from the same
environment are related. The set of tasks from the same environment was
assumed to be probabilistically distributed. He was able to find bounds of the
order $O(\sqrt{\log{n}/n})$ and $O(\sqrt{\log{T}/T})$, where $n$ is the number
of samples per task and $T$ is the number of related tasks.

Other convergence rates for MTL were reported in the literature. The bound in
[1] is based on Gaussian average and yields convergence rates of the order of
$O(\sqrt{1/n})$ and $O(\sqrt{1/T})$. However, in this work, task relatedness
is not explicitly explored. Liu et al. [12] explore the notion of algorithmic
stability rather than using complexity measures to derive algorithm-dependent
upper bounds. Their bounds achieve fast convergence rates of order $O(1/n)$
and of order $O(1/T)$. Another usual assumption of many theoretical studies in
MTL is to encode a priori that the relations between the related tasks are
linear. To overcome this limitation, Ciliberto et al. [13] study the
generalization properties of MTL with nonlinear task relations. They found
that nonlinear MTL for a set of unrelated tasks converges at a rate of
$\tilde{O}((T/n)^{1/4})$ and that for the set of related tasks, this rate is
$\tilde{O}((1/nT)^{1/4})$.

#### Quantifying Task Relatedness in MTL

The assumptions made about the level of relatedness between tasks is one of
the key aspects to obtain successful results in MTL. Naively learning
unrelated tasks together may be detrimental and outlier tasks may lead to
negative transfers that can damage performance. In non-neural MTL, some
preliminary research assumed that task relatedness was given a priori. Some
studies supposed that all tasks were uniformly related and were controlled
globally by imposing that the Frobenius norm of the differences between
parameters of all tasks is close to each other [14, 15]. Other studies adopted
a pairwise approach, evaluating the similarity between all pairs of tasks. A
pair of tasks are similar when the difference between the parameters of them
is upper bounded by some constant, or when they would use features in a
similar way [16]. The results from [16] were further explored, such that one
could model clusters [17] or hierarchies [18] among tasks. Methods for
inferring task relatedness from data were also proposed. Bonilla et al. [19]
explored learning task relatedness using Bayesian models and Zhang and Yeung
[20] used a regularization approach.

In the context of neural MTL for NLP, task relatedness is not as explored as
for non-neural methods. Some works assume relatedness is given a priori so
that relationships between tasks are linguistically motivated [6, 21]. Recent
works also adopt a brute-force approach by training all pairs of tasks to
realize which pairs lead to gain and losses when jointly trained [7, 8]. Other
methods are focused on neural architecture design. Changpinyo et al. [8]
include a task embedding to learn task similarity and Ruder et al. [22]
propose multitask architecture learning that lets the network itself learn
which layers to share, the amount of sharing per layer and how to best design
the output layer so as to maximize the benefit of mutually learning several
tasks. Finally, Clark et al. [23] rather focused on the multitask training
algorithm. They use concepts of knowledge distillation such that a “teacher”
model (single-task model) transfers knowledge to a “student” model (multitask
model) by training the student to mimic the teacher’s outputs.

#### Symbolic Regression

Symbolic Regression [24] is a learning method to find a formula that relates
some input-output pair. A common implementation of Symbolic Regression is done
using Genetic Programming (GP) [25]. The building blocks for GP-based methods
are constants, variables and algebraic operations. The terminal set is
composed by the input variables and constants. The primitive function set
consists of the relevant operators (e.g. op = {$+,-,$ $\div$, $\times$,
$inv,abs,cos,sin,log$}). In general, an expression $(log(x_{2})+0.5)x_{1}$ is
represented as syntax tree or an S-expression, such as
$y=\times(+(log(x_{2}),0.5),x_{1})$. The GP-based method for Symbolic
Regression consists of 5 phases: creation of an initial population of random
formulas, evaluation of the current population, selection of individuals,
crossover, and mutation. For more information about each of these phases, one
can refer to [25].

#### Sequence Labeling

Many NLP tasks can be cast as a sequence labeling problem, where the goal is
to assign a label to each token in a sequence. Sequence labeling tasks serve
also as an intermediate step to achieve more complex problems, which are
usually of more practical utility to understand language. As an example,
question answering systems commonly have Named Entity Recognition (NER) as an
embedded component. The goal of NER is to identify a word or a set of words as
named entities (e.g person, organization, location, time and quantity).

Traditionally, sequence labeling problems were solved by probabilistic
graphical models including Hidden Markov Models [26], Conditional Random
Fields [27] and Maximum Entropy Markov Models [28]. These models heavily
relied on hand-designed features carefully engineered for a certain task,
language and dataset. [29] tackled these limitations by proposing neural-based
methods. Instead of task-specific engineering, they let the network discover
internal representations that are useful for learning a given task. More
recently, Transformer based models [30, 31] have been used as one of the
state-of-the-art architecture for problems involving sequential data.

##  3 General Approach

To tackle the challenge of quantifying the benefit of MTL, we propose the
following methodology (visualized in Figure 1). First, we developed a task
simulator to create data for $T$ artificial tasks while controlling parameters
of interest, such as the number of tasks, the number of samples per task and
task relatedness. Then, we trained and evaluated more than 10000 multitask
models with varying architectures using the generated data from the task
simulator. We record the performance of these models and the values of the
controlled parameters into a new dataset. Using this new dataset, we apply
Symbolic Regression to experimentally find simple formulas relating the
empirical model performance and the controlled variables.

Figure 1: Empirical methodology

#### Task Simulator

We developed a simulator to generate tasks under controlled conditions. This
simulator allows a user to create a multitask environment while controlling:

  * •

the task relatedness measured by the adjusted mutual information (proxied by a
parameter $\alpha$)

  * •

the number of artificial tasks

  * •

the number of samples per task

Our task simulator uses as a starting point a real dataset (in our experiments
from the NER domain). Based on this dataset, related tasks are generated. The
user controls relatedness between pairs of the form {RealTask,
ArtificialTaskt} using a matrix $M$ that regulates the relationship between
the real task tag-set (with labels {$r_{1},r_{2},..,r_{k}$}, where $k$ is the
number of different labels in the real task), and a generated tag-set of
synthetic labels ({$s_{1},s_{2},..,s_{c}$}, where $c$ is the number of
different labels in the synthetic task). Each entry of $M$ represents
conditional probabilities of the form $P(ArtificialLabel=s_{i}\mid
RealLabel=r_{j})$, and is defined as:

$\displaystyle M[i,j]=\begin{cases}\alpha&\text{if }i=j\\\
\frac{1-\alpha}{c-1}&\text{otherwise}\end{cases}$

where $1\leq i\leq k$ and and $1\leq j\leq c$. $\alpha$ is a parameter that
controls how much probability mass is kept in the diagonal of the matrix.
Hence when $\alpha=1$ we have the maximum dependence between the tasks and
when $\alpha=1/c$, we can say that the tasks are independent. The parameter
$\alpha$ is a proxy for the Adjusted Mutual Information (AMI) as defined in
Equation 2.

$AMI(U,V)=\frac{I(U,V)-E\\{I(U,V)\\}}{max\\{H(U),H(V)\\}-E\\{I(U,V)\\}}$ (2)

In the equation above, $H(U)$ and $H(V)$ are the entropy of two clusterings U
and V, respectively. E{I(U,V)} is a factor to correct the original Mutual
Information (MI) metric for randomness. The Adjusted Mutual Information
quantifies how much information is shared by two random clusters. This
quantity is symmetric $AMI(U,V)=AMI(V,U)$ and non-negative, where $AMI(U,V)=0$
means that $U$ and $V$ are independent. This metric is always in the interval
[0,1]. We omit further details about this adapted measure and rather refer to
the work [32] for additional information.

We build several of these matrices in order to simulate artificial tasks with
different levels of task relatedness. For simplicity, we create artificial
tasks with the same number of labels in the tag-set of the real task.
Algorithm 1 further outlines the implementation details about the task
simulator.

Data: Array of Real Task R, Number of Tasks T, Number of Samples per Task
n,Array of alpha

Result: Matrix of $T$ Synthetic Tasks SYN, Array of Adjusted Mutual
Information AMI

1 $\texttt{LabelSet}\leftarrow GetLabelSet(\texttt{R})$ ;

2 for  _t = 1 …T iterations_ do

3 $\texttt{M}\leftarrow
GenerateLabelMatrix(\texttt{LabelSet},\texttt{alpha[t]})$ ;

4 for  _i = 1 … n iterations_ do

5 $\texttt{probabilities}\leftarrow\texttt{M[R[i]]}$ ;

6 $\texttt{SYN[i][t]}\leftarrow
RandomChoice(\texttt{LabelSet},\texttt{probabilities})$ ;

7 $\texttt{AMI[t]}\leftarrow
AdjustedMutualInformation(\texttt{R},\texttt{SYN[t]})$ ;

8

9 end for

10

11 end for

Algorithm 1 Task Simulator

#### Architecture of MTL Models

In an attempt to discover formulas that are not particular to one specific
model, we created different MTL architectures as a combination of a shared BI-
LSTM layer and fully connected layers (see Figure 2). The difference in each
architecture lies in the number of dense layers shared or task-specific in
each model. Initially, the input is fed to a shared BI-LSTM layer, which will
then learn generic parameters that are available for all tasks during
training. Then, the output of BI-LSTM flows into fully connected layers that
can be either shared or task-specific. In a multitask setting, the output
layers are always task-specific.

Figure 2: Model architecture

#### Applying Symbolic Regression

By training and evaluating the MTL models on the multitask datasets generated
by the task simulator, we obtain a model performance for each dataset (e.g.
$F_{1}score$). These can be used to obtain new, different datasets on which
the symbolic regressors can be fit. These datasets follow the format
[$I_{1},I_{2},...,I_{p},Out$], where {$I_{1},I_{2},...,I_{p}$} are the
controlled factors (e.g. {$n,T,AMI$}) and $Out$ is the measure of model
performance (e.g. $F_{1}score$). Using datasets of this format, we were able
fit formulas with different number of parameters, including $Out=f(I_{p})$ and
$Out=f(I_{1},I_{2},...,I_{p})$. We evaluated the obtained formulas with the
Mean Squared Error (MSE) by using held-out datasets.

##  4 Experimental Setup

For the purpose of this work, we adapted the simulator to output artificial
sequence labeling tasks that are related to the Named Entity Recognition
problem.

#### Training and Evaluating MTL models

All models were implemented and trained using Pytorch [33]. The input for all
models was the words of the CoNLL dataset [34] embedded with the pre-trained
FastText vectors available in [35]. The input words were converted to
sentences of the form [Left Context, Target Word, Right Context]. The size of
the context window was set to 3. The output has the format [NER label,
Synthetic Task1 label, … , Synthetic Taskt label] corresponding to the labels
for the center word on $t+1$ tasks. We converted the encoding of the original
CoNLL dataset from the IOB tagging scheme to IO encodings for the training set
and IOB2 for the development and test sets.

We adopted a similar version of the fully joint strategy in [4] to train our
models. For comparability, we trained all models using the Adam Optimizer with
default parameters over 100 epochs. The batch size is also fixed to 100. At
each step, we backpropagate the total cross-entropy loss defined as the sum
over the losses for all tasks in the related set. To evaluate the accuracy of
our models we used the official CoNLL evaluation script [36]. In this work, we
report only the overall $F_{1}score$ metric for the Named Entity Recognition
task. Each $F_{1}score$ reported is an average over 5 runs.

#### Fitting Formulas

To find the formulas we used the gplearn package [37]. This package performs
Symbolic Regression using a GP-based implementation.

The formulas were obtained from datasets created for fitting functions of the
form $F_{1}score=f(T,n,AMI)$. We used a 90/10 train/test split for each
experiment. We chose a set from $\\{+,-,\div,\times,\log,\sqrt{,}exp,neg\\}$
to represent the primitive set of operators. The exponential operator was not
available in the library, therefore, we had to add it manually. We train all
models for at least 20 generations while controlling the diversity of formulas
at each generation under a 70% probability of crossover and 10% probability of
point, hoist and subtree mutation. We set the parsimony coefficient to 0.2 in
most experiments. This last hyperparameter controls the tradeoff between
fitness and expression length. Finally, we fit the formulas using the actual
values of X, but a scaled y according to $c*y$, where c is a constant which
was set to 100 and 1000 in the different experiments. We have done this
because, in our experience with the gplearn module, when y values vary from
[0,1], the formulas always converged to a constant close to the mean value of
y.

##  5 Experimental Results

It is our intention to study the benefit of multitask learning by adopting an
empirical approach. Albeit the references mathematically exploring the
advantage of MTL, we take an empirical approach other than theoretical
methodologies and sound numerical proofs. We do that for two reasons. First,
we wanted to verify if a purely empirical tool, such as Symbolic Regression,
would find similar results already theoretically explored. Secondly, we wanted
to avoid being trapped by prior assumptions but rather letting a learning
algorithm discover explanatory patterns by itself.

We conducted experiments to inspect the patterns of model performance when
considering the influence of (1) sample size per task; (2) adjusted mutual
information; (3) number of tasks; and (4) all these parameters together. The
results of our experiments were used to induce simple expressions to quantify
the advantage of MTL models trained under these variable conditions. Two
scenarios were explored. First, we performed experiments with balanced
datasets, where all tasks have the same number of samples, which is a common
assumption in theoretical studies. Secondly, we explored the more realistic
scenario where tasks have datasets with a different number of samples.

###  5.1 Experiments with Balanced Datasets

To fit the functions we create the MTL Balanced Datasets which contains 1002
observations on 5 features (n, T, AMI, $\alpha$, $F_{1}score$ ). Each row in
our dataset represents the average $F_{1}score$ of 5 models trained with the
same variables configuration. Hence, we trained and evaluated 5010 models to
obtain this dataset. We vary $T\in\\{1,2,3,4,...,15\\}$,
$n\in\\{1000,1500,2000,2500,3000,....,20000\\}$,
$AMI\in\\{0,0.13,0.17,0.21,0.27,0.34,0.43,0.55,0.57,0.6,0.63,0.68,0.71,0.75,0.8,0.85,0.91,1\\}$
and
$\alpha\in\\{0.2,0.6,0.65,0.75,0.8,0.85,0.87,0.9,0.91,0.92,0.93,0.94,0.95,0.96,0.97,0.98,0.99,1\\}$.

(a)

(b)

(c)

Figure 3: Results for experiments with balanced datasets

#### Considering Number of Samples per Task

We study the effect of the number of samples per task on the $F_{1}$score of
the NER task for single task (STL) models and MTL models. We vary the number
of tokens per task from 1k up to 20k samples for both settings. Figure 3(a)
shows the performance of the STL and 4 MTL models ($T=2$), with different
amounts of adjusted mutual information. For comparability, we only plot MTL
models with the same number of tasks ($T=2$). The plot highlights the positive
impact of the increased sample size on both STL and MTL settings, where MTL
outperforms STL models when AMI between the related tasks is high.

We also observe that all models have a similar curve shape when varying this
parameter, which characterizes a good scenario for fitting symbolic
regressors. We fit a symbolic regressor to discover a function of the form
$F_{1}score=f(n)$, where $n$ is the number of tokens per task. The regressor
returns the following function with a test MSE of 0.027.

$F_{1}score=\sqrt{n}*0.69$ (3)

#### Considering Adjusted Mutual Information

In this experiment, we investigate the impact of the adjusted mutual
information on the performance of the NER task. We vary the AMI from 0.13 to
1. Each color in Figure 3(b) represents a MTL model with different amounts of
samples per task. From this plot, we note that for all multitask scenarios, we
obtained improvement in generalization as the AMI in the task set increases.
In this experiment, we obtained empirical evidence that the benefit of
multitask learning is almost certain to be observed in low resource scenarios
even for low task relatedness. This advantage is still observed for high
resource settings, but the task set needs a higher AMI value.

Once more we found a function of the form $F_{1}score=f(AMI)$. The formula
below has a test MSE of 0.034.

$F_{1}score=\sqrt{AMI}/0.015$ (4)

#### Considering Number of Tasks

We vary the number of tasks up to 15 and evaluate the impact in the
performance of the NER task. Each task compromised 2k tokens. Figure 3(c)
contains different MTL models with varying levels of adjusted mutual
information. Evaluating this plot, we can conclude that all MTL models
outperform their STL equivalent and that not always the multitask models
benefit from adding more tasks.

In this experiment, symbolic regression failed to converge to a meaningful
formula. By observing the shape of the curve and the results from the previous
experiment, we guessed factors of $\sqrt{T}$ as a fitting formula. We manually
tried a few constants and report the expression with the lowest test MSE of
0.072.

$F_{1}score=\sqrt{T}*20$ (5)

#### Considering All Parameters

In this experiment, we wanted to obtain an expression with all three factors
previously studied. Instead of searching for formulas of the form
$F_{1}score=f(n)$, where the input was a unique variable, our goal was
searching for an expression $F_{1}score=f(n,AMI,T)$. We report the simplest
formula with the best fit among all outcomes from the symbolic regressor. The
best candidate expression with the lowest test MSE of 0.009 is reported below.

$F_{1}score=\sqrt{(log(T)+0.45)*(AMI*n/(AMI+0.6)))}$ (6)

#### Discussion

When we evaluate the 4 formulas obtained in this set of experiments, we could
infer convergence rates of:

  * •

$O(\sqrt{n})$ in the number of samples based on the Equations 3 and 6,

  * •

$O(\sqrt{AMI})$ in the amount of adjusted mutual information from the
Equations 4 and 6,

  * •

$O(\sqrt{T})$ from Equation 5 and $\sqrt{(log(T)}$ from Equation 6 in the
number of tasks.

While contrasting with previous theoretical work, we conclude that our results
are not far off the outcomes from sound mathematical proofs adopted in [1, 9,
13, 12]. In fact, the convergence rates from Equations 3 and 5 are equivalent
to Maurer et al. [1]. Besides, we were able to derive novel results relating
model performance with the adjusted mutual information according to Equations
4 and 6.

###  5.2 Experiments with Unbalanced Datasets

In this second set of experiments, we explore the scenario where tasks
contribute with datasets of different sizes. We generated data where the NER
task has always less data than the synthetic tasks. We created the MTL
Unbalanced Datasets which consists of 1555 observations on 6 features
($n_{ner}$, $n_{syn}$, $T$, $AMI$, $\alpha$, $F_{1}score$). We vary
$T\in\\{1,2,3,...,30\\}$, $n_{ner}\in\\{1000,1500,2000,2500\\}$,
$n_{syn}\in\\{2000,3000,...,6000\\}$,
$AMI\in\\{0.13,0.22,0.43,0.55,0.61,0.71,0.8,0.91\\}$ and
$\alpha\in\\{0.6,0.7,0.85,0.9,0.92,0.95,0.97,0.99\\}$ . We performed the same
set of experiments as for the balanced datasets and the resulting formulas
were:

$F_{1}score=\sqrt{n_{ner}}$ (7) $F_{1}score=\sqrt{n_{syn}}*0.69$ (8)
$F_{1}score=\sqrt{2*AMI}/0.026$ (9) $F_{1}score=log(\sqrt{T}+log(T))/0.039$
(10) $F_{1}score=(\sqrt{log(n_{ner}/0.38)*(\sqrt{T}+T)})*\sqrt{n_{syn}*AMI}$
(11)

##  6 Conclusion

In this work, we developed a simulated environment of sequence labeling tasks
to explore the complexities of learning multiple tasks jointly. Using the data
generated from a task simulator, we were able to investigate the
generalization patterns of MTL models. Further, we resorted to a learning
algorithm to obtain symbolic expressions that can assist practitioners to
understand the factors which contribute to the generalization performance of
MTL models in the sequence labeling domain. Our findings were in accordance
with previous theoretical works, yet we were able to add additional insights
to the unexplored setting of unbalanced datasets. In addition to that, we
explicitly quantified task relatedness rather than making any assumption about
it or trying to learn it implicitly as a parameter of the model.

## Broader Impact

We bring awareness about the use of multitask learning and how to understand
its benefits using a novel methodology. Approaching theoretical questions via
simulation and the learning of formulas from data could ease the way towards
more explainable machine learning processes.

The possible negative impact of this work is the idea that our methodology
could be seen as a replacement instead of a support of theoretically grounded
work. Undoubtedly, rigorous proofs are an essential part of our domain.
Theoretical contributions are necessary since they make the use of machine
learning more reliable. Essentially, the only way one can safely say some
phenomena is true is by presenting a valid mathematical proof for it.

Building a better understanding of machine learning algorithms, and especially
deep learning methods, is an important requirement for applying them in a
responsible manner in real applications. More specifically, we believe that
knowledge sharing techniques such as MTL and Transfer Learning are important
components of many low-resource problems. Improving the use of these
techniques enables the use of machine learning for unsolvable problems due to
a lack of data. These problems include, for instance, the use of natural
language processing for languages in developing countries or general
classification tasks for small businesses.

## Acknowledgments

Gefördert durch die Deutsche Forschungsgemeinschaft (DFG) – Projektnummer
232722074 – SFB 1102 / Funded by the Deutsche Forschungsgemeinschaft (DFG,
German Research Foundation) – Project-ID 232722074 – SFB 1102 and the EU-
funded Horizon 2020 project ROXANNE under grant number 833635.

## References

  * Maurer et al. [2016] Andreas Maurer, Massimiliano Pontil, and Bernardino Romera-Paredes.  The benefit of multitask representation learning.  _Journal of Machine Learning Research_ , 17:2853–2884, Jan 2016.
  * Caruana [1997] Rich Caruana.  Multitask learning.  _Machine Learning - Special issue on inductive transfer_ , 28:41–75, Jul 1997.  doi: 10.1023/A:1007379606734.
  * Collobert and Weston [2008] Ronan Collobert and Jason Weston.  A unified architecture for natural language processing: Deep neural networks with multitask learning.  In _Proceedings of the 25th International Conference on Machine Learning_ , pages 160–167, Jun 2008.  doi: 10.1145/1390156.1390177.
  * McCann et al. [2018] Bryan McCann, Nitish Shirish Keskar, Caiming Xiong, and Richard Socher.  The natural language decathlon: Multitask learning as question answering.  2018\.
  * Crichton et al. [2017] Gamal Crichton, Sampo Pyysalo, Billy Chiu, and Anna Korhonen.  A neural network multi-task learning approach to biomedical named entity recognition.  _BMC Bioinformatics_ , 18, Aug 2017.  doi: 10.1186/s12859-017-1776-8.
  * Søgaard and Goldberg [2016] Anders Søgaard and Yoav Goldberg.  Deep multi-task learning with low level tasks supervised at lower layers.  In _Proceedings of the 54th Annual Meeting of the Association for Computational Linguistics: Volume 2, Short Papers_ , pages 231–235, Aug 2016\.  doi: 10.18653/v1/P16-2038.
  * Bingel and Søgaard [2017] Joachim Bingel and Anders Søgaard.  Identifying beneficial task relations for multi-task learning in deep neural networks.  In _Proceedings of the 15th Conference of the European Chapter of the Association for Computational Linguistics: Volume 2, Short Papers_ , pages 164–169, Apr 2017.  doi: 10.18653/v1/E17-2026.
  * Changpinyo et al. [2018] Soravit Changpinyo, Hexiang Hu, and Fei Sha.  Multi-task learning for sequence tagging: An empirical study.  In _Proceedings of the 27th International Conference on Computational Linguistics_ , pages 2965–2977, Aug 2018.
  * Baxter [2000] Jonathan Baxter.  A model of inductive bias learning.  _Journal of Artificial Intelligence Research_ , 12:149–198, Feb 2000.  doi: 10.1613/jair.731.
  * Maurer [2006] Andreas Maurer.  Bounds for linear multi-task learning.  _Journal of Machine Learning Research_ , 7:117–139, Jun 2006.
  * Lee et al. [2018] Haebeom Lee, Eunho Yang, and Sung Hwang.  Deep asymmetric multi-task feature learning.  In _Proceedings of the 35th International Conference on Machine Learning_ , pages 2956–2964, Jun 2018.
  * Liu et al. [2017] Tongliang Liu, Dacheng Tao, Mingli Song, and Stephen J. Maybank.  Algorithm-dependent generalization bounds for multi-task learning.  _IEEE Transactions on Pattern Analysis and Machine Intelligence_ , 39:227–241, Feb 2017.  doi: 10.1109/TPAMI.2016.2544314.
  * Ciliberto et al. [2017] Carlo Ciliberto, Alessandro Rudi, Lorenzo Rosasco, and Massimiliano Pontil.  Consistent multitask learning with nonlinear output relations.  In _Advances in Neural Information Processing Systems 30_ , pages 1986–1996. 2017.
  * Evgeniou and Pontil [2004] Theodoros Evgeniou and Massimiliano Pontil.  Regularized multi–task learning.  In _Proceedings of the 10th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining_ , pages 109–117, Aug 2004.  doi: 10.1145/1014052.1014067.
  * Parameswaran and Weinberger [2010] Shibin Parameswaran and Kilian Weinberger.  Large margin multi-task metric learning.  In _Proceedings of the 23rd International Conference on Neural Information Processing Systems - Volume 2_ , pages 1867–1875, Dec 2010.
  * Evgeniou et al. [2005] Theodoros Evgeniou, Charles A. Micchelli, and Massimiliano Pontil.  Learning multiple tasks with kernel methods.  _Journal of Machine Learning Research_ , 6:615–637, Apr 2005.
  * Kang et al. [2011] Zhuoliang Kang, Kristen Grauman, and Fei Sha.  Learning with whom to share in multi-task feature learning.  In _Proceedings of the 28th International Conference on International Conference on Machine Learning_ , pages 521–528, Jun 2011.
  * Görnitz et al. [2011] Nico Görnitz, Christian Widmer, Georg Zeller, Andre Kahles, Sören Sonnenburg, and Gunnar Rätsch.  Hierarchical multitask structured output learning for large-scale sequence segmentation.  In _Proceedings of the 25th International Conference on Neural Information Processing Systems_ , pages 2690–2698, Dec 2011.
  * Bonilla et al. [2007] Edwin V. Bonilla, Kian Ming A. Chai, and Christopher K. I. Williams.  Multi-task gaussian process prediction.  In _Proceedings of the 20th International Conference on Neural Information Processing Systems_ , pages 153–160, Dec 2007.
  * Zhang and Yeung [2014] Yu Zhang and Dit-Yan Yeung.  A regularization approach to learning task relationships in multitask learning.  _ACM TKDD_ , 8:12:1–12:31, Jun 2014.  doi: 10.1145/2538028.
  * Sanh et al. [2019] Victor Sanh, Thomas Wolf, and Sebastian Ruder.  A hierarchical multi-task approach for learning embeddings from semantic tasks.  _Proceedings of the AAAI Conference on Artificial Intelligence_ , 33:6949–6956, Jul 2019.  doi: 10.1609/aaai.v33i01.33016949.
  * Ruder et al. [2019] Sebastian Ruder, Joachim Bingel, Isabelle Augenstein, and Anders Søgaard.  Latent multi-task architecture learning.  _Proceedings of the AAAI Conference on Artificial Intelligence_ , 33:4822–4829, Jul 2019.  doi: 10.1609/aaai.v33i01.33014822.
  * Clark et al. [2019] Kevin Clark, Minh-Thang Luong, Urvashi Khandelwal, Christopher D. Manning, and Quoc V. Le.  Bam! born-again multi-task networks for natural language understanding.  In _Proceedings of the 57th Annual Meeting of the Association for Computational Linguistics_ , pages 5931–5937, Jul 2019.  doi: 10.18653/v1/p19-1595.
  * Orzechowski et al. [2018] Patryk Orzechowski, William La Cava, and Jason H. Moore.  Where are we now? a large benchmark study of recent symbolic regression methods.  In _Proceedings of the Genetic and Evolutionary Computation Conference_ , pages 1183–1190, Jul 2018.  doi: 10.1145/3205455.3205539.
  * Koza [1992] John R. Koza.  _Genetic Programming: On the Programming of Computers by Means of Natural Selection_.  MIT Press, Cambridge, MA, USA, 1992.  ISBN 0-262-11170-5.
  * Rabiner [1989] Lawrence R. Rabiner.  A tutorial on hidden markov models and selected applications in speech recognition.  pages 257–286, February 1989.  doi: 10.1109/5.18626.
  * Lafferty et al. [2001] John D. Lafferty, Andrew McCallum, and Fernando C. N. Pereira.  Conditional random fields: Probabilistic models for segmenting and labeling sequence data.  In _Proceedings of the 18th International Conference on Machine Learning_ , pages 282–289, Jun 2001.
  * Mccallum et al. [2000] Andrew Mccallum, Dayne Freitag, and Fernando Pereira.  Maximum entropy markov models for information extraction and segmentation.  In _Proceedings of the 17th International Conference on Machine Learning_ , pages 591–598, 2000.
  * Collobert et al. [2011] Ronan Collobert, Jason Weston, Leon Bottou, Michael Karlen, Koray Kavukcuoglu, and Pavel Kuksa.  Natural language processing (almost) from scratch.  2011\.
  * Vaswani et al. [2017] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Lukasz Kaiser, and Illia Polosukhin.  Attention is all you need.  pages 6000–6010, Dec 2017.
  * Devlin et al. [2019] Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova.  BERT: Pre-training of deep bidirectional transformers for language understanding.  In _Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long and Short Papers)_ , pages 4171–4186, Jun 2019.
  * Vinh et al. [2010] Nguyen Xuan Vinh, Julien Epps, and James Bailey.  Information theoretic measures for clusterings comparison: Variants, properties, normalization and correction for chance.  _Journal of Machine Learning Research_ , 11:2837–2854, Oct 2010.
  * Paszke et al. [2019] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zach DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.  Pytorch: An imperative style, high-performance deep learning library.  2019\.
  * Sang and Meulder [2003] Erik F. Tjong Kim Sang and Fien De Meulder.  Introduction to the conll-2003 shared task: Language-independent named entity recognition.  In _Proceedings of the Seventh Conference on Natural Language Learning at HLT-NAACL 2003_ , pages 142–147, Jun 2003.
  * Grave et al. [2020] Edouard Grave, Piotr Bojanowski, Prakhar Gupta, Armand Joulin, and Tomas Mikolov.  Learning word vectors for 157 languages, 2020.  URL https://fasttext.cc/docs/en/pretrained-vectors.html.
  * Sighsmile [2020] Sighsmile.  Conlleval-python, 2020.  URL https://github.com/sighsmile/conlleval.
  * Stephens [2019] Trevor Stephens.  Gp learn, 2019.  URL https://gplearn.readthedocs.io/en/stable/index.html.


else if $0.63\leq v_{k}\leq 0.94$ AND $0.15\leq\phi_{k}\leq 0.45$. then

Use Region $\mathcal{R}_{1}$ approximate maps (16)-(17):

else if $v_{k}>0.63-0.53\phi_{k}$ AND $v_{k}>0.55$ AND
$(v_{k},\phi_{k})\notin\mathcal{R}_{1}$, then

Use $\mathcal{R}_{2}$ approximate map (• ‣ iv).):

else if $v_{k}>0.63-0.53\phi_{k}$ AND $1.1<\phi_{k}<2.5$ and $v_{k}<0.55$,
then

Use $\mathcal{R}_{4}$ approximate map (A.5):

else if $2.5<\phi_{k}<\pi$ AND $v_{k}<0.55$, then

Use $\mathcal{R}_{5}$ approximate map (A.7):

else if $v_{k}<0.63-0.53\phi_{k}$, then

Use $\mathcal{R}_{3}$ approximate map (A.6):

end if