Radial velocities in the outermost disk toward the anticenter

TL;DR

This study measures stellar radial velocities in the outer Galactic disk toward the anticenter, revealing complex velocity patterns likely influenced by spiral arm overdensities rather than other dynamical factors.

Contribution

It provides new measurements of radial velocities in the outer disk and suggests spiral arm overdensities as a key factor influencing these motions.

Findings

Positive radial velocities at 9-13 kpc reaching ~6 km/s

Negative velocities beyond 17 kpc reaching ~-10 km/s

Velocity patterns consistent with spiral arm overdensities

Abstract

We measure the mean Galactocentric radial component of the velocity of stars ($v_R$) in the disk at 8 kpc$<R<28$ kpc in the direction of the anticenter. For this, we use the Apache Point Galactic Evolution Experiment (APOGEE). Furthermore, we compare the result with HI maps along the same line of sight. We find an increase in positive (expansion) $v_R$ at $R\approx 9-13$ kpc, reaching a maximum of $\approx 6$ km/s, and a decrease at large values of $R$ reaching a negative (contraction) value of $\approx -10$ km/s for $R>17$ kpc. Negative velocities are also observed in 21 cm HI maps, possibly dominated by local gas emission. Among the possible dynamical causes for these non-zero $v_R$, factors such as the effect of the Galactic bar, streams, or mergers do not seem appropriate to explain our observations. An explanation might be the gravitational attraction of overdensities in a spiral arm. As a matter of fact, we see a change of regime from positive to negative velocities around $R\approx 15$ kpc, in the position where we cross the Outer spiral arm in the anticenter. The mass in spiral arms necessary to produce these velocities would be about 3\% of the mass of the disk, consistent with our knowledge of the spiral arms. Another scenario that we explore is a simple class of out-of-equilibrium systems in which radial motions are generally created by the monolithic collapse of isolated self-gravitating overdensities.

Figures

Tables

Table 1. Table 1: Average radial Galactocentric velocity derived from Eq. ( 1 ) with radial heliocentric velocities from APOGEE within a region close to the Galactic Center–Sun line, including All of the stars or only with stars with [Fe/H] > > -1. Bins have N ≥ 6 𝑁 6 N\geq 6 and are separated by Δ ​ log ⁡ R ​ ( kpc ) = 0.03 Δ 𝑅 kpc 0.03 \Delta\log R({\rm kpc})=0.03 . The corresponding numbers of stars in each bin are indicated in the columns 4, 5; the r.m.s. of the radial velocities is given in columns 6, 7.

$R$ (kpc)	${⟨v_R⟩}_{\mathrm{all}}$ (km/s)	${⟨v_R⟩}_{\left[\mathrm{Fe}/H\right]>-1}$ (km/s)	$N_{\mathrm{all}}$	$N_{\left[\mathrm{Fe}/H\right]>-1}$	${\sigma }_{\mathrm{all}}$ (km/s)	${\sigma }_{\left[\mathrm{Fe}/H\right]>-1}$ (km/s)
7.9	-4.0$\pm$ 9.5(stat.)$\pm$ 1.0(syst.)	-4.0$\pm$ 9.5(stat.)$\pm$ 1.0(syst.)	16	16	38.1	39.0
8.5	-2.4$\pm$ 2.1(stat.)$\pm$ 1.0(syst.)	-2.4$\pm$ 2.1(stat.)$\pm$ 1.0(syst.)	289	289	36.5	36.5
9.1	-3.4$\pm$ 1.0(stat.)$\pm$ 1.0(syst.)	-3.5$\pm$ 1.0(stat.)$\pm$ 1.0(syst.)	1246	1244	33.9	34.0
9.8	0.4$\pm$ 0.8(stat.)$\pm$ 1.0(syst.)	0.6$\pm$ 0.8(stat.)$\pm$ 1.0(syst.)	2011	2007	34.6	34.7
10.5	2.1$\pm$ 0.7(stat.)$\pm$ 1.0(syst.)	2.2$\pm$ 0.7(stat.)$\pm$ 1.0(syst.)	2496	2493	32.7	33.0
11.2	3.3$\pm$ 0.6(stat.)$\pm$ 1.0(syst.)	3.3$\pm$ 0.6(stat.)$\pm$ 1.0(syst.)	2979	2979	31.4	31.5
12.0	4.7$\pm$ 0.6(stat.)$\pm$ 1.0(syst.)	4.7$\pm$ 0.6(stat.)$\pm$ 1.0(syst.)	3000	2998	30.5	30.5
12.9	6.1$\pm$ 0.6(stat.)$\pm$ 1.0(syst.)	6.1$\pm$ 0.6(stat.)$\pm$ 1.0(syst.)	2022	2021	27.1	27.1
13.8	5.6$\pm$ 0.9(stat.)$\pm$ 1.0(syst.)	5.4$\pm$ 0.9(stat.)$\pm$ 1.0(syst.)	1099	1094	30.2	28.9
14.8	2.3$\pm$ 1.5(stat.)$\pm$ 1.0(syst.)	1.6$\pm$ 1.4(stat.)$\pm$ 1.0(syst.)	451	449	32.6	30.7
15.8	-5.8$\pm$ 2.3(stat.)$\pm$ 1.1(syst.)	-6.6$\pm$ 2.1(stat.)$\pm$ 1.0(syst.)	196	191	32.1	30.0
17.0	-12.7$\pm$ 2.6(stat.)$\pm$ 1.0(syst.)	-12.0$\pm$ 2.6(stat.)$\pm$ 1.1(syst.)	109	106	27.6	26.6
18.2	-13.7$\pm$ 2.5(stat.)$\pm$ 1.8(syst.)	-13.0$\pm$ 2.3(stat.)$\pm$ 1.8(syst.)	87	85	23.4	21.3
19.5	-9.7$\pm$ 3.6(stat.)$\pm$ 1.2(syst.)	-8.8$\pm$ 3.6(stat.)$\pm$ 1.2(syst.)	54	53	26.4	26.7
20.9	-0.2$\pm$ 9.2(stat.)$\pm$ 2.2(syst.)	-12.9$\pm$ 5.0(stat.)$\pm$ 1.9(syst.)	32	28	52.2	26.5
22.4	-25.5$\pm$ 13.6(stat.)$\pm$ 1.6(syst.)	-13.0$\pm$ 8.9(stat.)$\pm$ 2.2(syst.)	9	8	40.8	24.8
23.9	–	–	0	0	–	–
25.7	–	–	1	1	–	–
27.5	-22.2$\pm$ 20.7(stat.)$\pm$ 2.5(syst.)	–	6	3	50.7	–

Equations

$$v_R=-\frac{v_r}{\cos (\varphi +\ell )\cos b}-\frac{\cos \ell }{\cos (\varphi +\ell )}{U}_{\odot }-\frac{\sin \ell }{\cos (\varphi +\ell )}{V}_{g,\odot }$$

$$+\frac{\tan b}{\cos (\varphi +\ell )}(W-W_{\odot })+\tan (\varphi +\ell )V_{\varphi },$$

$$V_{a,0}(R)=\frac{f_B(R){\sigma }_R(R{)}^2}{V_c(R)},$$

$$T(v)\propto {\int }_{R_{\odot }}^{\infty }{\rho }_{HI}(R)\exp \left(-\frac{(v-v_R(R){)}^2}{2{\sigma }_R(R{)}^2}\right)dR,$$

$$P(e)de\approx \frac{2^{1/2}de}{{\pi }^{3/2}\sigma V_{\varphi }{e}^2}{\int }_{-V_{\varphi }e}^{V_{\varphi }e}\frac{|x|}{\sqrt{1-\frac{x^2}{V_{\varphi }^2{e}^2}}}exp\left[-\frac{(x-v_R{)}^2}{2{\sigma }^2}\right]dx$$

$$\frac{\dot{M}}{M}(R)=\frac{-1}{R\sigma (R)}\frac{d[⟨v_R⟩(R)R\Sigma (R)]}{dR},$$

$$\frac{\dot{M}}{M}(R)=⟨v_R⟩(R)\left(\frac{1}{h_r}-\frac{1}{R}\right)-\frac{d⟨v_R⟩(R)}{dR}.$$

$$\Delta a_R=2G\sigma {\int }_{R_{\mathrm{sp}}-{\Delta }_{\mathrm{sp}}/2}^{R_{\mathrm{sp}}+{\Delta }_{\mathrm{sp}}/2}S\left[\frac{E(m)}{R(S-R)}-\frac{K(m)}{R(S+R)}\right]dS,$$

$$m=\frac{4RS}{(R+S{)}^2};\sigma =\frac{M_{\mathrm{sp}}}{2\pi R_{\mathrm{sp}}{\Delta }_{\mathrm{sp}}},$$

$$\text{mbox}Max[v_R]\sim 100\sqrt{\frac{M_{\mathrm{sp}}}{2\times 10^{11}{\mathrm{M}}_{\odot }}}\text{mbox}km/s,$$

$$v_r=\frac{\vec{v}\cdot \vec{r}}{|\vec{r}|}$$

$${\vec{v}}_t(r)=\frac{\vec{r}\times \vec{v}(r)}{|\vec{r}|},$$

$$v_a=v_c-\overline{v_{\varphi }}$$

$$v_c^2(r)=\frac{GM(r)}{r}.$$

$$\beta (r)=1-\frac{⟨v_t^2(r)⟩}{2⟨v_r^2(r)⟩},$$

Full text

Galaxy: kinematics and dynamics – Galaxy: disk

††journal: AJ

1 Introduction

Disks in spiral galaxies are usually considered to be in a quasi-equilibrium configuration, where stars and other emitters move on steady circular orbits around the center of the galaxy with almost constant velocity. However, most disks exhibit a wealth of non-axisymmetric structures (Rix & Zaritsky 1995; Laine et al. 2014); about one third of them are substantially lopsided at $>$2.5 times the scale length of the disk, although the spiral pattern couples significantly to the estimate of the intrinsic ellipticity and their measurement may represent an upper limit to the true potential triaxiality. Lopsidedness is quite typical in disk galaxies and this may be interpreted as a pattern of elliptical orbits (Baldwin et al. 1980; Song et al. 1983). Non-circular streaming motions have also been observed in the gas motions of other galaxies (Trachternach et al. 2008; Sellwood & Zanmar Sánchez 2010), but the stellar kinematics is usually more regular and symmetric than the gas kinematics (Pizzella et al. 2008). A new analysis using tracers of the stellar populations should be investigated. This can be more accurately done in the Milky Way, which is the goal of the present paper.

In our Galaxy, Siebert et al. (2011) and Williams et al. (2013) measured a significant average gradient in Galactocentric radial velocity ($\langle v_{R}\rangle$) outwards within 6 kpc$\lesssim R\lesssim 9$ kpc. Carrillo et al. (2018) reached $R=10$ kpc with Gaia-DR1 data. López-Corredoira & González-Fernández (2016; hereafter LG16) and Tian et al. (2017) measured the mean radial Galactocentric velocity in the anticenter direction, reaching distances up to $R=16$ kpc using red clump giants as standard candles, and they obtained significant positive values of $\lesssim+10$ km/s. Wang et al. (2018b) confirmed these results up to $R\approx 13$ kpc with K giant stars. Gaia Collaboration (2018, hereafter G18) provided a full kinematic information in all lines of sight for different populations with the recent release of Gaia-DR2 data that include parallaxes, proper motions, and radial velocities, but also are limited to $R\lesssim 13$ kpc; this range was extended up to $R\approx 20$ by López-Corredoira & Sylos Labini (2018, hereafter LS18) with a method to deconvolve parallax errors that allows the use of data with large errors in distance. G18 or Antoja et al. (2018) or LS18 could show that the Galactic disk as a time-independent axisymmetric component is definitively an incorrect image of our Galaxy. Harris et al. (2018) took A and F stars as standard candles and also tried to measure the radial velocities up to distances of 15 kpc, but the uncertainties of their average values of several km/s are too large to measure small values of $\langle v_{R}\rangle$. Using star counts that fit the different components of the Galaxy (López-Corredoira & Molgó 2014) or doing analyses of metallicity distributions in spectroscopic surveys (López-Corredoira et al. 2018), it was known that the outer disk extends up to Galactocentric distances much farther away than 20 kpc, so one might explore the kinematics in this very outer region, possibly dominated by non-equilibrium dynamics. In this paper, we aim to obtain measurements of the radial velocity up to these greater Galactocentric distances, using the stellar spectrophotometric parallaxes derived by several methods by the “Apache Point Galactic Evolution Experiment” (APOGEE). We give information about the data source and distance derivation in §2. The measurements of the mean radial velocity of these stars in the anticenter direction are given in §3, where we also include the comparison of HI gas in that direction.

These measurements may be interpreted in terms of ellipticity, net migration of stars, associated with causes such as streams, bars, spiral arms, mergers, or non-equilibrium states derived from the particular evolution of isolated disks. All of these kinematical and dynamical analyses are carried out in §4. Particular attention is paid to a scenario that has been proposed in a recent work by Benhaiem et al. (2017; hereafter B17), where it was shown that a far-from-equilibrium dynamics may give rise to a disk with transient spiral-like arms in a manner completely different to that envisaged by the perturbative mechanisms for generating such structure: the transient nature of the outer part of the disk and of the arms is related to the presence of coherent radial velocities, with both negative (i.e., contraction) and positive (i.e., expansion) signs. On the other hand, we show that galaxy-like structures formed by magnetohydrodynamical simulations (e.g., Grand et al. 2018) typically form axisymmetric structures without coherent radial velocities. Finally, we summarize our results and interpretation in a concluding section (§5).

2 Data

APOGEE (Majewski et al. 2017) is one of the surveys included in the Sloan Digital Sky Survey (SDSS) III (Eisenstein et al. 2011) and IV (Blanton et al. 2017). APOGEE employs a pair of custom-built high-resolution (R$\equiv\lambda/\delta\lambda\sim 22,500$) multi-object $H-$band (1.5–1.7 $\mu$m) spectrographs (Wilson et al. 2010, 2012) on 2.5-metre telescopes with a 3∘ FOV: the SDSS telescope at Apache Point Observatory (New Mexico, US; see Gunn et al. 2006) and the du Pont telescope at Las Campanas (Chile). The project is collecting hundreds of thousands of stellar spectra with high signal-to-noise ratios across the Milky Way, focusing on the regions where dust causes dramatic extinction at optical wavelengths, namely the disk and the central parts of the Galaxy. The APOGEE observations are providing a chemical map of our Galaxy with unprecedented quality (Badenes et al. 2018; Fernández-Alvar et al. 2018; Fu et al. 2018; García-Dias et al. 2018; Garcııa-Pérez et al. 2018; Hayes et al. 2018; López-Corredoira et al. 2018; Mackereth et al. 2019; Palicio et al. 2018; Souto et al. 2018; Weinberg et al. 2018; Wilson et al. 2018).

The 14th data release of the SDSS (DR14; Abolfathi et al. 2018; Holtzman et al. 2018; Jönsson et al. 2018) includes millions of APOGEE spectra for approximately 263,000 stars observed from Apache Point Observatory, which have been processed by the project’s pipelines (Nidever et al. 2015; García-Pérez et al. 2016) to derive radial velocities, atmospheric parameters, and abundances for more than 20 chemical elements, typically with a precision better than 0.1 dex. Distances to the stars have been estimated by four different methods and included in a value-added catalog released in conjunction with DR14: the StarHorse code (Santiago et al. 2016), the code described in Wang et al. (2016), the isochrone-matching technique described in Schultheis et al. (2014), and the distance code of Holtzman et al. (2018). The agreement among the four codes is fair, typically within 20%. We have adopted for our analysis the average values of the available estimates, so the error of this average distance is expected to be $\sim$10%.

The spatial distribution of the stars considered here spans a range of Galactic longitudes accessible from observatories in Northern-hemisphere. The range of Galactocentric distances explored in this paper is between $R=6$ and 28 kpc, always at $|z|<5$ kpc. In Fig. 1 we show the distribution in cartesian Galactocentric coordinates ($X_{\odot}=R_{\odot}$) of the 16 105 sources in the anticenter direction that we will use in this paper: a decreasing density of sources with increasing $X$ is evident. The number of sources for each value of $R$ is given in Table LABEL:Tab:anticenter.

The stars in this range are thin- and thick-disk stars; the disk reaches a distance of 25 kpc from the Galactic center or even farther (López-Corredoira et al. 2018), with some significant contamination by halo stars only for $R\gtrsim 20$ kpc (López-Corredoira & Molgó 2014; López-Corredoira et al. 2018). Although some authors attributed some excess stars in the plane at $R\approx 20$ kpc (the Monoceros Ring) to a putative tidal debris of a dwarf galaxy (Sollima et al. 2011), or a bending of the disk or a corrugation pattern linked to the interaction with a passing satellite (Gómez et al. 2016; Yanny & Newberg 2016), López-Corredoira & Molgó (2014) and Wang et al. (2018a) have shown that a flared disk+halo is the most likely explanation without the need for any extragalactic element. We assume this last position here.

3 Galactocentric radial velocities

3.1 Transformation of Heliocentric into Galactocentric Radial Velocities

The three components of the velocity derived from proper motions and radial velocities were already analyzed for Gaia data (G18; LS18). Here, however, we will only use the radial velocities. APOGEE provides independent information with a determination of distances different from parallax in Gaia that is more precise for the largest heliocentric distances, so we can reach greater distances with APOGEE than with Gaia. Moreover, given that our analyses will be restricted to the anticenter, proper motions make a negligible contribution to the radial component of the Galactocentric velocity. With only the radial heliocentric velocity $v_{r}$, we can obtain the Galactocentric radial velocity ($v_{R}$) through (LG16, Eq. 4):

[TABLE]

[TABLE]

where $(\ell,b)$ are the Galactic coordinates, $\phi$ is the Galactocentric azimuth (such that $x=R\cos\phi$, $y=R\sin\phi$, $\phi_{\odot}=0$); ($U_{\odot}$, $V_{\odot}$, $W_{\odot}$) is the velocity of the Sun with respect to the local standard of rest (LSR), $V_{g,\odot}=V_{\phi}(R_{\odot},z=0)+V_{\odot}$; $V_{\phi}(R,z)$ is the azimuthal velocity, and $W(R,z)$ is the vertical motion. In the calculation of the average or median of $v_{R}$, we neglect the possible asymmetries of random peculiar motions.

The disadvantage of Eq. (1) is that it is model-dependent, since we need to know the values of $V_{\phi}$ and $W$, but we can make an appropriate selection of regions in which the dependence on them is very small. In particular, as we will see next, toward the anticenter, the dependence on these two functions is negligible.

Here, we adopt the values $R_{\odot}=8$ kpc; $V_{g,\odot}=244\pm 10$ km/s, $U_{\odot}=10\pm 1$ km/s (Bovy et al. 2012), and $W_{\odot}=7.2\pm 0.4$ km/s (Schönrich et al. 2010). The azimuthal velocity $V_{\phi}(R,z=0)=V_{c,0}(R)-V_{a,0}(R)$ where $V_{c,0}$ is the rotation curve, taken from Sofue (2017, Fig. 4), and $V_{a,0}$ is the asymmetric drift at $z=0$. The contribution from the asymmetric drift $V_{a,0}$ is negligible, especially for the anticenter region; we include its small contribution as

[TABLE]

with $f_{B}(R)\approx 0.27\,R({\rm kpc})$ (Bovy et al. 2012, Fig. 4/top). $V_{a,0}$ has a value $\sim 10$ km/s for $R\lesssim 15$ kpc and somewhat larger for larger $R$. The generalized expression with dependence on $z$ that takes into account the lower rotation speed for $z\neq 0$ is taken as $V_{\phi}(R,z)=V_{\phi}(R,z=0)-19.2$ km/s$\left(\frac{|z|}{\rm kpc}\right)^{1.25}$ (Bond et al. 2010). As systematic error we take $\Delta V_{\phi}(R,z)=\Delta V_{g,\odot}=10$ km/s for $R<20$ kpc and $\Delta V_{\phi}(R,z)=(R({\rm kpc})-10)$ km/s for $R\geq 20$ kpc, in agreement with the errors given in Sofue (2017, Fig. 4). We set $W=0$, considering negligible the contribution of possible vertical motions; but we take into account the range of possible values in the calculation of systematic errors using $\Delta W\equiv W_{\rm best\ fit}$ of Eq. (15) of López-Corredoira et al. (2014).

3.2 Analysis in the Anticenter In-plane Region

First, we choose the lines of sight symmetrically around the anticenter, which minimizes the systematic errors. We choose $|b|<10^{\circ}$, which makes the contribution of $W$ or $W_{\odot}$ very low; and $|\ell-180^{\circ}|<20^{\circ}$, to make $\cos(\phi+\ell)$ close to one for all of the stars and $\langle\tan(\phi+\ell)\rangle$ close to zero, which reduces the impact of the error in $V_{\phi}$. The average Galactocentric velocities are shown in Fig. 2 for $R-R_{\odot}<20$ kpc, together with other results of previous studies for $R<18$ kpc in fair agreement with our data. Note that the selected regions in the different surveys are not exactly the same, although all are centered in the anticenter region, so no perfect agreement in the results is expected. The median of $V_{R}(R)$ is also plotted for the first sample, with results very similar to the average; thus, we will use the “average” hereafter. In Table 2 we give the numerical values. Bins are divided into $\Delta\log R({\rm kpc})=0.03$, so as to have $N\geq 6$ stars in each of them. We plot the average velocity for all the stars, and separately for those at [Fe/H]$>-1$. The second subsample, which is restricted in metallicity, is almost pure disk stars with a very small contamination of the halo population (López-Corredoira et al. 2018). Both of them show compatible results for $\langle V_{R}\rangle$ because, even for the first sample, disk stars dominate. The outcome is an increase in positive velocity (expansion of the disk) from $R\approx 9$ kpc to $R\approx 13$ kpc, reaching a maximum of $\approx 6$ km/s, and a decrease outwards, reaching $\approx-10$ kpc for $R>17$ kpc (negative means contraction of the disk). The signal is significantly larger than the statistical errors (dispersion in the measurements of $v_{r}$) and systematic errors (due to the uncertainties in the parameters used in Eq. (1)).

If we divide the sample with [Fe/H]$>-1$ into two subsamples corresponding to the northern and southern Galactic hemispheres, the results show some differences (Fig. 3) but within similar trends. The systematic errors are larger due to the $\Delta W$ term when $\langle z\rangle$ is significantly different from zero. If the division is carried out into $\ell<180^{\circ}$ and $\ell>180^{\circ}$, however, the results are significantly different (Fig. 4), indicating also that there is an important dependence of the radial velocity on the azimuth.

The errors in the determination of heliocentric distance for any given star (around 10% as stated; Fig. 5 shows the average error as a function of Galactocentric distance for this line of sight) might produce a smoothing in the function $\langle v_{R}(R)\rangle$, but, given that the width of the bins is of the order of these errors, the variation of the function should not be important. The positive or negative signal cannot be due to errors in distance estimates.

3.3 Radial velocity from the 21 cm line in the anticenter

We wonder whether this feature of non-zero Galactocentric radial velocity can also be observed in neutral hydrogen. The anticenter is a region that is usually avoided for the analysis of the 21 cm line, because distances should be derived from heliocentric velocities assuming circular orbits, and thus one would expect zero velocity (except for the corrections for the motion of the Sun with respect to the LSR) at any distance, so there is no way to disentangle the integrated flux. Nonetheless, there might be some net redshift/blueshift or some distortion that could make us suspect that the orbits are not perfectly circular in the gas in the outer disk.

A simple exercise can illustrate this. Assuming that gas and stars have a similar mean velocity, the integrated flux of the antenna temperature in HI in the anticenter should be (corrected to Galactocentric coordinates)

[TABLE]

where $\rho_{HI}(R)$ is the HI density at distance $R$ and $\sigma_{R}$ is the dispersion of $v_{R}$. We assume a constant velocity dispersion for the HI gas. The observational evidence for a constant dispersion is indirect. Dickey et al. (2009) have shown that in the outer part of the Milky Way up to the observational limit at radii of $\sim$25 kpc the ratio of the emission to the absorption, which determines the mean spin temperature of the gas, stays nearly constant. This implies on average a constant mixture between cold and warm neutral medium, hence a constant velocity dispersion.

The state of the HI is determined by gas-phase reactions, in equilibrium situations with isotropic velocity dispersions for individual gas clouds. The velocity dispersion of the gas is not related to the dynamical state of the stellar population, but depends on the star formation rate (SFR). Dib et al. (2006) investigated the relationship between the velocity dispersion of the gas, the supernova rate, and feedback efficiency with three-dimensional numerical simulations of supernovae-driven turbulence in the interstellar medium and explained the approximate constancy of observed velocity dispersion profiles in the outer parts of galactic disks. Thermal instabilities in the HI gas, driven by supernovae, determine spin temperatures, composition, and velocity dispersion of the HI.

Tamburro et al. (2009) studied a sample of galaxies with a resolved systematic radial decline in the HI line width, implying a radial decline in kinetic energy density of HI. They find a correlation between the kinetic energy of gas and SFR with a slope close to unity and a similar proportionality constant in all objects.

Observations of velocity dispersions in external galaxies are faced with difficulties in disentangling the radial HI surface density, the rotation curve, and the HI velocity dispersion. O’Brien et al. (2010) observed the HI in edge-on galaxies, where it is in principle possible to measure the force field in the halo vertically and radially from gas layer flaring and rotation curve decomposition respectively. They found (their Fig. 21) radial fluctuations of the velocity dispersions, but on average the dispersions are constant in the outer parts of the disks. In all cases a distinct flaring is observed (O’Brien et al. 2010, Fig. 25). In fact, for a constant velocity dispersion the gradient of the force field in the vertical direction implies exponential flaring for the gas layer. Such a dependence is also observed in the Milky Way (Kalberla et al. 2007).

In Fig. 6 we plot the expected profile of the HI line in observations toward the anticenter for $\sigma_{R}(R)=12$ km/s, $\rho_{HI}(R)\propto\exp(-R/h_{R,HI})$ with a gas scale length of $h_{R,HI}=3.75$ kpc (Kalberla & Dedes 2008), and mean radial velocities null or equal to the observed APOGEE data for all sources (second column in Table LABEL:Tab:anticenter). The result shows a small asymmetry with respect to perfect circular orbits with slightly larger dispersion for negative velocities: $T(-v)>T(v)$. The effect is not very large, but it might be significant. Note, however, that, with these parameters, 40% of the temperature comes from gas at $R<10$ kpc and 80% of the temperature from gas at $R<15$ kpc, so the outermost disk features are not dominant in the maps of gas flux.

As a matter of fact, the results of Liu (2008, Fig. 6) or those from the EBHIS survey (Winkel et al. 2016)111Data provided at https://www.astro.uni-bonn.de/hisurvey/AllSky$\_$profiles/index.php. See also HI4PI Collaboration (2016)., reproduced in our Fig. 6, show an even larger asymmetry, which would confirm this trend; the larger asymmetry might be due to a $h_{R,HI}>3.75$ kpc or to higher amplitude of negative $v_{R}$ for $R<10$ kpc (which is affected by very large errors in our APOGEE sample) or because the stars in our region ($160^{\circ}<\ell<200^{\circ}$, $|b|<10^{\circ}$) have on average lower radial velocities than stars in the region of the anticenter alone, or simply because gas presents higher asymmetries than the stars in the velocity distribution, as observed in other galaxies (Pizzella et al. 2008).

The velocity distribution of the 21 cm line emission is usually affected by turbulent motions (Kalberla & Kerp 2009). In Fig. 6, we compare observations with significantly different FWHM beam sizes of $\sim 6^{\circ}$ and $\sim 0.2^{\circ}$. It appears implausible that the asymmetries in both cases are caused by random turbulent motions, but we cannot exclude this possibility.

We extend our analysis using EBHIS observations for $\ell=180^{\circ}$ and latitude $-10^{\circ}<b<10^{\circ}$. Figure 7 shows a position-velocity diagram for this case. The observed HI emission is shown color-coded and overplotted with model data. The model is based on a global fit of the 3D HI distribution in the Milky Way, including flaring and warping of the HI disk (Kalberla et al. 2007 with updates from Kalberla & Dedes 2008 222Table 1 from Kalberla & Dedes (2008) contains erroneous densities for the cold and warm neutral media that are too high by a factor 1.4 due to an inappropriate correction for He. Here we corrected this error.). Clearly, from the equilibrium model of the HI distribution no asymmetries in velocity are expected for all latitudes at $\ell=180^{\circ}$. The observed HI emission shows, as usual, significant turbulent fluctuations in intensity and velocity that cannot be modeled. But for all latitudes there is a common trend that the emission is shifted to negative velocities. To quantify this velocity shift we determine velocity centroids ($v_{\mathrm{c}}=\frac{\sum Tv}{\sum T}$; Kerp et al. 2016), the first moment of the velocity field, shown by the black line in Fig. 7. The average centroid velocity for $|b|<10\arcdeg$ is -4 km/s, however, the HI emission for $|b|>5^{\circ}$ is affected by emission at negative velocities. For latitudes $b\gtrsim 5^{\circ}$ there are several clouds at $-25$ km/s$\lesssim v_{\mathrm{LSR}}\lesssim-10$ km/s with peak brightness temperatures $T_{\mathrm{B}}\sim$ 50 K. At $b\lesssim 5^{\circ}$ this emission at intermediate negative velocities is weaker, 20 K$\lesssim T_{\mathrm{B}}\lesssim 40$ K but more extended in velocity, -35 km/s$\lesssim v_{\mathrm{LSR}}\lesssim-10$ km/s. Emission at positive velocities $v_{\mathrm{LSR}}\gtrsim 15$ km/s is nearly completely absent. Thus the velocity centroid for $|b|>5^{\circ}$ is biased toward negative velocities. The observed emission at negative velocities is part of a population of intermediate-velocity clouds, which are most likely located in the transition region between the Galactic disk and halo (Stanimirović et al. 2006). Emission from similar clouds in the negative-velocity wings is insignificant close to the disk at $|b|<5^{\circ}$ (Fig. 7): we determine for this range a centroid velocity of -3.5 km/s, which we consider characteristic for the HI in the anticenter region.

Variations and systematic changes in the velocity dispersion of the HI along the line of sight and uncertainties in the density distribution imply uncertainties in the shape of calculated model profiles in Figs. 6 and 7, but they do not affect the velocity centroids. The center velocities of individual gas clouds along the line of sight at $\ell=180\arcdeg$ depend merely on the dynamics of the system.

3.4 Analysis in regions different from the anticenter

The same analysis may be carried out in regions away from the anticenter, but in these cases the systematic error is much larger because the uncertainties in the rotation curve have an important effect on the conversion of heliocentric to Galactocentric coordinates.

The results are plotted in Fig. 8, where we can see that the average velocities are of the order of the systematic errors, rendering these measurements unuseful. We need the proper motions to reduce these systematics in these directions, a work that was already done with Gaia-DR2 data in G18 up to $R\approx 13$ kpc or LS18 up to $R\approx 20$ kpc.

4 Interpretation of results

As shown by LG16, the detection of ${v}_{R}\neq 0$ in a certain sky region means that the average orbit of the stars is not perfectly circular, either because it is elliptical or because there is a component of secular expansion/contraction either of the whole disk or of a part of it.

4.1 Elliptical orbits

In the case of elliptical orbits, ${v}_{R}\neq 0$ would indicate an eccentricity, $e$, different from zero. We have no information on the dependence of $v_{R}$ on $\phi$, so we cannot derive the exact orbits, but we can statistically derive the most likely values of $e$ for the mean orbit333The mean orbit is that followed by the average motion of large groups of stars in the same position, which is not the same as the orbit for each of the individual stars. For instance, the motion of the LSR has a closed orbit around the Galactic center, which is different from the motion of the Sun.. In the approximation of low $e$, the probability of having a given value of the eccentricity $e$ when we measure a radial Galactocentric velocity 444Do not confuse this with the average of the radial velocity over the whole orbit: here we only have $v_{R}$ for a given azimuth. $v_{R}$ with rms $\sigma$, at a particular point of the orbit, is (LG16, Eq. (11)):

[TABLE]

where, as mentioned above, $V_{\phi}$ is the azimuthal velocity (as usual, we take the rotation curve data of Sofue 2017, included in his Fig. 4, and we subtract the asymmetric drift as above; we neglect here the variation with $z$). Using this expression, we can evaluate the most likely values of the eccentricity as a function of $R$ derived from our results in the previous section. We carry out the calculation for the anticenter data with [Fe/H]$>-1$, without taking into account the systematic errors (we take $\sigma$ equal to the statistical error), which are plotted in Fig. 9. We can observe that values of $e\approx 0.07$ are the most likely for $R\approx 18$ kpc, and lower values for other Galactocentric radii.

4.2 Secular Expansion/Contraction in the Anticenter

For a particular line of sight in the anticenter, the relative variation of stellar mass within radii between $R$ and $R+dR$ is (LG16, Eq. (13)):

[TABLE]

where $\Sigma(R)$ is the stellar surface density, assuming a constant average mass/luminosity ratio throughout the disk. With an exponential disk $\Sigma(R)\propto e^{-R/h_{R}}$, we obtain that the relative gain of stellar mass toward the anticenter at radius $R$ is

[TABLE]

Assuming $h_{R}=2.0\pm 0.4$ kpc for a thin disk (López-Corredoira & Molgó 2014), we get the relative mass gain/loss given in Fig. 10. It is indeed positive (gain) with a maximum amplitude at $R\approx 15$ kpc and negative (loss) with a maximum amplitude at $R\approx 19$ kpc.

The values of Fig. 10 indicate a very fast evolution of the mass distribution. Most of the mass will be concentrated in the central regions without much change, but the outer disk has a trend at present to lose stars at $R\lesssim 11$ kpc and $R\gtrsim 16$ kpc and gain them at $R$ between 11 and 16 kpc. Assuming a surface density of stars in the solar neighbourhood of 27 M⊙pc*-2* (McKee et al. 2015), the total variation of stellar mass within $\Delta\phi=40^{\circ}=\frac{2\pi}{9}$ in our line of sight toward the anticenter is: $-0.48\pm 0.13$ M⊙/yr between R=9 and 11 kpc; $+0.18\pm 0.06$ M⊙/yr between R=11 and 16 kpc; $-0.03\pm 0.02$ M⊙/yr between R=16 and 23 kpc. This is not possible if the Galaxy is in a steady state, since the life of the Galaxy is much longer than the time $\Delta T=\left(\frac{\dot{M}}{M}\right)^{-1}$ and the size of the Galaxy cannot change so fast. This means that either there is a different trend for other azimuths, or these velocities $v_{R}(R)$ are not constant with time, so we live now in a period of fast expansion/contraction but this period will be short. This mechanism may contribute to the variation of the disk size over cosmological times. Also some stars will escape the Galaxy and will be part of the halo once they abandon the disk.

This net migration of stars might be related to the hypothesis known as stellar radial migrations, which has indeed been useful to explain breaks in surface brightness (Sánchez-Blázquez et al. 2009), the formation of the thick disk (Sales et al. 2009), disk flaring (Vera-Ciro et al. 2014) or metallicity distributions (Casagrande et al. 2011; Grand et al. 2015) and age-metallicity distributions (Frankel et al. 2018).

4.3 Dynamical causes

Minor mergers may raise vertical waves (Gómez et al. 2013), but these perturbations may not intensively affect the in-plane velocity (Tian et al. 2017). A local stream cannot be the explanation for radial velocities along a wide range of $\sim 20$ kpc, but it could be a large-scale stream associated with the Galaxy in the Sun–Galactic center line. Since we do not have evidence of such a huge structure embedded in our Galaxy, we think this is not very likely. Moreover, there is a very small asymmetry between the northern and southern Galactic hemispheres (Fig. 3), thus we think it is very unlikely that the same stream independent of the Galaxy be so symmetric with respect to the plane.

In principle, the migrations of stars are theoretically expected as a consequence of the resonances of the bar and transient spiral arms in the disk (Sellwood & Binney 2002; Roskar & Debattista 2015; Halle et al. 2018; Hunt et al. 2018). No gradient in the radial velocity is expected from a bar effect for the observed distances (Monari et al. 2014), but it can be produced by spiral arms (Faure et al. 2014).

4.4 Spiral arms

There might be accelerations due to spiral arms, which pull the stars toward the arm. In such a case, we would see fluctuations around zero in the positions where we cross a spiral arm: around $R=15$ kpc (the value at which $\langle v_{R}\rangle$ crosses the value of zero in Fig. 2), and this is roughly coincident with the position of the Outer (Norma–Cygnus) spiral arm (Camargo et al. 2015; Vallée 2017). Therefore, the influence of the spiral arm looks like a plausible candidate for driving our asymmetry. This mechanism, with compression where the stars enter the spiral arm and expansion where they exit, was also proposed as a tentative explanation for the variations of more local radial velocities by Siebert et al. (2011), Monari et al. (2016) and G18. The association of spiral structure potential with some ripples and ridges in the kinematic data is also proposed (Faure et al. 2014; Hunt et al. 2018; Quillen et al. 2018) without the need to invoke an external perturbation, such as from the passage of the Sagittarius dwarf galaxy (Antoja et al. 2018). In this case, one thing is remarkable from our data in Fig. 2: the transition from the highest velocity to the lowest velocity is quite smooth, over 4-5 kpc (although we must also bear in mind that the errors in the distance determination of the order of 10% might also produce this smoothing), much larger than the typical observable width in spiral arms of $\sim 0.8$ kpc (Vallée 2014). Nonetheless, the transition from positive to negative velocities (or vice versa) is also smooth and short in the case of spiral arms and it is not directly linked to the inter-arm distance (Faure et al. 2014, Fig. 6). In the simulations by Faure et al., the value of $\langle v_{R}\rangle$ is expected to be quite stable in the inter-arm region.

Before we discuss this question, we want to comment briefly on large-scale deviations from circular motions as observed in the HI gas. The spiral density wave theory is based on the hypothesis of quasi-stationary spiral structures (for a recent review see Shu 2016). Accordingly the motion of the gas is assumed to be forced by the background gravitational potential of spiral arms. As the HI gas streams through the spiral pattern it gets compressed and shocked, leading to phase transitions from the warm to the cold neutral medium. Hydrodynamic instabilities may cause the growth of knots in the spiral arms, eventually leading to molecular clouds and star formation. Spiral structure in the HI distribution is observed predominantly from velocity crowding along the line of sight (Burton 1972).

These assumptions have been used by Roberts (1972) and Simonson (1976) (see also the review by Shu 2016), to model spiral structures in the gaseous Milky Way disk. From Figs. 5 and 6 of Roberts (1972), we estimate at $\ell=180^{\circ}$ an expected velocity shift of the HI line of -13 to -10 km/s for the Perseus arm, approximately three times the centroid velocity measured by us. The plots by Simonson (1976) show also some expected velocity shift, but it is hard to read off numerical values. In both cases the Outer arm is not part of the model. In any case, these numbers should depend on $R$.

Let us now consider this question: how much mass do we need in the spiral arm overdensity to produce the observed effects? In order to estimate it, let us assume a simple model in which a spiral arm is a ring with radius $R_{\rm sp}$, width $\Delta_{\rm sp}$ and total overdensity mass equal to $M_{\rm sp}$. Of course, we know arms are spiral-shaped but for a simple calculation of the local effect of the gravitational force this approximation is good enough. The radial component of the acceleration of a star at Galactocentric radius $R$ is (Feng & Gallo 2011)

[TABLE]

[TABLE]

$K(m)$ and $E(m)$ are respectively the first and second complete elliptical integrals. The radial velocity of a star follows from $v_{R}dv_{R}=a_{R}dR$. The distribution of mean velocities that we get in our Fig. 2 might be a combination of particles with different initial conditions, some of them coming from the inner parts of the Galaxy with initial $v_{R}=0$ and others coming from the fall-off of particles from very high values of $R$. A simple estimation of the order of magnitude integrating the previous motion equation numerically gives us that

[TABLE]

with a slight dependence on the width of the spiral arms $\Delta_{\rm sp}$. If we compare this number with the value of +6 km/s that we obtained in Fig. 2, we get that the mass of the whole spiral arm should of the order of $\sim 7\times 10^{8}$ M⊙.

If we have four spiral arms in our Galaxy with the same mass as this one, the total mass of the spiral arms would be $\sim 3\times 10^{9}$ M⊙, which is about 3% of the Galactic disk mass (Pérez-Villegas et al. 2015). This is precisely the mass that is attributed to the spiral arms in the models (Pérez-Villegas et al. 2015), in agreement with other observations.

Kalberla et al. (2007) also posited the existence of a dark ring at 13 kpc$<R<$ 18.5 kpc, in order to reproduce the observed flaring in the HI distribution; however, they attribute a mass to this dark ring of $(2.2-2.8)\times 10^{10}$ M⊙, which is much larger than our estimation of the necessary mass of a ring. The general assumption of Kalberla et al. (2007) was that the HI distribution is isothermal, which means velocity dispersions that cause the observed scale heights are constant with $R$; possibly this was a wrong assumption. Anyway, apart from the calculation of the mass, apparently something is taking place at $R\approx 15$ kpc that produces important kinematical effects.

4.5 Velocity Fields in out-of-equilibrium Structures: Radial Motions, Dispersion and Dynamical Equilibrium

In this section we discuss the kinematic properties of two different kinds of system that are caused by different dynamical processes and that can be compared, only from a qualitative point of view, with the observed velocity field. On the one hand we consider a mock galaxy catalog that was built from a cosmological $N$-body simulation, and on the other hand we discuss the properties of some systems created by the gravitational collapse of isolated overdensities. The signature of the dynamical processes at work in the two cases is imprinted in the amplitude of radial motions: let us start from the latter case.

In B17, it was shown that the purely self-gravitating systems evolving from quite simple initial configurations can give rise easily to produce a quasi-planar spiral structure surrounding a virialized core and qualitatively resembling a spiral galaxy. In particular, for a broad range of non-spherical and non-uniform rotating initial conditions, gravitational relaxation gives rise quite generically to a rich variety of structure characterized by spiral-like arms; the main characteristic of these structures is that they are features of the intrinsically out-of-equilibrium nature of the system’s collapse. That is, they represent long-lived non-stationary configurations characterized by predominantly radial motions in their outermost parts, but they also incorporate an extended flattened region that rotates coherently about a well virialized core of triaxial shape with an approximately isotropic velocity dispersion. Let us consider a few examples.

To characterize the velocity field, we consider the radial component of a particle’s velocity $\vec{v}$,

[TABLE]

and the vectorial “transverse velocity”

[TABLE]

i.e., the vector of which the magnitude is that of the non-radial component of the velocity, but oriented parallel to the particle’s angular momentum relative to the origin. We will refer hereafter to the average of a quantity in a spherical shell with a radius $r$: coherent rotation of all the particles in the shell about the same axis then corresponds to $\langle|\vec{v}_{t}|\rangle$= $|\langle\vec{v}_{t}\rangle|=v_{\phi}$. The asymmetric drift is defined as (Binney & Tremaine 2008),

[TABLE]

where $v_{c}$ is the circular speed that corresponds to perfectly circular orbits, i.e.,

[TABLE]

To characterize the kinematics further, we also consider the anisotropy parameter $\beta(r)$

[TABLE]

where $\langle v_{t}^{2}(r)\rangle$ and $\langle v_{r}^{2}(r)\rangle$ are respectively the average square values of the transverse and radial velocities.

Let us now consider three examples of this class of system (see Figs. 11-13): the initial conditions were represented by uniform and isolated ellipsoids (with different ratio between the semi-axes) and with normalized spin parameter (Peebles 1969) $\lambda=0.1-0.3$ (see B17 and Benhaiem et al 2018 for more details).

We note that at sufficiently small scales all simulations are characterized by a compact core with isotropic velocity dispersion, for which $\beta(r)\approx 0$. At larger distances from the system’s center, the simulation A1 is characterized by a regime where rotational motions dominate, i.e. $v\approx v_{\phi}\approx v_{c}$ and $\beta(r)\ll 0$. In this case, the radial velocity remains small even in the outermost region of the system. On the other hand, the simulation A2, in an intermediate range of distances larger than the core’s size but with $r/r_{0}\leq 1$, is dominated by rotational motions (i.e., $v_{\phi}>0$) and has negligible net radial motions, i.e. $\langle v_{r}\rangle\approx 0$: however, it shows a non-zero radial velocity dispersion so that $0<\beta(r)<1$ and thus the asymmetric drift becomes large. Then at larger scales, the system is completely dominated by radial motion so that $\beta(r)\approx 1$. Finally the simulation A3 shows an intermediate situation between A1 and A3: radial motions, in the very outermost region of the system, remain of the same order as the rotational ones.

By a closer inspection we may note that, in all cases, the behavior of the average radial velocity, when different from zero, is first positive but close to zero, then negative and then positive again. This is due to the fact that there is a stream of particles, spatially correlated, that are expanding outward and they produce an infall motion of the shells that lie in front of them (see Fig. 14). The amplitude of the negative radial velocities reduces with time. Thus, a region where the radial velocity is negative is thus a natural outcome of these out-of-equilibrium models.

The two cases where radial motions dominate in the outermost part of the systems, i.e. A2 and A3, are characterized by a noticeable break of axisymmetry, while in the case where radial motions are smaller (i.e, A1) the system remains close to axisymmetric. This is due to the fact that radial motions develop during the collapse because of the amplification of the initial anisotropy: indeed, the simulation A1 was started from an oblate ellipsoid and thus remains symmetric in the plane of the two largest semi-axes — while it contracts in the orthogonal direction parallel to the minor semi-axis. Thus, because of the correlation between the direction of the radial motion and the direction of the major semi-axis, one must take into account that a measurement of the radial velocities in the direction of the anticenter of a given observer (in the real case, the Sun) is not generally expected to give a strong radial signal in these simulations, unless such a direction is (by chance) aligned with the the direction of the largest semi-axis or if the departure from axisymmetry is small enough. On the other hand, the more the system is axisymmetric, the higher is the probability that a random observer sees at its anticenter a growing radial velocity but the smaller is its amplitude.

Finally, let us now briefly discuss the difference between the models that we have presented, where the formation of a disc-like flat structure is caused solely by a gravitational, and thus dissipationless, dynamics and models in which instead the formation of a disc-like structure is driven by dissipational effects. Fig. 15 shows the velocity components of a mock stellar catalog that was built from a cosmological simulation of the Auriga Project 555Data for the simulation called Au6 are available at https://auriga.h-its.org/gaiamock.html. This provides a large suite of high-resolution magnetohydrodynamical simulations of galaxies simulated in a fully cosmological environment by means of the ‘zoom-in’ technique (Grand et al. 2018). One may note in Fig. 15 that the average radial component of the velocity is almost zero at all distances, while the velocity of stars is dominated by the circular component. In addition, the radial velocity dispersion is of the same order as the tangential one at large radii, while $\beta(r)<0$ at smaller radii where the motion is predominantly circular.

Fig. 16 shows that the radial velocity field is not characterized by any coherent spatial structure at large enough distances from the center (i.e., $r>5$ kpc), where the system is very close to axisymmetric. As mentioned above, the departure from axisymmetry is related to the presence of a complex velocity field, where radial motions become predominant at large distance: systems of this kind can be generated by a dynamics that involves globally the whole structure, i.e. a collapse, rather than a bottom-up aggregation process of the type at work in cosmological simulations.

5 Conclusions

Slightly non-circular mean stellar orbits in the Galactic disk were already known for $R\lesssim 16$ kpc (Siebert et al. 2011; Williams et al. 2013; LG16; Tian et al. 2017; Wang et al. 2018b; Carrillo et al. 2018; G18; LS18), but the analysis carried out in the present paper for the anticenter direction spans a wider range of distances, 8 kpc$<R<28$ kpc, and better characterizes the Galaxy’s kinematics on large scales. The outcome is an increase in positive $v_{R}$ (expansion of the disk) from $R\approx 9$ kpc to $R\approx 13$ kpc, reaching a maximum of $\approx 6$ km/s, and a decrease outward, reaching $\approx-10$ kpc for $R>17$ kpc (negative means contraction of the disk). The analysis in other directions different from the anticenter is not possible, due to the large systematic errors resulting when we have only radial heliocentric velocities, 666Nonetheless, when accurate proper motions are available, it would be possible to obtain information for other lines of sight, like in G18 or LS18. associated with uncertainties in the rotation curve, which are generated in the conversion between heliocentric and Galactocentric radial velocities. This feature of higher negative velocities is also observed in 21 cm HI maps toward the anticenter, possibly dominated by local gas, thus corroborating for the gas the asymmetry of velocities that shows the departure of circularity in the mean orbits.

Following LG16, we have analyzed the results in terms of mean ellipticity of the orbits or secular expansion, possibly associated with stellar migration deduced by other data sets. We discuss some possibilities about the dynamical causes that create these non-zero mean $v_{R}$. One possible explanation might be the gravitational attraction of a spiral arm. As a matter of fact, we see a change of regime from positive to negative velocities in the position where we cross the Outer spiral arm (Camargo et al. 2015) in the anticenter, $R\approx 15$ kpc, as would be expected if an arm with total mass $\sim 7\times 10^{8}$ M⊙ (in agreement with the expected mass of a spiral arm) pushed the stars with $R<15$ kpc (thus, $v_{R}>0$) and pulled the stars with $R>15$ kpc (thus, $v_{R}<0$).

An unconventional theoretical model was explored using the B17 dynamical model: a simple class of out-of-equilibrium, rotating, and asymmetrical mass distributions that evolve under their own gravity to produce a quasi-planar spiral structure surrounding a virialized core. Non-circular orbits in the outer disk are precisely one of the predictions of this model. This scenario is also related to spiral arms and their formation, but it is not the direct gravitational attraction of the spiral arms that produces the mean radial velocities. Rather, it is the fact that orbits in the very outer disk are out of equilibrium so they have not reached circularity yet. Under some initial conditions, our simulations with the B17 model are able to reproduce the observed feature of $v_{R}$ in Fig. 2 and predict that at $R$ much larger than 30 kpc we should return again to the regime of positive $v_{R}$ with very large amplitude (of the order of the rotation speed). This is something that we cannot check with the present data, but it would be desirable to do in the future.

Certainly, further kinematic information at farther distances or along different lines of sight might constrain better the dynamical scenarios of our Galaxy. G18 already provided some $v_{R}(\phi)$ analysis, but constrained by $R<13$ kpc, and longer distances need to be explored. The extension of the kinematical maps with Gaia-DR2 made by LS18 gives a better insight to farther distances, and the future releases of Gaia data will extend those maps even more. And also in other galaxies, 2D spectroscopy [for instance, or radio data of THINGS (Sylos Labini et al. 2018), or using Multi Unit Spectroscopic Explore (MUSE) at VLT telescope] may allow us to carry out an analysis of non-circularity in mean orbits. At present, at least, we can say that non-zero radial velocities in the outer disk are a fact, and this should persuade us to abandon the general idea that spiral galaxies are in a quasi-equilibrium configuration where stars move on steady circular orbits around the center of the galaxy.

Thanks are given to the anonymous referee for helpful comments and to Sébastien Comerón for helpful comments on a draft of this paper. M.L.C. was supported by the grant AYA2015-66506-P of the Spanish Ministry of Economy and Competitiveness (MINECO). F.S.L. was granted access to the HPC resources of The Institute for scientific Computing and Simulation financed by Region Ile de France and the project Equip@Meso (reference ANR-10-EQPX-29-01) overseen by the French National Research Agency (ANR) as part of the Investissements d’Avenir program. Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy Office of Science, and the Participating Institutions. SDSS-IV acknowledges support and resources from the Center for High-Performance Computing at the University of Utah. The SDSS website is www.sdss.org. SDSS-IV is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS Collaboration including the Brazilian Participation Group, the Carnegie Institution for Science, Carnegie Mellon University, the Chilean Participation Group, the French Participation Group, Harvard-Smithsonian Center for Astrophysics, Instituto de Astrofísica de Canarias, The Johns Hopkins University, Kavli Institute for the Physics and Mathematics of the Universe (IPMU) / University of Tokyo, Lawrence Berkeley National Laboratory, Leibniz Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Astrophysik (MPA Garching), Max-Planck-Institut für Extraterrestrische Physik (MPE), National Astronomical Observatories of China, New Mexico State University, New York University, University of Notre Dame, Observatário Nacional / MCTI, The Ohio State University, Pennsylvania State University, Shanghai Astronomical Observatory, United Kingdom Participation Group, Universidad Nacional Autónoma de México, University of Arizona, University of Colorado Boulder, University of Oxford, University of Portsmouth, University of Utah, University of Virginia, University of Washington, University of Wisconsin, Vanderbilt University, and Yale University.