title: Identifying and Tracking Solar Magnetic Flux Elements with Deep Learning
authors: Haodi Jiang Jiasheng Wang Chang Liu Ju Jing Hao Liu Jason T. L. Wang Haimin Wang
text: 
  
  
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, 
  University Heights, Newark, NJ 07102-1982, USA
  	hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Department of Computer Science, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  	hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Center for Solar-Terrestrial Research, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  Big Bear Solar Observatory, New Jersey Institute of Technology, 40386 North Shore Lane, Big Bear City, CA 92314-9672, USA
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  		hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Center for Solar-Terrestrial Research, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  Big Bear Solar Observatory, New Jersey Institute of Technology, 40386 North Shore Lane, Big Bear City, CA 92314-9672, USA
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  	hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Center for Solar-Terrestrial Research, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  Big Bear Solar Observatory, New Jersey Institute of Technology, 40386 North Shore Lane, Big Bear City, CA 92314-9672, USA
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  	hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Department of Computer Science, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  	hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Department of Computer Science, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  
  
  Institute for Space Weather Sciences, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  	hj78@njit.edu, jw438@njit.edu,  chang.liu@njit.edu, ju.jing@njit.edu, hl422@njit.edu, wangj@njit.edu, haimin.wang@njit.edu
  Center for Solar-Terrestrial Research, New Jersey Institute of Technology, University Heights, Newark, NJ 07102-1982, USA
  Big Bear Solar Observatory, New Jersey Institute of Technology, 40386 North Shore Lane, Big Bear City, CA 92314-9672, USA
  
  
  Deep learning has drawn a lot of interest in recent years due to its effectiveness 
  in processing big and complex observational data gathered from diverse instruments. 
  Here we propose a new 
  deep learning method, called SolarUnet, 
  to identify and track solar magnetic flux elements or features in observed vector magnetograms
  based on the Southwest Automatic Magnetic Identification Suite (SWAMIS). 
  Our method consists of a data pre-processing component that prepares training data
  from the SWAMIS tool,
  a deep learning model implemented as a U-shaped convolutional neural network for fast and accurate image segmentation, 
  and a post-processing component that prepares tracking results.
  SolarUnet is applied to data from the 1.6 meter Goode Solar Telescope 
  at the Big Bear Solar Observatory.  
  When compared to the widely used SWAMIS tool,
  SolarUnet is faster while agreeing mostly with SWAMIS on feature size and flux distributions, 
  and complementing SWAMIS in tracking  long-lifetime features.
  Thus, the proposed physics-guided deep learning-based tool can be considered as an 
  alternative method for 
  solar magnetic tracking.
  
  
   
  
  
  
  § INTRODUCTION
   
  
  Tracking magnetic flux elements
  is an important subject in heliophysics research 
  <cit.>.[In the study presented here, we focus on tracking
  signed, including positive and negative, magnetic flux elements.]
  Identifying and tracking the surface
  magnetic elements is useful in deriving statistical parameters of the local and global solar dynamo, 
  allowing for sophisticated analyses of solar activity <cit.>. 
  It not only helps scientists understand the distribution of magnetic fluxes <cit.>, 
  but also helps estimate the amount of energy in acoustic waves, 
  which play an important part in the heating of the solar chromosphere and corona <cit.>. 
  In addition, magnetic tracking is useful in deriving boundary conditions of 
  magnetohydrodynamic (MHD) modeling of the solar corona and solar wind. 
  In the past, many researchers have studied the behaviors and patterns of magnetic flux elements.
  For instance, <cit.> developed a technique to detect and classify small-scale magnetic flux cancellations 
  and link them to chromospheric rapid blueshifted excursions. 
  <cit.> investigated the occurrence and persistence of magnetic elements in the quiet Sun 
  to understand the scales of organization at which turbulent convection operates. 
  <cit.> reported findings related to small-scale magnetic flux emergence in the quiet Sun. 
  
  In magnetic tracking, features are defined as a visually identifiable part of an image, 
  such as a clump of magnetic flux or a blob in a magnetogram. 
  One of the most popular software tools for magnetic feature tracking 
  across multiple images/frames
  is the Southwest Automatic Magnetic Identification Suite <cit.>. 
  SWAMIS takes five steps to track magnetic flux elements: 
  (1) feature discrimination for each frame; 
  (2) feature identification within a frame; 
  (3) feature association across frames; 
  (4) occasional noise filtering; 
  and (5) event detection <cit.>. 
  Magnetic events are broadly classified into two categories: death and birth <cit.>; 
  the former refers to the end of a magnetic feature’s existence while
  the latter refers to the start of a magnetic feature’s existence. 
  
  In this paper, we present a new tool, called SolarUnet, to track magnetic flux elements.
  Our tool is built using deep learning <cit.>.
  The tool can detect three different types of events in each category, 
  namely (i) disappearance and appearance,  
  (ii) merging and splitting, 
  and (iii) cancellation and emergence. 
  The event “disappearance" is defined as the end of a single unipolar magnetic feature 
  that “fades away" to nothing in the absence of nearby features 
  across two frames; 
  the opposite event “appearance" is defined as the origin of a single unipolar feature 
  where the unipolar feature does not exist in the previous frame.
  The event “merging" is defined as 
  the combination of two or more like-sign features into a single magnetic feature;
  the opposite event  “splitting" is defined as 
  the breakup of a single magnetic feature into at least two like-sign features, 
  where the total flux of all child features is roughly the same as that of the parent feature.
  The event “cancellation" is defined as the demise of a magnetic feature that collides with one or more opposite-sign features,
  resulting in the demise of these features or an alive feature carrying the remaining flux;
  the event “emergence" is defined as the appearance of  
  opposite-sign features with approximately the same magnitude 
  or a new feature adjacent to previously existing opposite-sign features
  in a nearly flux-conserving manner.
  
  Deep learning, which is a subfield of machine learning, has drawn a lot of interest in recent years <cit.>.
  Inspired by its success in computer vision, speech recognition and natural language processing, 
  researchers have started to use deep learning in astronomy and astrophysics 
  <cit.>.
  In contrast to the existing methods for magnetic tracking 
  <cit.>, 
  our SolarUnet tool is built using deep learning.
  Compared to the most closely related magnetic tracking tool, SWAMIS, 
  which uses hysteresis as the discrimination scheme and a 
  gradient-based “downhill" method to identify features in a frame,
  SolarUnet runs faster while producing similar or complementary results.
  
  The rest of this paper is organized as follows.
  Section 2 describes observations and data used in this study.
  Section 3 presents details of SolarUnet and tracking algorithms used by the tool.
  Section 4 reports experimental results.
  Section 5 concludes the paper.
  
  
  
  § OBSERVATIONS AND DATA PREPARATION
  
  
  
  We adopted two collections of observations in this study.
  The first collection was conducted by
  the Near InfraRed Imaging Spectropolarimeter <cit.>
  of the 1.6 m Goode Solar Telescope (GST) at the Big Bear Solar Observatory  
  <cit.>.
  This collection contained 
  observations of the magnetic polarity inversion region in NOAA AR 12665 (431, -131) 
  during 20:16–22:42 UT on 2017 July 13. 
  The obtained data included spectro-polarimetric observations of a full set of Stokes measurements at the Fe I 1564.8 nm line 
  (0.25 Å
  bandpass) by NIRIS with a FOV (field of view) of 80 at 0.24 resolution and 56 s cadence.
  Vector magnetic field products in local coordinates were constructed 
  after removing azimuth ambiguity <cit.>. 
  
  The second collection of observations was conducted with a clear seeing condition; 
  BBSO/GST achieved diffraction-limited imaging during 16:17–22:17 UT
  on 2018 June 07. 
  The obtained multi-wavelength observations revealed detailed structural and evolutionary properties 
  of small scale magnetic polarities in quiescent solar regions at north of the disk center  (-32, 294). 
  The essential data included in this collection were the images taken by the GST's NIRIS using a 2048×2048 pixels Teledyne camera 
  with a 80 FOV. 
  The spatial resolution (at diffraction limit θ= λ / D) of the NIRIS images was 02, 
  and the temporal cadence was 56 s. 
  The magnetograms were then aligned based on sunspot and plage features, 
  with an alignment accuracy within 0.3, which was the best accuracy by using interpolation.  
  
  We prepared our training and testing sets by using the magnetograms 
  taken from the two collections of observations described above. 
  Because the magnetograms taken on 2018 June 07 had higher quality than the observations conducted on 2017 July 13, 
  we used the higher-quality magnetograms to prepare our training data so as to obtain a better magnetic tracking model. 
  Specifically, we gathered all 202 frames from the second collection of observations, 
  and excluded 6 images with poor quality (these excluded images were very noisy). 
  The remaining 196 frames were used as training data for magnetic tracking.
  The testing set contained all 147 magnetograms from the first collection of observations. 
  Table <ref> summarizes the numbers of training and testing images used in this study.
  
  
  
  
  
  § METHODOLOGY
   
  
  
  
  
   §.§ Overview of SolarUnet
  
  
  
  
  
  Figure <ref> explains how SolarUnet works.
  Training magnetograms are pre-processed in steps 1 and 2, 
  and then used to train the deep learning model for image segmentation (step 3).
  The trained model takes a testing magnetogram (step 4) and produces a predicted mask (step 5).
  Through post-processing of the predicted mask, SolarUnet produces magnetic tracking results (step 6).
  
  Specifically, in step 1, we apply SWAMIS with the downhill option
  to the 196 training magnetograms to get 196 masks.
  These images, including the magnetograms and masks, 
  are converted to 8-bit grayscale images of 720×720  pixels, 
  suitable for our deep learning model.
  Pixels in the masks belong to three classes represented by three colors/labels respectively:  
  positive magnetic flux with a label of 1 (white), 
  negative magnetic flux with a label of -1 (black), and 
  non-significant flux with a label of 0 (gray). 
  During pre-processing, we convert the 196 three-class masks obtained from SWAMIS
  to 196 two-class (binary) masks by 
  (i) changing the label of the non-significant flux regions from 0 to 1;
  and (ii) changing both the positive magnetic flux regions and negative magnetic flux regions to
  significant flux regions with label -1 (step 2).
  
  The 196 magnetograms (images) and two-class (binary) masks 
  are then used to train the deep learning model, implemented in 
  TensorFlow <cit.> and Keras <cit.>,
  for image segmentation (step 3).
  Because our deep learning model needs a large amount of data in order to train successfully,
  the model invokes the ImageDataGenerator[<https://www.tensorflow.org/api_docs/python/tf/keras/preprocessing/image/ImageDataGenerator>] 
  in Keras to perform data augmentation, expanding the 
  training set by shifting, rotating, flipping and scaling the training images during the model training process.
  Shifting an image is to move all pixels of the image
  horizontally or vertically while keeping the dimensions of the image the same. 
  Rotating an image is to rotate the image clockwise by a given number of degrees from 0 to 360.
  Flipping an image is to reverse the rows or columns of pixels in the image.
  Scaling an image is to randomly zoom the image in and either add new pixel values around the image or interpolate pixel values in the image.
  We train the deep learning model using 1 epoch with 10,000 iterations/epoch.
  In each iteration, the model randomly selects one of the 196 training magnetograms and its binary mask,
  feeds them to the ImageDataGenerator to generate a synthetic magnetogram and binary mask, 
  and uses the synthetic magnetogram and binary mask to train the model.
  There are 10,000 iterations and hence 10,000 synthetic magnetograms and binary masks are generated through the data augmentation process,
  where the 10,000 generated magnetograms and binary masks 
  are used for model training.[Notice that SWAMIS is applied only to the 196 training magnetograms
  mentioned in Table <ref>; SWAMIS is never run on the 10,000 generated (synthetic) magnetograms.]
  We have chosen to use data augmentation as opposed to acquiring more training data because 
  the quality of ground-based observations is subject to many factors such as seeing conditions and observing time limits. 
  Obtaining large volumes of high-quality training data requires further observations. 
  Nevertheless, using the synthetic training images produces reasonably good results as shown in Section <ref>.
  
  When a testing magnetogram is submitted, it is converted to a 8-bit grayscale image of 720×720 pixels, 
  and fed to the trained deep learning model (step 4).
  The trained deep learning model predicts a two-class (binary) mask, 
  containing  non-significant flux regions with label 1 and
  significant flux regions with label -1 (step 5). 
  We convert the predicted two-class (binary) mask back to a three-class mask 
  via post-processing as follows.
  For the non-significant flux regions with label 1, we change their label from 1 to 0.
  For the significant flux regions with label -1, we use the information of radial components 
  in the vertical magnetic field image shown in Figure <ref>,
  where the radial components are perpendicular to the plane of the Sun, 
  to reconstruct positive and negative magnetic flux regions.  
  Specifically, for each pixel x in the significant flux regions in the predicted binary mask, 
  we check the magnetic strength of the pixel, y, at x's corresponding location in the vertical magnetic field image.
  If y's magnetic strength is greater than 150 G,
  we set x as positive magnetic flux and change the label of this pixel from -1 to 1. 
  If y's  magnetic strength is smaller than -150 G, 
  we set x as negative magnetic flux and the label of this pixel remains -1.
  If y's  magnetic strength is between -150 G and 150 G,
  we set x as non-significant flux and change the label of this pixel from -1 to 0.
  This yields a three-class mask with the polarity information.
  
  Finally, we apply our magnetic tracking algorithms described in Section <ref> 
  to the testing magnetogram and masks
  to get tracking results  (step 6).
  Magnetic tracking is often involved with more than one testing magnetogram, 
  and we output the tracking results in all of the testing magnetograms.
  
  
  
   §.§ Implementation of the Deep Learning Model in SolarUnet
  
  
  
  Figure <ref> illustrates the deep learning model used in SolarUnet, 
  which is a U-shaped convolutional neural network.
  We adapt U-Net <cit.> to our work, enhancing it to obtain our model.
  The model has an encoder, a bottleneck, a decoder, followed by a pixelwise binary classification layer.[Please see the Appendix for
  more detailed descriptions of the technical terms used here.]
  The encoder consists of 4 blocks: E1, E2, E3, E4. 
  Each block has two 3×3 convolution layers, represented by blue arrows, 
  followed by a 2×2 max pooling layer, represented by a red arrow. 
  In each convolution layer, we adopt batch normalization <cit.> after convolution, 
  followed by a rectified linear unit (ReLU) activation function.
  Furthermore, we add a dropout layer <cit.> after each max pooling layer.
  The four encoder blocks E1, E2, E3, E4 have 32, 64, 128 and 256 kernels, respectively.
  
  
  
  The bottleneck, denoted Bot, mediates between the encoder and the decoder. 
  It uses two 3×3 convolution layers followed by a 2×2 up-convolution layer, represented by a green arrow. 
  The bottleneck has 512 kernels.  
  Similar to the encoder, the decoder consists of 4 blocks: D1, D2, D3, D4. 
  Each block has two 3×3 convolution layers followed by a 2×2 up-convolution layer. 
  The four decoder blocks D1, D2, D3, D4 have 256, 128, 64 and 32 kernels, respectively. 
  The input of each decoder block is concatenated by the output of the corresponding encoder block
  where the concatenation is represented by a gray arrow. 
  A dropout layer is added after each concatenation.
  Finally, a 1×1 convolution layer, represented by a turquoise arrow, with 2 kernels 
  followed by a softmax activation function,
  is used to produce a segmentation mask. 
  During testing, the deep learning model
  takes as input a testing magnetogram and produces as output a two-class mask.
  
  The input resolution of the encoder block E1 is set to 720×720 pixels to match the size of the testing magnetogram. 
  Each max pooling layer reduces the size by a factor of 2.
  Hence, the input resolution of the encoder block E2 (E3, E4 respectively) is
  360×360 (180×180, 90×90 respectively) pixels.
  The input resolution of the bottleneck, Bot, is 45×45 pixels.
  Each up-convolution layer increases the size by a factor of 2.
  Thus, the input resolution of the decoder block D1 (D2, D3, D4 respectively) is 
  90×90 (180×180, 360×360, 720×720 respectively) pixels.
  
  The loss function, L, used by the deep learning model is 
  the binary cross-entropy function defined below:
  
      L = - ∑_xlog y_c(x, W).
  
  Here, W are the parameters of the convolutional neural network,
  y_c(x, W) is the output of the softmax layer of the convolutional neural network, 
  and c is the class label (1 vs. -1) of each pixel x.
  
  In training the deep learning model, we adopt  adaptive moment estimation (Adam) 
  to find the optimal parameters of the model.
  The learning rate of Adam is set to 0.0001. 
  Adam combines the advantages of two popular methods: AdaGrad and RMSProp <cit.>. 
  In most cases, Adam achieves better performance than other stochastic optimization methods 
  including the stochastic gradient descent (SGD) with momentum employed by U-Net <cit.>.
  
  Although both our deep learning model and U-Net <cit.> have the same U-shaped architecture, 
  they differ in several ways.
  First, U-Net used SGD with a momentum of 0.99 to train and optimize its model. 
  By contrast, we choose Adam because it achieves better performance in our case where 
  the training process would be trapped in a poor local minimum if SGD were used. 
  Second, U-Net focused on imbalanced datasets and used a weighted cross-entropy loss function 
  to tackle the imbalanced classification problem. 
  By contrast, because our training set is relatively balanced in the sense that 
  non-significant flux regions roughly have the same number of pixels as significant flux regions, 
  we use the binary cross-entropy loss function as defined in Eq. (<ref>). 
  Third, we adopt batch normalization and dropout layers, which were not used by U-Net. 
  Batch normalization improves model learning, 
  stabilizes the learning process, reduces the learning (training) time and improves prediction accuracy <cit.>. 
  Dropout prevents neural networks from overfitting <cit.>, where
  overfitting means that a trained model fits training data too well, and cannot generalize 
  to make predictions on unseen testing data.
  Finally, we reduce the numbers of kernels of the encoder, bottleneck and decoder blocks 
  by a factor of 2 compared to U-Net to speed up the training process and reduce GPU memory usage.
  
  
  
   §.§ Algorithms for Magnetic Tracking and Event Detection
  
  
  
  After describing the deep learning model used in SolarUnet,
  we now turn to the magnetic tracking algorithms employed by SolarUnet.
  Based on the positive magnetic flux regions and negative magnetic flux regions found in Section <ref>,
  we identify signed magnetic flux elements or features 
  in a magnetogram (image/frame)
  by utilizing a connected-component labeling algorithm <cit.> 
  to group all adjacent segments in the
  positive magnetic flux regions and negative magnetic flux regions respectively
   if their pixels in edges or corners touch each other.
  We filter out those magnetic features whose sizes are smaller than a user-determined threshold.
  The features eliminated from consideration are treated as noise.
  Then, we assign each of the remaining features a label number 
  and highlight the features with different bordering colors. 
  Finally, we consider the association of features (magnetic flux elements) across different frames 
  to perform event detection.
  
  Based on the observational data and instruments used, 
  we calculate the moving distance D (number of pixels) of a magnetic flux element X as follows:
  
      D = C ×725  km/arcsec ×Δ s
  
  where C is the transverse speed (km/s) on the photosphere 
  according to the observational environment and Sun's activity, and
  Δs is the pixel scale.
  In this study, C is set to 4 km/s.
  For the NIRIS magnetograms used here, Δs = 0.083/pixel.
  We then calculate the radius of the region of interest (ROI) with respect to the location or position of the magnetic flux element X, 
  denoted ROI_p(X), as follows:
  	
      radius(_p(X)) = 2 × D + r
  
  where r is the radius of the smallest region that covers the magnetic flux element.
  The ROI_p(X) defines the region which the magnetic flux element X can not move beyond between two contiguous frames.
  
  The magnetic flux (X) is calculated by the surface integral of 
  the normal component of magnetic field passing through X, as follows: 
  
      (X) = ∫_S B_z dS
  
  where B_z is the magnitude of the magnetic field from the vertical magnetic field image of the testing magnetogram and
  S is the area of the surface of X.
  	
  For any two features or magnetic flux elements X and Y in a frame, 
  we define the distance between X and Y, denoted Dist(X, Y), as follows:
  
      Dist(X, Y) = min_x ∈ X, y ∈ Y d(x, y)
  
  where x and y are pixels in X and Y respectively, 
  and d(x, y) is the Euclidean distance between x and y.
  For a given magnetic flux element X in a frame, 
  the adjacent features of X are defined as
  the k-nearest neighboring features of X in the frame.
  (In the study presented here, k is set to 10.)
  
  Let X_i be a magnetic flux element in the current frame F_1.
  Let Y_i be a magnetic flux element in the next frame F_2
  where Y_i occurs in the ROI_p(X_i) in F_2.
  We say X_i is approximately equal to Y_i, denoted X_i≈ Y_i, if 
  X_i and Y_i have the same sign, and
  
      |(X_i)-(Y_i)/(X_i)| ≤ϵ_1
  
  where ϵ_1 is a user-determined threshold based on the observation setting and tracking task requirement. 
  (In the study presented here, ϵ_1 is set to 0.33.)
  
  With the above terms and definitions, we are now ready to describe our algorithms for
  magnetic tracking and event detection. 
  For each magnetic flux element or feature X_i in the current frame F_1,
  the algorithms below determine and indicate whether X_i 
  exists in the next frame F_2, or X_i is involved in a merging or cancellation event, or
  X_i disappears. 
  
    
  
  [Main Algorithm] 
  
  (1) If there exists 
  a magnetic feature Y_i in the ROI_p(X_i) in the next frame F_2
  such that X_i ≈ Y_i, then indicate X_i exists in F_2
  (more precisely, X_i becomes Y_i in F_2)
  and go to (2); otherwise go to (3).
  
  (2) Check the sign of X_i, highlighting X_i by yellow bordering color if X_i is positive
  and by green bordering color if X_i is negative. Exit the Main Algorithm.
  
  (3) Find all magnetic features
  in the ROI_p(X_i) in the current frame F_1.
  Group those features in the ROI_p(X_i) whose signs are the same as the sign of X_i
  into G_same(X_i)
  and group those features in the ROI_p(X_i) whose signs are opposite to the sign of X_i
  into G_opposite(X_i). 
  
  (4) Go to the Merging Algorithm to check whether 
  X_i and the features in G_same(X_i) meet the merging criterion.
  If yes, perform the merging using the Merging Algorithm and then exit the Main Algorithm; 
  otherwise indicate X_i is not involved in a merging event.
  
  (5) Go to the Cancellation Algorithm to check whether
  X_i and the features in G_opposite(X_i) meet the cancellation criterion.
  If yes, perform the cancellation using the Cancellation Algorithm and then exit the Main Algorithm; 
  otherwise indicate X_i is not involved in a cancellation event.
  
  (6) If X_i is not involved in a merging event according to (4) and 
  X_i is not involved in a cancellation event based on (5), indicate X_i disappears and highlight 
  X_i by purple bordering color.[ For the events belonging to the death category, 
  namely disappearance, merging and cancellation, magnetic features involved in the events are highlighted by different bordering colors
  (purple for disappearance, amber for merging and pink for cancellation)
   in the current frame  F_1.
  For the events belonging to the birth category, 
  namely appearance, splitting and emergence, magnetic features involved in the events are highlighted by different bordering colors
  (blue for appearance, aqua for splitting and red for emergence)
   in the next frame  F_2.]
  Exit the Main Algorithm.
  
    
  
  [Merging Algorithm] 
  
  (1)  For each feature X_j
  in G_same(X_i),  check whether there exists a feature Y_j
  in the ROI_p(X_j) in the next frame F_2 such that X_j ≈ Y_j, and
  if yes, delete X_j from G_same(X_i).
  Call the remaining set, G'_same(X_i). 
  If there are too many features in G'_same(X_i), only keep those adjacent features of X_i in G'_same(X_i).
  
  (2) If there exist 
  a combination C_s of features in G'_same(X_i) and
  a magnetic feature Y_i in the ROI_p(X_i)
  in the next frame F_2 where Y_i and X_i have the same sign, 
  such that Eq. (<ref>) below is satisfied, then we
  say X_i and the features in C_s are merged into Y_i:
  
      |((X_i) + ∑_X ∈ C_s(X)) - (Y_i)/(Y_i)| ≤ϵ_2
  
  where (X_i) and (X) have the same sign,
  ϵ_2 is a user-determined threshold (which is set to 0.5).
  Indicate X_i and the features in C_s are merged into Y_i by
  highlighting X_i and the features in C_s using amber bordering color.
  Exit the Merging Algorithm. 
  
  (3)
  If there does not exist Y_i or a combination of features satisfying Eq. (<ref>)
  (i.e., the condition in (ii) is not satisfied),
  indicate X_i and the features in G_same(X_i) do not meet the merging criterion.
  Exit the Merging Algorithm. 
  
    
  
  [Cancellation Algorithm] 
  
  (1)  For each feature X_j
  in G_opposite(X_i),  check whether there exists a feature Y_j
  in the ROI_p(X_j) in the next frame F_2 such that X_j ≈ Y_j, and
  if yes, delete X_j from G_opposite(X_i).
  Call the remaining set, G'_opposite(X_i). 
  If there are too many features in G'_opposite(X_i), only keep those adjacent features of X_i in G'_opposite(X_i).
  
  (2)
  If there exists a combination C_o of features in G'_opposite(X_i) 
  such that Eq. (<ref>) below is satisfied, then we say X_i and the features in C_o cancel each other
  (referred to as balanced cancellation in <cit.>):
  
      |(X_i) + ∑_X ∈ C_o(X) /(X_i)| ≤ϵ_2
  
  where (X_i) and (X) have opposite signs.
  Indicate X_i and the features in C_o cancel each other by
  highlighting X_i and the features in C_o using pink bordering color. 
  Exit the Cancellation Algorithm. 
  
  (3) If there exist a combination C_o of features in G'_opposite(X_i) 
  and a magnetic feature Y_i in the ROI_p(X_i) in the next frame F_2, 
  such that Eq. (<ref>) below is satisfied, 
  then we say X_i and the features in C_o are canceled to yield a feature Y_i carrying the remaining flux
  (referred to as unbalanced cancellation in <cit.>):
  
      | |(X_i) + ∑_X ∈ C_o(X)| - |(Y_i)| |/|(Y_i)|≤ϵ_2
  
  where (X_i) and (X) have opposite signs.
  Indicate X_i and the features in C_o are canceled by
  highlighting X_i and the features in C_o using pink bordering color.
  Exit the Cancellation Algorithm. 
  
  (4)
  If neither the condition in (ii) nor the condition in (iii) is satisfied, 
  indicate X_i and the features in G_opposite(X_i) do not meet the cancellation criterion.
  Exit the Cancellation Algorithm. 
  
    
  
  To determine whether the magnetic feature X_i appears 
  or is involved in a splitting or emergence event, we use the same algorithms as described above
  except that we treat the next frame F_2 as the current frame and
  the current frame F_1 as the next frame.
  
  
  
  § RESULTS
   
  
  SolarUnet is implemented in Python.[<https://www.python.org/>]
  Our deep learning model is coded with TensorFlow[<https://www.tensorflow.org/>] <cit.> 
  and Keras[<https://keras.io/>] <cit.> libraries.
  The data processed by SolarUnet, with the aid of Astropy[<https://www.astropy.org/>] <cit.>,
   include FITS files containing vector magnetic fields, 
  PNG images of the observational data described in Section <ref>, 
  and image masks obtained from the SWAMIS tool presented in <cit.>.
  SWAMIS, written in PDL (Perl Data Language)[<http://pdl.perl.org/>] and
  available via SolarSoft[<https://sohowww.nascom.nasa.gov/solarsoft/>]  
  <cit.>, was run with the downhill option.
  Figures in this section were produced with the aid of matplotlib[<https://matplotlib.org/>] <cit.>. 
  Statistical tests were performed by SciPy[<https://www.scipy.org/>] <cit.>.
  All experiments were conducted on a Dell PC with i7-8700k CPU, 32GB RAM and a NVIDIA GeForce RTX 2080 GPU 
  for training and testing the deep learning model. 
  
  
  
   §.§ Magnetic Tracking and Event Detection Results
  
  
  In this series of experiments, we used the 196 magnetograms mentioned in Table <ref> 
  and the corresponding masks obtained from SWAMIS to train SolarUnet as described in Section <ref>,
  and then we performed testing on the set of 147 magnetograms mentioned in Table <ref>.
  The testing set contained observations in NOAA AR 12665 (431, -131) during 20:16–22:42 UT on 2017 July 13. 
  The filter size threshold was fixed at 10 pixels.
  Thus, in the experiments we considered features or magnetic flux elements having at least 10 pixels.
  
  We present figures to illustrate the six events studied here.
  Frames taken at 20:15:49 UT and 20:16:45 UT are used to illustrate a disappearance event.
  Frames taken at 21:00:48 UT and 21:01:45 UT are used to illustrate an appearance event. 
  Frames taken at 20:18:38 UT and 20:19:34 UT are used to illustrate a merging event.
  Frames taken at 20:19:34 UT and 20:20:30 UT are used to illustrate a splitting event and a cancellation event.
  Frames taken at 20:17:41 UT and 20:18:38 UT are used to illustrate an emergence event. 
  We present enlarged FOV results in these figures with a FOV of 7.5where each figure has two axes: E-W (x-axis) and S-N (y-axis).
  
  Figure <ref> shows a disappearance event.
  In Figure <ref>(A), a magnetic feature highlighted by purple bordering color exists at
   E-W = 421 and S-N = -191,
  which is pointed to by a red arrow
  in the frame from 20:15:49 UT.
  This feature disappears in the next frame from 20:16:45 UT as shown in Figure <ref>(B).
  Figure <ref> illustrates an appearance event.
  In Figure <ref>(A), there exists no feature at
  E-W = 420 and S-N = -192 in the frame from 21:00:48 UT.
  However, a new feature appears in the next frame from 21:01:45 UT,
  which is highlighted by blue bordering color and pointed to by a red arrow as shown in Figure <ref>(B).
  
  
  
  
  
  Figure  <ref> shows a merging event.
  In Figure <ref>(A), 
  two separate positive polarity magnetic features highlighted by amber bordering color and pointed to by a red arrow 
  in the frame from 20:18:38 UT
  are merged into a single positive polarity feature at
  E-W = 453 and S-N = -180 in the frame from 20:19:34 UT as shown in Figure  <ref>(B).
  Figure <ref> illustrates a splitting event.
  In Figure <ref>(A), a negative polarity feature exists at
  E-W = 448 and S-N = -179 in the frame from 20:19:34 UT.
  This negative polarity feature is split into two negative polarity features highlighted by aqua bordering color and 
  pointed to by a red arrow in the frame from 20:20:30 UT as shown in Figure <ref>(B).
  
  
  
  
  
  Figure <ref> shows an unbalanced cancellation event.
  In Figure <ref>(A), there exist two magnetic features with opposite signs around
  E-W = 437 and S-N = -158.5 in the frame from 20:19:34 UT.
  The two magnetic features with opposite signs, highlighted by pink bordering color,
  are canceled to yield a negative polarity magnetic feature carrying the remaining flux, pointed to by a red arrow, 
  in the frame from 20:20:30 UT as shown in Figure <ref>(B).
  Figure <ref> illustrates an unbalanced emergence event.
  A new negative polarity feature emerges, next to a pre-existing positive polarity feature, 
  in the frame from 20:18:38 UT as shown in Figure <ref>(B).
  The flux of the positive polarity feature pointed to by a red arrow in Figure <ref>(A) is approximately
  equal to the total flux of the two features with opposite signs, highlighted 
  by red bordering color and  pointed to by a red arrow, in Figure <ref>(B).
  
  
  
  
  
  
  
   §.§ Comparison with SWAMIS
  
  
  While both SolarUnet and SWAMIS <cit.> 
  aim to track magnetic features and detect magnetic events,
  they differ in two ways.
  
  
    * Their feature discrimination and identification algorithms are different.
  SWAMIS used hysteresis and a threshold-based method to separate non-significant flux regions, 
  positive magnetic flux regions and negative magnetic flux regions.
  Then it used direct clumping and a gradient based (“downhill”) method to identify magnetic features in these regions.
  By contrast, SolarUnet employs a U-shaped convolutional neural network to
  gain knowledge from training data, and then predicts a binary (two-class) mask
  containing non-significant flux regions and significant flux regions.
  Next, SolarUnet separates the significant flux regions
  into positive magnetic flux regions and negative magnetic flux regions
  through post-processing of the binary mask.
  Finally, SolarUnet uses a connected-component labeling algorithm <cit.> 
  to group all adjacent segments in the
  positive magnetic flux regions and negative magnetic flux regions respectively
  if their pixels in edges or corners touch each other to 
  identify positive and negative magnetic flux elements.
  
    * Their feature tracking and event detection algorithms are different.
  SWAMIS used a dual-maximum-overlap criterion to find persistent features across frames. 
  In contrast, SolarUnet defines the region of interest (ROI) of a magnetic feature and
  traces the flux changes of the magnetic features in the ROI to find the association
  of features across frames.
  
  
  It should be pointed out that, although the U-shaped network (i.e., the deep learning model) 
  in SolarUnet gains knowledge from the training data prepared by SWAMIS, 
  the model is able to generalize learned features to more generic forms.
  In our work, the model gains knowledge from the training images in quiescent solar regions collected on 2018 June 07 and
  uses the acquired knowledge to make predictions on unseen testing images from an active region (NOAA AR 12665) collected on 2017 July 13.
  With the generalization and inference capability, the model 
  may discover new magnetic flux elements not found by the SWAMIS method.
  For example, with the filter size threshold of SolarUnet fixed at 10 pixels, 
  SolarUnet detected two opposite-sign features not found by SWAMIS
  on the testing image (magnetogram) from
  AR 12665 collected on 2017 July 13 20:15:49 UT.
  Figure <ref>(A) highlights these two features;
  Figure <ref>(B) shows that the two features were not found by SWAMIS.
  
  
  
  Figure <ref> compares the feature size distributions of SWAMIS and SolarUnet
  on the testing image (magnetogram) where the features had at least 2 pixels (0.007242 Mm^2).[In
  this and subsequent experiments, features with 1 pixel were considered as noise and excluded.]
  The feature sizes of SWAMIS are represented by blue color and those of SolarUnet are represented by orange color.
  It can be seen from Figure <ref> that SolarUnet agrees mostly with SWAMIS on the feature size distributions.
  To quantify this finding, we conducted the Epps-Singleton two-sample test <cit.>.
  According to the test, the results of SolarUnet and SWAMIS have a significant difference when p ≤ 0.05. 
  In our case p = 0.858 > 0.05, and hence we conclude that the results of the two tools are similar.
  Table <ref> shows the minimum, maximum, median, mean and standard deviation (SD)
  of the feature sizes found by SWAMIS and SolarUnet, respectively.
  SWAMIS detected 548 features while SolarUnet identified 543 features.
  The largest magnetic feature,
  which was a negative feature,
  found by SWAMIS had 60213 pixels
  (218.03 Mm^2).
  This feature was also detected by SolarUnet, with a smaller size of 57662 pixels
  (208.80  Mm^2).
  This size difference occurs due to the different feature identification and tracking algorithms used by the two tools.
  
  
  
  
  
  Next, for each feature detected by the tools, we calculated its flux using the formula in Eq. (<ref>). 
  Figure <ref> compares the feature flux distributions of SWAMIS and SolarUnet.
  The results in Figure <ref> are consistent with those in Figure <ref>; 
  SolarUnet agrees mostly with SWAMIS on the feature flux distributions.
  According to the Epps-Singleton two-sample test,
  the feature flux distributions of SolarUnet and SWAMIS have a significant difference when p ≤ 0.05. 
  In our case p = 0.983 > 0.05, and consequently we conclude that the feature flux distributions of SolarUnet and SWAMIS are similar.
  Table <ref> shows the minimum, maximum, median, mean and standard deviation (SD)
  of the feature fluxes found by SWAMIS and SolarUnet, respectively.
  The feature fluxes detected by SWAMIS ranged from 0.009507×10^18 Mx to 1805.13×10^18 Mx.
  The feature fluxes detected by SolarUnet ranged from 0.011051×10^18 Mx to 1780.96×10^18 Mx.
  Some of the small fluxes could be noise whiles others might be involved in
  small-scale magnetic flux emergence <cit.> or small-scale magnetic flux cancellation <cit.>.
  Similar results on feature size and flux distributions were obtained from the other magnetograms in the testing set.
  
  
  
  To further understand the behavior of SolarUnet and compare it with SWAMIS, we performed additional
  experiments to examine the lifetimes of the features identified and tracked by the two tools.
  We applied SolarUnet and SWAMIS
  to all of the 147 testing magnetograms mentioned in Table <ref>.
  The lifetime of a feature X is defined as X's disappearance time minus X's appearance time.
  More precisely, assuming X appears in the mth frame and disappears after the nth frame 
  (i.e., X is not shown in the (n+1)th frame),
  the lifetime of X is defined to be n-m+1 frames.
  Feature lifetime is strongly dependent on the feature identification and tracking algorithms employed by a tool <cit.>,
  and can be used to measure flux turnover rate <cit.>.
  
  Figure <ref> compares the lifetimes of features 
  found by SWAMIS and SolarUnet.
  SWAMIS tracked 48145 features across the 147 testing magnetograms while
  SolarUnet tracked 42470 features.
  The lifetimes of features found by SWAMIS ranged from 1 frame (56 seconds) to 138 frames (128.8 minutes).
  The lifetimes of features detected by SolarUnet ranged from 1 frame to 147 frames (137.2 minutes).
  SWAMIS tracked more short-lifetime features than SolarUnet
  while SolarUnet tracked more long-lifetime features than SWAMIS.
  Specifically, among the 48145 features tracked by SWAMIS,
  37110 features had a lifetime of 1 frame while
  SolarUnet only identified and tracked 22657 such features.
  On the other hand, SolarUnet tracked 19813 features whose lifetimes lasted more than 1 frame
  while SWAMIS only identified and tracked 11035 such features.
  SolarUnet complements SWAMIS in tracking long-lifetime features.
  We note that the training data of SolarUnet are from SWAMIS.
  For those features with short lifetime in the training images,
  our deep learning model may not acquire enough knowledge about them, 
  and hence may miss similar features in the testing images.
  This may explain why SolarUnet detects fewer short-lifetime features than SWAMIS.
  
  
  
  
  
  § DISCUSSION AND CONCLUSIONS
   
  
  We develop a deep learning method, SolarUnet, for tracking signed magnetic flux elements (features) 
  and detecting magnetic events in observed vector magnetograms. 
  We apply the SolarUnet tool to data from the 1.6 meter Goode Solar Telescope (GST) 
  at the Big Bear Solar Observatory (BBSO). 
  The tool is able to identify the magnetic features and 
  detect three types of events, namely disappearance, merging and cancellation, in the death category and
  three types of events, namely appearance, splitting and emergence, in the birth category.
  We use the BBSO/GST images to illustrate how our tool works on feature identification and event detection,
  and compares with the widely used SWAMIS tool <cit.>.
  
  Our main results are summarized as follows:
  
          1. For the testing data considered, SolarUnet agrees mostly with SWAMIS on feature size (area) and flux distributions, 
            and complements SWAMIS in tracking long-lifetime features. 
  It is worth noting that because SolarUnet performs magnetic tracking through making predictions, it is
   faster than the current version of SWAMIS.
  In general, SolarUnet runs in seconds on a testing magnetogram while
  the current version of SWAMIS runs in minutes on the same testing magnetogram.
  
  	2. SolarUnet is a physics-guided tool in the sense that it incorporates physics knowledge into its model and algorithms in several ways.  
  First, the training data of SolarUnet are from the physics-based SWAMIS tool.
  Second, in designing the loss function for the deep learning model used by SolarUnet, based on the observation that 
  non-significant flux regions roughly have the same number of pixels as significant flux regions in the training set,
  we adopt a binary cross-entropy loss function as defined in Eq. (<ref>) instead of the weighted cross-entropy loss function
  used by the related U-Net model <cit.>.
  Third, in converting the binary (two-class) mask predicted by our deep learning model for a testing magnetogram
  to a three-class mask with polarity information, we use the information of radial components in the vertical magnetic field image
  of the testing magnetogram to reconstruct positive and negative magnetic flux regions in the predicted mask.
  Lastly, by exploiting physics knowledge and based on the observational data and instruments used,
  we introduce the moving distance as defined in Eq. (<ref>)
  and region of interest (ROI) as defined in Eq. (<ref>) of a magnetic flux element
  to find the association of features across frames 
  so as to track these features.
  
           3. Although SolarUnet gets training data from SWAMIS, 
           our tool may discover new features not found by the SWAMIS method.
  	For example, refer to Figure <ref>.
         SolarUnet may detect smaller opposite-polarity features, as shown and highlighted in Figure <ref>(A), 
        near larger magnetic flux elements.  
  Small-scale energy release phenomena, ranging from coronal jets down to
  spicules, may be responsible for providing the upward flux of energy and
  momentum for the observed heatings and flows in the corona, and may
  plausibly drive the small transients in the solar wind recently
  discovered by the Parker Solar Probe <cit.>. There is
  mounting evidence that these events are generated via small-scale
  magnetic reconnection <cit.>, the photospheric
  signature of which is flux cancellation involving opposite magnetic
  polarities <cit.>. The ability of SolarUnet in detecting smaller
  opposite-polarity features near larger magnetic flux elements in a 
  faster manner can result in an improved determination of
  magnetic reconnection rate, thus contributing to the understanding of
  the mechanisms of solar coronal heating and the acceleration of the solar wind.
   
       4. The deep learning model in SolarUnet 
  performs binary (two-class) classification, i.e., predicting a two-class mask, rather than three-class classification, 
  i.e., predicting a three-class mask, during image segmentation. 
  SolarUnet produces a three-class mask through post-processing of the predicted two-class mask
  as described in item 2 above and in Section <ref>. 
  As indicated in the machine learning literature, multiclass classification including three-class classification
  often adds more noise to the loss function 
  <cit.>, 
  and it is easier to devise algorithms for binary classification
  <cit.>.
  We conducted additional experiments to compare SolarUnet with a three-class classification method.
  This method trained its deep learning model using the three-class masks obtained directly from SWAMIS, and
  predicted three-class masks.
  Its model was the same as SolarUnet's model except that
  (i) its loss function was changed from the binary cross-entropy function defined in Eq. (<ref>)
  to a categorical cross-entropy loss with three class labels (1, 0, -1);
  (ii) its softmax activation function was modified to output three-class masks.
  The three-class classification method used the same tracking algorithms as described in Section <ref> for magnetic tracking and event detection.
  The results of the three-class classification method were not as good as those of SolarUnet.
  For example, the feature size distribution obtained from the three-class classification method was 
  significantly different from the feature size distribution obtained from SWAMIS with p =  0.025 ≤ 0.05 
  according to the Epps-Singleton two-sample test on the testing magnetogram from AR 12665 collected on 2017 July 13 20:15:49 UT.		
  			
  
  
  Based on our experimental results, we conclude that the proposed SolarUnet
  should be considered a novel and alternative
  method for identifying and tracking magnetic flux elements. 
  More testing of the method, using different training and test data, should be performed.
  With the advent of big and complex observational data gathered from diverse instruments
  such as BBSO/GST and the upcoming Daniel K. Inouye Solar Telescope (DKIST),
  it is expected that the physics-guided deep learning-based
  SolarUnet tool will be a useful utility for processing and analyzing the data.
  
  We thank the referee for very helpful and thoughtful comments.
  We also thank the BBSO/GST team for providing the data used in this study. 
  The BBSO operation is supported by the New Jersey Institute of Technology
  and U.S. NSF grant AGS-1821294. 
  The GST operation is partly supported by the Korea Astronomy and Space Science Institute, 
  the Seoul National University, and the Key Laboratory of Solar Activities of the Chinese Academy of Sciences (CAS) 
  and the Operation, Maintenance and Upgrading Fund of CAS for Astronomical Telescopes and Facility Instruments. 
  Portions of the analysis presented here made use of the Perl Data Language (PDL),
  which has been developed by K. Glazebrook, J. Brinchmann, J. Cerney, 
  C. DeForest, D. Hunt, T. Jenness, T. Lukka, R. Schwebel, and C. Soeller 
  and can be obtained from <http://pdl.perl.org>.
  This work was supported by U.S. NSF grants AGS-1927578 and AGS-1954737.
  C.L. and H.W. acknowledge the support of NASA under grants NNX16AF72G, 
  80NSSC18K0673, and 80NSSC18K1705.
  
  Big Bear Solar Observatory
  
  
  
  Here we explain the technical terms used in describing our deep learning model (i.e., the U-shaped convolutional neural network). 
  
  
  Encoder is a neural network, which takes an input image and generates a high-dimensional vector
  that is an abstract representation of the image <cit.>.
  By using the encoder, our model can better understand the content and context of the image.
  
  
  Decoder is a neural network, which takes a high-dimensional vector and
  generates a segmentation mask <cit.>.
  By using the decoder, our model can recover the spatial information in the input image.
  
  
  Bottleneck, also known as the “compressed code" <cit.>, 
   is a layer with less neurons than the layer below or above it <cit.>.
  In general, it can be used to obtain a representation of the input with reduced size (dimensionality).
  In our model, bottleneck mediates between the encoder and the decoder. 
  
  
  Convolution layer contains multiple kernels where 
  a kernel is a matrix whose elements (weights) need to be learned from training data
  <cit.>. 
  Each kernel is multiplied with an image vector X (via element-wise multiplications)
  to produce a new image vector that contains 
  only the important information in X <cit.>.
  
  
  Max pooling layer reduces the size of an image vector X while 
  retaining only the important information in X
  <cit.>.
  
  
  Up-convolution layer, containing learnable parameters (weights),
  increases the size of an image vector X.
  This layer, also called an upsampling <cit.> or deconvolution layer,
  can recover the spatial information in X  <cit.>. 
  
  
  Softmax activation function converts a vector of k real values to a vector of k real values that sum to 1
  <cit.>. 
  Softmax is useful because it converts the scores in the vector to a normalized probability distribution, which can be displayed to a user.
  In our model, softmax is used to output the class label (1 vs. -1 or non-significant flux vs. significant flux) of each pixel. 
  
  
  Rectified linear unit (ReLU) employs an activation function f(x), defined as f(x) = max(0, x),  
  where x is the input to a neuron, f(x) = x if x ≥ 0 and f(x) = 0 otherwise
  <cit.>.
  It is easy to train a model that uses ReLUs, which often achieves good performance.
  
  aasjournal
  
  
  