Fitting pseudo-S${\rm\acute{e}}$rsic (Spergel) light profiles to galaxies in interferometric data: the excellence of the $uv$-plane

As a generalization of the $r^{1/4}$ law, the $r^{1/n}$ profile first proposed by Sersic (1968) is one of the most common functions used to describe how the intensity of a galaxy varies with distance from its center. The surface density (or equivalently surface brightness) of the S${\rm\acute{e}}$rsic profile can be written as

where $R$ is the projected distance to the source center, $R_{\rm e}$ is the effective radius containing half of the total luminosity, $\Sigma_{\rm e}$ is the surface brightness at $R_{\rm e}$, and $n$ is the S${\rm\acute{e}}$rsic index which determines the shape of the light profile (see Fig. 1). The parameter $\kappa$ is a function of S${\rm\acute{e}}$rsic index and is such that $\Gamma(2n)=2\gamma(2n,\kappa)$, where $\Gamma$ and $\gamma$ represent the complete and incomplete gamma functions (Ciotti & Bertin, 1999), respectively. We will use the terms half-light radius or effective radius interchangeably to refer to the radius within which half of a galaxy’s luminosity is contained.

The S${\rm\acute{e}}$rsic index, $n$, determines the degree of curvature of the profile, with $n=0.5$ giving a Gaussian profile, $n=1$ an exponential disk profile and $n=4$ generally associated with galaxy bulges. As the index $n$ increases, the core steepens more rapidly for $R<R_{\rm e}$, and the intensity of the outer wing at $R>R_{\rm e}$ is significantly extended.

However, as mentioned in the introduction, the general S${\rm\acute{e}}$rsic profile is not analytically transformable in Fourier space for most values of the parameter $n$, as the parameter $\kappa$ cannot be solved in closed form. Various techniques have been developed to address this issue when calculations of the Fourier transform are needed, including numerical integration methods, approximations, and asymptotic expressions (e.g., Ciotti & Bertin, 1999; Mazure & Capelato, 2002; Baes & Gentile, 2011). The inability to solve the S${\rm\acute{e}}$rsic profile analytically in Fourier space creates challenges in certain scenarios. For example, when performing convolutions, such as corrections for seeing, or when directly working in Fourier space, like with interferometric data.

Spergel (2010) introduced an alternative to the S${\rm\acute{e}}$rsic model for galactic luminosity profiles with functional form

where $f_{\nu}(x)=\left(\frac{x}{2}\right)^{\nu}\frac{K_{\nu}(x)}{\Gamma(\nu+1)}$, $\Gamma$ is the Gamma function, $K_{\nu}$ is a modified spherical Bessel function of the third kind, $c_{\nu}$ is a constant, $R_{\rm e}$ is the half-light radius, and $\nu$ is known as Spergel index that controls the relative peakiness of the core and the relative prominence of the wings (similar to S${\rm\acute{e}}$rsic  $n$), with a theoretical limit of $\nu>-1$.

This family of functions is found to provide a good fit for galaxy light profiles and resembles the S${\rm\acute{e}}$rsic function over a range of indices. The Spergel profile at $\nu=0.5$ is identical to an exponential profile, which is equivalent to a S${\rm\acute{e}}$rsic profile with $n=1$. However, the two functions do not exactly coincide for different Spergel $\nu$ (see Fig. 1).

The Spergel profile has a significant advantage over the S${\rm\acute{e}}$rsic profile because it is analytic in both real space and Fourier space (Spergel, 2010), which means that it can be described mathematically using equations, making it easier to work with and analyze in the $uv$-plane.

In the left panel of Fig. 1, a comparison is shown between S${\rm\acute{e}}$rsic and Spergel profiles. This comparison covers a range of $n$ and $\nu$ indices. All the profiles are normalized at $R_{\rm e}$. Within a certain range of $\nu$, the Spergel profiles resemble S${\rm\acute{e}}$rsic profiles in shape such that a relation can be conceived between $\nu$ and $n$. For example, the Spergel profile at $\nu=-0.6$ and in the radial range near the effective radius, exhibits similarities to a de Vaucouleurs $n=4$ profile (see also Spergel, 2010). However, when elongating far from the normalization point, these profiles start to differ at both the innermost (i.e., $R<0.1\ R_{\rm e}$) and outermost (i.e., $R>5\ R_{\rm e}$) regions, indicating that converting one index into the other depends on the observed scales. This is generally the case for all S${\rm\acute{e}}$rsic models with $n>1$, respect to their first-order matching Spergel profiles with $-1<\nu<0.5$: they display steeper inner profiles and drop faster at large radii. This is discussed more further down in this Section.

We also see from Fig. 1 that Spergel models are not meant to reproduce S${\rm\acute{e}}$rsic models with a flatter shape than $n=1$, i.e., $n<1$. For example, even choosing $\nu=0.85$ the Spergel model only slightly flattens compared to the $\nu=0.5$ ($n=1$) case, and it remains far from resembling a Gaussian model $n=0.5$. Given that fitting Gaussian models in the $uv$-plane is a straightforward approach, we suggest using this method instead of attempting Spergel fits with high values of $\nu>1$.

In the right panel of Fig. 1, we show cumulative distributions for the S${\rm\acute{e}}$rsic and Spergel profiles. When truncated at $10\ R_{\rm e}$, an exponential profile (S${\rm\acute{e}}$rsic $n=1$ and Spergel $\nu=0.5$) retains almost all of its flux ($\sim 100\%$). In contrast, a S${\rm\acute{e}}$rsic profile with $n=4$ retains 96.1$\%$ of its flux, while 99.7$\%$ of the flux is contained within this radius for a Spergel profile with $\nu=-0.6$. At the small radius, within the inner $0.05\ R_{\rm e}$, an exponential profile contains 0.3$\%$ of the flux, while $n=4$ contains 3.2$\%$ of the flux. For a Spergel profile with $\nu=-0.6$, 5.4$\%$ of the flux is contained within this radius. As the index $n$ increases (or decreases for $\nu$), the differences between the fluxes contained within the radius of the small and large ends also somewhat increase, while remaining overall contained. A comparison of the S${\rm\acute{e}}$rsic profile with the Spergel profile using mathematical simulations is presented in Appendix A (see Fig. 16).

To empirically derive the conversion (spatial-scale dependent) between Spergel $\nu$ and the S${\rm\acute{e}}$rsic $n$, we create noise-free²²2 The model sources were created with an extremely high S/N to ensure a robust galfit measurement images of Spergel two-dimensional models, as described in the next section, and use galfit to measure the corresponding S${\rm\acute{e}}$rsic index. In Fig. 2, the cross symbols represent the galfit measurements, which we refer to the intrinsic best value and are labelled as $n_{\rm galfit}$, compared to the input $\nu$ for different source sizes of $R_{\rm e}/\theta_{\rm b}$ ranging from 0.1 to 2.0. Here $\theta_{\rm b}$ is defined as the synthesized circularized beam size, given by $\sqrt{ab}$, where $a$ and $b$ represent the FWHMs of the major and minor axes of the synthesized beam, respectively.

As a zero-order check of the procedure, Fig. 2 shows that for a Spergel profile with $\nu=0.5$, sources with different sizes all return a S${\rm\acute{e}}$rsic index of $n=1$ in galfit, as expected. This result is accurate within $3\%$, which represents the inherent maximal precision of this empirical calibration. We find that $\nu=-0.6$ corresponds to $n=4$ for R${}_{\rm e}/\theta_{\rm b}\sim 0.9$–1.0, which is close to when the FWHM scale of the galaxy and of the beam are identical, as expected (Spergel, 2010). However, when R${}_{\rm e}/\theta_{\rm b}<0.5$ an input $\nu=-0.6$ converts to $n\sim 2$–3, while models with $\nu=-0.6$ and R${}_{\rm e}/\theta_{\rm b}>1$ correspond to steeper $n>4$. The dependence of the conversion on R${}_{\rm e}/\theta_{\rm b}$ systematically decrease with increasing the Spergel indices, and nearly vanishes for $\nu>0$ (Fig. 2).

We find that the whole set of measurements can be well described by the form:

for sources with $R_{\rm e}/\theta_{\rm b}$ ranging from 0.1 to 2.0, respectively. The solid lines in Fig. 2 represent the best-fit relation with coefficients $p_{1}=0.0249$, $p_{2}=-7.72$, $p_{3}=0.191$, $p_{4}=-0.721$, and $p_{5}=1.32$. The residuals of the fit indicate that the uncertainties of S${\rm\acute{e}}$rsic $n$ are largely within 10% for model sources with $R_{\rm e}/\theta_{\rm b}$ of 0.1–2.0 at $\nu>-0.7$ (see Fig. 2). We confirm this trend when comparing the Spergel index with the S${\rm\acute{e}}$rsic index analytically through a mathematical matching between the S${\rm\acute{e}}$rsic and Spergel functions (see Fig. 17).

This analysis demonstrates that converting a Spergel index to a S${\rm\acute{e}}$rsic index depends on the ratio between the angular size of the galaxy and that of the beam of the observations being examined. We have mapped this relationship over a reasonably large range of this ratio, i.e., $R_{\rm e}/\theta_{\rm b}$ of 0.1–2.0. Anticipating results from forthcoming paper, which focuses on analyzing the morphologies of an ALMA archival sample of about 100 distant star-forming galaxies in the submillimeter bands (Q. Tan et al., in preparation), we find that around 82% (93%) sources fall within the range of $R_{\rm e}/\theta_{\rm b}=0.1-1.0\ (0.1-2.0)$, while approximately 88% (99%) of sources exhibit best-fitting $\nu>-0.7$ ($\nu<1$). The median $R_{\rm e}/\theta_{\rm b}$ is 0.32 and the semi-interquartal range is 0.2–0.6. This suggests that the calibration presented in Fig. 2 is representative of general observations of distant galaxies with ALMA and NOEMA.

The differences in the profiles imply that also for sizes (half-light radii) and total fluxes, a conversion might be required when comparing results based on Spergel to those from S${\rm\acute{e}}$rsic. Based on our simulations, we find that both the size and flux density estimated by galfit using a S${\rm\acute{e}}$rsic profile tend to be larger than when using Spergel, as the S${\rm\acute{e}}$rsic index becomes larger. Instead, axis ratios and position angles are unaffected. Fig. 3 shows the ratios of $R_{\rm e}$ and flux densities measured from fitting S${\rm\acute{e}}$rsic models to simulated galaxies created as Spergel models, as a function of Spergel index. Both the ratio of $R_{\rm e}$ and flux density exhibit a similar increasing trend as the profile get steeper. The ratio becomes larger for better resolved sources (i.e., larger $R_{\rm e}/\theta_{\rm b}$). For example, in the case of a large-sized ($R_{\rm e}/\theta_{\rm b}=1$) galaxy with a de Vacouleurs-like profile of $\nu=-0.6$, the $R_{\rm e}$ estimated by the S${\rm\acute{e}}$rsic model can be larger by about 30–40%, while the excess flux density is around 15–20%. By comparison, the relative increase in $R_{\rm e}$ and flux density are less significant for sources with flatter profile or smaller sizes compared to the beam.

To correct for these systematical biases and thus to enable comparison of measurements using Spergel and S${\rm\acute{e}}$rsic profiles, we fit the distribution of the ratio of half-light radius and total flux measured from S${\rm\acute{e}}$rsic to the Spergel-based input values. To distinguish the fitted parameters between those obtained from Spergel versus S${\rm\acute{e}}$rsic profile modeling, we label $R_{\rm e}$ and $S_{\rm sp}$ the half-light radius and total flux measured from Spergel profile fitting, respectively. For S${\rm\acute{e}}$rsic-based measurements, we use instead $R_{\rm e,se}$ and $S_{\rm se}$. We find that both ratios of $R_{\rm e,se}/R_{\rm e}$ and $S_{\rm se}/S_{\rm sp}$ can be accurately described by a similar form:

for sources with $R_{\rm e}/\theta_{\rm b}$ ranging from 0.1 to 2.0. For the size ratio of $R_{\rm e,se}/R_{\rm e}$, the best-fit gives coefficients $p_{1}=0.00138$, $p_{2}=-8.96$, $p_{3}=0.260$, $p_{4}=-0.0260$, and $p_{5}=0.996$, while for the flux ratio of $S_{\rm se}/S_{\rm sp}$, the best-fit coefficients are $p_{1}=0.00217$, $p_{2}=-7.43$, $p_{3}=0.149$, $p_{4}=0.00942$, and $p_{5}=1.00$ (see the solid lines in Fig. 3). Analytical calculations, matching the Spergel profile with S${\rm\acute{e}}$rsic profiles numerically (see Appendix A), fully confirm the trends encoded in the Eq.(4) above and in Fig. 3.

We caution that while we believe that our methodology captures the bulk of the systematic effects in the conversion as encoded in the $R_{\rm e}/\theta_{\rm b}$ ratio, some further systematics might be expected depending on higher order terms describing the actual shape of the beam. We derived best fitting parameters for Eq.(3) and Eq.(4) averaging over three different ALMA array configurations. By comparing to results from single ALMA array configurations, we estimate that further systematic uncertainties are small. The details of the three ALMA array configuration are summarized in Table 1 (see Appendix B).

Each simulated galaxy was created by inserting a model source signal into an empty dataset with realistic noise. The noise was obtained from real data using ALMA band 7 observed visibilities from galaxies in a recent survey (e.g., Puglisi et al., 2019, 2021; Valentino et al., 2020).

To analyze the data, the calibrated ALMA visibilities for several targets were exported from CASA using the exportuvfits. The exported data was then converted to $uv$-tables using the GILDAS fits$\_$to$\_$uvt task. The four spectral windows were combined using the uv$\_$continuum and uv$\_$merge tasks. Prior to introducing the model source into the $uv$-plane, any detected sources were subtracted from the visibility data using the best-fitting Spergel model to the visibilities. This produced a residual data set which, upon inspection, did not reveal any further source.

Fig. 4 (left) shows the map of the residual $uv$ data for a case in our simulations: the primary beam field of view (FOV) is 18^′′, and the FWHM of the synthesized beam using natural weighting is $\theta_{\rm b}=1.0^{\prime\prime}$. The right panel of Fig. 4 shows the pixel distribution, which can be fitted well with a Gaussian profile. This indicates that the data is mostly noise without any other significant features. The noise level measured from the Gaussian fit is 40 $\mu$Jy beam^-1. We refer to the rms noise derived in this way as $\sigma_{\rm b}$ in the following sections. This is an objective characterization of the noise in the data, independent on source properties.

We emphasize that each empty galaxy map can be used for a large number of independent simulations. This can be achieved by placing simulated signal at random positions within the primary beam. Considering that half-synthesized-beam offsets produce an independent noise realization, the number of such independent realizations are of order 1000 for each empty ALMA band-7 dataset.

The model sources were created using GILDAS through the MAPPING task uv$\_$fit. We generate elliptical Spergel model sources by fixing their seven free parameters: centroid position, flux density, effective radius $R_{\rm e}$, minor-to-major axis ratio, position angle (PA), and Spergel index $\nu$³³3In practice this is done by fitting to the empty data a model with all parameters fixed and negative flux. The residual image will then have the desired model added (and positive).. The flux density varies over a range with step sizes of factors of two. This range corresponds to a total flux density that is normalized by the noise, varying from 25 to 400, which is represented by the ratio of integrated flux density to the pixel rms noise (S_tot/$\sigma_{\rm b}$). We define S/N as the ratio of the integrated flux density to the noise per beam S_tot/$\sigma_{\rm b}$: the advantage of this choice lays in its model-independence and reliance on very basic properties of the source and of the noise. Note that for extended sources, the S/N defined in this way is obviously higher than the S/N eventually recovered for the integrated flux density coming from a full Spergel/S${\rm\acute{e}}$rsic profile.

The simulation process involves setting the size of the sources in units of $\theta_{\rm b}$, i.e., $R_{\rm e}/\theta_{\rm b}$ = 0.1, 0.2, 0.4, and 0.7, to represent very compact, small, intermediate, and large-sized (relative to the beam) sources, respectively. For each set of Monte Carlo (MC) sources that share a fixed effective radius and axis ratio $q$ ($\equiv b/a$), the Spergel model sources vary in both their flux density (thus, S/N) and Spergel index with values of $\nu=$ 0.5, $-0.3$, $-0.5$, $-0.6$, and $-0.7$. The position angle remain constant. To mimic real observations, we add the model sources to the realistic noise data that is derived from observed visibilities, to produce a simulated data set (see Fig. 4 for an example).

Fitting elliptical Spergel models to the simulated dataset is performed using the task uv$\_$fit. All the seven fitted parameters are allowed to vary. As typically done within MAPPING’s uv$\_$fit in GILDAS, we generate a range of initial guesses for each fitted parameter, which are built into a N-dimensional list of combinations of these guesses. This approach helps to explore a wider range of possible solutions and thus identify the best fit model and return a well-sampled range for uncertainties. We test the fit by setting the starting range parameters of initial guesses within a factor of 2 centered on the input mock values and the number of start parameters (e.g., 3 guesses for each parameter). We found that for simulated sources (where model parameters are known a priori), the results using single initial guesses identical to the real parameters are not significantly different from when the code is run using multiple initial guesses for each parameter. However, using uv$\_$fit with a large range of initial guesses is critical for real observations, where the source’s properties are not known beforehand.

Each $uv$-plane simulation is imaged and then fitted with a single-component S${\rm\acute{e}}$rsic profile using galfit (Peng et al., 2002). The galfit run is performed on the dirty maps, which are created by Fourier transforming visibilities without cleaning. The (full) dirty beam is used as the galfit PSF. The known parameters of the mock sources are used as initial guesses for galfit. All parameters are left free without constraints in the fitting. The initial guess of the S${\rm\acute{e}}$rsic $n$ is calculated by converting the input $\nu$ using Eq. (3). The S${\rm\acute{e}}$rsic fits also provide measurements of seven free parameters: central position, total magnitude, effective radius, S${\rm\acute{e}}$rsic index, axis ratio, and PA. In some cases, galfit may fail to provide accurate measurements. To ensure the validity of our comparisons, we remove measurements in both the image- and $uv-$plane for sources with output parameters marked as problematic in galfit. All this procedure is extremely favorably biased towards positively amplifying the performances of galfit.

Fig. 5 illustrates a typical case of a simulated source generated with a Spergel profile and with $S_{\rm tot}/\sigma_{\rm b}$ ratio of 50. Additionally, it shows the best-fit models derived from uv$\_$fit and galfit, respectively. Both methods provide good constraints to the simulated source at this S/N ratio, as no significant component is visible in the residual map.

The uv_fit command uses the SLATEC/DNLS1E implementation of the Levenberg-Marquardt algorithm (see Press et al., 1992, for an intuitive presentation) to minimize the reduced $\chi^{2}$ of this non-linear least-square problem. This algorithm only requires to deliver a routine that computes the complex function and its partial derivatives with respect to the different fitted parameters. Appendix C delivers the equations for the elliptical Spergel profile and its partial derivatives. Once the minimum of the least-square problem is found, the routine SLATEC/DCOV is called to compute the covariance matrix at this minimum. The diagonal elements of this covariance matrix are the $\pm 1\sigma$ uncertainties on each fitted parameters.

In estimation theory, the Fisher matrix, $\boldsymbol{I}_{F}$, quantifies the amount of information in the least-square problem. When the noise on the measurements (the visibilities) is well modeled by an uncorrelated centered white Gaussian random variable of standard deviation $\sigma$, the computation of the Fisher matrix reduces to (Stoica & Moses, 2005)

where $[\boldsymbol{I}_{F}]_{ij}$ stands for the term $(i,j)$ of the Fisher matrix, $\sigma_{k}$ and $V_{k}$ are the noise and the fitted visibility function for visibility $k$, and $(\varphi_{i})$ is the vector of fitted parameters. In our case, $(\varphi_{i})=(x_{0},y_{0},L_{0},R_{\mathrm{maj}},R_{\mathrm{min}},\mathrm{\phi},\nu)$, i.e., the central position of the Spergel profile as an offset with respect to the phase center, its luminosity, its major and minor half-light radius, its position angle, and its index. The Cramer-Rao Bound (CRB) for each fitted parameter, ${\cal B}(\varphi_{i})$, is defined as the $i$th diagonal element of the inverse of the Fisher matrix

The CRB is the reference precision of the least-square problem. Indeed, the variance of any unbiased estimator of the parameter $i$ will always be larger than the associated CRB (Garthwaite et al., 1995), or

In other words, an efficient fitting algorithm will deliver variances for the estimated parameters, which reach the associated CRB values. In practice, sufficiently large signal-to-noise ratios are required to ensure that the $\chi^{2}$ minimization converges towards the actual solution. Additional explanations and an example of application to the fit of CO(1-0) profiles in the local inter-stellar medium can be found in Roueff et al. (2021).

In Fig. 6, we present the results of our simulations, using only two models of galaxies as typical examples, for clarity. The input mock source was chosen with $R_{\rm e}/\theta_{\rm b}$ of 0.4 (0.2), axis ratio $q$ of 0.75 (0.6), and Spergel index $\nu$ of $-0.6$ (0.5). The distribution of recovered fitted parameter is shown in the top-two (bottom-two) panels. The $\nu=0.5$ case is particularly useful as it provides perfect coincidence with the Spergel and Sersic ($n=1$) models.

We extracted the distribution of the structural parameters obtained by fitting general Spergel profiles (panels labelled with uvfit) and S${\rm\acute{e}}$rsic profiles (panels labelled with galfit) with all parameters free, respectively (Fig. 6). The distribution of recovered best fitting values is shown for flux density, size, axis ratio, and S${\rm\acute{e}}$rsic index, where the true values are shown by vertical lines. The dispersion of recovered values can be used as a gauge of the measurement uncertainty. The average difference between recovered and true values probes any measurement biases and accuracy. In all cases, the scatter of the distributions decreases as the S/N of the simulated sources increases, as expected. Similarly, any systematic bias decreases with S/N, both in the image plane and in the $uv$-plane.

These two examples demonstrate that, when comparing fits in the $uv$-plane to measurements of each structural parameter obtained from profile fitting in the image-plane with galfit, it is clear that the latter exhibit larger scatter and have larger uncertainties. Also, there is evidence for systematic biases that are more pronounced for image-plane fitting with galfit, especially at low S/N and clearly visible for both the low and high S${\rm\acute{e}}$rsic cases (albeit larger for the latter).

The systematic relative deviations of all key parameters of the fit (biases) for a larger variety of input parameters are shown in Fig. 7, which again shows how the S${\rm\acute{e}}$rsic fits gets increasingly biased at low S/N much more rapidly than uv-plane fits. As the biases vanishes at the highest S/N even for the S${\rm\acute{e}}$rsic case, we conclude that these are not due to the previously discussed systematic differences in the profiles.

Fig. 8 compares the measurement uncertainty between $uv$-plane and image-plane. To enhance the readability of the figure, we plotted only the measurements obtained from fitting model sources with a size of $R_{\rm e}/\theta_{\rm b}$ = 0.2 and 0.7, and flux density of $S_{\rm tot}/\sigma_{\rm b}$ = 25, 100, and 400 in Figs. 7 and 8. The following sections provide detailed results on the recovery bias of individual key parameters. Then, the focus shifts to the comparative uncertainty in parameters estimates. Again, we emphasize how the case $\nu=0.5$/$n=1$ is included, where the S${\rm\acute{e}}$rsic and Spergel models are identical, and it behaves fully similar to the other $\nu$/$n$ cases, demonstrating that the results are not driven by small differences between models at higher $n$.

Apart from systematic biases, it is also important to evaluate the scatter of the measurements, both in the uv-plane and in the image plane, to verify if any of the two returns better measurements with lower scatter.

In order to evaluate such relative measurement uncertainty, we use the scatter of the distribution of recovered values, evaluated as the median absolute deviation (MAD) of the data around the true value (i.e., input mock value) converted to $\sigma$ using $\sigma=$1.48$\times$MAD, which is what expected for a Normal distribution. The use of MAD is preferred in order to to be less dependent on outliers while capturing the bulk of the spread in the sample. By defining the deviation with respect to the true value, the measurement bias is also taken into account when evaluating the overall uncertainty in the measurements. Fig. 8 compares the scatter of measurements obtained from $uv$-fitting with those obtained in the image-plane for each key structural parameter.

Again, we find that the random uncertainties in all structure parameter recoveries are systematically larger for image-based estimates than those obtained by visibility analysis. The median scatter ratio of measurements obtained from $uv$-fitting to that obtained from image-plane for the recovered flux densities is 0.5, while for the recovered effective radius, axis ratio, and S${\rm\acute{e}}$rsic $n$, the median scatter ratios are 0.5, 0.6, and 0.5, respectively. We did not find a significant correlation between the scatter ratio and the S/N of data. These findings are in line with a study that compared the performance of stacking data in the $uv$-plane and image-plane, which found that $uv$-stacking resulted in significantly improved accuracy of size estimates, with typical errors less than half compared to image-stacking (Lindroos et al., 2015). We emphasize that the small residual systematics shown in Figs. 2 and 3, inherent in the conversion of Spergel-based parameters of our simulations to S${\rm\acute{e}}$rsic parameters, only account for a negligible amount of the excess noise coming from galfit (and have no effect at all for the $n=1$ case).

In conclusion, we find that measurements in the image-plane are not only more subject to bias, but also returning parameters that are more substantially affected by noise, compared to those in the $uv$-plane.

We check the achievable accuracy for the parameter estimates from Spergel fits in the uv-plane by computing the ratio of the uncertainties of the parameter estimates to the input mock value, given by $\sigma$(para)/para. The uncertainties of the parameter estimates, $\sigma$(para), can be evaluated as 1.48$\times$MAD, as discussed in Section 4.2.

In Fig. 9, the accuracy of parameter estimates is shown to vary with the S/N of the data. The parameters being estimated are flux density, size, axis ratio, and S${\rm\acute{e}}$rsic index. To enhance readability of the figure, only Spergel models with an exponential ($\nu=0.5$) profile and steep profile with $\nu=-0.7$ (close to a de Vaucouleurs profile) were considered. These models are shown in the panels between the second and fifth rows. It is worth noting that for other cases where the Spergel index is assumed to be between $\nu=0.5$ and $-0.7$, the measurements are found to be either within or close to the values obtained from the above two cases. Therefore, the results presented in Fig. 9 can be considered representative of the Spergel model.

Fig. 9 shows that the uncertainties in measuring galaxy shape parameters, such as size and Spergel index, are significantly larger than those for flux density. The most difficult parameter to estimate is the Spergel index. The uncertainty in estimating the Spergel index can be as high as $70\%$ for model sources with a S_tot/$\sigma_{\rm b}$ of 25 or lower, suggesting that the measurement of Spergel index can be highly uncertain. As the S_tot/$\sigma_{\rm b}$ increases to 50, the accuracy of Spergel index estimates significantly improves, with a median value of $\sigma$($n$)/$n$ of 0.36, which we consider as the bare minimum to define a meaningful estimate. At S_tot/$\sigma_{\rm b}$ of 50, the uncertainties of the flux density, size, and axis ratio are significantly smaller, with median values of 9%, 18%, and 22% of estimates, respectively. To obtain a meaningful and reliable profile fitting result with the Spergel model in $uv$-plane, our simulations strongly suggest that a S_tot/$\sigma_{\rm b}$ of at least 50 should be required.

In addition, we find that the accuracy of both the flux density and size is lower for simulated data with a steep profile compared to the model data with a flat profile. On average, the accuracy of flux densities and sizes for the galaxy with an exponential profile ($\nu=0.5$) is 50$\%$ higher than those of the galaxy with a de Vaucouleurs-like profile ($\nu=-0.7$). However, it is important to note that the measurements of Spergel $\nu$ tend to be less accurate for smaller sources ($R_{\rm e}/\theta_{\rm b}\leqslant$0.2) compared with larger sources in the data. We have also found that both the flux densities and sizes of the small-sized sources are more accurately measured, with a typical factor of about 40%. The differences in the accuracy of parameter estimates imply that all the fitted parameters are interrelated.

In Fig. 10, we compare the average value of parameter errors in simulated galaxies to the posterior scatter of the recovered distributions to explore their reliability. We find that in most cases, the uncertainties of parameters measured from Spergel profile fitting are underestimated, although typically within a factor of 2. The underestimation is generally less important for steep profiles ($\nu=-0.7$) than for disks ($\nu=0.5$) and for extended versus compact galaxies. At each flux S/N bin, the ratio between the error obtained from simulations and the median of the error measured by uv$\_$fit for the whole set of simulated sources is in the range of 1.4–1.9, 2.0–3.2, 1.9–2.5, and 1.4–1.6 (see the black bars in Fig. 10), with a mean value of 1.6, 2.5, 2.1, and 1.5 for the parameter estimates of flux density, effective radius, axis ratio, and Spergel index, respectively.

In this section, we examine the covariance between fitted parameters estimated from uv$\_$fit using a Spergel model. The top panels in Fig. 11 show the correlations between the fitted parameters, which are flux density $S_{\rm tot}/\sigma_{\rm b}$, effective radius $R_{\rm e}/\theta_{\rm b}$, axis ratio $q$, and Spergel index $\nu$, for a simulated dataset with an input $S_{\rm tot}/\sigma_{\rm b}$ of 100, $R_{\rm e}/\theta_{\rm b}$ of 0.2, $q$ of 0.6, and $\nu$ of $-0.5$, respectively. In this case, the Spearman correlation coefficient, which is calculated by dividing the covariance by the intrinsic scatter of each parameter, show weak correlations ($|\rho|\lesssim 0.3$) between fitted parameters of size, axis ratio, and $\nu$, while moderate correlations ($|\rho|\sim 0.5$) are found between the flux density and both $\nu$ and size.

The bottom part of Fig. 11 summarizes the pairwise correlation coefficients calculated for all the data sets in our simulation. Generally, we find that the correlation between variables becomes more prominent as the S/N increases, except for the correlation between flux density and size. For sources that are significantly more compact than the beam ($R_{\rm e}/\theta_{\rm b}<0.2$), the correlation between the flux density and measured source size is weaker when the source is detected with a high S/N. We did not find any significant correlation between q and other parameters. This indicates that the measured axis ratios are almost independent of these parameters.

There is a positive correlation between flux density and S${\rm\acute{e}}$rsic $n$ (anti-correlated with Spergel $\nu$), indicating that sources with higher measured S${\rm\acute{e}}$rsic $n$ tend to have a larger flux density. In addition, a strong positive correlation is seen between flux density and size, except for very compact sources with $R_{\rm e}/\theta_{\rm b}$ of 0.1, where the correlation is relatively weak. This means that sources for which sizes are overestmated have a tendency to be also boosted in flux density.

A comparison between the Spergel and S${\rm\acute{e}}$rsic profiles shows that the Spergel model has a steeper core and faster declining wings compared to the S${\rm\acute{e}}$rsic model (see Fig. 1 and Fig. 16), except for the case of Spergel $\nu=$0.5 which is equivalent to an exponential profile (S${\rm\acute{e}}$rsic $n$=1). The question of whether Spergel or S${\rm\acute{e}}$rsic models provide better fits to actual galaxies remains open and requires future investigation. The current available data quality may not be sufficient to determine a clear preference between the two functional forms, in absolute terms. For the time being and with the typical data available, we deem the two functional forms as equivalent.

Converting a Spergel index to a S${\rm\acute{e}}$rsic one (as well as $R_{\rm e}$ and total flux) requires knowledge of the intrinsic size of the galaxy and the angular resolution of the observations. In the case of this study, the angular resolution is determined by the synthesized beam of the interferometric data. We believe this requirement also applies implicitly to optical observations. When a S${\rm\acute{e}}$rsic index is derived from optical data, it applies to the range of scales actually observed in the data. It may not apply beyond the observed scales by definition. It is possible that a somewhat different Spergel/S${\rm\acute{e}}$rsic index might be recovered when re-observing the same galaxy with a much different surface brightness sensitivity.

It is also relevant to question whether the differences in the profiles at the outer and inner ranges could affect the systematic biases of the parameter estimates, particular for a steep profile. Along these lines, considering that we have simulated Spergel models and then fitted them with galfit in the image-plane, one might wonder if the under-performance of galfit in the image-plane could be simply related to the discussed differences between Spergel and S${\rm\acute{e}}$rsic models.

As already mentioned through-out the description of results in the previous section, we believe that the argument can be dismissed, but it’s worth recalling the key evidences here. Based on the bottom panels of Fig.7, the simulations with the highest S/N (400 in this case) show that any average bias in S${\rm\acute{e}}$rsic modeling is vanishing or strongly reduced, as expected by construction (Eq.(3)), showing that any residual systematics beyond our conversions do not have a measurable impact. The systematic differences in the recovered profile parameters ($R_{\rm e}$ and S${\rm\acute{e}}$rsic $n$, crucially) are in fact strongly S/N-dependent, and therefore mostly an effect of noise, rather than arising from structural differences. Additionally, Fig.8 shows that the excess noise in parameter recovery from image-plane fitting is not a function of S/N, while the systematic deviations between the profiles are expected to become more evident at the highest S/N. Also importantly, Fig.8 shows that an entirely consistent situation is seen for the $n=1$ case, which is fully identical in shape to a Spergel $\nu$ of 0.5.

Fig.12 extends further this point by showing the performance comparison of uv$\_$fit versus galfit in the case of Gaussian sources and PSF models. In this case, we directly simulated these shapes in the $uv$-plane, which are not Spergel models, and there is no shape difference with respect to the model fitted by galfit. However, the overall behaviour is qualitatively the same as in Figs.7 and 8. There are stronger biases in the image plane (top panels), and the effective noise is also much higher there (lower panels).

It has been already shown that ignoring noise correlation in interferometric images can lead to a significant underestimation of the statistical uncertainty of the results (e.g., Tsukui et al., 2023). To test this idea further, we carried out aperture photometry as an alternative to galfit for point-source simulations. Aperture photometry is one of techniques used to measure flux density in astronomical observations, which is also used in interferometric images (e.g., Lang et al., 2019; Gómez-Guijarro et al., 2022). It should be noted that, for extended sources, there is an added complication of correcting for flux losses, and this technique has the limitation of not allowing to estimate morphological parameters. However, it is well-defined for flux density measurements of point sources.

Figure 13 shows the systematic bias and the measurement uncertainty on the recovery of the flux density for point sources using different fitting methodologies, i.e., uv$\_$fit with a Point source model, galfit with a PSF profile that is identical to the dirty beam of the simulated data, and aperture photometry. The latter is performed with an aperture radius equal to the circularized radius of beam size and at the position returned by galfit PSF fitting. It is clear that the flux density estimates obtained through aperture photometry are fully consistent to those measured in the $uv$-plane, with a similar level of accuracy in the estimates. Both of these methods exhibit smaller systematic errors and scatters when compared to fitting with a PSF model using galfit. In a way, aperture photometry measurements being fairly basic and raw just are not fooled by false coherence in the signal induced by correlation, and return the full S/N performances in the limiting case of point-like emission (Fig.13-right).

It is common practice to perform deconvolution of low-frequency radio interferometric data to remove strong sidelobes from bright sources in very large fields, a process known as cleaning. We emphasize that in ALMA and NOEMA observations the field of view is tiny compared to the radio and the source density is basically always manageable, so that cleaning can generally be safely avoided. Cleaning, as any deconvolution approach, is a model-dependent process which in most implementations arbitrarily interpolates the observed source with a superposition of point sources, whose side-lobes are then subtracted from the data using the dirty beam. Due to its arbitrariness and its approximation of real sources as multiple point-sources, it is clearly not ideal for analysis of galaxies morphological profiles. We finally emphasize that, obviously, the cleaning process would not remove correlations in the signal in the image plane, which are due to Fourier transforming the visibilities.

To investigate the morphologies of galaxies observed in interferometric images, Spergel modeling directly in the $uv$-plane should be the preferred approach. However, a high S/N ratio of at least 50, and ideally much more, must be required to perform a Spergel fit that meaningfully constrain the Spergel index, as we have previously commented based on Fig.9. This requirement is similar to optical observations of galaxies and their galfit modeling, where a high S/N of order 100 is needed to attempt derivation of a S${\rm\acute{e}}$rsic index (van der Wel et al., 2014; Magnelli et al., 2023).

These considerations should be extended further to include other parameters of interest, such as the size or even the flux density. These parameters are easier to derive and require less S/N compared to a Spergel/S${\rm\acute{e}}$rsic index. However, there is a tension between extracting a meaningful measure from the data and ensuring an unbiased measure. It is important to determine the best approach in balancing these factors.

The top-panel of Fig.9 provides information on the uncertainties in flux density, size, and axis ratio for Gaussian models, as a function of their S/N ratios. By comparing this information to the bottom panels, one can make a decision between a less complex (Gaussian) fit and a more complex (Spergel) fit to a given dataset. For example, Gaussian fits can provide size uncertainties of about 20% down to an S/N ratio ($S_{\rm tot}/\sigma_{\rm b}$) of 20, whereas Spergel fits do not offer the same level of accuracy at lower S/N ratios. Similar considerations apply to flux density and axis ratio. In addition, there are no significant differences in the accuracy and uncertainty estimates between circular and elliptical Gaussian models (see Fig. 9 and Fig. 12). This demonstrates that the increase in degrees of freedom in Spergel models leads to greater uncertainty in the measured parameters, including flux density and size.

Figure 14 shows the distribution of the ratio of the size obtained from a Spergel fit to that from a Gaussian fit, measured for a sample of about 100 high-$z$ star-forming galaxies observed with the ALMA (Q. Tan et al., in preparation). The median ratio between the $R_{\rm e}$ size obtained from a Spergel fit and the FWHM size obtained from a Gaussian fit is $0.47_{-0.08}^{+0.10}$, where uncertainties correspond to the 16th and 84th percentiles of the distribution. This implies that in most cases, the FWHM size measured from a Gaussian model can be used as a good approximate of the effective radius using the relation $R_{\rm e}$=FWHM/2. However, it is worth noting that almost all of the outliers with $R_{\rm e}$/FWHM ratio far from the median value in Fig. 14 are sources measured with large S${\rm\acute{e}}$rsic index $n$ ($n>2$; Q. Tan et al., in preparation). This suggests that the difference in size measured from Gaussian fit and the Spergel fit could be significant when the source has a large S${\rm\acute{e}}$rsic index.

Finally, to achieve a global optimal solution for the Spergel fit, we recommend utilizing the fitting results obtained from Gaussian or exponential fits as prior knowledge for the initial guesses of each structure parameter. This provides good starting guesses, which helps the solver in achieving convergence and in providing a reasonable level of uncertainty. Such preliminary information will also inform the user about the actual merit of proceeding to a more demanding full-Spergel fit.

Fitting Spergel models to interferometric data can be complex and requires deep high S/N data. However, it can also potentially brings powerful insights and enable investigations into new scientific questions.

We aim to demonstrate the potential of our approach by re-examining two cases of published ALMA observations of distant sources that were claimed to contain giant halos surrounding the central galaxies (PACS-787 from Silverman et al. (2018), including several coauthors of this paper, and HZ7 from Lambert et al. (2023)). These two specific examples are taken from a growing body of results that report the existence of large halos, often observed in [CII]$\lambda 158\mu$m and other tracers (e.g., Ginolfi et al., 2017; Fujimoto et al., 2019, 2020; Pizzati et al., 2020; Cicone et al., 2021; Herrera-Camus et al., 2021; Jones et al., 2023; Li et al., 2023; Posses et al., 2023; Scholtz et al., 2023, etc). Extended halos around galaxies are generally interpreted as evidence for accretion, outflows, tidal stripping, or other phenomena affecting distant galaxies. This presents a relevant opportunity for further constraining these processes. However, it is worth considering whether these halos are genuine different structures from the galaxies or simply the outer scale extension of high S${\rm\acute{e}}$rsic-index profiles. In many cases, simple Gaussian fitting was attempted for the central galaxies, and high S${\rm\acute{e}}$rsic index profiles are well known to display large halos (e.g., Mancini et al., 2010).

We have downloaded the data for both datasets as described below. In ALMA Cycle 6, HZ7 was observed with 80 minutes of on-source integration time in the C43-4 configuration with 47/45 12m antennas in band 7 (Project 2018.1.01359.S; PI: M. Aravena). The [CII] emission at 303.93GHz falls into one of four SPWs with a native channel width of 15.625 MHz. PACS-787 was observed with high (C40-6, 43 12m antennas with a maximum baseline of 1.1 km) and low (C40-1, 42 12m antennas with a maximum baseline of 278.9 m) resolution configuration in Band 6 with 19.7 and 10.2 minutes of on-source integration time in ALMA Cycle 4 (Project 2016.1.01426.S; PI: J. Silverman). The CO(5-4) emission at 228.17GHz falls into a native channel width of 3.91 MHz in high-resolution observation and 15.62 MHz in low-resolution observation. The calibration targets for the three observations above are J1058+0133 for bandpass, pointing, and flux, and J0948+0022 for phase.

The left panels of Fig. 15 show the dirty images of PACS-787 and HZ7, where emission on large scales of several arcsec (tens of kpc) can be readily seen. We modeled the emission in the $uv$-space with simple Spergel profiles. For PACS-787, which contains two galaxies in the process of merging, we used two Spergel components (one for each galaxy), while for HZ7, we used a single Spergel. The Spergel fitting in both cases is able to fully account for the emission from the galaxies, simultaneously reproducing the inner emission and the outer halos. The residuals are clean, and no further components are needed, as shown in the right panels of Fig.15. For the case of PACS-787, the two Spergel component have $\nu=-0.36,-0.14$ and $R_{\rm e}/\theta_{\rm b}\sim 0.7$ for both. For HZ7, we find $\nu=-0.42$ and $R_{\rm e}/\theta_{\rm b}\sim 1.3$. Based on Eq.(3) and Fig.2, this corresponds to S${\rm\acute{e}}$rsic $n\sim 1.5$–3, which is well above the Gaussian approximation ($n=0.5$) and also above the exponential case ($n=1$), although not as steep as a de Vacouleurs profile ($n=4$).

We have shown that the full emission in these two systems can be fit with a single component model, which raises questions about the interpretation of the outer emission as a halo. Although high S${\rm\acute{e}}$rsic values indicate the presence of both a central component and a profile extending further out than a disk model with both smoothly connected, there is no apparent solution of continuity between the inner parts and the outer halos. This is similar to the inner versus outer parts of elliptical galaxies, and therefore suggestions of different physical origins for them are less substantiated.

We emphasize that we have only re-evaluated two cases of halos from the literature (out of many more existing) to exemplify the power of fitting more complex Spergel models to ALMA data. It is beyond the scope of this work to re-analyse all similar observations from the ALMA archive. However, we obviously anticipate that in other cases, the halos might disappear once fitted with a Spergel model. This is not only because of the analysis presented here but also based on one of the key results from our forthcoming companion paper (Q. Tan et al. in preparation), which suggests that $n>1$ models are required for most ALMA observations of distant galaxies.