Proceedings of the CTD 2023 PROC-CTD2023-43

An Application of HEP Track Seeding to Astrophysical Data

Mi̇ne Gökçen

Department of Physics Engineering
Istanbul Technical University, Turkey

Maurice Garcia-Sciveres, Xiangyang Ju

Lawrence Berkeley National Lab.
Cyclotron Rd., Berkeley, CA 94720, USA

> ABSTRACT
>
>
>
>
> We apply methods of particle track reconstruction in High Energy Physics
> (HEP) to the search for distinct stellar populations in the Milky Way, using
> the Gaia EDR3 data set. This was motivated by analogies between the 3D space
> points in HEP detectors and the positions of stars (which are also points in
> a coordinate space) and the way collections of space points correspond to
> particle trajectories in the HEP, while collections of stars from distinct
> populations (such as stellar streams) can resemble tracks. Track
> reconstruction consists of multiple steps, the first one being seeding. In
> this note we describe our implementation and results of the seeding step to
> the search for distinct stellar populations, and we indicate how the next
> step will proceed. Our seeding method uses machine learning tools from the
> FAISS library, such as the kNN nearest neighbor search.

> PRESENTED AT
>
>
>
>
> Connecting the Dots Workshop (CTD 2023)
>
> October 10-13, 2023

##  1 Introduction

The Milky Way (MW) was formed by merging of many smaller objects and continues
to accrete matter from satellite systems as it ages. Studying the progenitors
of our galaxy is key to understanding its formation and evolution. Recent
merger events are particularly useful in revealing insights about the the
gravitational structure our galaxy, including its dark matter halo. As stars
from satellite systems tidally get disrupted by our galaxy’s gravitational
potential, stellar streams form around the galactic disk, inside its halo.
Alongside with star clusters, they have been extensively studied by
astronomers as relics of recent merger events.

There have been numerous research efforts on identifying stellar streams and
star clusters using statistical and computational techniques [1, 2]. Although
traditional methods of constraining stars by their color, age and kinematics
are still in use [3], simulations [4, 5] are widely utilized to determine star
candidates for discoveries of these astrophysical structures. More novel
methods of anomaly detection using machine learning (ML) [6] and deep learning
[7] are also being actively investigated with the common purpose of automating
searches for distinct stellar populations.

Orbits of stellar streams are reminiscent of particle tracks in high energy
physics (HEP) detectors. In this work, we explore applying (HEP) track
reconstruction methods to astrophysical data, inspired by that resemblance,
and more generally to identify stars from distinct populations in an automated
way.

In HEP experiments, the trajectory a particle follows inside the detector
needs to be reconstructed—a process called particle tracking—in order to study
the particle interactions and events that have occurred inside. In essence,
tracking consists of three major steps: forming seeds, extrapolating the seeds
to form track candidates, and fitting tracks using these seeds. Particle
detectors have multiple coaxial layers that record the location a particle
passes by. These locations that are marked when a particle hits the detector
layer are called space-points (SP). Seeds are clustered sets of SPs that
potentially belong to the same particle track, which are then combined and
extended to arrive at final track candidates.

We aim to establish an analogous method to track stars—tracking meaning
identifying stars with a common origin in this context. Stars belonging to a
certain stellar population could be characterized by their intrinsic
properties, such as their ages and chemical compositions, as well as their
kinematics. Here, to construct seeds of stars with shared origins, we will
cluster them based on their kinematics in the sky. In a later step we would
extend and combine the seeds to arrive at final distinct populations.

This paper describes our ML-based seeding method as the first stage of
tracking stars. In section 2, we introduce the astrophysical data set used for
our analysis. Section 3 goes over the seeding method. In section 4, we share
our current seed results and give examples of known objects that collections
of seeds clearly correspond to. Finally, in section 5, we conclude with a
discussion on the performance of our seeding method and the steps for
implementing the second tracking stage to automatically extend and combine
seeds.

##  2 Gaia Dataset

We use the EDR3 catalog [8] of the Gaia mission [9]. Gaia EDR3 has
observations of roughly two billion stars and their properties. In particular,
we use the following five positional and kinematic parameters: galactic
longitude and latitude, proper motion (PM) in right ascension and declination,
and parallax. It is important to note that parallax is a parameter with high
uncertainties [10] and could yield unreliable results if not treated
carefully. Therefore, we mainly use the first four parameters in our seeding
method. The distribution of catalog stars in the 2D position and kinematics
space is illustrated in fig. 1.

(a)

(b)

Figure 1: All Gaia EDR3 stars with parallax (total of 1,467,744,818 stars)
binned in (a) position- and (b) proper motion-space, with respective bin sizes
of 0.25 deg and 0.1 mas/yr in each axis. The two small red crosses at (l, b)
$\sim$ (350, 20) degs in (a) mark the patches of sky that will be shown in
fig. 3 and discussed in section 3.

##  3 Seeding Methods

###  3.1 Subsetting the Data

The seeding will look for differences in PM as a signature of distinct
populations. The basic assumption of our seeding method is that most stars
have a PM given by the galaxy average, which we call background stars, while a
fraction of stars that originated form recent mergers or accretion (signal
stars) will move differently. Additionally and critically, we assume each
signal population is localized in one part of the sky. Therefore, if the whole
sky PM distribution were considered at once (fig. 1(b)), the correlation
between PM and position would be hidden.

On the other hand, in a small patch of sky (from here on referred to as a
chunk), signal stars have a chance to stand out in PM space. Thus, we carry
out PM space clustering of stars chunk-by-chunk, rather than globally.

Chunking, or subsetting stars in position-space is not as trivial as using a
fixed grid. The number of stars in a chunk can have a considerable effect on
the clustering results, depending on clustering parameters. We pick the
minimum and maximum number of stars to be contained in a chunk based on a
sample size threshold that the clustering algorithm we use is optimized for.
The chunk sizes vary within these limits and the distribution of number of
stars in chunks is shown in fig. 2(b) along with centers of the chunks in fig.
2(a). Note that we make no attempt to match the chunk size to the expected
physical extent of some distinct population or another, because we are only
searching for seeds at this stage, which can be thought of as hints of where
distinct populations may reside in both PM and position. Eventually, when
populations are fully reconstructed, there will be no memory of the seeds and
chunks, just as there is no memory of the seeds in the HEP final reconstructed
tracks.

Because the density of stars in the galactic disk region ($\mid$b$\mid\leq$ 20
degrees) is much higher than that of the halo region ($\mid$b$\mid>$ 20
degrees), the chunk dimensions in these two regions differ significantly.
Analysis of seed results also requires subtraction of extra-galactic
backgrounds in the halo region; so disk and halo seeds must be analyzed
separately. Due to space limitations, we focus the rest of this paper on
galactic disk region.

(a)

(b)

Figure 2: Gaia stars subsetted in the position space. (a) shows position-space
centers of the chunks and (b) is the binned number of stars in chunks. Orange
represents galactic disk region ($\sim$ 30,500 chunks) and blue represents the
galactic halo region ($\sim$ 2,500 chunks).

###  3.2 Clustering Stars

Clustering is done in each chunk individually based on the two dimensional PM
of stars. Fig. 3 shows the PM distributions in two selected neighboring
chunks, which have been singled out in fig. 1(a). A very clear signal can be
seen in fig. 3(b). Note that there is no signal in the chunk just next door
(a), which illustrates the point that signals can have high PM-position
correlation. The goal of clustering is to find such signal peaks in an
automated way; not only the obvious ones ones like in fig. 3(b), but also ones
that may be hard to spot visually.

We cluster stars using the open source FAISS library [11], which can conduct
the clustering on GPU. More specifically, we use the FAISS k-means algorithm
to obtain suggested centroids for a given number of clusters, $k$, and the
k-NN algorithm to uniquely assign each star in the chunk to one of the $k$
centroids. Fig.3(c) and (d) show the results of clustering with $k=250$ on (a)
and (b), respectively.

If all stars in a chunk have a background distribution, which has a broad
shape as in fig. 3(a), then the number of stars in every centroid will be
approximately $N/k$ plus or minus statistical fluctuations, where $N$ is the
number of stars in the chunk. If a signal is present, which is more peaked
than the background, the centroids containing the signal will have a number of
stars determined by the signal density, and not by $N/k$. This makes it clear
that $1/k$ plays the role of a signal detection threshold.

(a)

(b)

(c)

(d)

Figure 3: The stars from two exemplary neighboring chunks marked by red
crosses in fig. 1(a) are shown in PM space in (a) and (b). Each chunk spans
(dl, db) $\sim$ (4, 1) degrees in position space. Chunk (a) contains only
background stars, whereas (b) has a secondary peak of signal stars. Panels (c)
and (d) show the same stars with the 250 centroids found by k-means clustering
with k=250 overplotted.

###  3.3  Testing for Seeds

After we perform the clustering over each chunk, we automatically test if a
chunk contains a seed by checking for outliers in the binned distribution of
number of stars per centroid, as shown in fig. 4. We look for high-side gaps
between the bins, and cut on the counts within the outlier bins. We
empirically settled on using 20 bins, having checked that varying between 15
and 25 bins made little difference. We investigated more formal statistical
methods and found that this simple test was the most robust against
differences in the shape of the main distribution. As in HEP track
reconstruction, seeding does not have to be perfect, as seeds are intended to
be redundant hints to be followed up with extension and merging later. In fig.
4, example (a) fails the seed test (no signal found) while example (b) passes.

(a)  (b)

Figure 4: Distribution of number of stars per centroid in k-means clusters for
k=250, with the inserts showing the cluster centroids in PM space (crosses).
Plot (a) shows the chunk of fig. 3(a), which failed the seed test (no seed
detected in this chunk), while (b) is that of fig. 3(b), which passed seed
test. The seed stars can be clearly seen in the insert and as well separated
bins on the right side of the histogram.

The simple gap selection method works well for small signals (where most of
the stars in the chunk are not signal), but can fail when the number of signal
stars in a chunk approaches or exceeds the number of background stars, in
which case the histogram in number of stars per centroid can have two peaks
with no gap in between. We have not added a test to detect those extreme
seeds, because such a large signal will spread over many chunks and will
already have many redundant seeds. As in HEP track reconstruction, one needs a
sufficient number of seeds to cover the objects to be reconstructed, and
adding more redundant seeds will not change the outcome.

Once a chunk passes the seed test, we save the stars within the top three
k-means clusters with highest number of stars as seed stars. One reason for
saving multiple cluster stars is the possibility of having more than one
signal peak in a chunk. Another is that there are often more than one centroid
per signal peak.

Fig. 5 shows how seed stars identified by our method in PM space also cluster
in the position space, as expected. A second example in which the signal peak
overlaps the background distribution in PM space is illustrated in fig. 6.

(a)

(b)

(c)

Figure 5: Example of a seeded chunk. All chunk stars were plotted in fig. 3
(reproduced here in the insert of (a)). The seed stars are shown in PM space
in (a) and position space in (c), while (b) shows all chunk stars in position
space.

(a)

(b)

(c)

Figure 6: Another example of a seeded chunk. The insert in (a) shows all chunk
stars in PM space, while (a) shows only the seed stars. (b) Shows all chunk
stars in position space while (c) shows only the seed stars.

##  4 Seed Results

This section shows the seeding results for the galactic disk region. Fig. 7
and fig. 8 show the results in position and PM space respectively, with
statistics given in table 1. The halo region still requires the second
chunking step described below and the removal of extra-galactic backgrounds
that will cluster near zero proper motion (i.e. fixed in the sky), because
they are very distant objects.

| Total Chunks | Seeded Chunks | Percentage of Seeded Chunks | Seeded Stars
---|---|---|---|---
Galactic disk | 30,264 | 972 | 3.21 | 1,013,408
Galactic halo | 2,375 | 618 | 26.0 | 354,557
Table 1: Overview of galactic disk (top row) and preliminary galactic halo
(bottom row) seed results. Galactic halo results are only preliminary: they
are missing the addition of seeds from chunks with offset centers. Seeds
corresponding to extra-galactic backgrounds with zero proper motion (i.e.
fixed in the sky because they are very distant) are removed.

(a)

(b)

Figure 7: Seed results in position space, where gray points are all chunk
centers and black represents the seeded ones. (a) is that of galactic disk
region and (b) is that of galactic halo region. In (a), all seed results are
shown, including the seeds from chunks with offset centers, with the
corresponding original chunk centers marked for visualization purposes.

(a)

(b)

Figure 8: Binned galactic disk seed stars in PM space: (a) shows all seed
stars and (b) zoomed in central region showing only the seed stars from the
most populated k-means clusters in seeded chunks (as opposed to the top
three).

For the galactic disk region, we performed chunking twice. Stellar populations
residing at the corners of chunks may not have enough signal stars in any
chunk and thus, tend to lack seeds. Shifting chunk centers diagonally such
that the corners of original chunks become the centers for new ones accounts
for this vulnerability. For each set of chunks, we ran our seeding algorithm
for $k=250$ and $k=600$, which amounts to four runs of seeding in total. At
the end, seed stars from all runs are combined and get checked for redundancy.

Extending (meaning adding non-seed stars) and combining seed stars into
candidate objects is the next step in the HEP-inspired reconstruction chain
and the automation of this step is work in progress. Such object candidates
can be compared to astrophysical object catalogs. We give two examples of
known objects that match seed star collections, where the collection was done
manually.

The first example is Sagittarius dwarf galaxy [12], which stands out due to
the many neighboring seeded chunks circled in fig. 9(a). We compare the PM of
seed stars in these chunks, plotted in fig. 9(c), along with their galactic
coordinates and parallax, to the SIMBAD data base [13] objects and find they
are well matched.

The second example is the NGC 6397 globular cluster [14], which matches the
distinguishable PM-space distribution of its associated seed stars circled in
fig. 10(a). The position-space center of these seed stars, plotted in fig.
10(c), along with their PM and parallax center, indicate that these are
potentially NGC 6397 stars. The gap in the position space region of NGC 6397
seed stars is due to the extreme signal strength in that chunk, as explained
in sec. 3.3.

(a)

(b)

(c)

Figure 9: Seed stars and seeded chunks matching the Sagittarius dwarf galaxy
and stream. In (a), associated seeded chunk centers are circled with red. In
(b) position space and in (c) PM space distributions are shown for the seed
stars within the selection in (a).

(a)

(b)

(c)

Figure 10: Seed stars matching the known NGC 6397 globular cluster. In (a),
associated seed stars are circled with red among all the galactic disk seed
stars. In (b), seed stars are shown in position space. The central chunk in
(b), whose PM space distribution is shown in (c), failed the seed test due to
an extremely high signal content, as explained in the text.

##  5 Conclusion and Future Work

We established an ML-based, automated seeding method to identify distinct
stellar populations. Seeding is the first step of a HEP track reconstruction
inspired object search. We tested our method on Gaia EDR3 stars and presented
our seed results for the galactic disk region.

The next stage after seeding, analogous to HEP track reconstruction, is under
development. For this purpose, we have begun exploring methods to combine and
extend seeds as well as to refine them. Re-clustering seed stars in position
and PM space stands out as a promising start. However, we face challenges due
to the varying structures and sizes of stellar populations in these parameter
spaces.

This work suggests that it is possible to develop an automated reconstruction
flow that, starting from survey data, can potentially identify in an instant
all known stellar populations that have been manually cataloged over centuries
of varied observations. If that is the case, previously unknown populations
could also be identified. Significant further work is needed to reach this
point.

##  6 Acknowledgements

We thank Jessica Lu for helpful discussion. This work was supported in part by
U.S. Department of Energy, Office of Science under contract DE-AC02-05CH11231.

This work has made use of data from the European Space Agency (ESA) mission
Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing
and Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the Gaia Multilateral Agreement.

## References

  * [1] Yuan Zhen et al, “StarGO: A New Method to Identify the Galactic Origins of Halo Stars,” ApJ 863 26 (2018) [arXiv:1806.06341 [astro-ph.GA]]
  * [2] Malhan K. and Ibata R., “STREAMFINDER – I. A new algorithm for detecting stellar streams,” MNRAS 477 4063 (2018)
  * [3] Vedant Chandra et al, “Discovery of the Magellanic Stellar Stream Out to 100 Kiloparsecs,” ApJ 956 110 (2023) [arXiv:2306.15719 [astro-ph.GA]].
  * [4] Nora Shipp et al, “Streams on FIRE: Populations of Detectable Stellar Streams in the Milky Way and FIRE,” ApJ 949 44 (2023) [arXiv:2208.02255 [astro-ph.GA]].
  * [5] Ana Bonaca et al, “Milky Way Mass and Potential Recovery Using Tidal Streams in Realistic Halo,” ApJ 795, 94 (2014) [arXiv:1406.6063 [astro-ph.GA]].
  * [6] Mariel Pettee et al, “Weakly-Supervised Anomaly Detection in the Milky Way,” [arXiv:2305.03761 [astro-ph.GA]].
  * [7] David Shih et al, “Via Machinae 2.0: Full-Sky, Model-Agnostic Search for Stellar Streams in Gaia DR2,” [arXiv:2303.01529 [astro-ph.GA]].
  * [8] Gaia Collaboration, “Gaia Early Data Release 3: Summary of the contents and survey properties,” A&A 650 C3 (2021) [arXiv:2012.01533 [astro-ph.GA]].
  * [9] Gaia Collaboration, “The Gaia Mission,” A&A 595 A1 (2016) [arXiv:1609.04153 [astro-ph.IM]].
  * [10] L. Lindegren et al, “Gaia Early Data Release 3: Parallax bias versus magnitude, colour, and position,” A&A 649 A4 (2021) [arXiv:2012.01742 [astro-ph.IM]].
  * [11] Jeff Johnson et al, “Billion-scale similarity search with GPUs,” IEEE Transactions on Big Data 7 535 (2019) [arXiv:1702.08734 [cs.CV]].
  * [12] C. Alard, “Evidence for the Sagittarius Dwarf Galaxy at Low Galactic Latitudes,” ApJ 458 L17 (1996).
  * [13] Marc Wenger et al, “The SIMBAD astronomical database,” Astron. Astrophys. Suppl. Ser. 143 9 (2000) [arXiv:astro-ph/0002110].
  * [14] H. Baumgardt et al, “Mean proper motions, space orbits, and velocity dispersion profiles of Galactic globular clusters derived from Gaia DR2 data,” MNRAS 482 5138 (2019) [arXiv:1811.01507 [astro-ph.GA]].