Inferring Plasma Flows at Granular and Supergranular Scales With a New Architecture for the DeepVel Neural Network

The wealth of observational data available has been instrumental in investigating physical features relevant to solar granulation, supergranulation and Active Regions. Meanwhile, numerical models have attempted to bridge the gap between the physics of the solar interior and such observations. However, there are relevant physical quantities that can be modeled but that cannot be directly measured and must be inferred. For example, direct measurements of plasma motions at the photosphere are limited to the line-of-sight component. Methods have consequently been developed to infer the transverse plasma motions from continuum images in the case of the Quiet Sun and magnetograms in the case of Active Regions. Correlation-based tracking methods calculate the optical flows by correlating series of images locally while other methods like “Coherent Structure Tracking” or “Balltracking” exploit the coherency of photospheric granules to track them and use the group motions of the granules as a proxy of the average plasma flows advecting them. Recently, neural network computing has been used in conjunction with numerical models of the Sun to be able to recover the full velocity vector in photospheric plasma from continuum images. We experiment with a new architecture for the DeepVel neural network which takes inspiration from the U-Net architecture. Simulation data of the Quiet Sun and Active Regions are then used to evaluate the response at granular and supergranular scales of the aforementioned method.


INTRODUCTION
The Quiet Sun (i.e., in the absence of significant magnetic activity), hereafter QS, is filled with patterns of flows at multiple spatial and temporal scales. Granules, which typically have a diameter of 1 Mm and a lifespan of the order of 6 min, are associated with hot plasma upwellings whereas the intergranular lanes surrounding them, smaller in comparison, are associated with cold plasma sinking back in the interior. Supergranular flows are found at greater scale, i.e., above 20 Mm in diameter and with lifespans that range from hours to nearly 2 full days. This supercell-like pattern in the Quiet Sun is only revealed indirectly, for example by analyzing Dopplergrams and magnetograms or by tracking the average motion of granulation. It is believed that supergranulation originates from the deeper layers in the convection zone (see , and references therein). More advanced techniques can provide actual images of supergranular cells (Potts and Diver, 2008;Attie et al., 2009Attie et al., , 2016. Sunspots, which appear as dark spots in intensitygrams with diameters ranging from 1 to 50 Mm and lifespans of the order of days, are signatures of the magnetic activity of the Sun (i.e., Active Regions, hereafter ARs), as revealed by magnetograms. The strong magnetic field in sunspots inhibits the convective motions of the plasma. Centered around sunspots are moat flows which start from the penumbra and expand horizontally and radially outward (i.e., away from the sunspot). This flow was originally revealed through the analysis of spectroheliograms and moving magnetic features in magnetograms (Sheeley, 1969;Hagenaar, 2005). When the moat flow encounters the neighboring supergranular cells of the quiet sun, it forms a ring-shaped boundary around the sunspot that looks like a rather dynamic, yet coherent supergranule-like cell structure (Attie et al., 2018).
Dopplergrams measure only the line-of-sight component of the aforementioned flows (Welsch et al., 2013). Methods were thus developed to recover the missing transverse component of the velocity vector by either inversion of the velocity vector as is done here or by finding the electric fields instead and then applying the ideal MHD Ohm's law to find velocity (Kazachenko et al., 2014). Subsequent applications include the estimation of energy fluxes , the inferral of boundary conditions from which to drive simulations of the Sun  or the inferral of synthetic observations such as velocity field reconstructions for data assimilation in a model of the solar photosphere, as suggested in Abbett and Fisher (2010). Reconstruction algorithms for ARs solve the magnetic induction equation with the vertical velocity and magnetic field vector being assigned Dopplergrams and vector magnetograms (e.g., Longcope, 2004;Schuck, 2005Schuck, , 2006Schuck, , 2008. Intensitygrams are used to track optical flows in the Quiet Sun, i.e., the displacement that needs to be applied to one image subfield to recover the image subfield at the following timestep. Optical flows were shown to be highly correlated with actual plasma motions averaged over spatial scales of the order of a megameter or larger and timescales of the order of 30 min (Rieutord et al., 2001). Correlation-based tracking methods such as "Local Correlation Tracking" (LCT: November and Simon, 1988) and "Fourierbased Local Correlation Tracking" (FLCT: Fisher and Welsch, 2008) calculate optical flows by correlating series of images locally while other methods like "Coherent Structure Tracking" (CST: Rieutord et al., 2007) or "Balltracking" (Potts et al., 2004) exploit the coherency of photospheric granules to track them and use the group motions of the granules as a proxy of the average plasma flows advecting them. Recently, a deeplearning algorithm was trained with computations performed by a radiative magnetohydrodynamics (MHD) simulation of the solar photosphere to emulate the physics that relate the continuum intensity to the velocity vector that appears in the model equations (i.e., not an optical flow). DeepVel 1 (Asensio 1 DeepVel is an open-source neural network: https://github.com/aasensio/deepvel. Ramos et al., 2017) is a fully-convolutional neural network that infers instantaneous depth-dependent transverse plasma motions from pairs of intensitygrams (i.e., the same inputs as tracking methods). Outputs are dependent on the spatial resolution, cadence and physics of the model presented during training. Versions of the neural network have been trained for the QS (e.g., Asensio Ramos et al., 2017;Tremblay et al., 2018) and ARs (e.g., Tremblay et al., submitted) and have been used to generate synthetic SDO/HMI observations, i.e., estimates of the plasma motions that reflect a numerical simulation but appear as though they were derived from the Helioseismic Magnetic Imager (HMI: Schou et al., 2012) onboard the Solar Dynamics Observatory (SDO: Hoeksema et al., 2014).
A comparison between a sample of intensity-based methods as a function of spatial scales identified DeepVel as best capturing the physics of the Quiet Sun at granular scales whereas it appeared less effective at supergranular scales (Tremblay et al., 2018). In this paper, we use simulation data at the HMI instrument spatial resolution to test whether adapting the architecture of the DeepVel neural network to that of a U-Net (Ronneberger et al., 2015) could further improve reconstructions at supergranular scales in the QS and ARs. Additionally, the architecture is modified to accept a combination of intensitygrams, magnetograms, and Dopplergrams as input to account for the physics and spatial features they encompass and their impact on the inferred flows.

METHODOLOGY
U-nets are widely used for image segmentation, e.g., the segmentation of coronal holes in solar data (Illarionov and Tlatov, 2018). This neural network inherits its name from the shape of its architecture which features a contracting branch, a bottleneck and an expansive branch (Ronneberger et al., 2015). The contracting branch identifies the dominant features at a given spatial scale through convolutional layers and downsampling operations to halve the resolution. The number of channels is doubled at each level. The expansive branch mirrors the contracting branch, upsampling from the low spatial resolution output of the bottleneck to higher resolution to provide context. Skip-connections concatenate the outputs of the contracting branch with inputs of each level of the expansive branch to localize features (Ronneberger et al., 2015). Moreover, U-nets have the ability to train for specific spatial scales by freezing (i.e., stop training) the weights and biases of other layers.
Photospheric flows range from subgranular scales (< 1 Mm) to supergranular scales (> 10 Mm). We adapt the architecture of DeepVel to that of a U-net to probe spatial scales that range from the pixel-size to the size of the sub-images presented to the neural network during training. We refer to Asensio Ramos et al. (2017) for a detailed description of the original DeepVel architecture. Each level of the contracting branch features in succession a 2D convolutional layer with a kernel of 3 by 3 pixels 2 , batch normalization (Ioffe and Szegedy, 2015), a ReLU activation function and a Dropout layer of 50 % to avoid overtraining. Downsampling is then performed using 2D convolutional layers with strides of two pixels in all directions. The expanding branch mirrors the contracting branch, with the output of each level being concatenated with the output from the contracting branch at the same spatial scale. Upsampling in the expansion branch is performed by repeating each row and each column twice. Unlike the classic U-net architecture, the number of channels is maintained throughout all layers (i.e., 128 channels), with the exception of the bottleneck (i.e., 256 channels), to diminish the number of free parameters to adjust.
Versions of the U-net architecture (hereafter DeepVelU) are trained, validated and tested using the Tremblay et al. (2018) dataset for the Quiet Sun (hereafter the QS dataset) and the Tremblay et al. (submitted) dataset which features a mixture of ARs and QS (hereafter the ARs dataset). Training of the DeepVel and DeepVelU neural networks relies on a similar approach to that of measuring optical flows: two consecutive images (i.e., intensitygrams, Dopplergrams, and/or magnetograms) are presented to the network and their differences and dominant features are used to reconstruct a velocity vector. The ARs dataset is required in addition to the QS dataset as a neural network trained exclusively with examples of the QS will fail to extrapolate the flows in Sunspots whose physical mechanisms governing plasma motions in relation to the observed light intensity is very different from those in the QS Tremblay et al. (submitted). On the other hand, although it is presented a dataset that differs from the one used for training, the neural network is capable of generalizing the behavior for granulation (Tremblay et al., 2018, Tremblay et al., submited).
The QS dataset was derived from the STAGGER magnetoconvection simulation 2 of solar granulation (Stein, 2012;Stein and Nordlund, 2012). It features maps of the continuum intensity I c at 500 nm and the velocity vector v at optical depths τ = {1, 0.1, 0.01} with a field of view of dimensions 96.768 by 96.768 Mm 2 , spatial resolution x = y = 96 km pixel −1 and time step t = 60 s for a total duration of 6 h (i.e., many turnover times). The ARs dataset for I c , v and B was generated by the MURaM simulation of a Sunspot (Rempel and Cheung, 2014). The Sunspot has a diameter of ≈ 25 Mm and is featured at the center of a field of view of 98.304 by 98.304 Mm 2 . The spatial resolution of the dataset is the same as the QS dataset (i.e., x = y = 96 km pixel −1 ), but the cadence t is higher at 45 s which coincides with the cadence of SDO/HMI level-2 products (Hoeksema et al., 2014).
The simulation data was convolved with the SDO/HMI point spread function (PSF; Wachter et al., 2012) before being resampled to the SDO/HMI spatial resolution near disk center ( x = y ≈ 368 km pixel −1 ) using nearest-neighbor sampling. Patches of 48 by 48 pixels 2 were then extracted at random positions in planes of constant τ and at consecutive times Computations were performed on a NVIDIA-GTX-960 GPU using the Keras library with the Tensorflow backend. Weights and biases were only updated when the mean square error for the validation set improved.

RESULTS
Full field-of-view maps in sequences of 30 to 80 time steps from the QS and ARs datasets are used as test sets to evaluate the performance of DeepVelU at scales equal to and greater than those analyzed in Tremblay et al. (2018).

Metrics
Agreement between the simulation vector field [ v ref ] (i.e., the reference case) and the reconstructed vector field [ v D ] is quantified by the mean absolute errors and mean relative errors where · is the spatial average operator. Additionally, the Pearson linear correlation coefficient between v ref and v D , denoted [C], is introduced as a measure of similarity (or discrepancy) that takes into account different spatial scales of the vector amplitude. The averaged normalized dot product is used to assess the global orientation of the inferred velocity vectors with respect to the simulation, with values A = ±1 for parallel/anti-parallel vectors and A = 0 for perpendicular vectors. Response as function of the spatial scales is evaluated through the power spectrum of the total kinetic energy of the transverse plasma motions. Energy densities E(k) are computed following the definition of  for a square dataset (i.e., n x = n y and x = y): wherev x andv y are the discrete Fourier transform of v x and v y , and the wavenumber k is an inverse measure of the spatial scale. Finally, for each velocity field, we compute the unsigned "shear" component of the vertical Poynting flux [S z ] and then integrate it over the field-of-view (Liu and Schuck, 2012;Welsch, 2015): where B ref,t and B ref,z are the transverse component and the vertical component of the magnetic field, respectively, computed from the reference simulation. Note that the total Poynting flux also includes the "emergence" term, which contains Doppler velocity. However the division into "shearing" and "emergence" terms is done conceptually, since both terms involve emergence of magnetized plasma across the photosphere (Welsch, 2015).

Quiet Sun Test Set
Plasma motions v D,t (τ ≈ 1) inferred by DeepVelU are consistent with v ref,t , producing divergent velocity vectors at the center of granules and converging vectors in the intergranular lanes ( Figure 1A). Furthermore, the new architecture results in reduced errors for v D,t (τ ≈ 1), ( ∇ × v D,t ) z and ∇ t · v D,t in comparison to DeepVel ( Table 1). Both methods underestimate flow amplitudes, as suggested by the scatterplot in Figure 1B, but the similarity with the reference velocity amplitudes increases from C = 0.841 to 0.947 when transitioning from DeepVel to the DeepVelU architecture (Table 1). Furthermore, the RMSE, E abs and E rel at τ ≈ 1 decrease from 0.777 km s −1 , 0.694 km s −1 and 70.6% to 0.501 km s −1 , 0.442 km s −1 and 43.7 %, respectively. Similar improvements are noted for the divergence and the curl of the flow fields ( Table 1). The global orientation of the velocity vectors is also improved, with A increasing from 0.786 when using DeepVel to 0.914 when using DeepVelU. Further analysis identified regions of downflows as the largest sources of errors (not shown). These are typically associated with intergranular lanes which feature more complex flow structures confined in small areas. Similar conclusions are drawn for the inversion of flows at optical depths τ ≈ {0.1, 0.01} (i.e., at higher geometrical heights above the surface) from I c (τ ≈ 1) (Figure 1 and Table 1). More specifically, τ ≈ 0.1 is just beyond the height where the reversal of the granulation pattern occurs (i.e., a few hundred kilometers above the surface: Cheung et al., 2007), with the center of granules becoming colder than the intergranular network. At τ ≈ 0.01, this pattern is more diffused but remains well-correlated with the surface granulation pattern (not shown). By extension, the structures in v t (τ ≈ {0.1, 0.01}) correlate with I c (τ ≈ 1).
The power spectrum of the kinetic energy (Equation 4) as a function of the spatial scales shows an improved response for DeepVelU at supergranular scales (k −1 ≈ 100 Mm) with respect to the test set (green and blue curves in Figure 2A). The transition from a loss to a gain in signal with respect to DeepVel (orange curve) occurs close to the spatial scale that is achieved through downsampling at the bottleneck of the U-net architecture (k −1 ≈ 2.944 Mm). This change is interpreted as the new architecture favoring signal at supergranular scales over the pixel-size and granular scales to further optimize the cost function during the training process.
In addition, Figure 2C shows a significant improvement for DeepVelU in the correlation between the inferred velocities and the simulation at all time and spatial scales for the QS. DeepVel's flows have the surprising disadvantage of losing correlation at increasing time averages by a few percents, whereas the correlation increases for DeepVelU. Similarly, DeepVel's correlations significantly decreases at greater spatial averages, e.g., from (resp.) C ≈ 0.85 to C ≈ 0.62 between 368 km and 5 Mm (resp.), which is not the case for DeepVelU which consistently correlates very well with the simulation at C > 0.9. The latter, however, plateaus near 3 Mm which coincides with the pixel-size in the bottleneck of the U-net architecture. Future work will test if increasing the field-of-view of the training images affects favorably the ability to improve the correlation further at greater spatial scales.

Active Region Test Set
Figures 3A,C,E show a subset of velocity field inversions generated by DeepVelU from single-quantity inputs that relate to SDO/HMI level-2 data products i.e., intensitygrams, line-ofsight magnetograms,and Dopplergrams, respectively. Only a subfield of 50 by 50 pixels 2 is shown for clarity. The position of this patch in the field-of-view was selected to highlight distinct flow structures in the presence of ARs: the center of the sunspot (upper right corner), the flows in the penumbra, the moat-like flows around the latter (close to the diagonal connecting the upper left corner to the lower right corner), and finally, the granulation like in the QS dataset (lower left corner). Scatterplots in Figure 3 were computed over the sub-field whereas the metrics compiled in Table 1 were computed over the entire field-ofview and include all tested combinations of inputs. An arbitrary threshold of |B ref,z | ≥ 100 G is set to compute metrics that are specific to the AR and the magnetic field network. Regions where |B ref,z | < 100 G are interpreted as QS.
DeepVel and DeepVelU generate very similar velocity fields from consecutive intensitygrams, with DeepVel performing slightly better where |B ref,z | < 100 G and DeepVelU (Figure 3A) improving absolute and relative errors slightly where |B ref,z | ≥ 100 G (Table 1) but underestimating low amplitude velocities ( Figure 3B). In fact, the performances of the two neural networks are comparable to that of DeepVel with the QS dataset ( Table 1). For both architectures, the signal for the continuum intensity inside the Sunspot results in less effective flow inversions where |B ref,z | ≥ 100 G. More specifically, the mean relative error almost doubles when transitioning from the QS to strong field regions and the metric A describing the global orientation of the vectors decreases to a value of about 0.55 (Table 1). For this reason, different combinations of physical inputs were tested to measure their impact on photospheric flows. Although these tests were performed using the DeepVelU architecture, we expect that similar conclusions can be drawn for the DeepVel neural network.
The use of line-of-sight magnetograms as inputs instead of intensitygrams slightly improves the errors where |B ref,z | ≥ 100 G (Table 1), i.e., where there is the most signal in the input data. Despite weak magnetic fields in the QS, the neural network is capable of distinguishing individual granules and their flow patterns (Figure 3C), but the errors are significantly larger with greatly underestimated flow amplitudes ( Figure 3D) and misaligned vectors ( Table 1). The overall performance is thus worse due to a larger fraction of the field-of-view being covered by QS (Table 1).
Dopplergrams, such as the ones depicted in the background of Figures 3A,C,E, are the inputs that best capture the behavior   of flows in both QS and ARs (Table 1 and Figure 3F). The granulation pattern from intensitygrams is clearly outlined in Dopplergrams by the cold sinking plasma in the intergranular lanes and the hot rising plasma at the center of granules. Meanwhile convective motions are inhibited inside Sunspots. However, flows in the penumbra are predicted by DeepVelU to be almost purely radial with respect to the center of the Sunspot (Figure 3E), whereas the simulation penumbral flows resemble more closely those seen in Figures 3A,C. Although the metric A was slightly improved in strong field regions, vector orientations remain much less accurately reproduced than in the QS (Table 1). Combining B ref,z and v ref,z , which yielded the best reconstructions for ARs and QS, respectively, significantly improves the performance for ARs and the integrated Poynting flux ( Table 1). Magnetograms provide to Dopplergrams additional signal inside Sunspots. Both quantities are coupled physically to v t through the magnetic induction equation, which could further explain the increase in performance. The addition of intensitygrams provides more context for the QS and improves the metrics where |B ref,z | < 100 G, with very little to no decrease in performance where |B ref,z | ≥ 100 G ( Table 1).
The use of Dopplergrams or combinations of inputs in DeepVelU improves the response in the power spectrum of the kinetic energy for spatial scales larger than k −1 ≈ 3 Mm (Figure 2B), which could again be related to the spatial scales probed in the U-net's bottleneck. Although the power is generally underestimated, its variations as a function of k are matched more consistently by DeepVelU, with DeepVel generating more signal than the simulation at supergranular scales.

CONCLUSION
We trained DeepVelU, a U-net-inspired architecture for the DeepVel neural network, using simulations of the QS and ARs and evaluated the method's response as a function of spatial and temporal scales. DeepVelU shows significant improvement over DeepVel for the QS test set. The correlations for the latter falls close to 0.6 at spatial scales of 5 Mm whereas it stays consistently above 0.9 for DeepVelU; increasing at greater time averages and plateauing at spatial scales above 3 Mm. Thus DeepVelU's QS model appears more effective than the other tracking methods tested in Tremblay et al. (2018) against the same dataset, with increased correlations and lower errors being achieved over DeepVel at granular scales and over FLCT at supergranular scales. The results for the ARs dataset are not as conclusive, but may be further improved by training for granulation, penumbra and sunspots separately. This approach could, however, introduce discontinuities at the edges of the different structures. Frontiers in Astronomy and Space Sciences | www.frontiersin.org Current efforts are meant to be a proof of concept. Limitation of the method include the input image dimensions in each direction which must be a factor of 2 n where n is the number of downsampling or upsampling layers in the network architecture (here n = 3). For example, the dimensions of the sub-images (48×48 px 2 ) presented to the network during training limited the number of downsampling/upsampling layers in the architecture to n = 3, corresponding to a spatial of about 3 Mm. This could explain the plateauing of the correlation above this spatial scale. In this context, development is underway for a deeper version of DeepVelU (i.e., with more downsampling/upsampling layers) that will be trained on (almost) full field-of-view images (256 × 256 px 2 or 94.192 × 94.192 Mm 2 ) of the QS dataset and that will also double (resp. halve) the number of filters after each downsampling (resp. upsampling) operation. In addition, the simulations do not model actual supergranulation which is known to advect granules over spatial scales greater than 3 Mm. Therefore it will be worth exploring further inferences on actual photospheric observations. For example we plan to compare the inferred supergranular (QS) and moat flow patterns (ARs) with those of Attie et al. (2018) which uses a new implementation of the "Balltracking" method that is more accurate than the one of Potts et al. (2004) and that had not been tested by Tremblay et al. (2018). Similarly, the MURaM simulation includes a deep seated flow system with velocities in the 200-500 m/s range which extends about 10 Mm past the sunspot boundary, which seems in line with many observations of moat flows (Rempel, 2015). Using the data of both experiments (observations and simulation), we will compare the moat flow patterns revealed by DeepVel, DeepVelU and Balltracking.
The inclusion of Dopplergrams and line-of-sight magnetograms as inputs in the neural network architecture, both of which provide more signal than intensitygrams inside sunspots and are tied to the transverse plasma motions through the magnetic induction equation, have improved reconstructions in the context of ARs and response at supergranular scales. Evaluating the impact of additional inputs such as the transverse magnetic field vector B t which also appears in the magnetic induction equation or its strength | B t | which is not subjected to ambiguities is left as future work.
Despite the changes in architecture and the success of the neural network at capturing the spatial distribution of flows, velocity amplitudes are generally underestimated in the QS and overestimated in ARs by both DeepVel and DeepVelU. Future efforts will be dedicated to improving the inferral of amplitudes.
All versions of the neural networks were trained to generate synthetic data that is consistent with a given simulation and the SDO/HMI cadence (if the velocity is multiplied by a factor of 45/60 for the QS dataset) and spatial resolution near disk center. SDO/HMI level-2 products can thus be used as input, however one should first assess if the physics and preprocessing of the training set is consistent with the observations presented as input. For instance, SDO/HMI Dopplergrams measure the superposition of the line-of-sight components of the satellite motion, meridional flows, differential rotation, p-mode oscillations, and plasma motions. Simulation Dopplergrams only feature the latter two. Additional preprocessing steps are thus required to correct SDO/HMI Dopplergrams for the aforementioned effects (e.g., Welsch et al., 2013) or to project the simulation data in the observations space prior to training. The resulting velocity fields may then serve as synthetic observations or first estimates when performing data assimilation in an MHD model of the photosphere, or as boundary conditions driving a simulation. The method may also be used to estimate and evolve a Poynting flux vector that is representative of a given epoch of the Sun.
The velocity vector that the neural networks are trying to recover is the same vector as physics-based velocity inversion methods, i.e., v t such that the magnetic induction equation is satisfied. Furthermore, DeepVelU best performed when using Dopplergrams and magnetograms as input, with both quantities appearing in the magnetic induction equation alongside v. The training process may be revisited in the future to incorporate more effectively the physics when estimating the plasma motions, e.g., through the loss function or a physics-informed network (Raissi et al., 2019).

DATA AVAILABILITY STATEMENT
The QS and ARs weights and biases as well as the DeepVel and DeepVelU outputs generated and analyzed for this study can be found in the DeepVelU_Frontiers repository (https://github. com/tremblaybenoit/DeepVelU_Frontiers).

AUTHOR CONTRIBUTIONS
BT generated the QS and ARs datasets from simulation data. BT also trained and tested all versions of the DeepVel and DeepVelU neural networks. RA acquired the ARs dataset prior to BT downsampling the data to the SDO/HMI spatial resolution. RA also provided the analysis of the QS dataset.