propagate to the random motions through the BS effect, but on small scales, and it is more likely to find phase angles of the *k* = 2 mode orientated along the minor axis. The same effect applies to the major axis in the presence of a tangential large-scale anisotropy as caused by stronger streamings in the tangential dimension. When both $V_R$ and $V_\theta$ streaming motions are present, the phase angles of the bisymmetry of velocity dispersion reflect the competition between the two directions. Therefore, producing distributions of *k* = 2 phase angles similar to the observations would probably require specific mixing between the radial and tangential perturbations in discs, such as having different regions of dominance of streamings in one direction over the other. Also, no *k* = 4 phase angles centred around $\pi/4$, as in the observations, could be produced within the toy models. Our results suggest that the fourth order asymmetry is likely a combination between BS residuals and harmonics of the second order asymmetry in the velocity dispersion.

This analysis thus shows that it is possible to constrain which asymmetry between $V_R$ or $V_\phi$ dominates the streamings by studying the systematic orientations of the *k* = 1 and *k* = 2 phase angles of $\sigma_{los}$. In particular, to explain the systematic orientation of *k* = 1 phase angles along the major axis observed in the samples, this work suggests that lopsidedness in $V_{los}$ must be dominated by tangential perturbations, while the systematic orientations of *k* = 2 phase angles outside the major axis require elliptical streamings in $V_{los}$ affected by a stronger contribution from radial perturbations. By consequence, the observation of velocity asymmetries more significant in one direction in the plane than the other implies that H I velocities are highly anisotropic, at least at the angular scales probed by the Fourier analysis and the resolution of the observations. Yet it is possible that the velocity ellipsoid of gas remains isotropic at angular scales much smaller than those probed by the data, thus very locally, perhaps at the scale of individual clouds. Moreover, we predict that the signature of velocity anisotropy on large scales is also present in velocity fields of galaxies, again in the shape of systematic orientations of asymmetries near the principal axes of the discs. This prediction should greatly benefit from analyses of residual velocity fields of galaxies.

We further show that an anisotropic velocity ellipsoid with correlations between velocity components can generate asymmetries in the velocity dispersion with any orientation. Within this formalism, we predict that the probability to find coefficients of a correlation between $V_R$ and $V_\theta$ must be larger towards $|\rho_{R\theta}| \sim 0.6$, in order to explain the larger occurrence of gaseous $\phi_2$ at $\sim 0.35\pi$ and $\sim 0.65\pi$, assuming $\sigma_\theta/\sigma_R \sim 0.7$ as for the youngest populations of stars in the Milky Way. It would mean that the correlations between velocities in nearby gaseous discs are not random values, and thus that perturbations of the gravitational potential share common properties among the galaxies. This prediction could be compared with the velocity correlations measured for the youngest stellar populations of the Milky Way, to which the kinematics of interstellar gas should be comparable in many aspects.

This pioneering work has highlighted the importance of velocity anisotropy in shaping the asymmetries of H I velocity dispersions, through velocity streamings in the H I gas which are dominant in particular directions in the plane and/or tilted velocity ellipsoids. The present study is intended to serve as a stepping stone towards the design and the analysis of future large-scale surveys of galaxies, observed at higher sensitivity and resolution

than those samples studied here, such as those upcoming with the Square Kilometre Array and its precursors for the H I gas, or those already observed in the molecular interstellar medium of galaxies by ALMA, as well as in ionised gas and stellar velocity fields. In that aspect, comparisons of the shapes of asymmetries in $\sigma_{los}$ from both the stellar and gaseous components will be crucial to support the proposal that gas velocities are anisotropic. Harmonic decompositions of $V_{los}$ will also be helpful. This will be the subject of future papers following this work.

**Acknowledgements.** We are grateful to F. Walter for sharing the maps of the dirty beams of THINGS observations. P. Adamczyk and L. Chemin warmly acknowledge the Comité Mixto European Southern Observatory-Gobierno de Chile. The research of L. Chemin is funded by the Chilean Agencia Nacional de Investigación y Desarrollo (ANID) through the Fondo Nacional de Desarrollo Científico y Tecnológico (FONDECYT) Regular Project 1210992.

## References

* Adamczyk, P. 2021, Theses, Aix-Marseille Université, France

* Baldwin, J. E., Lynden-Bell, D., & Sancisi, R. 1980, MNRAS, 193, 313

* Bershady, M. A., Verheijen, M. A. W., Swaters, R. A., et al. 2010, ApJ, 716, 198

* Binney, J., & Merrifield, M. 1998, Galactic Astronomy (Princeton: Princeton University Press)

* Bosma, A. 1978, PhD Thesis, Groningen University, The Netherlands

* Chemin, L., Huré, J.-M., Soubiran, C., et al. 2016, A&A, 588, A48

* Chemin, L., Braine, J., Combes, F., Kam, Z. S., & Carignan, C. 2020, A&A, 639, A145

* Combes, F., & Becquaert, J. F. 1997, A&A, 326, 554

* Daigle, O., Carignan, C., Amram, P., et al. 2006, MNRAS, 367, 469

* Debattista, V. P., Gonzalez, O. A., Sanderson, R. E., et al. 2019, MNRAS, 485, 5073

* de Blok, W. J. G., Walter, F., Brinks, E., et al. 2008, AJ, 136, 2648

* Dicaire, I., Carignan, C., Amram, P., et al. 2008, MNRAS, 385, 553

* Di Teodoro, E. M., & Fraternali, F. 2015, MNRAS, 451, 3021

* Epinat, B., Amram, P., Balkowski, C., & Marcelin, M. 2010, MNRAS, 401, 2113

* Epinat, B., Tasca, L., Amram, P., et al. 2012, A&A, 539, A92

* Gaia Collaboration (Prusti, T., et al.) 2016, A&A, 595, A1

* Gaia Collaboration (Brown, A. G. A., et al.) 2018a, A&A, 616, A1

* Gaia Collaboration (Helmi, A., et al.) 2018b, A&A, 616, A12

* Gaia Collaboration (Brown, A. G. A., et al.) 2021a, A&A, 649, A1

* Gaia Collaboration (Luri, X., et al.) 2021b, A&A, 649, A7

* Gaia Collaboration (Drimmel, R., et al.) 2023, A&A, 674, A37

* Hagen, J. H. J., Helmi, A., de Zeeuw, P. T., & Posti, L. 2019, A&A, 629, A70

* Hayashi, E., & Navarro, J. F. 2006, MNRAS, 373, 1117

* Hunter, D. A., Elmegreen, B. G., & van Woerden, H. 2001, ApJ, 556, 773

* Kennicutt, R. C., Jr., Armus, L., Bendo, G., et al. 2003, PASP, 115, 928

* Koyama, H., & Ostriker, E. C. 2009, ApJ, 693, 1346

* Krajnović, D., Cappellari, M., de Zeeuw, P. T., & Copin, Y. 2006, MNRAS, 366, 787

* Kuijken, K., & Tremaine, S. 1991, in Dynamics of Disc Galaxies, ed. B. Sundelius, 71

* Kuijken, K., & Tremaine, S. 1994, ApJ, 421, 178

* Oh, S.-H., Hunter, D. A., Brinks, E., et al. 2015, AJ, 149, 180

* Oman, K. A., Marasco, A., Navarro, J. F., et al. 2019, MNRAS, 482, 821

* Schoenmakers, R. H. M., Franx, M., & de Zeeuw, P. T. 1997, MNRAS, 292, 349

* Siebert, A., Bienaymé, O., Binney, J., et al. 2008, MNRAS, 391, 793

* Smith, M. C., Evans, N. W., & An, J. H. 2009, ApJ, 698, 1110

* Spekkens, K., & Sellwood, J. A. 2007, ApJ, 664, 204

* Swaters, R. A., Schoenmakers, R. H. M., Sancisi, R., & van Albada, T. S. 1999, MNRAS, 304, 330

* Trachternach, C., de Blok, W. J. G., Walter, F., Brinks, E., & Kennicutt, R. C., Jr. 2008, AJ, 136, 2720

* van der Hulst, J. M., van Albada, T. S., & Sancisi, R. 2001, in Gas and Galaxy Evolution, eds. J. E. Hibbard, M. Rupen, & J. H. van Gorkom, ASP Conf. Ser., 240, 451

* Visser, H. C. D. 1980, A&A, 88, 159

* Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ, 136, 2563

* Wegg, C., Gerhard, O., & Bieth, M. 2019, MNRAS, 485, 3296

* Wilcots, E. M., & Thurow, J. C. 2001, ApJ, 555, 758

* Zhou, Z., Wu, H., Zhou, X., & Ma, J. 2018, PASP, 130, 094101
propagate to the random motions through the BS effect, but on small scales, and it is more likely to find phase angles of the *k* = 2 mode orientated along the minor axis. The same effect applies to the major axis in the presence of a tangential large-scale anisotropy as caused by stronger streamings in the tangential dimension. When both $V_R$ and $V_\theta$ streaming motions are present, the phase angles of the bisymmetry of velocity dispersion reflect the competition between the two directions. Therefore, producing distributions of *k* = 2 phase angles similar to the observations would probably require specific mixing between the radial and tangential perturbations in discs, such as having different regions of dominance of streamings in one direction over the other. Also, no *k* = 4 phase angles centred around $\pi/4$, as in the observations, could be produced within the toy models. Our results suggest that the fourth order asymmetry is likely a combination between BS residuals and harmonics of the second order asymmetry in the velocity dispersion.

This analysis thus shows that it is possible to constrain which asymmetry between $V_R$ or $V_\theta$ dominates the streamings by studying the systematic orientations of the *k* = 1 and *k* = 2 phase angles of $\sigma_{los}$. In particular, to explain the systematic orientation of *k* = 1 phase angles along the major axis observed in the samples, this work suggests that lopsidedness in $V_{los}$ must be dominated by tangential perturbations, while the systematic orientations of *k* = 2 phase angles outside the major axis require elliptical streamings in $V_{los}$ affected by a stronger contribution from radial perturbations. By consequence, the observation of velocity asymmetries more significant in one direction in the plane than the other implies that H I velocities are highly anisotropic, at least at the angular scales probed by the Fourier analysis and the resolution of the observations. Yet it is possible that the velocity ellipsoid of gas remains isotropic at angular scales much smaller than those probed by the data, thus very locally, perhaps at the scale of individual clouds. Moreover, we predict that the signature of velocity anisotropy on large scales is also present in velocity fields of galaxies, again in the shape of systematic orientations of asymmetries near the principal axes of the discs. This prediction should greatly benefit from analyses of residual velocity fields of galaxies.

We further show that an anisotropic velocity ellipsoid with correlations between velocity components can generate asymmetries in the velocity dispersion with any orientation. Within this formalism, we predict that the probability to find coefficients of a correlation between $V_R$ and $V_\theta$ must be larger towards $|\rho_{R\theta}| \sim 0.6$, in order to explain the larger occurrence of gaseous $\phi_2$ at $\sim 0.35\pi$ and $\sim 0.65\pi$, assuming $\sigma_\theta/\sigma_R \sim 0.7$ as for the youngest populations of stars in the Milky Way. It would mean that the correlations between velocities in nearby gaseous discs are not random values, and thus that perturbations of the gravitational potential share common properties among the galaxies. This prediction could be compared with the velocity correlations measured for the youngest stellar populations of the Milky Way, to which the kinematics of interstellar gas should be comparable in many aspects.

This pioneering work has highlighted the importance of velocity anisotropy in shaping the asymmetries of H I velocity dispersions, through velocity streamings in the H I gas which are dominant in particular directions in the plane and/or tilted velocity ellipsoids. The present study is intended to serve as a stepping stone towards the design and the analysis of future large-scale surveys of galaxies, observed at higher sensitivity and resolution

than those samples studied here, such as those upcoming with the Square Kilometre Array and its precursors for the H I gas, or those already observed in the molecular interstellar medium of galaxies by ALMA, as well as in ionised gas and stellar velocity fields. In that aspect, comparisons of the shapes of asymmetries in $\sigma_{los}$ from both the stellar and gaseous components will be crucial to support the proposal that gas velocities are anisotropic. Harmonic decompositions of $V_{los}$ will also be helpful. This will be the subject of future papers following this work.

**Acknowledgements.** We are grateful to F. Walter for sharing the maps of the dirty beams of THINGS observations. P. Adamczyk and L. Chemin warmly acknowledge the Comité Mixto European Southern Observatory-Gobierno de Chile. The research of L. Chemin is funded by the Chilean Agencia Nacional de Investigación y Desarrollo (ANID) through the Fondo Nacional de Desarrollo Científico y Tecnológico (FONDECYT) Regular Project 1210992.

## References

* Adamczyk, P. 2021, Theses, Aix-Marseille Université, France

* Baldwin, J. E., Lynden-Bell, D., & Sancisi, R. 1980, MNRAS, 193, 313

* Bershady, M. A., Verheijen, M. A. W., Swaters, R. A., et al. 2010, ApJ, 716, 198

* Binney, J., & Merrifield, M. 1998, Galactic Astronomy (Princeton: Princeton University Press)

* Bosma, A. 1978, PhD Thesis, Groningen University, The Netherlands

* Chemin, L., Huré, J.-M., Soubiran, C., et al. 2016, A&A, 588, A48

* Chemin, L., Braine, J., Combes, F., Kam, Z. S., & Carignan, C. 2020, A&A, 639, A145

* Combes, F., & Becquaert, J. F. 1997, A&A, 326, 554

* Daigle, O., Carignan, C., Amram, P., et al. 2006, MNRAS, 367, 469

* Debattista, V. P., Gonzalez, O. A., Sanderson, R. E., et al. 2019, MNRAS, 485, 5073

* de Blok, W. J. G., Walter, F., Brinks, E., et al. 2008, AJ, 136, 2648

* Dicaire, I., Carignan, C., Amram, P., et al. 2008, MNRAS, 385, 553

* Di Teodoro, E. M., & Fraternali, F. 2015, MNRAS, 451, 3021

* Epinat, B., Amram, P., Balkowski, C., & Marcelin, M. 2010, MNRAS, 401, 2113

* Epinat, B., Tasca, L., Amram, P., et al. 2012, A&A, 539, A92

* Gaia Collaboration (Prusti, T., et al.) 2016, A&A, 595, A1

* Gaia Collaboration (Brown, A. G. A., et al.) 2018a, A&A, 616, A1

* Gaia Collaboration (Helmi, A., et al.) 2018b, A&A, 616, A12

* Gaia Collaboration (Brown, A. G. A., et al.) 2021a, A&A, 649, A1

* Gaia Collaboration (Luri, X., et al.) 2021b, A&A, 649, A7

* Gaia Collaboration (Drimmel, R., et al.) 2023, A&A, 674, A37

* Hagen, J. H. J., Helmi, A., de Zeeuw, P. T., & Posti, L. 2019, A&A, 629, A70

* Hayashi, E., & Navarro, J. F. 2006, MNRAS, 373, 1117

* Hunter, D. A., Elmegreen, B. G., & van Woerden, H. 2001, ApJ, 556, 773

* Kennicutt, R. C., Jr., Armus, L., Bendo, G., et al. 2003, PASP, 115, 928

* Koyama, H., & Ostriker, E. C. 2009, ApJ, 693, 1346

* Krajnović, D., Cappellari, M., de Zeeuw, P. T., & Copin, Y. 2006, MNRAS, 366, 787

* Kuijken, K., & Tremaine, S. 1991, in Dynamics of Disc Galaxies, ed. B. Sundelius, 71

* Kuijken, K., & Tremaine, S. 1994, ApJ, 421, 178

* Oh, S.-H., Hunter, D. A., Brinks, E., et al. 2015, AJ, 149, 180

* Oman, K. A., Marasco, A., Navarro, J. F., et al. 2019, MNRAS, 482, 821

* Schoenmakers, R. H. M., Franx, M., & de Zeeuw, P. T. 1997, MNRAS, 292, 349

* Siebert, A., Bienaymé, O., Binney, J., et al. 2008, MNRAS, 391, 793

* Smith, M. C., Evans, N. W., & An, J. H. 2009, ApJ, 698, 1110

* Spekkens, K., & Sellwood, J. A. 2007, ApJ, 664, 204

* Swaters, R. A., Schoenmakers, R. H. M., Sancisi, R., & van Albada, T. S. 1999, MNRAS, 304, 330

* Trachternach, C., de Blok, W. J. G., Walter, F., Brinks, E., & Kennicutt, R. C., Jr. 2008, AJ, 136, 2720

* van der Hulst, J. M., van Albada, T. S., & Sancisi, R. 2001, in Gas and Galaxy Evolution, eds. J. E. Hibbard, M. Rupen, & J. H. van Gorkom, ASP Conf. Ser., 240, 451

* Visser, H. C. D. 1980, A&A, 88, 159

* Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ, 136, 2563

* Wegg, C., Gerhard, O., & Bieth, M. 2019, MNRAS, 485, 3296

* Wilcots, E. M., & Thurow, J. C. 2001, ApJ, 555, 758

* Zhou, Z., Wu, H., Zhou, X., & Ma, J. 2018, PASP, 130, 094101
# Appendix A: Significance of FFT results with noisy data

In inner regions of galaxies, rings are smaller and contain less pixels than in the outermost parts. Moreover, the average velocity dispersion tends to decrease with radius, and can thus vary significantly from one ring to another, so that the properties of rings are not necessarily uniform as a function of radius. The results of a FFT, and in particular the precision that can be obtained on Fourier amplitudes, is impacted these parameters. It is essential to understand how the noise in the data affects the results of the FFT and to estimate to which limits we can recover the real signal, or in other words the reliability of the recovered signal. In our study, we apply the FFT implemented in the Python package FFTPACK from the Scipy library on the squared velocity dispersion $\sigma^2$. To have a realistic estimate of the noise $\sigma^2$, we express the uncertainties on $\sigma^2$ using the differential of the squared velocity dispersion as:

$$d(\sigma^2) = 2\sigma_0 d\sigma, \qquad (A.1)$$

where $\sigma_0^2 = \langle\sigma^2\rangle$ is the value of the zeroth order of the FFT, because all higher orders cancel out on average over a ring, and where $d\sigma$ is the uncertainty on the observed velocity dispersion and that roughly corresponds to the uncertainty on the instrumental resolution at a given signal-to-noise ratio. Since we do not have uncertainty maps due to the processing techniques (moments method) on THINGS data, $d\sigma$ is not known precisely. Since the spectral resolution of both THINGS and WHISP data corresponds to about $\sigma_{LSF} \sim 2-5$ km s$^{-1}$, and because the actual uncertainty on the line width depends on the S/N, we assume that the uncertainty on the velocity dispersion is of the order $d\sigma \sim 1$ km s$^{-1}$.

The way the uncertainties propagate from the data to the FFT coefficients depends on how the FFT is implemented. A discrete Fourier transform is used in FFTPACK to compute Fourier coefficients:

$$C_k = \sum_{m=0}^{N-1} c_m \cos\left(-\frac{2\pi km}{N}\right) + j \sum_{m=0}^{N-1} c_m \sin\left(-\frac{2\pi km}{N}\right), \qquad (A.2)$$

where j is the unit imaginary number, $C_k$ is the Fourier coefficient of order k as computed by FFTPACK, and $c_m$ are the N elements in the discrete array. In our study, $c_m$ are the squared velocity dispersion values along rings, $\sigma^2(m)$, and we deduce the amplitude of each order $k \neq 0$ as $\sigma_k^2 = 2|C_k|/N$, and $\sigma_0^2 = |C_0|/N$ as the average of the squared velocity dispersion values along the ring. From this formalism, it can be shown that for $k \neq 0$:

$$\mathrm{Var}[\sigma_k^2] = \frac{4}{N^2} \sum_{m=0}^{N-1} \mathrm{Var}[\sigma(m)^2] + \frac{8}{N^2} \sum_{\Delta m=1}^{N-1} \sum_{m=0}^{N-1-\Delta m} \mathrm{Cov}[\sigma(m)^2, \sigma(m+\Delta m)^2] \cos\left(\frac{2\pi k \Delta m}{N}\right), \qquad (A.3)$$

with Var and Cov being respectively the variance and covariance. Assuming that the variance is constant across the ring ($\mathrm{Var}[\sigma(m)^2] = \mathrm{Var}[\sigma^2]$) and that covariance is null between different pixels, we can get at first approximation that:

$$\mathrm{Var}[\sigma_k^2] = \frac{4}{N} \mathrm{Var}[\sigma^2], \qquad (A.4)$$

<figure><img src="image_17.png" /><figcaption>Fig. A.1. Square root of the noise on the squared velocity dispersion amplitude for orders $k = 1$ to 4 of the FFT as a function of $\sigma_0$ for two distinct values of $d\sigma$ and various number of points $N$ in the FFT. The curves correspond to Eq. A.5 with $d\sigma = 1$ km s$^{-1}$ (dashed) and $d\sigma = 4$ km s$^{-1}$ (dotted), each corresponding to a different value of $N$ (from top to bottom: 1000, 2000, 4000 and 8000). Coloured symbols correspond to different orders measured on toy models. The similar behaviour of orders leads to superposed points.</figcaption></figure>

leading to the uncertainty on coefficients:

$$d(\sigma_k^2) = \frac{4\sigma_0 d\sigma}{\sqrt{N}}. \qquad (A.5)$$

Interestingly, this equation shows that the uncertainty is the same for all orders, that it is proportional to the zeroth order, and that it decreases with radius since the number of pixels increases linearly with radius and because $\sigma_0$ rarely tends to decrease with radius.

We have checked numerically the consistency between this expression and the actual results from FFT using a toy model. We define the noise $\delta$ on $\sigma_{obs}^2$ as a Gaussian distribution centred on 0 and with a standard deviation of $\sigma_{Gauss} = d(\sigma^2)$. We create a toy model of the observed velocity dispersion with a constant velocity dispersion to which we added noise:

$$\sigma_{obs}^2 = \sigma_0^2 + \delta, \qquad (A.6)$$

with $\delta$ a random value in the noise distribution. Doing so, the distribution of each observed point is taken independently. We used several values for $d\sigma$, $\sigma_0$, and $N$ the number of points in the FFT, since $N$ varies from one ring to another and from one galaxy to another. We do 1000 iterations of the FFT and analyse the amplitudes for orders 1 to 4, which are null in our model, as shown in Eq. A.6. Fig. A.1 shows the measured uncertainties of the toy models and confirms that Eq. A.5 is a good description of the uncertainties. Equation A.5 can therefore be used to infer at first order the significance of Fourier coefficients, $s_k = \sigma_k^2/d(\sigma_k^2)$, depending on the number of rings and on the mean velocity dispersion value.
## Appendix A: Significance of FFT results with noisy data

In inner regions of galaxies, rings are smaller and contain less pixels than in the outermost parts. Moreover, the average velocity dispersion tends to decrease with radius, and can thus vary significantly from one ring to another, so that the properties of rings are not necessarily uniform as a function of radius. The results of a FFT, and in particular the precision that can be obtained on Fourier amplitudes, is impacted these parameters. It is essential to understand how the noise in the data affects the results of the FFT and to estimate to which limits we can recover the real signal, or in other words the reliability of the recovered signal. In our study, we apply the FFT implemented in the Python package FFTPACK from the Scipy library on the squared velocity dispersion $\sigma^2$. To have a realistic estimate of the noise $\sigma^2$, we express the uncertainties on $\sigma^2$ using the differential of the squared velocity dispersion as:

$$d(\sigma^2) = 2\sigma_0 d\sigma, \qquad (A.1)$$

where $\sigma_0^2 = \langle\sigma^2\rangle$ is the value of the zeroth order of the FFT, because all higher orders cancel out on average over a ring, and where $d\sigma$ is the uncertainty on the observed velocity dispersion and that roughly corresponds to the uncertainty on the instrumental resolution at a given signal-to-noise ratio. Since we do not have uncertainty maps due to the processing techniques (moments method) on THINGS data, $d\sigma$ is not known precisely. Since the spectral resolution of both THINGS and WHISP data corresponds to about $\sigma_{LSF} \sim 2-5$ km s$^{-1}$, and because the actual uncertainty on the line width depends on the S/N, we assume that the uncertainty on the velocity dispersion is of the order $d\sigma \sim 1$ km s$^{-1}$.

The way the uncertainties propagate from the data to the FFT coefficients depends on how the FFT is implemented. A discrete Fourier transform is used in FFTPACK to compute Fourier coefficients:

$$C_k = \sum_{m=0}^{N-1} c_m \cos\left(-\frac{2\pi km}{N}\right) + j \sum_{m=0}^{N-1} c_m \sin\left(-\frac{2\pi km}{N}\right), \qquad (A.2)$$

where j is the unit imaginary number, $C_k$ is the Fourier coefficient of order k as computed by FFTPACK, and $c_m$ are the N elements in the discrete array. In our study, $c_m$ are the squared velocity dispersion values along rings, $\sigma^2(m)$, and we deduce the amplitude of each order $k \neq 0$ as $\sigma_k^2 = 2|C_k|/N$, and $\sigma_0^2 = |C_0|/N$ as the average of the squared velocity dispersion values along the ring. From this formalism, it can be shown that for $k \neq 0$:

$$\mathrm{Var}[\sigma_k^2] = \frac{4}{N^2} \sum_{m=0}^{N-1} \mathrm{Var}[\sigma(m)^2] + \frac{8}{N^2} \sum_{\Delta m=1}^{N-1} \sum_{m=0}^{N-1-\Delta m} \mathrm{Cov}[\sigma(m)^2, \sigma(m+\Delta m)^2] \cos\left(\frac{2\pi k \Delta m}{N}\right), \qquad (A.3)$$

with Var and Cov being respectively the variance and covariance. Assuming that the variance is constant across the ring ($\mathrm{Var}[\sigma(m)^2] = \mathrm{Var}[\sigma^2]$) and that covariance is null between different pixels, we can get at first approximation that:

$$\mathrm{Var}[\sigma_k^2] = \frac{4}{N} \mathrm{Var}[\sigma^2], \qquad (A.4)$$

<figure><img src="image_17.png" /><figcaption>Fig. A.1. Square root of the noise on the squared velocity dispersion amplitude for orders $k = 1$ to 4 of the FFT as a function of $\sigma_0$ for two distinct values of $d\sigma$ and various number of points $N$ in the FFT. The curves correspond to Eq. A.5 with $d\sigma = 1$ km s$^{-1}$ (dashed) and $d\sigma = 4$ km s$^{-1}$ (dotted), each corresponding to a different value of $N$ (from top to bottom: 1000, 2000, 4000 and 8000). Coloured symbols correspond to different orders measured on toy models. The similar behaviour of orders leads to superposed points.</figcaption></figure>

leading to the uncertainty on coefficients:

$$d(\sigma_k^2) = \frac{4\sigma_0 d\sigma}{\sqrt{N}}. \qquad (A.5)$$

Interestingly, this equation shows that the uncertainty is the same for all orders, that it is proportional to the zeroth order, and that it decreases with radius since the number of pixels increases linearly with radius and because $\sigma_0$ rarely tends to decrease with radius.

We have checked numerically the consistency between this expression and the actual results from FFT using a toy model. We define the noise $\delta$ on $\sigma_{obs}^2$ as a Gaussian distribution centred on 0 and with a standard deviation of $\sigma_{Gauss} = d(\sigma^2)$. We create a toy model of the observed velocity dispersion with a constant velocity dispersion to which we added noise:

$$\sigma_{obs}^2 = \sigma_0^2 + \delta, \qquad (A.6)$$

with $\delta$ a random value in the noise distribution. Doing so, the distribution of each observed point is taken independently. We used several values for $d\sigma$, $\sigma_0$, and $N$ the number of points in the FFT, since $N$ varies from one ring to another and from one galaxy to another. We do 1000 iterations of the FFT and analyse the amplitudes for orders 1 to 4, which are null in our model, as shown in Eq. A.6. Fig. A.1 shows the measured uncertainties of the toy models and confirms that Eq. A.5 is a good description of the uncertainties. Equation A.5 can therefore be used to infer at first order the significance of Fourier coefficients, $s_k = \sigma_k^2/d(\sigma_k^2)$, depending on the number of rings and on the mean velocity dispersion value.
## Appendix B: Toy model description to study projection and resolution effects on random motions

Here, we describe the detailed generation of the toy models presented in Sect. 3.3 and used in Sections 5.2.1 and 5.3.1 to investigate the impact of projection effects on the study of random motions in both axisymmetric and asymmetric cases as well as with and without local anisotropy. We first built disc-like kinematics with $5 \times 10^6$ points representing individual particles or gas clouds uniformly distributed over a disc of 200'' radius. The adopted axisymmetric contribution of azimuthal velocities $V_\theta$ follows $V_\theta(R) = v_0/(1 + (r_s/R)^\gamma)^{1/\gamma}$ for two different sets of parameters to describe either (i) a rotation curve slowly rising to a moderate maximum rotation velocity with $v_0 = 200 \text{ km s}^{-1}$, $r_s = 100''$ and $\gamma = 1$, or (ii) a rotation curve with a steep inner gradient and a sharp transition to a high velocity plateau reached at small radius with $v_0 = 250 \text{ km s}^{-1}$, $r_s = 20''$ and $\gamma = 2$. We considered those two rotation curves to have either a low or a high impact of large-scale rotation on BS. The disc was not assumed contracting or expanding, that is, $V_R(R) = 0$, and a negligible average vertical component, $V_z(R) = 0$. A few solutions exist to describe the departure of $V_\theta$ and $V_R$ from axisymmetry at a given $R, \theta$ position in the plane. For example, planar asymmetries can follow Fourier harmonics of second order to mimic elliptical orbits in a bar or spiral potential (e.g. Spekkens & Sellwood 2007; Gaia Collaboration 2023). This is useful to describe large-scale perturbations of $V_\theta$ and $V_R$, but not necessarily appropriate to mock sharp velocity gradients across and along spiral arms (Chemin et al. 2016). Therefore, we adopted a different prescription than the simple cylindrical harmonics, and defined a spiral pattern by:

$$ \theta = \phi_{sp} + 2\pi R/r_0, \qquad (B.1) $$

where $R$ and $\theta$ are radial and azimuthal coordinates in the plane of the galaxy, $\phi_{sp}$ is the phase of the pattern at the centre, and $r_0$ is the radius interval necessary for the spiral to complete $2\pi$. We used $r_0/2\pi = 50''$ in all models, and since the interesting parameter to probe here is the difference between the orientation of the velocity perturbation and a reference axis in the disc, we vary the initial angle of the spiral perturbation with a sampling of $15^\circ$, from $0^\circ$ to $165^\circ$ with respect to the major axis of the galaxy, leading to 12 possible orientations for the pattern. For each $R, \theta$ position, the closest point along the spiral in the galaxy frame was found to infer the distance $D_{sp}$ to the spiral and the unitary vector $u$ between those points. A density perturbation following the spiral pattern was introduced as a constant over-density for all points having their distance $D_{sp}$ lower than a maximum distance $\Delta D$ where the spiral has an impact on the axisymmetric model. We note that the distribution of points is not important in this exercise, as it does not affect the kinematics. It can only impact the derivation of the BS model, though with negligible overall effect. We did not model the case of asymmetric vertical motions, and defined the amplitude of the velocity perturbation by $\Delta V = \Delta V_0 \times (1 - D_{sp}/\Delta D)$ for $D_{sp} \le \Delta D$, $\Delta V_0$ being the maximum amplitude of the velocity perturbation. We used $\Delta V_0 = 25 \text{ km s}^{-1}$ and $\Delta D = 25''$. The value of $25 \text{ km s}^{-1}$ was chosen to match the maximum variation of velocity through the spiral arms of the grand-design spiral Messier 99 with respect to the axisymmetric circular velocity (Fig. 8 of Chemin et al. 2016), or the maximum strength of the elliptical motions of RGB stars in the Galactic bar (Fig. 19 of Gaia Collaboration 2023). We assumed the six following configurations:

<figure><img src="image_18.png" /><figcaption>Fig. B.1. Velocity ellipsoid of the mock disc with local velocity isotropy for the cases without large-scale anisotropy (upper panel), and with a large-scale anisotropy induced by the spiral $V_R$ streamings (bottom panel), assuming radial motions are null on average. The white contours represent densities of 100, 1000, 5000, and 10000 toy-model particles.</figcaption></figure>

- (1) Axisymmetric models with no perturbation ($\Delta V = \Delta V_R = \Delta V_\theta = 0$). These models are used as references from which comparisons can be made, with respect to previous asymmetric models.

- Asymmetric velocity models with a spiral velocity perturbation along either (2) the azimuthal direction ($\Delta V_\theta = \Delta V \times u \cdot u_\theta / |u \cdot u_\theta|$, $\Delta V_R = 0$), $u_\theta$ being the unitary azimuthal vector, or (3) the radial direction ($\Delta V_R = \Delta V \times u \cdot u_R / |u \cdot u_R|$, $\Delta V_\theta = 0$), $u_R$ being the unitary radial vector. These cases are necessary to identify the impact of velocity gradients along one planar direction, independently from the other direction. This thus represents a departure from $V_\theta(R)$ of at least ~13%, and significantly more from $V_R(R)$.

- Asymmetric velocity models with a spiral velocity perturbation along both azimuthal and radial directions for the
## Appendix B: Toy model description to study projection and resolution effects on random motions

Here, we describe the detailed generation of the toy models presented in Sect. 3.3 and used in Sections 5.2.1 and 5.3.1 to investigate the impact of projection effects on the study of random motions in both axisymmetric and asymmetric cases as well as with and without local anisotropy. We first built disc-like kinematics with $5 \times 10^6$ points representing individual particles or gas clouds uniformly distributed over a disc of 200'' radius. The adopted axisymmetric contribution of azimuthal velocities $V_\theta$ follows $V_\theta(R) = v_0/(1 + (r_s/R)^\gamma)^{1/\gamma}$ for two different sets of parameters to describe either (i) a rotation curve slowly rising to a moderate maximum rotation velocity with $v_0 = 200 \text{ km s}^{-1}$, $r_s = 100''$ and $\gamma = 1$, or (ii) a rotation curve with a steep inner gradient and a sharp transition to a high velocity plateau reached at small radius with $v_0 = 250 \text{ km s}^{-1}$, $r_s = 20''$ and $\gamma = 2$. We considered those two rotation curves to have either a low or a high impact of large-scale rotation on BS. The disc was not assumed contracting or expanding, that is, $V_R(R) = 0$, and a negligible average vertical component, $V_z(R) = 0$. A few solutions exist to describe the departure of $V_\theta$ and $V_R$ from axisymmetry at a given $R, \theta$ position in the plane. For example, planar asymmetries can follow Fourier harmonics of second order to mimic elliptical orbits in a bar or spiral potential (e.g. Spekkens & Sellwood 2007; Gaia Collaboration 2023). This is useful to describe large-scale perturbations of $V_\theta$ and $V_R$, but not necessarily appropriate to mock sharp velocity gradients across and along spiral arms (Chemin et al. 2016). Therefore, we adopted a different prescription than the simple cylindrical harmonics, and defined a spiral pattern by:

$$ \theta = \phi_{sp} + 2\pi R/r_0, \qquad (B.1) $$

where $R$ and $\theta$ are radial and azimuthal coordinates in the plane of the galaxy, $\phi_{sp}$ is the phase of the pattern at the centre, and $r_0$ is the radius interval necessary for the spiral to complete $2\pi$. We used $r_0/2\pi = 50''$ in all models, and since the interesting parameter to probe here is the difference between the orientation of the velocity perturbation and a reference axis in the disc, we vary the initial angle of the spiral perturbation with a sampling of $15^\circ$, from $0^\circ$ to $165^\circ$ with respect to the major axis of the galaxy, leading to 12 possible orientations for the pattern. For each $R, \theta$ position, the closest point along the spiral in the galaxy frame was found to infer the distance $D_{sp}$ to the spiral and the unitary vector $u$ between those points. A density perturbation following the spiral pattern was introduced as a constant over-density for all points having their distance $D_{sp}$ lower than a maximum distance $\Delta D$ where the spiral has an impact on the axisymmetric model. We note that the distribution of points is not important in this exercise, as it does not affect the kinematics. It can only impact the derivation of the BS model, though with negligible overall effect. We did not model the case of asymmetric vertical motions, and defined the amplitude of the velocity perturbation by $\Delta V = \Delta V_0 \times (1 - D_{sp}/\Delta D)$ for $D_{sp} \le \Delta D$, $\Delta V_0$ being the maximum amplitude of the velocity perturbation. We used $\Delta V_0 = 25 \text{ km s}^{-1}$ and $\Delta D = 25''$. The value of $25 \text{ km s}^{-1}$ was chosen to match the maximum variation of velocity through the spiral arms of the grand-design spiral Messier 99 with respect to the axisymmetric circular velocity (Fig. 8 of Chemin et al. 2016), or the maximum strength of the elliptical motions of RGB stars in the Galactic bar (Fig. 19 of Gaia Collaboration 2023). We assumed the six following configurations:

<figure><img src="image_18.png" /><figcaption>Fig. B.1. Velocity ellipsoid of the mock disc with local velocity isotropy for the cases without large-scale anisotropy (upper panel), and with a large-scale anisotropy induced by the spiral $V_R$ streamings (bottom panel), assuming radial motions are null on average. The white contours represent densities of 100, 1000, 5000, and 10000 toy-model particles.</figcaption></figure>

- (1) Axisymmetric models with no perturbation ($\Delta V = \Delta V_R = \Delta V_\theta = 0$). These models are used as references from which comparisons can be made, with respect to previous asymmetric models.

- Asymmetric velocity models with a spiral velocity perturbation along either (2) the azimuthal direction ($\Delta V_\theta = \Delta V \times u \cdot u_\theta / |u \cdot u_\theta|$, $\Delta V_R = 0$), $u_\theta$ being the unitary azimuthal vector, or (3) the radial direction ($\Delta V_R = \Delta V \times u \cdot u_R / |u \cdot u_R|$, $\Delta V_\theta = 0$), $u_R$ being the unitary radial vector. These cases are necessary to identify the impact of velocity gradients along one planar direction, independently from the other direction. This thus represents a departure from $V_\theta(R)$ of at least ~13%, and significantly more from $V_R(R)$.

- Asymmetric velocity models with a spiral velocity perturbation along both azimuthal and radial directions for the
<figure><img src="image_19.png" /><figcaption>**Fig. B.2.** Kinematics maps and velocity dispersion Fourier analysis histograms for the toy models with a weak velocity gradient rotation curve ($v_0 = 200$ km s$^{-1}$, $r_s = 100''$ and $\gamma = 1$) without any perturbation (top), and with a spiral velocity perturbation towards either radial (middle) or azimuthal (bottom) components. We show, from left to right the density map, the velocity field and the BS corrected dispersion maps in the isotropic and anisotropic cases, obtained with an inclination of 60° and with $\phi_{sp} = 135°$. We show from left to right amplitude histograms of orders $k=2$ and $k=4$ and normalised phase angle histograms of these orders, obtained from the FFT analysis of toy models dispersion maps. Blue and green histograms correspond to isotropic and anisotropic cases respectively.</figcaption></figure>
<figure><img src="image_19.png" /><figcaption>**Fig. B.2.** Kinematics maps and velocity dispersion Fourier analysis histograms for the toy models with a weak velocity gradient rotation curve ($v_0 = 200$ km s$^{-1}$, $r_s = 100''$ and $\gamma = 1$) without any perturbation (top), and with a spiral velocity perturbation towards either radial (middle) or azimuthal (bottom) components. We show, from left to right the density map, the velocity field and the BS corrected dispersion maps in the isotropic and anisotropic cases, obtained with an inclination of 60° and with $\phi_{sp} = 135°$. We show from left to right amplitude histograms of orders $k=2$ and $k=4$ and normalised phase angle histograms of these orders, obtained from the FFT analysis of toy models dispersion maps. Blue and green histograms correspond to isotropic and anisotropic cases respectively.</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>**Fig. B.3.** Same as Fig. B.2, but for the toy models with a steep velocity gradient rotation curve ($v_0 = 250$ km s$^{-1}$, $r_s = 20''$ and $\gamma = 2$).</figcaption></figure>

cases (4) $\Delta V_{\theta} = 2\Delta V_{R} = 2\Delta V \times u \cdot u_{\theta} / |u \cdot u_{\theta}|$, (5) $\Delta V_{R} = 2\Delta V_{\theta} = 2\Delta V \times u \cdot u_{R} / |u \cdot u_{R}|$, and (6) with the spiral velocity perturbation $\Delta V$ along $u$, that is $\Delta V_{R} = \Delta V \times u \cdot u_{R}$ and $\Delta V_{\theta} = \Delta V \times u \cdot u_{\theta}$. In the latter case, the orientation of

the velocity perturbation changes with radius, being more azimuthal in the inner parts and more radial in the outer parts. These three models are useful to identify the combined effect of planar asymmetries with different strengths.
<figure><img src="image_20.png" /><figcaption>**Fig. B.3.** Same as Fig. B.2, but for the toy models with a steep velocity gradient rotation curve ($v_0 = 250$ km s$^{-1}$, $r_s = 20''$ and $\gamma = 2$).</figcaption></figure>

cases (4) $\Delta V_{\theta} = 2\Delta V_{R} = 2\Delta V \times u \cdot u_{\theta} / |u \cdot u_{\theta}|$, (5) $\Delta V_{R} = 2\Delta V_{\theta} = 2\Delta V \times u \cdot u_{R} / |u \cdot u_{R}|$, and (6) with the spiral velocity perturbation $\Delta V$ along $u$, that is $\Delta V_{R} = \Delta V \times u \cdot u_{R}$ and $\Delta V_{\theta} = \Delta V \times u \cdot u_{\theta}$. In the latter case, the orientation of

the velocity perturbation changes with radius, being more azimuthal in the inner parts and more radial in the outer parts. These three models are useful to identify the combined effect of planar asymmetries with different strengths.
The azimuthal velocity at the coordinate R, θ was thus obtained by $V_{\theta}(R, \theta) = V_{\theta}(R) \pm |\Delta V_{\theta}(R, \theta)|$, and similarly for the radial component. This creates a sharp kinematic discontinuity through the spiral feature, and the asymmetries in $V_{\theta}$ or $V_{R}$ make the velocity ellipsoid intrinsically anisotropic on scales larger than that of single gas clouds. To illustrate this, Fig. B.1 shows the velocity ellipsoid in the case of local velocity isotropy without large-scale anisotropy due to a spiral streaming (top panel), and in the case of local velocity isotropy but with anisotropy induced at large scale by the $V_{R}$ streaming described above (bottom panel). The signature of the radial bias of the anisotropy is clearly seen as contours of density elongated in the radial dimension of the velocity tensor.

Each particle was then assigned $V_{\theta}$, $V_{R}$ and $V_{z}$ components, which were drawn randomly following a Gaussian probability distribution centred on the local values described above, $V_{\theta}(R, \theta)$, $V_{R}(R, \theta)$ and $V_{z}(R)$, and with standard deviation $\sigma_{\theta}$, $\sigma_{R}$ and $\sigma_{z}$, respectively. The two following cases have been explored: $\sigma_{\theta} = \sigma_{R} = \sigma_{z} = 8$ km s$^{-1}$, corresponding to traditional isotropic velocity ellipsoids in the absence of asymmetries, and $\sigma_{R} = 8$ km s$^{-1}$, $\sigma_{\theta} = 0.7\sigma_{R} = 5.6$ km s$^{-1}$, and $\sigma_{z} = 5$ km s$^{-1}$, which is a direct way to also generate a velocity anisotropy, uniformly, and in addition to the velocity streamings. To keep the modelling simple, we thus did not allow direct variations of $\sigma_{\theta}$ and $\sigma_{R}$ within the disc. Combined to the six previously enumerated configurations, this leads to a total of 12 cases for each rotation curve. We did not explore the possibility $\sigma_{\theta} > \sigma_{R}$, as this situation only occasionally occurs in the stellar populations of the Milky Way or the Large Magellanic Cloud (Gaia Collaboration 2021b, 2023).

Then, the mock cubes of data were built, giving each particle a $V_{\text{los}}$ following Eq. C.1. Each data cube contains 200 spectral channels of 3 km s$^{-1}$ width, and 400 × 400 squared pixels, with a pixel scale of 1". The assigned channel map of a particle corresponds to the one closest to $V_{\text{los}}$, and the line profile in a spatial pixel was obtained by summing the counts of all individual particles from this pixel. Constant position angle of the disc major axis and disc inclination were assumed as a function radius. The position angle was fixed at 45° in all models, and inclinations of 45°, 60°, and 75° were used. For each inclination, 12 mock data cubes were generated, corresponding to the 12 orientations of the spiral perturbation described above. To mimic an observational BS effect, the data cubes were then smoothed by a Gaussian function with a FWHM of 8 pixels (corresponding to ≈465 pc, assuming a distance of 12 Mpc). These are idealised data cubes to which no further smoothing by synthesised light and point spread functions were applied, and no further random noise was added. The first and second moments of the data cubes were then derived, the BS correction applied, using the low-resolution flux and velocity maps, to remain consistent with the methodology applied to the observed data.

In Figures B.2 (weak velocity gradient, moderate maximum rotation velocity) and B.3 (steep inner gradient, high maximum rotation velocity reached at small radius), we show both the effects of models of uniform isotropy and anisotropy for (1) the axisymmetric configuration, and for the spiral velocity perturbation along either (2) radial or (3) azimuthal directions. An example with the velocity perturbation along **u** is shown in Fig. 3 for the rotation curve with the steep gradient. We show for each case the maps obtained with an inclination of 60° and with $\phi_{sp} = 135°$. Whereas no obvious BS residuals due to rotation are visible in Fig. B.2, stronger residuals are observed in Fig. B.3, especially in the inner parts. This is due to both the shape and amplitude of the rotation curve. Nevertheless, we checked

that using the moment maps of the high-resolution mock data cubes, that is those obtained before the spatial smoothing, leads, as expected, to an accurate correction of the BS effect (see Sect. 3.1), except at the very centre for the steeper rotation curve. We also point out that when the spiral perturbation is set along the radial (azimuthal) direction, the impact is minimum along the major (minor) axis in both the velocity field and the velocity dispersion map. Fig. B.2 also shows the histograms of amplitudes and phases for orders two and four resulting from the Fourier analysis of the corresponding models, described and discussed in Sections 5.2.1 and 5.3.1.

## Appendix C: Line-of-sight velocity dispersion for a disc in cylindrical coordinates

In this Appendix, we introduce and discuss the geometrical and mathematical framework used to study asymmetries and anisotropies in the neutral gas spatially resolved velocity dispersion fields. Velocity dispersion is a second order moment, and therefore is a more complex quantity to study than flux and radial velocity which are respectively the zero and first moment order.

Due to its hydrodynamical properties, the gas likely lies in a plane. The velocity vector in the frame of the galactic plane is described by two components both lying in this plane. Thus, in the cylindrical frame, the observed velocity along the line-of-sight is expressed as:

$$V_{\text{los}} = V_{\text{sys}} + V_{\theta} \cos \theta \sin i + V_{R} \sin \theta \sin i + V_{z} \cos i, \quad (C.1)$$

where $V_{\theta}$ is the azimuthal velocity, that is the rotation motions; $V_{R}$ is the radial velocity, that is the inward or outward motions; and $V_{z}$ is the vertical velocity, that is motions perpendicular to the galactic plane.

The observed velocity $V_{\text{los}}$ (los for line-of-sight) is linked to the projection of $V_{\theta}$, $V_{R}$ and $V_{z}$ along the line-of-sight through five additional parameters: PA, the position angle of the major axis of the galaxy (measured anticlockwise from the north to the direction of receding side of the galaxy); i, the inclination of the galactic disc with respect to the sky plane; $V_{\text{sys}}$, the systemic velocity of the galaxy; and $x_{c}$ and $y_{c}$, the coordinates of the rotation centre in Cartesian coordinates (sky projection).

Both the radial, azimuthal, and vertical components can vary with R and θ, which are the polar coordinates in the plane of the galaxy with respect to the centre, choosing the major axis as reference θ = 0 (receding side). The azimuth in the plane of the galaxy, θ, is linked to the position angle PA, the inclination i, the position x (east-west), y (north-south) and centre $x_{c}$, $y_{c}$ in the sky by the set of equations C.2 to C.7:

$$R \cos \theta = r \cos \psi, \quad (C.2)$$

$$R \sin \theta = r \frac{\sin \psi}{\cos i}, \quad (C.3)$$

$$\cos \psi = \frac{(y - y_c) \cos PA - (x - x_c) \sin PA}{r}, \quad (C.4)$$

$$\sin \psi = - \frac{(x - x_c) \cos PA + (y - y_c) \sin PA}{r}, \quad (C.5)$$

$$r = \sqrt{(x - x_c)^2 + (y - y_c)^2}, \quad (C.6)$$

$$R = r \sqrt{\cos^2 \psi + \frac{\sin^2 \psi}{\cos^2 i}}, \quad (C.7)$$

with ψ being the anticlockwise angle from the major axis, and r being the distance to the centre, both in the plane of the sky.
The azimuthal velocity at the coordinate R, θ was thus obtained by $V_{\theta}(R, \theta) = V_{\theta}(R) \pm |\Delta V_{\theta}(R, \theta)|$, and similarly for the radial component. This creates a sharp kinematic discontinuity through the spiral feature, and the asymmetries in $V_{\theta}$ or $V_{R}$ make the velocity ellipsoid intrinsically anisotropic on scales larger than that of single gas clouds. To illustrate this, Fig. B.1 shows the velocity ellipsoid in the case of local velocity isotropy without large-scale anisotropy due to a spiral streaming (top panel), and in the case of local velocity isotropy but with anisotropy induced at large scale by the $V_{R}$ streaming described above (bottom panel). The signature of the radial bias of the anisotropy is clearly seen as contours of density elongated in the radial dimension of the velocity tensor.

Each particle was then assigned $V_{\theta}$, $V_{R}$ and $V_{z}$ components, which were drawn randomly following a Gaussian probability distribution centred on the local values described above, $V_{\theta}(R, \theta)$, $V_{R}(R, \theta)$ and $V_{z}(R)$, and with standard deviation $\sigma_{\theta}$, $\sigma_{R}$ and $\sigma_{z}$, respectively. The two following cases have been explored: $\sigma_{\theta} = \sigma_{R} = \sigma_{z} = 8$ km s$^{-1}$, corresponding to traditional isotropic velocity ellipsoids in the absence of asymmetries, and $\sigma_{R} = 8$ km s$^{-1}$, $\sigma_{\theta} = 0.7\sigma_{R} = 5.6$ km s$^{-1}$, and $\sigma_{z} = 5$ km s$^{-1}$, which is a direct way to also generate a velocity anisotropy, uniformly, and in addition to the velocity streamings. To keep the modelling simple, we thus did not allow direct variations of $\sigma_{\theta}$ and $\sigma_{R}$ within the disc. Combined to the six previously enumerated configurations, this leads to a total of 12 cases for each rotation curve. We did not explore the possibility $\sigma_{\theta} > \sigma_{R}$, as this situation only occasionally occurs in the stellar populations of the Milky Way or the Large Magellanic Cloud (Gaia Collaboration 2021b, 2023).

Then, the mock cubes of data were built, giving each particle a $V_{\text{los}}$ following Eq. C.1. Each data cube contains 200 spectral channels of 3 km s$^{-1}$ width, and 400 × 400 squared pixels, with a pixel scale of 1". The assigned channel map of a particle corresponds to the one closest to $V_{\text{los}}$, and the line profile in a spatial pixel was obtained by summing the counts of all individual particles from this pixel. Constant position angle of the disc major axis and disc inclination were assumed as a function radius. The position angle was fixed at 45° in all models, and inclinations of 45°, 60°, and 75° were used. For each inclination, 12 mock data cubes were generated, corresponding to the 12 orientations of the spiral perturbation described above. To mimic an observational BS effect, the data cubes were then smoothed by a Gaussian function with a FWHM of 8 pixels (corresponding to ≈465 pc, assuming a distance of 12 Mpc). These are idealised data cubes to which no further smoothing by synthesised light and point spread functions were applied, and no further random noise was added. The first and second moments of the data cubes were then derived, the BS correction applied, using the low-resolution flux and velocity maps, to remain consistent with the methodology applied to the observed data.

In Figures B.2 (weak velocity gradient, moderate maximum rotation velocity) and B.3 (steep inner gradient, high maximum rotation velocity reached at small radius), we show both the effects of models of uniform isotropy and anisotropy for (1) the axisymmetric configuration, and for the spiral velocity perturbation along either (2) radial or (3) azimuthal directions. An example with the velocity perturbation along **u** is shown in Fig. 3 for the rotation curve with the steep gradient. We show for each case the maps obtained with an inclination of 60° and with $\phi_{sp} = 135°$. Whereas no obvious BS residuals due to rotation are visible in Fig. B.2, stronger residuals are observed in Fig. B.3, especially in the inner parts. This is due to both the shape and amplitude of the rotation curve. Nevertheless, we checked

that using the moment maps of the high-resolution mock data cubes, that is those obtained before the spatial smoothing, leads, as expected, to an accurate correction of the BS effect (see Sect. 3.1), except at the very centre for the steeper rotation curve. We also point out that when the spiral perturbation is set along the radial (azimuthal) direction, the impact is minimum along the major (minor) axis in both the velocity field and the velocity dispersion map. Fig. B.2 also shows the histograms of amplitudes and phases for orders two and four resulting from the Fourier analysis of the corresponding models, described and discussed in Sections 5.2.1 and 5.3.1.

## Appendix C: Line-of-sight velocity dispersion for a disc in cylindrical coordinates

In this Appendix, we introduce and discuss the geometrical and mathematical framework used to study asymmetries and anisotropies in the neutral gas spatially resolved velocity dispersion fields. Velocity dispersion is a second order moment, and therefore is a more complex quantity to study than flux and radial velocity which are respectively the zero and first moment order.

Due to its hydrodynamical properties, the gas likely lies in a plane. The velocity vector in the frame of the galactic plane is described by two components both lying in this plane. Thus, in the cylindrical frame, the observed velocity along the line-of-sight is expressed as:

$$V_{\text{los}} = V_{\text{sys}} + V_{\theta} \cos \theta \sin i + V_{R} \sin \theta \sin i + V_{z} \cos i, \quad (C.1)$$

where $V_{\theta}$ is the azimuthal velocity, that is the rotation motions; $V_{R}$ is the radial velocity, that is the inward or outward motions; and $V_{z}$ is the vertical velocity, that is motions perpendicular to the galactic plane.

The observed velocity $V_{\text{los}}$ (los for line-of-sight) is linked to the projection of $V_{\theta}$, $V_{R}$ and $V_{z}$ along the line-of-sight through five additional parameters: PA, the position angle of the major axis of the galaxy (measured anticlockwise from the north to the direction of receding side of the galaxy); i, the inclination of the galactic disc with respect to the sky plane; $V_{\text{sys}}$, the systemic velocity of the galaxy; and $x_{c}$ and $y_{c}$, the coordinates of the rotation centre in Cartesian coordinates (sky projection).

Both the radial, azimuthal, and vertical components can vary with R and θ, which are the polar coordinates in the plane of the galaxy with respect to the centre, choosing the major axis as reference θ = 0 (receding side). The azimuth in the plane of the galaxy, θ, is linked to the position angle PA, the inclination i, the position x (east-west), y (north-south) and centre $x_{c}$, $y_{c}$ in the sky by the set of equations C.2 to C.7:

$$R \cos \theta = r \cos \psi, \quad (C.2)$$

$$R \sin \theta = r \frac{\sin \psi}{\cos i}, \quad (C.3)$$

$$\cos \psi = \frac{(y - y_c) \cos PA - (x - x_c) \sin PA}{r}, \quad (C.4)$$

$$\sin \psi = - \frac{(x - x_c) \cos PA + (y - y_c) \sin PA}{r}, \quad (C.5)$$

$$r = \sqrt{(x - x_c)^2 + (y - y_c)^2}, \quad (C.6)$$

$$R = r \sqrt{\cos^2 \psi + \frac{\sin^2 \psi}{\cos^2 i}}, \quad (C.7)$$

with ψ being the anticlockwise angle from the major axis, and r being the distance to the centre, both in the plane of the sky.
It is possible to define velocity dispersion components:

$$ \sigma_{ij}^2 = \overline{V_i V_j} - \overline{V_i} \overline{V_j}, \qquad (C.8) $$

where the subscripts i and j denote the different coordinate directions<fnref n="6" />: radial (R), azimuthal (θ), and vertical (z), whereas $V_i$ and $V_j$ give the corresponding velocities. When $i \neq j$, the term is the covariance between $V_i$ and $V_j$, whereas in the specific case of $i = j$, the term is the variance of $V_i$ and we note $\sigma_{ii}^2 = \sigma_i^2$. By definition, the line-of-sight velocity dispersion is thus defined as

$$ \sigma_{\text{los}}^2 = \overline{V_{\text{los}}^2} - \overline{V_{\text{los}}}^2, \qquad (C.9) $$

which leads, by replacing the terms by their corresponding expression using Eq. C.1, by using Eq. C.8, and because $V_{\text{sys}}$ is constant, to:

$$ \sigma_{\text{los}}^2 = \sigma_{\theta}^2 \cos^2 \theta \sin^2 i + \sigma_R^2 \sin^2 \theta \sin^2 i + \sigma_z^2 \cos^2 i \\ + 2(\sigma_{R\theta}^2 \cos \theta \sin \theta \sin^2 i + \sigma_{\theta z}^2 \sin \theta \cos i \sin i) \\ + \sigma_{Rz}^2 \cos \theta \cos i \sin i, \qquad (\text{C.10}) $$

where $\sigma_{R\theta}^2$, $\sigma_{Rz}^2$ and $\sigma_{\theta z}^2$ are the squared cross terms of the velocity dispersion tensor. Those terms can be expressed as $\sigma_{R\theta}^2 = \rho_{R\theta}\sigma_R\sigma_\theta$, $\sigma_{Rz}^2 = \rho_{Rz}\sigma_R\sigma_z$, $\sigma_{\theta z}^2 = \rho_{\theta z}\sigma_\theta\sigma_z$, owing to the definition of the correlation coefficients between the components $\rho_{R\theta}$, $\rho_{Rz}$ and $\rho_{\theta z}$ that are bounded between -1 and 1. The term $\sigma_{R\theta}^2$ can thus be negative, depending on the sign of the correlation coefficient.

Equation C.10 can be recast in a sum of trigonometric polynomials of second degree:

$$ \sigma_{\text{los}}^2 = \left( \frac{\sigma_{\theta}^2 - \sigma_R^2}{2} \cos 2\theta + \sigma_{R\theta}^2 \sin 2\theta \right) \sin^2 i \\ + \left( \sigma_{Rz}^2 \cos \theta + \sigma_{\theta z}^2 \sin \theta \right) \sin 2i \\ + \sigma_z^2 \cos^2 i + \frac{\sigma_{\theta z}^2 + \sigma_R^2}{2} \sin^2 i. \qquad (\text{C.11}) $$

We can formulate differently this equation to match better the formalism introduced with Fourier Transforms (see Sect. 3.2):

$$ \sigma_{\text{los}}^2 = k_0 + k_1 \cos(\theta - \alpha_1) + k_2 \cos(2(\theta - \alpha_2)), \qquad (\text{C.12}) $$

with:

$$ k_0 = \sigma_z^2 \cos^2 i + \frac{\sigma_{\theta}^2 + \sigma_R^2}{2} \sin^2 i, \qquad (\text{C.13}) $$

$$ k_1 = \sqrt{\sigma_{Rz}^4 + \sigma_{\theta z}^4} \times \sin 2i, \qquad (\text{C.14}) $$

$$ \cos \alpha_1 = \frac{\sigma_{Rz}^2}{k_1} \sin 2i, \qquad (\text{C.15}) $$

$$ \sin \alpha_1 = \frac{\sigma_{\theta z}^2}{k_1} \sin 2i, \qquad (\text{C.16}) $$

$$ k_2 = \sqrt{\left(\frac{\sigma_{\theta}^2 - \sigma_R^2}{2}\right)^2 + \sigma_{R\theta}^4} \times \sin^2 i, \qquad (\text{C.17}) $$

$$ \cos 2\alpha_2 = \frac{(\sigma_{\theta}^2 - \sigma_R^2)}{2k_2} \sin^2 i, \qquad (\text{C.18}) $$

$$ \sin 2\alpha_2 = \frac{\sigma_{R\theta}^2}{k_2} \sin^2 i. \qquad (\text{C.19}) $$

This formalism further motivated our choice of performing Fourier Transforms on squared velocity dispersion. With these definitions, we find

$$ \tan 2\alpha_2 = 2\sigma_{R\theta}^2/(\sigma_{\theta}^2 - \sigma_R^2), \qquad (\text{C.20}) $$

so that $\alpha_2$ is the tilt angle of the velocity ellipsoid in $(R, \theta)$ (see e.g. Smith et al. 2009), that is, the vertex deviation (Kuijken & Tremaine 1991, 1994; Binney & Merrifield 1998), analogue to the tilt angle of the velocity ellipsoid in $(R, z)$ plane used, for example, in Hagen et al. (2019).

In case there are asymmetries, those are embedded in the terms $k_i$ and $\alpha_i$. For the sake of simplicity, let's assume that cross terms are null, which leads to

$$ k_1 = 0, k_2 = \frac{\sigma_{\theta}^2 - \sigma_R^2}{2} \times \sin^2 i, \text{ and } \alpha_2 = 0. $$

If we only have an asymmetry of order $k$ (sinusoidal) in both $\sigma_R^2$ and $\sigma_{\theta z}^2$ with similar amplitudes $\sigma_{R,k} = \sigma_{\theta,k} = \sigma_k$ and with phases $\phi_{R,k}$ and $\phi_{\theta,k} - \phi_{R,k} = \Delta\phi_k$, we can therefore express $k_0$ and $k_2$ as:

$$ k_0 = \sigma_z^2 \cos^2 i + \sigma_k^2 \cos(\Delta\phi_k/2) \cos(k\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i, \qquad (\text{C.21}) $$

$$ k_2 = \sigma_k^2 \sin(\Delta\phi_k/2) \sin(k\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i, \qquad (\text{C.22}) $$

and Eq. C.12 becomes:

$$ \sigma_{\text{los}}^2 = \sigma_z^2 \cos^2 i + \sigma_k^2 \cos(\Delta\phi_k/2) \cos(k\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i \\ + \sigma_k^2/2 \sin(\Delta\phi_k/2) \sin((k+2)\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i \\ + \sigma_k^2/2 \sin(\Delta\phi_k/2) \sin((k-2)\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i. \qquad (\text{C.23}) $$

Orders $k, k-2$ and $k+2$ appear. Nevertheless, even if $\Delta\phi_k$ is constant, the phase $\phi_{R,k}$ is supposed to remain random, so no peculiar direction should be favoured.

The general case without the assumption that $\sigma_{R,k} = \sigma_{\theta,k}$ leads to:

$$ k_2 = c_k \cos(k\theta - \phi_{R,k}) + s_k \sin(k\theta - \phi_{R,k}), \qquad (\text{C.24}) $$

with

$$ c_k = \frac{\sigma_{\theta,k}^2 \cos \Delta\phi_k - \sigma_{R,k}^2}{2}, \quad \text{and} \qquad (\text{C.25}) $$

$$ s_k = \frac{\sigma_{\theta,k}^2 \sin \Delta\phi_k}{2}, \qquad (\text{C.26}) $$

which can also be written as

$$ k_2 = x_k \cos(k\theta - \phi_{R,k} - \phi_k), \qquad (\text{C.27}) $$

with

$$ x_k = \sqrt{c_k^2 + s_k^2}, \quad \cos \phi_k = c_k/x_k, \quad \text{and} \quad \sin \phi_k = s_k/x_k. $$

The general case would therefore also lead to orders $k-2$ and $k+2$.

If $\sigma_{R,k} = \sigma_{\theta,k} = \sigma_k$, we obtain that $c_k = -\sigma_k^2 \sin^2(\Delta\phi_k/2)$ and $s_k = \sigma_k^2 \sin(\Delta\phi_k/2) \cos(\Delta\phi_k/2)$, and therefore $x_k = \sigma_k^2 \sin(\Delta\phi_k/2)$ and $\phi_k = \Delta\phi_k/2 + \pi/2$, which enables us to recover Eq. C.22.

<footnote n="6">⁶ by construction, $\sigma_{ij}^2$ is not necessarily positive, this notation was nevertheless used for the sake of homogeneity and because it was used by other authors (e.g. Smith et al. 2009).</footnote>
It is possible to define velocity dispersion components:

$$ \sigma_{ij}^2 = \overline{V_i V_j} - \overline{V_i} \overline{V_j}, \qquad (C.8) $$

where the subscripts i and j denote the different coordinate directions<fnref n="6" />: radial (R), azimuthal (θ), and vertical (z), whereas Vᵢ and Vⱼ give the corresponding velocities. When i ≠ j, the term is the covariance between Vᵢ and Vⱼ, whereas in the specific case of i = j, the term is the variance of Vᵢ and we note σᵢᵢ² = σᵢᵢ². By definition, the line-of-sight velocity dispersion is thus defined as

$$ \sigma_{\text{los}}^2 = \overline{V_{\text{los}}^2} - \overline{V_{\text{los}}}^2, \qquad (C.9) $$

which leads, by replacing the terms by their corresponding expression using Eq. C.1, by using Eq. C.8, and because Vsys is constant, to:

$$ \sigma_{\text{los}}^2 = \sigma_{\theta}^2 \cos^2 \theta \sin^2 i + \sigma_R^2 \sin^2 \theta \sin^2 i + \sigma_z^2 \cos^2 i + 2(\sigma_{R\theta}^2 \cos \theta \sin \theta \sin^2 i + \sigma_{\theta z}^2 \sin \theta \cos i \sin i + \sigma_{Rz}^2 \cos \theta \cos i \sin i), \qquad (\text{C.10}) $$

where σ<sub>Rθ</sub><sup>2</sup>, σ<sub>Rz</sub><sup>2</sup> and σ<sub>θz</sub><sup>2</sup> are the squared cross terms of the velocity dispersion tensor. Those terms can be expressed as σ<sub>Rθ</sub><sup>2</sup> = ρ<sub>Rθ</sub>σ<sub>R</sub>σ<sub>θ</sub>, σ<sub>Rz</sub><sup>2</sup> = ρ<sub>Rz</sub>σ<sub>R</sub>σ<sub>z</sub>, σ<sub>θz</sub><sup>2</sup> = ρ<sub>θz</sub>σ<sub>θ</sub>σ<sub>z</sub>, owing to the definition of the correlation coefficients between the components ρ<sub>Rθ</sub>, ρ<sub>Rz</sub> and ρ<sub>θz</sub> that are bounded between -1 and 1. The term σ<sub>Rθ</sub><sup>2</sup> can thus be negative, depending on the sign of the correlation coefficient.

Equation C.10 can be recast in a sum of trigonometric polynomials of second degree:

$$ \sigma_{\text{los}}^2 = \left( \frac{\sigma_{\theta}^2 - \sigma_R^2}{2} \cos 2\theta + \sigma_{R\theta}^2 \sin 2\theta \right) \sin^2 i + \left( \sigma_{Rz}^2 \cos \theta + \sigma_{\theta z}^2 \sin \theta \right) \sin 2i + \sigma_z^2 \cos^2 i + \frac{\sigma_{\theta}^2 + \sigma_R^2}{2} \sin^2 i. \qquad (\text{C.11}) $$

We can formulate differently this equation to match better the formalism introduced with Fourier Transforms (see Sect. 3.2):

$$ \sigma_{\text{los}}^2 = k_0 + k_1 \cos(\theta - \alpha_1) + k_2 \cos(2(\theta - \alpha_2)), \qquad (\text{C.12}) $$

with:

$$ k_0 = \sigma_z^2 \cos^2 i + \frac{\sigma_{\theta}^2 + \sigma_R^2}{2} \sin^2 i, \qquad (\text{C.13}) $$

$$ k_1 = \sqrt{\sigma_{Rz}^4 + \sigma_{\theta z}^4} \times \sin 2i, \qquad (\text{C.14}) $$

$$ \cos \alpha_1 = \frac{\sigma_{Rz}^2}{k_1} \sin 2i, \qquad (\text{C.15}) $$

$$ \sin \alpha_1 = \frac{\sigma_{\theta z}^2}{k_1} \sin 2i, \qquad (\text{C.16}) $$

$$ k_2 = \sqrt{\left(\frac{\sigma_{\theta}^2 - \sigma_R^2}{2}\right)^2 + \sigma_{R\theta}^4} \times \sin^2 i, \qquad (\text{C.17}) $$

$$ \cos 2\alpha_2 = \frac{(\sigma_{\theta}^2 - \sigma_R^2)}{2k_2} \sin^2 i, \qquad (\text{C.18}) $$

$$ \sin 2\alpha_2 = \frac{\sigma_{R\theta}^2}{k_2} \sin^2 i. \qquad (\text{C.19}) $$

This formalism further motivated our choice of performing Fourier Transforms on squared velocity dispersion. With these definitions, we find

$$ \tan 2\alpha_2 = 2\sigma_{R\theta}^2/(\sigma_{\theta}^2 - \sigma_R^2), \qquad (\text{C.20}) $$

so that α₂ is the tilt angle of the velocity ellipsoid in (R, θ) (see e.g. Smith et al. 2009), that is, the vertex deviation (Kuijken & Tremaine 1991, 1994; Binney & Merrifield 1998), analogue to the tilt angle of the velocity ellipsoid in (R, z) plane used, for example, in Hagen et al. (2019).

In case there are asymmetries, those are embedded in the terms kᵢ and αᵢ. For the sake of simplicity, let's assume that cross terms are null, which leads to

$$ k_1 = 0, k_2 = \frac{\sigma_{\theta}^2 - \sigma_R^2}{2} \times \sin^2 i, \text{ and } \alpha_2 = 0. $$

If we only have an asymmetry of order k (sinusoidal) in both σ<sub>R</sub><sup>2</sup> and σ<sub>θ</sub><sup>2</sup> with similar amplitudes σ<sub>R,k</sub> = σ<sub>θ,k</sub> = σ<sub>k</sub> and with phases φ<sub>R,k</sub> and φ<sub>θ,k</sub> - φ<sub>R,k</sub> = Δφ<sub>k</sub>, we can therefore express k₀ and k₂ as:

$$ k_0 = \sigma_z^2 \cos^2 i + \sigma_k^2 \cos(\Delta\phi_k/2) \cos(k\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i, \qquad (\text{C.21}) $$

$$ k_2 = \sigma_k^2 \sin(\Delta\phi_k/2) \sin(k\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i, \qquad (\text{C.22}) $$

and Eq. C.12 becomes:

$$ \sigma_{\text{los}}^2 = \sigma_z^2 \cos^2 i + \sigma_k^2 \cos(\Delta\phi_k/2) \cos(k\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i + \sigma_k^2/2 \sin(\Delta\phi_k/2) \sin((k+2)\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i + \sigma_k^2/2 \sin(\Delta\phi_k/2) \sin((k-2)\theta - \phi_{R,k} - \Delta\phi_k/2) \sin^2 i. \qquad (\text{C.23}) $$

Orders k, k - 2 and k + 2 appear. Nevertheless, even if Δφₖ is constant, the phase φₖ,ₖ is supposed to remain random, so no peculiar direction should be favoured.

The general case without the assumption that σ<sub>R,k</sub> = σ<sub>θ,k</sub> leads to:

$$ k_2 = c_k \cos(k\theta - \phi_{R,k}) + s_k \sin(k\theta - \phi_{R,k}), \qquad (\text{C.24}) $$

with

$$ c_k = \frac{\sigma_{\theta,k}^2 \cos \Delta\phi_k - \sigma_{R,k}^2}{2}, \text{ and} \qquad (\text{C.25}) $$

$$ s_k = \frac{\sigma_{\theta,k}^2 \sin \Delta\phi_k}{2}, \qquad (\text{C.26}) $$

which can also be written as

$$ k_2 = x_k \cos(k\theta - \phi_{R,k} - \phi_k), \qquad (\text{C.27}) $$

with

$$ x_k = \sqrt{c_k^2 + s_k^2}, \cos \phi_k = c_k/x_k, \text{ and } \sin \phi_k = s_k/x_k. $$

The general case would therefore also lead to orders k - 2 and k + 2.

If σ<sub>R,k</sub> = σ<sub>θ,k</sub> = σ<sub>k</sub>, we obtain that cₖ = -σₖ² sin²(Δφₖ/2) and sₖ = σₖ² sin(Δφₖ/2) cos(Δφₖ/2), and therefore xₖ = σₖ² sin(Δφₖ/2) and φₖ = Δφₖ/2 + π/2, which enables us to recover Eq. C.22.

<footnote n="6">⁶ by construction, σ<sub>ij</sub><sup>2</sup> is not necessarily positive, this notation was nevertheless used for the sake of homogeneity and because it was used by other authors (e.g. Smith et al. 2009).</footnote>
**Appendix D: Results of the FFTs for individual galaxies in THINGS**

In this section we show the various kinematic product for each THINGS galaxy as presented in Figure 1 and Figure 4. The upper row plots represent from left to right the observed flux density map, the observed velocity field, the observed velocity

dispersion, the BS model, and the velocity dispersion map corrected from the BS effect. The second row shows the observed squared velocity dispersion and its modelling through FFT up to $k = 4$, and their residuals. The third row shows the individual orders of the FFT in the galaxy. The bottom row presents the amplitudes (left) and phase angles (right) of the FFT coefficients as a function of radius.

<figure><img src="image_21.png" /></figure>

**Fig. D.1.** H I density and velocity maps and FFT results of NGC925. (Top panel, from left to right) Observed flux density map, observed velocity field, observed velocity dispersion, beam smearing model, and velocity dispersion map corrected from the beam smearing effect. (Second panel, from left to right) Squared observed velocity dispersion and its corresponding modelling through FFT up to the order 4, and the residuals between those two maps. (Third panel) Individual squared orders of the FFT modelling projected in the plane of the galaxy. (Bottom panel) Radial variation of the FFT amplitudes (left) and phase angles (right) for orders $k = 0$ to $k = 4$.
**Appendix D: Results of the FFTs for individual galaxies in THINGS**

In this section we show the various kinematic product for each THINGS galaxy as presented in Figure 1 and Figure 4. The upper row plots represent from left to right the observed flux density map, the observed velocity field, the observed velocity

dispersion, the BS model, and the velocity dispersion map corrected from the BS effect. The second row shows the observed squared velocity dispersion and its modelling through FFT up to $k = 4$, and their residuals. The third row shows the individual orders of the FFT in the galaxy. The bottom row presents the amplitudes (left) and phase angles (right) of the FFT coefficients as a function of radius.

<figure><img src="image_21.png" /><figcaption>**Fig. D.1.** H I density and velocity maps and FFT results of NGC925. (Top panel, from left to right) Observed flux density map, observed velocity field, observed velocity dispersion, beam smearing model, and velocity dispersion map corrected from the beam smearing effect. (Second panel, from left to right) Squared observed velocity dispersion and its corresponding modelling through FFT up to the order 4, and the residuals between those two maps. (Third panel) Individual squared orders of the FFT modelling projected in the plane of the galaxy. (Bottom panel) Radial variation of the FFT amplitudes (left) and phase angles (right) for orders $k = 0$ to $k = 4$.</figcaption></figure>
<figure><img src="image_22.png" /><figcaption>Fig. D.2. Same as Fig. D.1 for NGC2403.</figcaption></figure>
<figure><img src="image_22.png" /><figcaption>Fig. D.2. Same as Fig. D.1 for NGC2403.</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>Fig. D.3. Same as Fig. D.1 for NGC2903.</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>Fig. D.3. Same as Fig. D.1 for NGC2903.</figcaption></figure>
<figure><img src="image_24.png" /><figcaption>Fig. D.4. Same as Fig. D.1 for NGC2976.</figcaption></figure>
<figure><img src="image_24.png" /><figcaption>Fig. D.4. Same as Fig. D.1 for NGC2976.</figcaption></figure>
<figure><img src="image_25.png" /><figcaption>Fig. D.5. Same as Fig. D.1 for NGC3031.</figcaption></figure>
<figure><img src="image_25.png" /><figcaption>**Fig. D.5.** Same as Fig. D.1 for NGC3031.</figcaption></figure>
<figure><img src="image_26.png" /><figcaption>Fig. D.6. Same as Fig. D.1 for NGC3198.</figcaption></figure>
<figure><img src="image_26.png" /><figcaption>Fig. D.6. Same as Fig. D.1 for NGC3198.</figcaption></figure>
<figure><img src="image_27.png" /><figcaption>Fig. D.7. Same as Fig. D.1 for NGC3521.</figcaption></figure>
<figure><img src="image_27.png" /><figcaption>Fig. D.7. Same as Fig. D.1 for NGC3521.</figcaption></figure>
<figure><img src="image_28.png" /><figcaption>Fig. D.8. Same as Fig. D.1 for NGC3621.</figcaption></figure>
<figure><img src="image_28.png" /><figcaption>Fig. D.8. Same as Fig. D.1 for NGC3621.</figcaption></figure>
<figure><img src="image_29.png" /><figcaption>Fig. D.9. Same as Fig. D.1 for NGC3627.</figcaption></figure>
<figure><img src="image_29.png" /><figcaption>Fig. D.9. Same as Fig. D.1 for NGC3627.</figcaption></figure>
<figure><img src="image_30.png" /><figcaption>Fig. D.10. Same as Fig. D.1 for NGC4736.</figcaption></figure>
<figure><img src="image_30.png" /><figcaption>Fig. D.10. Same as Fig. D.1 for NGC4736.</figcaption></figure>
<figure><img src="image_31.png" /><figcaption>Fig. D.11. Same as Fig. D.1 for NGC5055.</figcaption></figure>
<figure><img src="image_31.png" /><figcaption>Fig. D.11. Same as Fig. D.1 for NGC5055.</figcaption></figure>
<figure><img src="image_32.png" /><figcaption>Fig. D.12. Same as Fig. D.1 for NGC6946.</figcaption></figure>
<figure><img src="image_32.png" /><figcaption>Fig. D.12. Same as Fig. D.1 for NGC6946.</figcaption></figure>
<figure><img src="image_33.png" /><figcaption>Fig. D.13. Same as Fig. D.1 for NGC7331.</figcaption></figure>
<figure><img src="image_33.png" /><figcaption>Fig. D.13. Same as Fig. D.1 for NGC7331.</figcaption></figure>
<figure><img src="image_34.png" /><figcaption>Fig. D.14. Same as Fig. D.1 for NGC7793.</figcaption></figure>
<figure><img src="image_34.png" /><figcaption>Fig. D.14. Same as Fig. D.1 for NGC7793.</figcaption></figure>