Combining Machine Learning with Physics: A Framework for Tracking and Sorting Multiple Dark Solitons

In ultracold atom experiments, data often comes in the form of images which suffer information loss inherent in the techniques used to prepare and measure the system. This is particularly problematic when the processes of interest are complicated, such as interactions among excitations in Bose-Einstein condensates (BECs). In this paper, we describe a framework combining machine learning (ML) models with physics-based traditional analyses to identify and track multiple solitonic excitations in images of BECs. We use an ML-based object detector to locate the solitonic excitations and develop a physics-informed classifier to sort solitonic excitations into physically motivated sub-categories. Lastly, we introduce a quality metric quantifying the likelihood that a specific feature is a kink soliton. Our trained implementation of this framework – SolDet – is publicly available as an open-source python package. SolDet is broadly applicable to feature identification in cold atom images when trained on a suitable user-provided dataset.



There are no comments yet.


page 1

page 2

page 5

page 7

page 8

page 9

page 10

page 11


Enhancing predictive skills in physically-consistent way: Physics Informed Machine Learning for Hydrological Processes

Current modeling approaches for hydrological modeling often rely on eith...

Jet-Images -- Deep Learning Edition

Building on the notion of a particle physics detector as a camera and th...

Using Machine Learning for Model Physics: an Overview

In the overview, a generic mathematical object (mapping) is introduced, ...

Vamsa: Tracking Provenance in Data Science Scripts

Machine learning (ML) which was initially adopted for search ranking and...

Automated Dissipation Control for Turbulence Simulation with Shell Models

The application of machine learning (ML) techniques, especially neural n...

QRMine: A python package for triangulation in Grounded Theory

Grounded theory (GT) is a qualitative research method for building theor...

Highway Traffic State Estimation Using Physics Regularized Gaussian Process: Discretized Formulation

Despite the success of classical traffic flow (e.g., second-order macros...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

I Introduction

Machine learning (ML) techniques promise improved data analysis and enhanced performance for today’s quantum devices and technologies. Recent developments include: quantum noise characterization harper_efficient_2020 ; quantum state detection carrasquilla_machine_2017 ; zhang_quantum_2017 ; torlai_neural-network_2018 ; venderley_machine_2018 ; rem_identifying_2019 ; miles_correlator_2020 ; cha_attention-based_2020 ; metz_deep_2020 ; venderley_harnessing_2021 ; maskara_learning_2021 ; guo_machine-learning_2021 ; parameter space exploration and optimization wigley_fast_2016 ; tranter_multiparameter_2018 ; kalantre_machine_2019 ; zwolak_autotuning_2020 ; barker_applying_2020 ; saito_creation_2020 ; and quantum control baum_experimental_2021 ; ai_experimentally_2021 . Together, these results show that ML techniques can extract information from ambiguous data, efficiently search large parameter spaces, and optimally control quantum systems.

Traditional statistical analysis using physics-based models, such as least-square fitting and hypotheses testing, have been go-to techniques for data analysis since the 1800’s Legendre1805 and remain widely applied in quantum physics research today ketterle1999making ; fritsch_creating_2020 ; purdy2017quantum

. The outcome of physics-model-based fits are intuitive, physically meaningful, and can help identify patterns present in the data; even fits based on more heuristic functions can have coefficients that are derived in obvious ways from the data. By contrast, ML methods work as “black boxes,” making their operation difficult to interpret. Conventional statistical methods use fixed algorithms in conjunction with preconceived models for data reduction. Overfitting occurs when the number of fit parameters is comparable or larger than the number of independent data points. In this context, the process of training an ML tool essentially co-designs the fitting algorithm and the data model, as encoded by a large number of internal parameters. Training ML models is itself a fitting process that can be susceptible to overfitting, for example when the training dataset has too little variability or the ML model has too many internal parameters. ML involves a class of data-driven techniques that do not rely on preexisting models, but also add additional opportunities for overfitting that can make them less reliable on new data than conventional techniques.

Here, we describe the hybrid two-module feature identification framework, shown in Fig. 1, that combines the flexibility of ML techniques with the intuition and robustness of conventional fitting methods. Furthermore the separate outputs of these two very different modules allows us to assess data quality by cross-validation. Hybrid approaches have been employed in other settings, for example for landslide prediction huang2020comparisons , medical image processing ghosh2012new , and cyber attack detection sakhnini2019smart .

The framework begins with a labeled dataset that is used to train the ML module and initialize the physics-based module. Before trusting either module, we independently validate each module on an unused subset of the labeled data that was not used for training. Model-redesign may be needed until satisfactory performance of each module is reached. We then combine both modules into an integrated system able to analyze new data.

Figure 1: Overview of the framework. The colored arrows link the preparation (Secs. II.1, II.2, and II.3), validation (Secs. III.1. III.2), and application (Sec. III.4) phases of the framework. The red path represents the preparation and implementation of the physics-based-approximation module of the framework. The blue path represents the ML modules.

We demonstrate the performance of our framework using data from atomic Bose-Einstein condensates (BECs), quintessential quantum systems. Quantum research with BECs, and cold-atom quantum gases more broadly, is multifaceted with examples ranging from realizing collective many-body physics Greiner2002 to creating today’s most accurate atomic clocks Bothwell2019

. In the vast majority of these experiments, data is acquired in the form of noisy images that typically have undergone evolution—such as a time-of-flight—before measurement, which obfuscates the computation of the quantities of interest. Cold quantum gases therefore make an ideal test-bed for our methodology that combines physically motivated, but heuristic, fitting functions with established computer vision techniques.

We focus on the specific problem of locating dark solitons—spatially compact excitations that manifest as reductions in the atomic density—as they move in BECs burger_dark_1999 ; Denschlag2000 ; fritsch_creating_2020 . This allows us to leverage our established soliton dataset guo_machine-learning_2021 ; solitons-data to train and validate our framework; representative elements of the dataset are shown in Fig. 2. These data consist of elliptical atom clouds (first row) where solitons appear as vertically aligned density depletions (second row). Not all vertically aligned density depletions are created equal (bottom row): while deep depletions mark the location of kink-solitons, shallow and asymmetric depletions can result from solitonic vortices Mateo2015 , and chains of stripes can result from highly excited phonon modes. Our framework is the first tool that can distinguish between these cases as well as identify the location of all the solitons and solitonic vortices present in each image.

Our ML module leverages and extends established computer vision techniques. Computer vision is a broad field with applications ranging from image classification to semantic segmentation and object detection voulodimos2018deep . Object detection refers to the capability of software systems to locate and identify objects in an image. Convolutional neutral networks (CNNs) underlie solutions to all of these tasks, and unsurprisingly were employed in our previous work classifying soliton image data into three categories: no soliton, one soliton, and other guo_machine-learning_2021 . Our ML module goes beyond simple classification and uses a CNN based object detector (OD) to provide the location of all candidate excitations in a given image.

By contrast our physics-based module employs a least-squares fit of an inverted and skewed Mexican-hat function to 1D background-subtracted projections of soliton-candidates (shown in Fig. 

2). We initialized this module using our previously labeled single soliton data and employ a Yeo-Johnson transformation yeo_new_2000

to produce a multivariate normal distribution yielding the likelihood that an unknown feature is a soliton.

This approach yielded three immediate benefits. First, a careful analysis of the coefficients from the physics based module identified previously overlooked correlations that allow us to distinguish between some solitonic excitations (kink and solitonic vortex burger_dark_1999 ; Denschlag2000 ; donadello_observation_2014 ; Mateo2015 ). Second, combining the results of the ML and fitting modules allowed us to automatically create a larger, more reliable dataset that includes fine-grained information such as the soliton position and type of excitation. This dataset is described in Ref. Fritsch21-DSD and published in the NIST data repository solitons-data . Third, our hybrid framework was prepared solely from a training dataset whose images contain either zero or one solitonic excitation, however, it is performant on complex data containing multiple excitations.

The remainder of this manuscript is structured as follows: Sec. II introduces both modules and describes their training and initialization. Sec. III describes the validation of both modules and their performance on new data that include multiple solitonic excitations. In Sec. III.5, we describe an open-source Python reference implementation of our framework: SolDet SolDet . Lastly, in Sec. IV we conclude and discuss possible future directions.

Figure 2: Representative data. The top panels plot pre-processed images from our dataset and the bottom panels plot profiles: (green) profile of full image, (black) TF fits, blue (blue) density fluctuations. The red lines mark the location of the deepest depletion in the density fluctuations, while the orange lines mark the solitons locations found from our OD. (a) An element of the no-excitation class. (b) Three elements of the single excitation class: (i) a single kink solition, (ii) an off-center kink soliton, and (iii) a solitonic vortex. (c) Two representative elements of the other excitations class.

Ii Data and modules

In addition to recent success of ML methods metz_deep_2020 ; guo_machine-learning_2021 ; leykam2021dark , solitonic excitations have also be located and characterized using traditional fitting techniques. For example, Ref. fritsch_creating_2020 began with the background-removed atom density profiles (blue curves in Fig. 2) described in Sect. II.1, then identified the deepest depletion (orange dashed line), and fit to a Gaussian function (a physically motivated, but heuristic choice) centered near the deepest depletion. This yielded physical information including soliton width, depth, and position. Unfortunately, this simple approach is failure prone, for example in Fig. 2(b-ii) the deepest depletion is far from the actual soliton. Moreover, it detects only single solitonic features, making human intervention necessary when many excitations are present. Rather than finding the deepest minimum, our framework first uses an OD (described in Sect. II.2

) to provide an initial estimate of the positions all solitonic excitations, and then uses a skewed Mexican hat fit function (Sect. 

II.3) that accurately describes their density profiles. The resulting fit coefficients serve two purposes: qualitative likelihood assessment and fine-grained categorization.

ii.1 Data

We trained and initialized our framework using a dataset consisting of about manually labeled experimental images of BECs with and without solitonic excitations solitons-data . The experimental setup and pre-processing techniques are described in fritsch_creating_2020 .

Figure 2 shows six selected sample images from the labeled dataset. The dataset includes labels for three classes: “no excitation”, images that do not contain any excitations; “single excitation” images containing one solitonic excitation; and “other excitations,” images not in the preceding classes (including those with multiple solitonic excitations, high degrees of noise, and those human labelers could not agree upon). In addition, the dataset includes the horizontal position of excitations in the single excitation class.

Figure 2(a) displays an image from the no excitation class, which lacks the pronounced stripes present in the remaining examples. In (b), we show three elements of the single excitation class, each containing a single dark vertical fringe: (b-i) an iconic kink soliton; (b-ii) an off-center single kink soliton; and (b-iii) a solitonic vortex (see Sect. II.3.2). In (c), we show two elements of the other excitations class containing more than one dark vertical fringe.

Horizontal 1D profiles (bottom row of Fig. 2) also have features associated with our vertically aligned solitonic excitations and are amenable to least squares fitting. We obtain these profiles by first summing the pixel values vertically to compress 2D images to 1D; this sum can be over all (green curves) or part (see Sec. II.3.2) of the vertical extent of the image. We then fit a 1D Thomas-Fermi (TF) model


to each summed 1D profile, where is the horizontal pixel index, and , , , and are fitting parameters representing peak density, center position, TF radius, and an overall offset, respectively. This fit (black curves) serves as an overall background that we subtract from the 1D profiles, leaving behind the 1D density fluctuations (blue curves). The orange dashed lines represent the location of deepest depletion in the 1D fluctuations.

ii.2 ML Modules

Our previous dark soliton detection and positioning system used a CNN classifier guo_machine-learning_2021 that returned one of the three classes: no excitation, single excitation, or other excitations. However,this detector did not locate the excitations. Rather, the positioning step located the deepest depletion and fit to a Gaussian, as described above. This algorithm has two limitations: (1) The soliton may not be the deepest depletion [as in Fig. 2(b-ii)]; and (2) multiple solitons cannot be located [as in Fig. 2(c)]. Here we retain the CNN classifier to globally organize the data, but inspired by a highly successful recent result using an OD to locate vortices in numerically simulated 2D BECs metz_deep_2020 , we employ an OD to locate solitonic excitations in experimental images of highly elongated BECs.

The OD complements the CNN classifier in two ways: (1) it identifies soliton positions rather than classifying; and (2) even thought it is trained with single-soliton data, it can locate multiple excitations in the same image. We employ a neural network based OD with six convolution layers and four max-pooling layers but no fully connected layers (see App. 

A for more detail). The OD has an order of magnitude fewer trainable parameters than our previous CNN ( versus parameters), accelerating the training process and making it lightweight to deploy. Because the OD simply requires a dataset with many representative instances of the object to be detected, it requires far less training data than the CNN classifier (which by design required substantial data from all the considered classes).

In our data the solitonic excitations are roughly 4 pixels in width. Since our images are pixels wide we designed our OD to aggregate the image into spatial cells, each with outputs in the range ; the OD therefore returns a array . For our dataset this aggregation guarantees that each output cell can describe the state of at most one soliton.

is a probability estimate that cell

contains a soliton, and is the fractional position of the soliton center within that cell, where or correspond to the left or right edge of the cell, respectively. The OD considers any cell with as containing an excitation, and then obtains its position from .

When comparing to the training data set with labels denoted by , we use the cost function metz_deep_2020


for each training image, where the label identifying the presence of an excitation in a cell is fully confident, i.e., either 0 or 1. The coefficients , are hyper-parameters controlling the relative importance of each term. The terms increase the cost function when the the OD misidentifies solitons, while the quadratic term contributes when a soliton is mislocated within a cell. Our training set uses images with at most one soliton, so cells with are much less frequent than those with , as a result we expect that to give similar overall weight to the three terms in Eq. 2. We train the OD by minimizing the cost function summed over all training images, updating the predicted OD values in each iteration. Because the cell size is comparable to the soliton size, a single soliton can span two cells. To prevent this double counting, we merge detections occurring in adjacent cells and take the position to be their average.

We deem the OD’s detection successful if our training data contains a labeled soliton close to the detected one (within pixels in our implementation). The two failure modes are failing to detect a solitonic excitation and reporting an excitation that is not present.

ii.3 Physics-based modules

In this subsection, we introduce our physics-based module that uses constrained least-squares fitting to estimate soliton parameters, and following a Yeo-Johnson transformation yeo_new_2000 , produces a quality estimate giving the likelihood of a given feature being solitonic.

We fit the Ricker wavelet Ricker43 , i.e., a “Mexican hat” function


to the 1D density fluctuations described Sect. II.1, where is evaluated with . The function takes six parameters: normalized logarithmic amplitude , center position , width , logarithmic symmetrical shoulder height , asymmetrical shoulder height , and an offset . When and are zero this function is a simple Gaussian, making non-zero adds symmetric shoulders to the distribution and introduces an asymmetry. Our solitonic features are well described by this function as the second term in Eqn. 3 is negative, since our excitations manifest as density depletions.

Our constrained least squares fit requires initial guesses for all of these parameters. The guess for the center position also provides the initial guess for by setting it equal to the 1D density fluctuations evaluated at . We found the initial values , , , and to lead to convergent fits across the whole dataset. In order to produce reliable fits we apply the following constraints: must remain within three pixels from the initial guess, and , to prevent numerical fitting errors.

ii.3.1 Physics-informed excitation classifier

Many candidate solitonic excitations are not vertically symmetric as might be expected, see e.g., Fig. 2(b-iii). The location of the largest “shoulder” in top half of the excitation is reversed with respect to the bottom half; in addition, the location of the minimum is slightly displaced going from the top half to the bottom. Inspired by these differences, we bisect each image into top and bottom halves (labeled by and

, respectively) and separately apply the Mexican hat fit to fluctuations in these data, giving vectors

. Using this observation, we develop a physics-informed excitation (PIE) classifier based on the single-soliton data set and discover that correlations between these vectors allow for a more fine-grained excitation classification.

Figure 3 shows the distribution of parameters from a single soliton dataset that were useful for classifying excitations. No meaningful correlations were found for parameters and , thus these did not assist in classification. The markers in the top panel show the amplitude ratio versus the top-bottom position difference , and show that they are not correlated. By contrast, the bottom panel shows that the asymmetric shoulder height difference is clearly anti-correlated with . Both panels are colored based on the cut-off points discussed in Sec. III.2 (see also Fig. 5).

This distribution and its correlation guide the classification rules described in Sec. III.2, yielding a PIE classifier based on cutoffs defined by human examination of the data.

Figure 3: Correlations between parameters implemented in PIE classifier. The top panel shows the distribution of center position difference versus the amplitude ratio (on a logarithmic scale). The bottom panel shows the correlation between the center position difference and the asymmetrical shoulder height difference for the gray points from top panel. Both panels are colored based on the cut-off points discussed in Sec. III.2.

ii.3.2 Quality estimation

Here we describe a quality estimate that a candidate excitation in an image is solitonic. We derive the likelihood that a vector of fit outcomes is drawn from a dimensional prior distribution spanning the set of representative solitonic excitations 111We found that was strongly correlated with the remaining five parameters and did not improve the the quality estimate performance. Ideally this distribution would be an uncorrelated multivariate normal distribution, but it is not, and we developed the following procedure to bring the distribution into this desired form.

We first fit a Yeo-Johnson power transformation yeo_new_2000

to each separate parameter distribution (having summed the 5D distribution along the remaining parameters) to transform them into independent zero-mean 1D Gaussian distributions with unit variance. Note that this treatment cannot transform the parameter distributions into perfect Gaussians, nevertheless each resulting distribution is balanced, contains a single-peak, and has long tails. The covariance matrix

is uncorrelated after this treatment and the distribution is qualitatively Gaussian in shape.

To calculate the quality estimate for a candidate excitation detected in an image, we:

  1. Fit the subtracted background 1D profile to Mexican hat function 3 giving .

  2. Use the established power transformation on to obtain .

  3. Return our quality estimate: , the likelihood between 0 and 1 that the excitation is solitonic.

The chi-squared cumulative distribution function

relates the Mahalanobis distance mahalanobis1936generalized to the likelihood that an outcome was drawn from the specified distribution 222This argument assumes no prior knowledge about the distribution of fit outcomes for structures that are not solitonic excitations.. is unbounded above and decreases to zero as approaches , the average over the prior distribution.

Iii Results

iii.1 ML modules

Figure 4: OD performance compared to ground truth (top), and the CNN classifier prediction (bottom), For ground truth and CNN classifier, the ticks “0”, “1”, “other” represent no, single, and other excitation classes. For OD, ticks represent the total number of positive excitations within an image.

We train both the CNN classifier and the OD using refined dataset with added soliton position labels (see Ref. Fritsch21-DSD ). The CNN classifier is trained using the full dataset while the OD training uses only the no excitation and single excitation classes. We assess the performance of both modules using 5-fold cross validation, that is using 80 % of the data to train a given module and the remaining 20 % to test it, and repeating the process 5 times to fully cover the dataset (see App. A for training details).

The results are summarized in the three cumulative confusion matrices plotted in Fig. 4. The top panel compares the outcome of the OD to the initial labels, showing near perfect delineation between no excitations and single excitations classes. However, the OD further subdivides the excitations class, counting anywhere from 0 to 4 candidate solitonic excitations within it. This results from the existence of excitations in this class that are not solitonic, as well as the possibility of having multiple solitons in the same image. The analogous comparison to CNN classification labels in the bottom panel is nearly indistinguishable from the one presented in the top panel, evidencing the quality of the CNN predictions.

Together, these ML tools effectively classify these data and locating excitations, however, it does not provide any fine-grained information on the nature nor the quality of the identified excitations. This is addressed in the following subsections.

iii.2 PIE classifier

The PIE classifier operates by applying a sequence of “cuts” driven by different combinations of the top-bottom fit outcomes . The exact parameter values described below are arrived at manually by exploring the data accepted and rejected by the cut to minimize the number of false positive kink soliton identifications.

The following cuts are applied sequentially, and the PIE classifier stops as soon as a classification is assigned.

Figure 5: The flow of the PIE classifier with example images for classification categories. Flow pathways and nodes are square-root scaled.

The amplitude parameters , and their ratio allow us to identify excitations that do not span the whole cloud. Cases with are classified as no excitation as they are associated with a density maximum, not a depletion. Cases with and are classified as “top partial excitation”, while those with and are classified as “bottom partial excitation”. Data with are also classified as “top partial excitation” and those with are classified as “bottom partial excitation”. This latter threshold identifies large fractional jumps in depth between the top and bottom that likely are off-axis solitonic vortices. Applying these cuts first is important because partial excitations interfere with the subsequent steps.


Figure 2(b-iii) illustrates a case with large shoulder height difference ; Ref. donadello_observation_2014 showed that such data result from solitonic vortices. As a result, we classify data with as “counterclockwise solitonic vortex” and as “clockwise solitonic vortex”.


Since kink solitons have a vertically aligned density depletion 333Our images are slightly rotated with respect to the vertical and horizontal axes leading to a small angle in our data. classify data with or as “kink soliton.”

Weaker cuts

Figure 3 shows that differences and are anti-correlated, indicating that asymmetries in position and shoulder height are related. A closer look at Fig. 2(b-iii) indicates that it is such a case, with and . We therefore add images with and to the “counterclockwise solitonic vortex” class and those with and to the “clockwise solitonic vortex” class.

Other data

The remaining images have but as “canted excitations”, likely kink solitons in the process of decay.

Figure 6: Quality estimate performance on no excitation and single excitation classes. In all cases we use the following color scheme: kink solitons (orange), all other solitonic excitations (green), and all non-solitonic local minima (blue). (a) Power transformed fit coefficient distributions, with untransformed variables labeled on the top axis. (b) Distribution of quality estimate. (c) Performance of quality estimate quantified by F1 score, the stars indicate the optimal F1 value and the circles mark the threshold we use for classification. The inset shows performance / recall curves.

The flow chart in Fig. 5 shows the application of this classifier to a single-soliton dataset; as a result the no excitation class is empty. We found that of the initial 3197 images about 1/3 failed a cut and were rejected as kink soliton candidates.

This classification was also used in the preparation of Ref. Fritsch21-DSD in which we present a refined soliton dataset, which includes improved single kink soliton labels. The cuts above are fairly aggressive to avoid false positives in the kink soliton classification. This implies possible misclassification in the other categories in order to ensure a high quality kink soliton dataset.

iii.3 Quality estimator

We initialized our quality estimator on the subset of the single excitation class identified as kink soliton using the PIE classifier. Figure 6(a) shows the power-transformed distribution of Mexican hat fit coefficients , with non-transformed coordinates marked on the top axis for reference. As would be expected, the data from the initialization dataset (orange) are nearly-normally distributed; interestingly, the remaining elements of the single excitation class (partial solitons, canted excitations, and solitonic vortices, as labeled by the PIE filter) follow very similar distributions (green). By contrast, the coefficients from every local minimum 444We require that the local minima be at least 7 pixels wide, i.e., a minimum at pixel must obey for . in the initialization set except solitonic excitations (blue curve) obey a qualitatively different distribution. These distributions demonstrate the ability of the quality estimator to discriminate between solitonic excitations and other features in the data, reinforcing the importance of the PIE filter for fine-grained classification.

Figure 7: Performance of quality estimate on other excitation classes. In all cases we use the same color scheme: purple plots data from all excitations identified as kink solitons by the PIE classifier and blue plots data from the set of local minima not identified as solitonic by the OD. (a) Power transformed fit coefficient distributions, with untransformed variables labeled on the top axis. (b) Distribution of quality estimate. (c) Representative images with (i-iii) OD+PIE identified one kink soliton (orange arrows) and one top partial soliton (green arrows), with the quality estimates of kink solitons equal ; and (iv) the OD+PIE found three kink solitons.

Using this initialization, we compare quality estimates obtained from the single excitation class in Fig. 6(b). The orange data show for kink solitons, and as intended the majority of this data is associated with larger values of

. The green data for the remaining solitonic excitations are nearly uniformly distributed, and the non-soliton minima (blue) are highly concentrated at small

. We note that the small peak in kink soliton distribution near-zero contains a negligible fraction of the kink soliton dataset (about ). However, this peak is more pronounced for the remaining excitations, which is not surprising because the power transform was initialized using kink soliton data.

We quantify the performance of the quality estimator in terms of the F1 scores plotted in Fig. 6(c), for kink solitons (orange) and all other solitonic excitations (green). The F1 score for kink solitons is maximized with a threshold of just (stars), however, in practice we minimized false positives and assign a feature to be solitonic when (circles). This choice gives only small change in the F1 score, however, it gives a marked increase in precision with only a small reduction in recall, as shown in the inset. The performance of the quality estimate on the other solitonic excitations, while far better than random, is subpar; this reemphasizes the importance of the PIE classifier in our framework.

iii.4 Application to other excitation class

Here we apply our SolDet framework to the other excitations class. This class includes images such as in Fig. 2(c) with multiple solitonic excitations, as well as confusing structures that made human labeling difficult. As such, this is an ideal test dataset since it defeated previous labeling attempts.

As a reminder, after the CNN classification step, the framework first uses the OD to locate soliton candidates that are then sorted by the PIE classifier, and we focus on only those features thus identified as kink solitons. Figure 7(a) plots the frequency of transformed Mexican hat fit outcomes , giving distributions that are qualitatively the same as those in Fig. 6(a) for the labeled single solitons. Furthermore, panel (b) histograms the quality estimate for these kink solitons, again exhibiting the same behavior seen in Fig. 6. Together these observations demonstrate that the efficacy of the OD/PIE combination in generating kink soliton candidates is comparable to manual identification.

Figure 7(c i-iii) shows three images that share the same labels. All three images were classified as other excitations by the CNN classifier and have two soliton candidates identified by the OD that were then classified as one kink soliton and one top partial by the PIE classifier.

Qualitative differences between these images are quantified by the quality estimate. Image (i), with is reminiscent of kink soliton but was assigned a low quality because it is off-center and shallow. By contrast image (ii) had modest value , penalized because it was shallow for a soliton near the center of the cloud. Image (iii) image contains two objects, but only one was identified as a kink soliton, and as expected by its high quality factor, is a clear kink soliton. Finally, (iv) shows an image with three relatively high quality kink solitons with from left to right. The left most one shows a clear kink soliton, while the other two have a lower quality estimate either because it is wide (middle), or it is shallow and has asymmetric shoulders (right).

iii.5 SolDet: Open-source Python package for solitonic excitation detection

Figure 8: The SolDet flow chart. The black line follows the SolDet dataflow and contains the labels added by each module (rectangles). Blue blocks represents ML modules, red blocks represent physics-based modules.

In this section, we describe our software package SolDet that integrates both the ML modules (CNN classifier and OD) with the fitting physics-based modules (PIE classifier and quality estimator), as we described in previous sections. The above discussion showed that the ML modules classify images effectively and can accurately locate one or many candidate solitons. The physics-based modules can sort these candidates into subclasses and provide a quality estimate for kink soliton candidates. Therefore, the ML and physics-based modules contribute to the task of soliton detection in different ways, and the SolDet infrastructure leverages their complementing strengths. We emphasize that soliton detection is one of a larger class of feature identification in quantum gases and that SolDet was designed to be broadly applicable.

The SolDet distribution includes a CNN classifier, OD, PIE classifier and quality estimator trained and initialized using the soliton dataset solitons-data . In addition, we provide training scripts to enable the ready application to user-defined data with custom preprocessors, ML models, fitting functions, and even the overall process flow.

Figure 8 illustrates a single use of SolDet for the specific example of kink soliton detection, where the individual blocks operate as follows:

Data processing

Preprocess raw data into image format that just enclose the elliptical atom clouds guo_machine-learning_2021 . The preprocessing particulars are not generic and instead are specific to both our task as well as the experimental parameters.

CNN classifier

Apply a trained CNN classifier to processed data, yield labels no excitation, single excitation, or other excitations.

Object detection

Apply trained OD to processed data, yield a list of positions of solitonic excitations.


If either the CNN classifier or OD finds no soliton, SolDet terminates.

PIE classifier

The PIE classifier is applied to each solitonic excitation.

Quality estimator

The quality estimator is applied to each excitation identified as “kink soliton” by the PIE classifier.

This algorithm is designed to be usable in a laboratory environment where one needs real-time identification, as well as for automated labeling of large datasets, as in Ref. Fritsch21-DSD .

Iv Discussion & Conclusion

In this work, we introduced a high level framework that combines ML methods with physics-based analysis, providing an integrated platform for studying experimental data. Our implementation of this framework, SolDet, currently targets the identification, classification, and tracking of features in image data generated by cold atom experiments. We demonstrated its initialization and performance using a publicly available dark soliton dataset Fritsch21-DSD . This investigation focused only on properties of individual images, however, the dataset also includes a label giving time elapsed since the excitation’s were created. This opens the door for studies correlating system control parameters and the SolDet labels.

While our initialization used only the no excitation and single excitation classes, SolDet’s feature detection successfully generalizes the learned patterns. This is confirmed by its performance on the other excitations class that was not part of training, where the CNN classifier gave ambiguous results and human classifiers often disagreed. Going beyond simple classification tasks, SolDet allowed us to identify unexpected structure in the data, for example dividing the single excitation into physically relevant subclasses, including solitonic vortices and partial solitons. This illustrates the power of our combined framework as a data analysis tool for discovery. Because the OD presented here is trained with the full single excitation, it identifies a range of excitations. An interesting extension of this work would be to train an OD on a dataset containing a single subclass found by the PIE classifier, e.g., kink solitons, or solitonic vortices.

From the ML perspective, adding modules based on unsupervised celebi2016unsupervised

, active learning 

sun2010survey , and synthetic data generation with generative models gui2020review may further enhance the performance of the SolDet framework.


This work was partially supported by NIST and NSF through the Physics Frontier Center at the JQI. We acknowledge the careful reading of the manuscript by D. J. D’Amato and Z. Grey.

Figure 9: Illustration of (a) OD and (b) CNN classifier neural network architectures. Yellow-orange boxes show convolutional layers while orange-red boxes show max-pooling layers. The horizontal lengths of boxes represent number of filters and the other two dimensions represent the image sizes. The horizontal blue and purple rectangles in (a) denote output vectors. Each cell of the blue vector describes the probability that it contains a soliton and and the purple vector contains the position of a soliton within the cell. And the vertical blue-green rectangles in (b) are three fully connected layers and the output layer. The lengths of edges are logarithmically scaled.
Layer 1 2 3 4 5 output
Filter 8 16 32 64 128 2
Kernel 55 55 55 15 15 15
Padding same same same same same same
Activation ReLu ReLu ReLu ReLu ReLu sigmoid
Pool size 42 42 41 21 N/A N/A
Strides 42 42 41 21 N/A N/A
Padding valid valid same same N/A N/A
Table 1: The OD architecture parameters. Middle rows are for for the convolutional 2D layers and bottom rows are for maxpooling 2D layers.

Appendix A Parameters of Machine Learning Models

Both machine learning modules are built and trained using the TensorFlow (v.2.5.0) Keras Python API 

tensorflow2015-whitepaper . Fig. 9(a) and (b) show the visualization of the network architecture for the OD and the CNN classifier, respectively. The model parameters of OD are presented in Tab. 1. The model parameters for the CNN classifier are presented in the Appendix of Ref. guo_machine-learning_2021 .

As can be seen in Fig. 9, there are three main differences between the two architectures: (1) the OD outputs 41 local probabilities and positions while the CNN classifier only outputs 1 of 3 possible classes; (2) the CNN classifier contains three fully-connected layers, which dramatically increase the number of trainable parameters, while OD does not; (3) the OD has asymmetric pool size and strides for vertical and horizontal directions, which are customized to the features in our dataset; the pool size and strides are symmetric for the CNN classifier. As a result, the OD has more than an order of magnitude fewer trainable parameters () than the CNN classifier ().


  • [1] R. Harper, S. T. Flammia, and J. J. Wallman. Efficient learning of quantum noise. Nat. Phys., 16:1184–1188, 2020.
  • [2] J. Carrasquilla and R. G. Melko. Machine learning phases of matter. Nat. Phys., 13(5):431–434, 2017.
  • [3] Y. Zhang and E.-A. Kim. Quantum Loop Topography for Machine Learning. Phys. Rev. Lett., 118(21):216401, 2017.
  • [4] G. Torlai, G. Mazzola, J. Carrasquilla, M. Troyer, R. Melko, and G. Carleo. Neural-network quantum state tomography. Nat. Phys., 14(5):447–450, 2018.
  • [5] J. Venderley, V. Khemani, and E.-A. Kim. Machine learning out-of-equilibrium phases of matter. Phys. Rev. Lett., 120(25):257204, 2018.
  • [6] B. S. Rem, N. Käming, M. Tarnowski, L. Asteria, N. Fläschner, C. Becker, K. Sengstock, and C. Weitenberg.

    Identifying quantum phase transitions using artificial neural networks on experimental data.

    Nat. Phys., 15(9):917–920, 2019.
  • [7] C. Miles, A. Bohrdt, R. Wu, C. Chiu, M. Xu, G. Ji, M. Greiner, K. Q. Weinberger, E. Demler, and E.-A. Kim.

    Correlator Convolutional Neural Networks: An Interpretable Architecture for Image-like Quantum Matter Data.

    Nat. Commun., 12:3905, 2020.
  • [8] P. Cha, P. Ginsparg, F. Wu, J. Carrasquilla, P. L. McMahon, and E.-A. Kim. Attention-based Quantum Tomography. arXiv:2006.12469, 2020.
  • [9] F. Metz, J. Polo, N. Weber, and T. Busch. Deep learning based quantum vortex detection in atomic Bose-Einstein condensates. Mach. Learn.: Sci. Technol., 2:035019, 2021.
  • [10] J. Venderley, M. Matty, K. Mallayya, M. Krogstad, J. Ruff, G. Pleiss, V. Kishore, D. Mandrus, D. Phelan, L. Poudel, A. G. Wilson, K. Weinberger, P. Upreti, M. R. Norman, S. Rosenkranz, R. Osborn, and E.-A. Kim. Harnessing Interpretable and Unsupervised Machine Learning to Address Big Data from Modern X-ray Diffraction. arXiv:2008.03275, 2021.
  • [11] N. Maskara, M. Buchhold, M. Endres, and E. van Nieuwenburg. A learning algorithm with emergent scaling behavior for classifying phase transitions. arXiv:2103.15855, 2021.
  • [12] S. Guo, A. R. Fritsch, C. Greenberg, I. B. Spielman, and J. P. Zwolak. Machine-learning enhanced dark soliton detection in Bose-Einstein condensates. Mach. Learn.: Sci. Technol., 2:035020, 2021.
  • [13] P. B. Wigley, P. J. Everitt, A. van den Hengel, J. W. Bastian, M. A. Sooriyabandara, G. D. McDonald, K. S. Hardman, C. D. Quinlivan, P. Manju, C. C. N. Kuhn, I. R. Petersen, A. N. Luiten, J. J. Hope, N. P. Robins, and M. R. Hush. Fast machine-learning online optimization of ultra-cold-atom experiments. Sci. Rep., 6(1):25890, 2016.
  • [14] A. D. Tranter, H. J. Slatyer, M. R. Hush, A. C. Leung, J. L. Everett, K. V. Paul, P. Vernaz-Gris, P. K. Lam, B. C. Buchler, and G. T. Campbell. Multiparameter optimisation of a magneto-optical trap using deep learning. Nat. Commun., 9(1):4360, 2018.
  • [15] S. S. Kalantre, J. P. Zwolak, S. Ragole, X. Wu, N. M. Zimmerman, M. D. Stewart, and J. M. Taylor. Machine learning techniques for state recognition and auto-tuning in quantum dots. npj Quantum Inf., 5(1):1–10, 2019.
  • [16] J. P. Zwolak, T. McJunkin, S. S. Kalantre, J.P. Dodson, E.R. MacQuarrie, D.E. Savage, M.G. Lagally, S.N. Coppersmith, M. A. Eriksson, and J. M. Taylor. Autotuning of double-dot devices in situ with machine learning. Phys. Rev. Appl., 13:034075, 2020.
  • [17] A. J. Barker, H. Style, K. Luksch, S. Sunami, D. Garrick, F. Hill, C. J. Foot, and E. Bentine. Applying machine learning optimization methods to the production of a quantum gas. Mach. Learn.: Sci. Technol., 1(1):015007, 2020.
  • [18] H. Saito.

    Creation and Manipulation of Quantized Vortices in Bose–Einstein Condensates Using Reinforcement Learning.

    J. Phys. Soc. Jpn., 89(7):074006, 2020.
  • [19] Y. Baum, M. Amico, S. Howell, M. Hush, M. Liuzzi, P. Mundada, T. Merkh, A. R. R. Carvalho, and M. J. Biercuk. Experimental Deep Reinforcement Learning for Error-Robust Gateset Design on a Superconducting Quantum Computer. arXiv:2105.01079, 2021.
  • [20] Mi.-Z. Ai, Y. Ding, Y. Ban, J. D. Martín-Guerrero, J. Casanova, J.-M. Cui, Y.-F. Huang, X. Chen, C.-F. Li, and G.-C. Guo. Experimentally Realizing Efficient Quantum Control with Reinforcement Learning. arXiv:2101.09020, 2021.
  • [21] A. M. Legendre. Nouvelles méthodes pour la détermination des orbites des cometes. F. Didot, 1805.
  • [22] W. Ketterle, D. S. Durfee, and D. M. Stamper-Kurn. Making, probing and understanding bose-einstein condensates. arXiv:cond-mat/9904034, 1999.
  • [23] A. R. Fritsch, M. Lu, G. H. Reid, A. M. Piñeiro, and I. B. Spielman. Creating solitons with controllable and near-zero velocity in Bose-Einstein condensates. Phys. Rev. A, 101(5):053629, 2020.
  • [24] T. P. Purdy, K. E. Grutter, K. Srinivasan, and J. M. Taylor. Quantum correlations from a room-temperature optomechanical cavity. Science, 356(6344):1265–1268, 2017.
  • [25] F. Huang, Z. Cao, J. Guo, S.-H. Jiang, S. Li, and Z. Guo. Comparisons of heuristic, general statistical and machine learning models for landslide susceptibility prediction and mapping. Catena, 191:104580, 2020.
  • [26] S. Ghosh, M. R Malgireddy, V. Chaudhary, and G. Dhillon. A new approach to automatic disc localization in clinical lumbar mri: combining machine learning with heuristics. In 2012 9th IEEE International Symposium on Biomedical Imaging (ISBI), pages 114–117. IEEE, 2012.
  • [27] J. Sakhnini, H. Karimipour, and A. Dehghantanha.

    Smart grid cyber attacks detection using supervised learning and heuristic feature selection.

    In 2019 IEEE 7th International Conference on Smart Energy Grid Engineering (SEGE), pages 108–112. IEEE, 2019.
  • [28] M. Greiner, O. Mandel, T. Esslinger, T.W. Hänsch, and I. Bloch. Quantum phase transition from a superfluid to a mott insulator in a gas of ultracold atoms. Nature, 415:39–44, 2002.
  • [29] T. Bothwell, D. Kedar, E. Oelker, J. M. Robinson, S. L. Bromley, W. L. Tew, J. Ye, and C. J. Kennedy. JILA SrI optical lattice clock with uncertainty of . Metrologia, 56(6):065004, 2019.
  • [30] S. Burger, K. Bongs, S. Dettmer, W. Ertmer, K. Sengstock, A Sanpera, G. V. Shlyapnikov, and M. Lewenstein. Dark Solitons in Bose-Einstein Condensates. Phys. Rev. Lett., 83(25):5198–5201, 1999.
  • [31] J. Denschlag, J. E. Simsarian, D. L. Feder, C. Clark, L. A. Collins, J. Cubizolles, L. Deng, E. Hagley, K. Helmerson, W. P. Reinhardt, S. L. Rolston, B. I. Schneider, and W. D. Phillips. Generating Solitons by Phase Engineering of a Bose-Einstein Condensate. Science, 287(5450):97–101, 2000.
  • [32] J. Zwolak. Dark solitons in BECs dataset. National Institute of Standards and Technology, 2020.
  • [33] A. M. Mateo and J Brand. Stability and dispersion relations of three-dimensional solitary waves in trapped bose–einstein condensates. New J. Phys., 17(12):125013, 2015.
  • [34] A. Voulodimos, N. Doulamis, A. Doulamis, and E. Protopapadakis. Deep learning for computer vision: A brief review. Comput. Intell. Neurosci., 2018:7068349, 2018.
  • [35] I.-K. Yeo and R. A. Johnson. A New Family of Power Transformations to Improve Normality or Symmetry. Biometrika, 87(4):954–959, 2000.
  • [36] S. Donadello, S. Serafini, M. Tylutki, L. P. Pitaevskii, F. Dalfovo, G. Lamporesi, and G. Ferrari. Observation of Solitonic Vortices in Bose-Einstein Condensates. Phys. Rev. Lett., 113(6):065302, 2014.
  • [37] A. R. Fritsch, S. M. Koh, S. Guo, I. B. Spielman, and J. P. Zwolak. Dark Solitons in Bose-Einstein Condensation: Real Dataset for Many-Body Physics Research. (in preparation).
  • [38] SolDet Team, 2021.
  • [39] D. Leykam, I. Rondon, and D. G. Angelakis. Dark soliton detection using persistent homology. arXiv:2107.14594, 2021.
  • [40] N. Ricker. Further developments in the wavelet theory of seismogram structure. Seismol. Soc. Am., Bull., 33(3):197–228, 1943.
  • [41] We found that was strongly correlated with the remaining five parameters and did not improve the the quality estimate performance.
  • [42] P. C. Mahalanobis. On the generalized distance in statistics. National Institute of Science of India, 1936.
  • [43] This argument assumes no prior knowledge about the distribution of fit outcomes for structures that are not solitonic excitations.
  • [44] Our images are slightly rotated with respect to the vertical and horizontal axes leading to a small angle in our data.
  • [45] We require that the local minima be at least 7 pixels wide, i.e., a minimum at pixel must obey for .
  • [46] M Emre Celebi and Kemal Aydin. Unsupervised learning algorithms. Springer, 2016.
  • [47] Li-Li Sun and Xi-Zhao Wang. A survey on active learning strategy. In 2010 International Conference on Machine Learning and Cybernetics, volume 1, pages 161–166. IEEE, 2010.
  • [48] J. Gui, Z. Sun, Y. Wen, D. Tao, and J. Ye. A review on generative adversarial networks: Algorithms, theory, and applications. arXiv:2001.06937, 2020.
  • [49] M. Abadi, A. Agarwal, P. Barham, E. Brevdo, Z. Chen, C. Citro, G. S. Corrado, A. Davis, J. Dean, M. Devin, S. Ghemawat, I. Goodfellow, A. Harp, G. Irving, M. Isard, Y. Jia, R. Jozefowicz, L. Kaiser, M. Kudlur, J. Levenberg, D. Mané, R. Monga, S. Moore, D. Murray, C. Olah, M. Schuster, J. Shlens, B. Steiner, I. Sutskever, K. Talwar, P. Tucker, V. Vanhoucke, V. Vasudevan, F. Viégas, O. Vinyals, P. Warden, M. Wattenberg, M. Wicke, Y. Yu, and X. Zheng. TensorFlow: Large-scale machine learning on heterogeneous systems, 2015. Software available from