Real-time Flare Prediction Based on Distinctions between Flaring and Non-flaring Active Region Spectra

Selecting an informative basis to represent the data, in this case spectra, is an integral part of machine learning. If we pass the entire profile to our algorithms, uninformative features such as the noisy small-scale structure of the continuum will hamper the algorithm's ability to learn, since it first has to learn what information is useful and what information should be neglected. Although there are many deep-learning models that automate this feature selection process, we find it more informative to describe each line profile in terms of a set of 10 parameters with known associations to observables. We list here the 10 selected features that provide us with a physical basis for describing each profile. Many of the listed features below (features, 1, 2, 3, 9, and 10), were shown by Leenaarts et al. (2013) in a forward modeling paper with a 10-level-plus-continuum model atom, to be descriptive of temperature, velocity, and density gradients within the chromosphere. To extract robust measures of the line center, line width, and asymmetry, we calculate each profile's normalized cumulative distribution function (NCDF) as illustrated in Figure 2. All features except the continuum emission and k/h ratio are taken exclusively with respect to the k-line, since both lines are highly symmetric.

• 1.  
  Intensity: for temperatures in excess of 6 kK, the peak intensity is strongly correlated to the gas temperature, with a correlation that starts to decrease for lower temperatures on account of the optical depth unity forming higher up in the chromosphere where the Planck and source function decouple. We divided each profile by its exposure time and selected the intensity as the maximum DN/s value.
• 2.  
  Line center: the Doppler shift of the line core from its rest frame is an excellent diagnostic of the line-of-sight velocity in the upper chromosphere. This was selected as the second quartile (q2) of each profile's cumulative distribution function (CDF).
• 3.  
  Line width: because of a small thermal width, line broadening is associated with nonthermal velocities. The line width was taken as (q3 − q1), where q1 and q3 are the first and third quartiles, respectively, of each profile's CDF.
• 4.  
  Line asymmetry: line asymmetry is an expression of nonthermal velocity flows and can take on values in the range±1, with positive and negative values being associated with up and downflows, respectively. The line asymmetry expressed in terms of quartiles is given by $[(q3-q2)-(q2-q1)]/(q3-q1)$.
• 5.  
  Total continuum: since the local continuum radiation between the k and h line cores is formed in LTE, the total intensity is dictated by the plasma temperature through the relation $I\simeq B\left[T\left(\tau =2/3\right)\right]$, where B is the Planck function and τ is the optical depth. For solar flares, small-scale heating events, and observations over sunspots, the continuum often appears flatter because of a different temperature structure. We define this feature as the sum of the hashed region in Figure 1 of a normalized profile. Therefore, larger total continuum values should be interpreted as intensity gains relative to the line core emission and not taken in the absolute sense.
• 6.  
  Triplet emission: the triplet emission on the red wing of the k-line core has been associated with lower atmospheric heating in the quiet Sun (Pereira et al. 2015). To distinguish triplet emission from continuum emission, we defined it as log(I_trip/I_wing), where I_trip is the height at the red vertical line in Figure 1 and I_wing is the height of the profile at the blue vertical line. Zhu et al. (2019) showed that the triplet lines at flare ribbons have a formation height similar to that of the h and k cores and, therefore, may not always be descriptive of the lower atmosphere.
• 7.  
  k/h ratio: the k/h ratio can be used as a measure of opacity, with ratios of 1:1 indicating a more opaque atmosphere and ratios of 2:1 indicating optically thin line formation (see Panos et al. 2018). We calculated the k/h ratio as the ratio of the integrated unhashed areas in Figure 1.
• 8.  
  k3-height: Rubio da Costa & Kleint (2017) showed that the central reversal can go into emission if large densities or temperatures are introduced at the core formation height. In the former case, the increased densities recouple the Planck and source functions, while in the later case, the temperature increase compensates for the decoupling allowing the source function to increase with height. Additionally, introducing layered velocity fields could fill in the reversal to create a single peak. The k3-height is taken to be the height of the k3 central minimum of a normalized profile. In the case of single-peaked profiles, k3 = 1.
• 9.  
  Peak ratios: the peak intensity ratios are correlated with upper chromospheric velocities, with large peak ratios indicating the existence of large velocity flows, and moderate ratios being linked to the difference between the average upper chromospheric velocity and the velocity at the peak formation height. We measure this as k2v/k2r. Carlsson & Stein (1997) showed relative intensity gains in the violet peaks of Ca ii to be correlated with downwards moving atmospheres in a region just above the peak formation height, an attribute shared by the Mg ii line as shown by Leenaarts et al. (2013). Therefore, values in excess of one correspond to larger 2v peaks and, consequently, indicate downflows that shift the peak opacity to longer wavelengths. Conversely, values smaller than one are associated with upflows above the peak formation height.
• 10.  
  Peak separation: the separation in wavelength between the 2v and 2r peaks are well correlated with the difference between the minimum and maximum velocities in the region of atmosphere just above the peak formation height. We define this separation as the number of wavelength grid points between k2v and k2r.

An evolution of these features leading up to and during an X1.6-class flare can be seen in Figure 3. Before implementing any machine-learning techniques, all features x_i have been standardized by

$$\begin{eqnarray}&&{x}_{i}^{{\prime} }=\displaystyle \frac{{x}_{i}-{\overline{x}}_{i}}{{\sigma }_{i}},\end{eqnarray} \tag{ 1 }$$

where ${\overline{x}}_{i}$ is the mean, and σ_i is the standard deviation of said feature. This standardization allows the algorithms to converge faster and removes biases based on the specific ranges of each feature.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

This study is concerned with the distinguishability of spectra collected from different solar regions. The spectral data were therefore divided into four subclasses collected from 25-minute full field of view (FOV) observation windows over four solar regions: QS, sunspot (SS), non-flaring AR and ARs resulting in flares, referred to as pre-flare regions (PF). The total number of instances analyzed in the conjoined data sets comes to 9 million spectral profiles, with the QS, AR, PF, and SS data sets comprising 32%, 28%, 25%, and 16% of the total data, respectively. For the details of each observation, we refer the reader to Tables 1 and 2. Each spectrum was then labeled according to the region they were extracted from. For instance, all spectra from an SS observation carry a label of 0, while all PF spectra carry a label of 1. These four different labeled data sets are referred to as classes in the standard machine-learning vernacular. The PF 25 minute observation windows were positioned such that the end of each window coincided precisely with the 1 minute mark to flare onset, as defined by the Geostationary Operational Environmental Satellite (GOES). Only observations leading to X- or M-class flares were considered for the PF data set, while the AR data set is a composite of ARs that never resulted in a major flare within the full (unparsed) IRIS observation time window, which often exceeded 25 minutes. Ordinarily, the duration of each observation varies, with typical values ranging from a few minutes to several hours, depending on the data rates and several other limiting factors unique to the specific research goals. The small time interval of 25 minutes for each observation is a limitation imposed by the PF data set, as IRIS on average only observed 25 minutes before each major flare. As a consequence of using the full FOV, the data sets are expected to carry an intrinsic number of shared profiles. Additionally, the data sets are naturally coupled to some degree, since our selected ARs contain sunspots, and PFRs ordinarily emanate from ARs. However, this does not preclude the possibility of these regions being remarkably distinct when observed in the near-ultraviolet.

The IRIS observations were chosen in a way that minimized the possibility of introducing unwanted biases. If the observations from the four target regions are not evenly distributed in time and on the solar disk, machine-learning algorithms could artificially distinguish between them based on some gross underlying temporal and geometric features not inherent to the physics of the Sun. For instance, if all AR observations are collected at the limb and all PF observations at disk center, then one might distinguish AR profiles from PF profiles based on the sum of the continuum emission, since center to limb variations result in a systematic decrease of the continuum (at least for plane parallel models). We note that this limb-darkening effect was not directly observed by us, possibly because the contribution to the Mg ii line shapes from local atmospheric conditions far exceeds the effects introduced by the observation angle. Nevertheless, to avoid introducing "artifacts," we attempted to match the distributions of each region's observations in time and on the solar disk as closely as possible. Additionally, we ensured that AR and PF region activity was similarly represented by carefully monitoring SJI animations and keeping only those IRIS observations for which the GOES soft X-ray flux remained below the 10⁻⁵ W m⁻² mark. The distributions of each of the four data sets in terms of the unstandardized 10 descriptive features can be seen in Figure 4. The data has been fitted with kernel density estimators (KDEs). A KDE can be viewed as the continuous counterpart of a histogram. For the case of a histogram, the height of each bar is determined by the number of points falling within the sharp boundaries of each bin. Histograms therefore do not make any assumptions about the underlying distribution of the data and have a single free parameter (number of bins), which can only be used to degrade the resolution of the data. In contrast, a KDE uses a subtler weighting function with smooth edges, called a kernel. The height or KDE ${\hat{f}}_{h}$, at a particular point x is determined by the number of contributing points falling within the kernel:

$$\begin{eqnarray}&&{\widehat{f}}_{h}(x)=\displaystyle \frac{1}{{nh}}\sum _{i=1}^{n}K\left(\displaystyle \frac{x-{x}_{i}}{h}\right),\end{eqnarray} \tag{ 2 }$$

where K is the kernel function (in our case Gaussian) and h the so-called bandwidth. The bandwidth determines the shape of the kernel function and, therefore, affects how each point within the kernel is weighted. We have set h = σ n^−1/5, where n is the number of data points for a specific region and σ, the variance for a particular feature of the data set. An alternative but equivalent interpretation is to imagine that kernels are centered at each data point and then summed together to produce a smooth continuous function. The denominator in Equation (2) divides out the area of the kernel, ensuring that the integral of the KDE equals unity. Note that the density function of a continuous variable at a single point need not be less than 1. Unlike histograms, KDEs do not degrade the data but generate missing data in a controlled manner. Therefore, they are used to estimate probability density functions of unknown distributions.

Zoom In Zoom Out Reset image size

The separation between each distribution in Figure 4 can be quantified with the use of a univariate F-score, or Fisher ranking score, which assumes uncoupled Gaussian distributed features. We have displayed the F-scores of each feature in Figure 5 for the case of PF/AR distributions. Higher F-scores imply larger dissimilarities between a feature's classes. In this first approximation, one can guess which feature is most or least important in distinguishing PF from AR spectra.

Zoom In Zoom Out Reset image size

In the PCA scheme (Hotelling 1933), high-dimensional data is projected onto a surface that best preserves the spread or variance of the data. As way of explanation, Figure 6 shows a bivariate Gaussian distribution with points living in a two-dimensional space. One-dimensional representations can be obtained by projecting the data orthogonally onto lines defined by unit vectors. Two such projections are shown in the inset at the bottom of Figure 6. The low-dimensional representation resulting from projecting the data onto the cyan line is superior to that of the red representation, since it captures more of the data's variance and is therefore more descriptive of the original data set. The eigenvectors of the covariance matrix, commonly referred to as principal components, form a convenient linearly independent basis for describing the data's variance, and eigenvectors with larger eigenvalues describe more of this variance. Instead of projecting our data onto a line, we constructed surfaces from the first two eigenvectors of the correlation matrix derived from the conjoined data sets (PF + QS), (PF + SS), and (PF + AR) separately. We then projected the spectra from 10-dimensions into two-dimensions, as seen in Figure 7, and compared the sample distributions of QS, SS, and AR spectra with spectra derived from the PF data set (orange). In all three cases, the first two principal components captured roughly 54% of the total variation within the data. The first panel on the left is for illustrative purposes, with the abstract data points being replaced by the spectra that they represent. We used QS spectra (cyan) and spectra from a flare (orange). Note that the two principal components appear to be descriptive of continuum emission and red wing enhancements, with large downflow profiles being found in the upper extremities of the plot. This clean correlation with single physical attributes is in no way guaranteed with PCA, and often the basis is composed of a complicated mixture of features. All remaining low-dimensional plots show only the abstract data points, since generating plots with 20,000 plus subplotted spectra is not feasible.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The plots indicate that PF and QS spectra are clearly different; however, the distinction between the first two principal components of PF and AR spectra is less obvious, with a large degree of overlap between their distributions. This does not mean that the PF and AR data sets are inseparable from one another, it may simply be that the PCA technique is not sophisticated enough to visualize a distinction that exists.

In contrast to PCA, t-SNE is a nonlinear map that adapts to the underlying structure of the data (van der Maaten & Hinton 2008). For a mathematical description of the algorithm, we refer the reader to the Appendix. What sets t-SNE apart from other dimensionality reduction algorithms is its ability to probe data at different length scales. A parameter called the perplexity, α, balances the importance between local and global structure representation. For instance, a large α penalizes the algorithm for representing the relationship between distant points poorly, whereas a small α only incurs a cost when nearby points are incorrectly represented. t-SNE is a complex dimensionality reduction algorithm that can lead to some unintuitive results. To better understand the output from t-SNE, we refer the reader to the excellent real-time simulations of Wattenberg et al. (2016). The separation of two classes can be judged with reference to Figure 8. The left panel shows a t-SNE embedding of flare and QS spectra, in contrast to the embedding on the right, based on spectra from the same data set. Spectra generated in solar flares are dramatically different in shape than those generated in the QS. As a result, the left panel shows an embedding between classes that are highly separable. Analogously, the right hand panel shows a t-SNE embedding between classes that are maximally inseparable, since they come from the same data set. Unlike many dimensionality reduction techniques, the distance between classes (blue and gray data points), is a poor measure of separability, and in the t-SNE framework, separability is measured by the degree to which the classes are mixed, with cleaner unmixed groups indicating higher separability between classes. The reader can judge the degree of separation between two classes with reference to these baselines. To gain insight into the topology of the data set, multiple runs at a variety of perplexity settings are required. We took the same combination of data sets used to generate Figure 7, and constructed low-dimensional representations at a variety of perplexity settings. The results complement that of PCA, where one can clearly separate QS from PF profiles, with a distinction that becomes more faded for SS and especially AR profiles. It is important to note that the sense of distance between clusters as well as relative cluster size lose all meaning in the t-SNE framework. Topological concepts such as containment, however, are preserved. Additionally, the algorithm is not guaranteed to converge to the same results, with multiple t-SNE runs at a fixed perplexity leading to different representations. Although different runs produce slightly different representations, we found no additional insight between these differences and, therefore, only display the results from a single t-SNE run in Figure 9.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The magnitudes of each profile's features have been plotted in Figure 10 for the (AR+PF) data set at a perplexity setting of 450. Each panel corresponds to one of the 10 standardized descriptive features, with the upper left panel indicating which points belong to PF (orange) and AR (green) spectra. The figures indicate that triplet emission is the most descriptive feature when trying to separate AR from PF spectra.

Zoom In Zoom Out Reset image size

A technique called Uniform Manifold Approximation and Projection, introduced by McInnes et al. (2018), is a scalable nonlinear manifold learning technique that uses concepts from Riemannian geometry and algebraic topology to produces results competitive to those of t-SNE, while arguably retaining more of the global structure. However, we found little discrepancy between the two methods on our data set and preferenced the t-SNE algorithm based on its simpler and more comprehendible mathematical formalism.

Neural networks are what we call supervised learning techniques that require a labeled data set {X, y}. In our case, xⁱ of X is a single Mg ii spectrum, with a corresponding label yⁱ of y indicating which solar region the spectrum originated from. From now on, we are only concerned with distinguishing PF spectra from non-PF spectra; therefore, our problem becomes a binary classification problem where yⁱ can take on the value 0 or 1. The goal of supervised machine learning is to use a fraction of the labeled data, called the training set {X_train, y_train}, to learn what is referred to as a model. This model can then be used to predict whether a profile is a PF profile 1 or not 0, from a never-before-seen test set {X_test, y_test}. The model that the algorithm learns is analogous to the human idea of a concept. In much the same way as a human can identify an apple with slightly different features from any apple they have seen before, so too can a trained NN make informed decisions about an unseen data set. This "mental agility" is what separates machine learning from other forms of artificial intelligence.

The upper image in Figure 11 shows a schematic of the simplest type of NN called a perceptron. This NN takes the 10 descriptive features associated with a single profile and forms a linear combination assigning relative importance to each feature through a set of pre-factors called weights w_i. This linear combination is then passed to an activation function, in this case, a sigmoid function $\sigma (z)=1/(1-{e}^{-z})$, which outputs a ${\hat{y}}^{i}$ value between 0 and 1. The linear combination plus activation function is called a node, and in NN diagrams is commonly represented by a circle. If a profile xⁱ belonging to the PF data set is fed into the NN, and the NN predicts that the profile comes from one of the other three data sets, i.e., produces an incorrect output close to 0, then the NN makes corrections to the weights w_i through a process known as back propagation, which allows the network to perform gradient descent across some cost function ${ \mathcal L }({\hat{y}}^{i},{y}^{i})$ sensitive to the severity of the mismatch between the guessed label ${\hat{y}}^{i}$ and the actual label yⁱ. Profiles from the training set are fed through the network, and the weights are continually adjusted until the model starts understanding the difference between PF/non-PF profiles. If there is no difference between these two classes, then the NN will converge to a model that always outputs a value of 0.5.

Zoom In Zoom Out Reset image size

The perceptron forms the basis of more complicated NNs in the same way that a biological neuron represents the building block for an organic brain. For more complicated NNs, the input features could connect to multiple nodes, each with their own associated system of weights such as those in the lower two panels of Figure 11. NNs that have two or more hidden layers (layers not including the input and output layer) are commonly referred to as Deep Neural Networks, with each layer learning a more abstract internal representation of the data. The number of layers in a network, number of nodes per layer and particular activation functions define what is called the network's architecture. It is difficult to know what architecture best suits a particular problem, and therefore, the architecture is chosen based on trial and error. For a more comprehensive review about the general mathematics behind NNs, we direct the reader to the detailed book of Murphy (2012).

We tested three NN architectures: architectures A, a simple perceptron with zero hidden layers, architecture B with two hidden layers, and architecture C with four hidden layers, as seen in Figure 11. We also included an SVM with a nonlinear kernel. SVMs are linear classifiers used for binary classification problems, which can be made to fit nonlinear decision boundaries by mapping the data into a higher dimensional space via the use of an appropriate kernel function. An SVM of this type should be comparable with a two-hidden-layer NN, where the first layer projects the data into a more abstract space, while the second layer classifies the data. For a detailed description of an SVM applied to solar flare prediction, see, for example, Bobra & Couvidat (2015). The three models in Figure 11 are what we call fully connected NNs, on account that every node between adjacent layers is connected.

Each model is trained on a fraction of the data, and the model's performance is measured on the test set. Several different performance metrics can be obtained by comparing the actual set of labels y_test = (0, 1, 1, ⋯, 1, 0, 1) to the model's output guesses $\hat{y}=(.1,.8,.6,\cdots ,\; .\,7,.2,.9)$. The results of a binary classification problem can be characterized by a confusion matrix, where the number of true positives TP (${y}^{i}=1\,| \,{\hat{y}}^{i}=1$), false negatives FN (${y}^{i}=1\,| \,{\hat{y}}^{i}=0$), true negatives TN (${y}^{i}=0\,| \,{\hat{y}}^{i}=0$), and false positives FP (${y}^{i}=0\,| \,{\hat{y}}^{i}=1$) are used in different combinations to express different aspects of the model's performance. We list all of the metrics used in this study, starting with

$$\begin{eqnarray}&&\mathrm{precision}=\displaystyle \frac{\mathrm{TP}}{\mathrm{TP}+\mathrm{FP}},\end{eqnarray} \tag{ 3 }$$

which is a measure of reliability when the model predicts that a profile comes from the PF data set.

$$\begin{eqnarray}&&\mathrm{recall}=\displaystyle \frac{\mathrm{TP}}{\mathrm{TP}+\mathrm{FN}},\end{eqnarray} \tag{ 4 }$$

on the other hand, is a measure of how many PF profiles were correctly labeled. Looking at just one of the above two metrics could be misleading, since it is possible to have a high precision and a low recall at the same time. The harmonic mean of both measures addresses this ambiguity by rolling both metrics into the so-called F1 score, given by

$$\begin{eqnarray}&&{\rm{F}}1=2\left(\displaystyle \frac{\mathrm{precision}\times \mathrm{recall}}{\mathrm{precision}+\mathrm{recall}}\right).\end{eqnarray} \tag{ 5 }$$

Notice that this formulation gives an interpretation of performance based on the coupled behavior of both quantities, i.e., favorable performance is assigned to classifiers/models who have both high precision and recall scores. Another useful metric is the ratio of correct predictions over the total number of predictions, called the accuracy:

$$\begin{eqnarray}&&\mathrm{accuracy}=\displaystyle \frac{\mathrm{TP}+\mathrm{TN}}{P+N},\end{eqnarray} \tag{ 6 }$$

where P and N are the number of observations in the positive and negative class, respectively. These measures are susceptible to misinterpretation due to the ambiguity introduced by the relative size of each class, called the class imbalance. For instance, say we have a high class imbalance with only a single PF profile in comparison to 99 non-PF profiles. Labeling every profile as non-PF results in a deceptively positive accuracy of 99%. To address the problem of class imbalance, Bloomfield et al. (2012) suggest the use of a ratio invariant skill score (Woodcock 1976), known as the true skill statistic (TSS). This measure uses all of the confusion matrix elements and promotes a standard for the comparison of studies with different flare/no-flare ratios. The TSS is defined as the difference between the recall and false alarm rate:

$$\begin{eqnarray}&&\mathrm{TSS}=\displaystyle \frac{\mathrm{TP}}{\mathrm{TP}+\mathrm{FN}}-\displaystyle \frac{\mathrm{FP}}{\mathrm{FP}+\mathrm{TN}},\end{eqnarray} \tag{ 7 }$$

with a TSS score of 1, 0, and −1 indicating a model with perfect, random, and adverse labeling, respectively. The problem of class imbalance is not directly applicable to our case, since we ensured that the number of spectral profiles in the PF/non-PF classes were held equal at all times; however, for completeness, we have included as many measures as possible.

The performance of each model over the binary classification problem PF/AR can be seen in Figure 12, with each column corresponding to the performance of a specific model on an observationally split data set. The term "observationally split" means that each IRIS observation has either been isolated to the training or test set, with none of that observation's spectra being shared between the two data sets. If the training and test sets were divided on a spectral basis only, without regard for the division of observations, the models could gain leverage on the shared data and learn spectral peculiarities unique to each particular observation. Splitting our data sets observationally ensures that our results are not artificial and remain as close to reality as possible. Figure 12 shows that architecture B, with two hidden layers, has the strongest performance across all metrics. The splitting of data into a training and test set injects a level of uncertainty into the results, with different splittings resulting in slightly different model performances. In machine learning, an estimation of this uncertainty can be captured by performing multiple runs with the same model, while continually swapping out the observations that appear in the test and training sets. Therefore, to construct the error bars in Figure 12, we performed the computationally expensive exercise of (1) randomly splitting the AR and PF observations seen in Tables 1 and 2, such that each random splitting resulted in the delegation of 4 AR and 4 PF observations in the test set, while the spectra from the remaining observations formed the training set; (2) training each model architecture over 15 epochs (the number of times the network is allowed to see the entire data set to update its weights) before extracting the model at an epoch that scored the highest accuracy on the test set; and (3) calculating the model's performance in terms of all of the above mentioned metrics over the data from the current unseen test set. This three-step procedure was iterated 100 times for each model, resulting in a distribution of each model's metric scores. The mean of each metric score was then plotted in Figure 12, and the standard deviations were used as error bars.

Zoom In Zoom Out Reset image size

Each of the metrics mentioned thus far are sensitive to the so-called threshold value δ, which defines the "watershed" between a PF/AR decision. In general, the threshold is set to δ = .5, so that output values $\hat{y}\lt .5$ are classified as non-PF, while values $\hat{y}\geqslant .5$ are classified as PF profiles. Choosing a threshold value turns a model into a binary classifier. Adjusting the threshold has a direct effect on the confusion matrix elements, which in turn cascade into all of the performance metrics. For this reason, an additional metric that integrates over the threshold value is preferable. This is commonly achieved by monitoring the relationship between the true and false positive rates given by

$$\begin{eqnarray}\begin{array}{rcl}\mathrm{TPR}(\delta ) & = & {\int }_{\delta }^{\infty }{f}_{1}(y)\,{dy},\\ \mathrm{FPR}(\delta ) & = & {\int }_{\delta }^{\infty }{f}_{0}(y)\,{dy},\end{array}\end{eqnarray} \tag{ 8 }$$

where f₁ and f₀ are the probability distributions assigned by the model to the positive and negative classes, respectively. Both of these distributions can be seen in the left panel of Figure 13, with f₀ (gray) corresponding to the probabilities associated with AR spectra, and f₁ (orange) to the PF spectra of a never-before-seen observationally split test set. As the threshold value changes, so too does the true positive and false negative rate, giving rise to the so-called receiver operating characteristic curve, or (ROC) curve seen in the right hand panel of Figure 13. Each point on the ROC curve represents a trade off between a high number of false positives and a high number of false negatives, with the dashed meridian line representing a random guess. A perfect score is represented by a point in the upper left hand corner of the parameter space, and the area under the curve

$$\begin{eqnarray}&&\mathrm{AUC}={\int }_{\delta =0}^{1}\mathrm{TPR}\left({\mathrm{FPR}}^{-1}(\delta \right)\,d\delta ,\end{eqnarray} \tag{ 9 }$$

is commonly used as a threshold-integrated meta statistic to judge the performance of a model on a binary classification problem. In this way, the AUC can be used to quantify the overlap between the two probability distributions, with larger areas corresponding to less overlap and therefore stronger model performance. Because we integrate out the threshold, the AUC score can be interpreted as the probability that the model will rank a randomly chosen PF spectrum higher than a randomly chosen non-PF spectrum. Therefore, the AUC is a metric that is both scale and threshold invariant. For this reason, the AUC scores have been computed and included in Figure 12.

Zoom In Zoom Out Reset image size

The first clue to support the notion of separability can be found in Figure 4. Many of the feature distributions show distinct differences between the four data sets. The QS profiles for instance have, on average, higher continuum emission, lower triplet emission, deeper central reversals, and appear to be formed under comparatively optically thicker conditions as suggested by their reduced k/h ratios. On the other hand, the distinguishing features of PF profiles appear to be enhanced intensities, triplet emission, and decreased continuum emission. We caution the reader when making literal inferences from these distributions, as they are subject to both the performance of our feature-finding algorithm, as well as the particular definition of each feature, which, in the case of triplet emission, can be quite obscure.

Low-dimensional representations of the data provide an alternative but complementary approach to viewing the possible separability of the four data sets, with both PCA and t-SNE supporting the premise that PF profiles can be distinguished from QS, SS, and AR profiles. The t-SNE representations in Figure 9 appear to display more dissimilarity than their PCA counterparts, while at the same time maintaining richer structures that may be indicative of the original data set. There is a clear distinction between QS and PF profiles for perplexities in excess of 2. The distinction between SS, AR, and PF profiles is slightly more strained, with a greater overlap between profiles from different data sets; however, profiles from different regions tend, on average, to occupy distinct subspace, especially for higher perplexities.

Increasing the perplexity from left to right results in an evolution from noisy circular representations with the two relevant classes being largely indistinguishably distributed (first column), to more ordered and complex structures. One can understand this evolution as follows: the perplexity defines the volume of space that a particular data point is aware of. The larger the perplexity, the more neighbors are included in each point's reference frame. t-SNE then embeds the data into a low-dimensional space, such that each point's frame of reference is minimally augmented. In this way, the data points can be distorted on a global scale, while still maintaining the same local structure. The first column in Figure 9 has the lowest perplexity, meaning each point is surrounded by a small volume containing only its nearest neighbors. Consequently, the global structures, such as the data's border, are of no significance to the map, and the representation only captures local clusters. In contrast, large perplexities have to preserve large-scale structures, leading to more elaborate borders.

Beyond their immediate aesthetic appeal, low-dimensional embeddings often serve no purpose other than to display the separability of classes or speed up the training of certain models by removing superfluous features. Projecting feature values onto the final embeddings, like we have done in Figure 10, allow us to gain traction back into the world of physics. We find it best to analyze Figure 10 in terms of symmetries. For instance, triplet emission has a rough symmetry about the horizontal axis, corresponding to the maximum direction of change of the triplet emission variable. In contrast, the line center divides the plot in such a way that it has an approximate symmetry about the vertical axis. Since the reference image in the upper left plot indicates a symmetry about the horizontal axis, any feature with horizontal symmetry is more useful for class prediction, because the direction of change of said feature overlaps maximally with how the PF and AR profiles are distributed. Following this line of reasoning, in order of predictive capacity, we have: triplet emission, followed by intensity, total continuum, line width, k/h ratio, peak separation, k3-height, peak ratios, line asymmetry, and finally, line center position. Conversely, those features displaying symmetry about the vertical axis are least valuable as class predictors, since profiles from either the PF or AR class can be found anywhere along their range. This is the case for all velocity features, which clearly divide the upper and lower hemispheres of the plot into down and upflows. For instance, assuming that a mean of 0 (gray) corresponds to an unshifted spectral profile, the feature associated with line center, (bottom right feature plot), divides the upper and lower hemispheres of the maps into down and upflows, respectively. A similar division is made in the line asymmetry map, with red negative values associated with downflows and blue positive values with upflows. The peak ratios divide the map in a similar way, but indicate a reversal of velocities, with the lower region being associated with downflows and the upper region with upflows. We are not sure how to resolve this discrepancy, and are hesitant to dismiss it as an oddity arising from our particular methods of measurement. What we can conclude however, is that, on average, velocity features seem indistinguishably distributed between PF and AR profiles. All of these observations are consistent pictorial realization of the F-scores displayed in Figure 5.

The t-SNE feature projection plots allow us to examine the coupled behavior of features, something that cannot be directly obtained from the univariate Fisher scores. Features corresponding to opposite symmetry groups represent loosely uncoupled variables. For instance, a profile can be highly symmetric or antisymmetric regardless of how much triplet emission it has. The uncoupled sets of variables are horizontal: (triplet emission, intensity, total continuum, line width, k/h ratio), middle: (peak separation, k3-height), and vertical: (line center, line asymmetry, peak ratios). Each feature has been ordered such that the leading features of each group have the highest degree of freedom. In contrast, features that fall into the same symmetry groups are coupled to some degree. An interpretation of the coupling from the horizontal group is as follows: more intense profiles are often accompanied by triplet emission and decreased continuum emission, and they have a tendency to form under optically thinner conditions (higher k/h ratios), with a weak propensity toward closer possibly single-peaked spectra. Additionally, peak separation and k3-height are strongly anticorrelated, with deeper reversals being associated with broader profiles.

The strong coupling between total continuum and line width should be dismissed as artificial, since it is clear that some continuum is included within the line windows of Figure 1. Therefore, an increase in total continuum results in an increase in line width. Once again, these interpretations are based on the performance of our feature-finding algorithm, and the statements we have made are statistical in nature and only capture general trends.

From Figure 12, it is clear that model B is superior on all fronts when it comes to distinguishing PF/AR profiles, with the four-hidden-layer NN performing below expectations. This may be due to overfitting, since the network has many more free parameters that could be used to form extremely complex decision boundaries that do not generalize well. However, no such overfitting was apparent from the learning curves, and subsequent trials involving the insertion of dropout layers returned little to no improvement. We can conclude that the division between PF/AR profiles is relatively simple, requiring only a light-weight NN to produce satisfying results. The success of model B may also indicate the relevance of our choice of basis, since additional feature-learning layers do not seem to be required. We expect that the intrinsic number of dimensions necessary to represent the data is no more than 10. Our basis was selected on the grounds that each feature can be used as a diagnostic for some observable. If a feature cannot be linked unambiguously to some statement about the solar atmosphere, then it is by definition random and nondescriptive.

The high metric scores warrant the conclusion that there is a great deal of continuity between IRIS observations of the same type. That is, different IRIS PF observations produce similar-looking spectra.

It appears that there is a persistent overlap in Figure 13 between the classification distributions of PF/AR scores. We postulate that this overlap has two distinct sources. The overlap in the middle of the distribution is likely due to the inclusion of QS profiles within both the PF and AR data sets. This is understandable since the IRIS spectrogram can cover an FOV of 130 × 175 arcsec², while an ordinary active region only occupies at most half that FOV, so that all three data sets (AR, SS, PF) are not entirely clean, but rather contain a certain number of mislabeled profiles that cannot be untangled. The overlap toward the tails of the distributions however, indicate the degree to which a non-PF region can produce a PF like profile, and visa versa. For the PF/AR case in Figure 13, we see that the probability of producing successively more convincing profiles from the opposite class steadily decreases.

The metric scores and ROC curves for the spatial experiment indicate that the model has learned to correlate PF spectra with higher levels of activity, as seen in the corresponding 1400 SJIs. This activity can be seen toward the top of the SJI in panel (B) of Figure 14. At 80'' away from these visible brightenings, the model does not know if the spectra are from PF or AR data sets, however, the scores and AUCs drastically increase when we present the model with spectra from successively closer regions. This equates to an f₁ distribution that starts out centered at 0, and migrates toward 1. Since the metric scores and curves go from random guess performance to perfect classification, fears related to possible biases introduced by yearly varying baseline UV emissions can be dispelled. Furthermore, these results fortify the assumption that mislabeled QS profiles are responsible for the central overlap seen in the probability distribution plot in Figure 13. This conclusion comes from the fact that at 80 arcseconds off of an active region (where one would expect QS spectra), the model produces metric scores that are entirely consistent with a random guess. This is precisely the score given in the center of the probability distribution plot. With this in mind, the entire enterprise of the spectral flare prediction would greatly benefit from a careful pruning of QS profiles, in a non-full FOV/targeted implementation.

The scores for the spatial experiment do not predict the future location of the flare ribbon that propagates directly over the slit at flare maximum, as seen in Panel D of Figure 14. This leads us to believe that the model is simply keeping track of the number of energetic events per time sample.

The metric scores for the temporal experiment do not vary as extensively, however, the scores and ROC curves indicate that the f₁ distribution shifts closer to 1 at about 35 minutes before flare onset, and continues to do so for successive future time steps. One can see a significant increase in metric scores in Figure 15 between t = −60 and t = −50. This increase occurs in conjunction with the first signs of intensity enhancements over the slit, displayed in panel (A) of Figure 14. Each successive step forward in time after the 30 minute mark produces higher metric scores and larger areas under the curve. In other words, the spectra become more PF like in nature the closer we move to flare onset, and thus easier for the model to correctly label. There are a number of reasons why long-range temporal tests like these are rare. First, there are a limited number of large flare observations with the IRIS slit positioned directly over the active region that is about to flare. Second, this set of observations contains only a small subset where IRIS was recording well in advance of flare onset. Third, the remaining long-range PF observations have to be split into training and test sets, which further diminishes this rarefied pool. Lastly, energetic events occurring sufficiently before a flare spoil certain observations. This is indeed the case for this particular flare, although IRIS was recording for a long time prior to flare onset, a C1.0-class flare spoils the possibility of extending the test further back in time. This may also account for the high metric scores observed even 80 minutes before flare onset.

For consistency, we performed a temporal test on a non-flaring active region. The results in the bottom panels of Figures 15 and 16 demonstrate how the model consistently labels the region as non-flaring, with metric scores below 0.5 indicating the prediction of AR profiles. ROC curves forming on the opposite side of the dividing random choice line indicate the prediction of non-flare profiles, with larger areas under the curve being associated with larger confidences. The scores deviate only slightly with time, in spite of visible small-scale energy releases.