A catalog of polychromatic bulge-disk decompositions of ~ 17.600 galaxies in CANDELS
Paola Dimauro Marc Huertas-Company Emanuele Daddi Pablo G. Pérez-González Mariangela Bernardi Guillermo Barro Fernando Buitrago Fernando Caro Andrea Cattaneo Helena Dominguez-Sánchez Sandra M. Faber Boris Häußler Dale D. Kocevski Anton M. Koekemoer David C. Koo Christoph T. Lee Simona Mei Berta Margalef-Bentabol Joel Primack Aldo Rodriguez-Puebla Mara Salvato Francesco Shankar Diego Tuccillo


firstpage–lastpage
On the importance of the nonequilibrium ionization of  and  and 
the line-of-sight in solar surges
    [
    December 30, 2023
==================================================================================================











Understanding how bulges grow in galaxies is critical step towards unveiling the link between galaxy morphology and star-formation. To do so, it is necessary to decompose large sample of galaxies at different epochs into their main components (bulges and disks). This is particularly challenging, especially at high redshifts, where galaxies are poorly resolved. 
This work presents a catalog of bulge-disk decompositions of the surface brightness profiles of ∼17.600 H-band selected galaxies in the CANDELS fields (F160W<23, 0<z<2) in 4 to 7 filters covering a spectral range of 430-1600nm. This is the largest available catalog of this kind up to z = 2. By using a novel approach based on deep-learning to select the best model to fit, we manage to control systematics arising from wrong model selection and obtain less contaminated samples than previous works. We show that the derived structural properties are within ∼10-20% of random uncertainties.  We then fit stellar population models to the decomposed SEDs (Spectral Energy Distribution) of bulges and disks and derive stellar masses (and stellar mass bulge-to-total ratios) as well as rest-frame colors (U,V,J) for bulges and disks separately. All data products are publicly released with this paper and through the web page  <https://lerma.obspm.fr/huertas/form_CANDELS> and will be used for scientific analysis in forthcoming works.








galaxies: fundamental parameters, galaxies: high-redshift, galaxies: bulges








§ INTRODUCTION

Most galaxies have two main components: disks and bulges. These two components are believed to have very different formation mechanisms. Disks are generally rotationally supported and confined to a thin plane. They are believed to be the consequence of gas infall into halos, which transfers angular momentum to the baryons. Bulges generally have a 3D shape and larger stellar velocity dispersions. Their formation requires dissipative processes and a loss of angular momentum. Mergers of two disks is the classical channel to grow bulges (e.g.). However numerical models show that disks, especially at high redshift when they are more unstable and gas rich, can also self generate a bulge through instabilities (e.g. ) and/or inflow of cold gas towards the center (e.g. ). Properly understanding how all these different processes come together to assemble galaxies into their main components requires identifying bulges and disks in galaxies and studying their evolution across cosmic time. Since disks and bulges have different projected surface brightness distributions, the decomposition of the light by fitting analytic models  to the 1D or 2D light profiles has been widely used in the literature. Extending this approach to large datasets arising from deep surveys, where objects cannot be checked individually, is particularly challenging. Not only because of the computing time, but also due to the large amount of systematics that need to be controlled. 
Many works have obtained bulge-disk decompositions of several hundreds of thousands of galaxies at low redshift, where galaxies are reasonably well resolved <cit.>. A significant amount of post processing is required to assess the quality of the fits and eventually identify unphysical solutions. One key issue, for instance, is deciding whether two components are really needed to model the light profile or if one unique component is better suited. This is usually addressed by performing a-posteriori statistical tests to measure if the addition of an extra component improves the fit (e.g. ). 

At high redshift, the situation is even more dramatic, both because of lower S/N and because galaxies start to be less well resolved even with space based imaging. This is why most of the works involving surface brightness fitting of large samples of distant galaxies tend to use one single Sérsic component, reducing the amount of free parameters (e.g. ). Two component fitting is generally done on smaller datasets (e.g. ). Even there, degeneracies are reduced by adding more constraints on the parameters. For example, <cit.> forced the index of the bulge to be 4. However, many works have shown that bulges have a wider distribution of the index (e.g ) so this might not be the ideal solution. 

An additional issue of bulge-disk decompositions is that they are performed on the light profiles, while models predict stellar mass distributions. Deriving stellar masses from light distributions requires assuming an mass/light (thereafter M/L) ratio, which can be different for bulges and disks and also from galaxy to galaxy. Assuming that the light profile directly trace the stellar mass (a unique M/L for all galaxies/components) is clearly an oversimplification which can introduce additional systematics. This is especially true for high redshift studies, where very different cosmic epochs are probed. 
 
In this paper we present a catalog of bulge-disk decompositions of ∼ 17.600 galaxies in the CANDELS fields. This is the largest catalog of this kind for objects spanning such a wide redshift range (z < 2).  This work also introduces several novel techniques to address some of the issues discussed above. Firstly, we develop a method based on deep-learning to estimate the optimal model that should be used to fit the light profile (namely one or two components). While other techniques existing in the literature which rely on the fitting residuals, our method acts before the modeling phase, at the pixel level. Additionally, our fits are done simultaneously in 4 to 7 (depending on the fields) high resolution filters using the modified version of Galfit, GalfitM  <cit.>. This multi-wavelength fitting allows us to increase the Signal/Noise (S/N) and to reduce the random uncertainties but also to estimate Spectral Energy Distributions (SEDs) of bulges and disks and a M/L for every component by fitting stellar populations models with the FAST code <cit.>. We  thus provide stellar population properties (stellar masses, SFRs) and rest-frame colors (U,V,J) for bulges and disks. This should enable a less biased comparison with predictions of galaxy formation models. The catalog will be made public upon publication of this work[<https://lerma.obspm.fr/huertas/form_CANDELS>]. 

The paper proceeds as follows. We describe the dataset in section <ref>. The methodology used for profile fitting is discussed in sections <ref> and sections <ref>. The accuracy of the catalog is quantified in section <ref>. The  stellar population properties are described in section <ref>. All magnitudes are measured in the AB system. 



§ DATA


Our starting point for the selection are the official CANDELS <cit.> H-band (F160W) selected catalogs ( for UDS,  for GOODS-S, Barro et al. (2017)  for GOODS-N and  for COSMOS and AEGIS). For this study, we only consider galaxies brighter than F_160W=23. This magnitude selection is applied to ensure reliable two component decompositions as detailed in section <ref>. In addition to the three NIR images (F105, F125, F160), observed as part of the CANDELS survey, we use ancillary data in four additional bands for GOODS-N and GOODS-S (F435W, F606W, F775W, F850L) and two in  the AEGIS, UDS and COSMOS fields (F606W, F814W). All images are resampled to a common pixel scale of 0.06 arcsec/pixel. This is required to perform simultaneous multi-wavelength fits to the surface brightness profiles as described in section <ref>. 

We also use the 2D single fits published in <cit.> in three NIR filters (F105W, F125W, F160W) and the deep-learning based visual morphologies published in <cit.>. The official CANDELS redshifts are used. More details can be found in <cit.>. Spectroscopic redshifts are used when available. If not, we use photometric redshifts derived through SED fitting by combining different available codes.  Although we derive stellar masses of bulges of disks (described in section <ref>) we also use total stellar masses from previous works for comparison <cit.>. Namely, the best available redshift is used to estimate stellar masses based on the PEGASE01 stellar population models (Fioc & Rocca-Volmerange 1999). We assume solar metallicity, exponentially declining star formation histories, Salpeter (1995) IMF and a Calzetti et al. (2000) extinction law. The M_*/L ratios are converted to a Chabrier (2000) IMF by applying a constant 0.22 dex offset. Rest-frame magnitudes (U,V,J) are also derived as part of the SED fitting procedure and are used to divide galaxies between star-forming and quiescent systems with the now standard way (see ).

Our final sample consists of ∼17.600 galaxies out of ∼19.000 with F160W<23 in the CANDELS survey. The difference comes from the fact that there is not a perfect overlap between the fields in the different filters. We restricted our analysis to galaxies observed with at least 4 filters. Figure <ref> shows the distribution of the selected galaxies in the M_*-z plane. We estimate the stellar mass completeness using the method described in <cit.>. The lower limit stellar mass (M^*_lim) is computed as: log(M^*_lim)=log(M_*)+0.4×(m_f160W-23). m_f160W is the apparent magnitude in the F160W filter where the selection is done. The completeness is then estimated as the 90th percentile of the distribution of M^*_lim. We repeat the same procedure for all galaxies as well as for passive and star-forming galaxies (selected according to their UVJ colors). The obtained values are plotted in figure <ref> and summarized in table <ref>. 
We estimate a mass completeness limit of ∼ 10.7 at z∼2. We notice that for the passive population, very few galaxies lie below the completeness limit, despite the bright magnitude cut. This is expected given the evolution of the stellar mass function of passive galaxies (e.g. ), which at z>1.5, is dominated by massive galaxies. As a sanity check, since the CANDELS catalogs are several magnitudes deeper than our selected catalog, the stellar mass completeness is also estimated by computing the stellar mass above which at least 80% of galaxies from the deep catalog at a given mass and redshift are selected (see also ). We obtain similar results. We also provide in table <ref> an estimate of the completeness for embedded components (bulges and disks). 



The key measurement added by this work is a one/two component (bulge/disk) decomposition of the surface brightness profiles in different filters. In the following we describe the details of the methodology as well as the procedure used to estimate the uncertainties. 

 



§ BEST-MODEL SELECTION WITH DEEP LEARNING



A fundamental problem with surface brightness profile fitting is to decide how many components are needed to model the galaxy light. By using two component models, we force the fitting algorithm to find a solution with 2 components within a given set of constraints even if the galaxy might be better fitted with a single profile or with another combination of profiles. This can lead to unphysical solutions, introducing a systematic error in our subsequent analysis of bulges and disks properties. The reason is that some light might be associated with a bulge and/or a disk even if there is not such a component in the galaxy. This systematic uncertainty can potentially dominate over random uncertainties when performing a scientific analysis (e.g. ).

Several works have used a statistical approach to tackle this problem. By looking at the residuals of the resulting fits it is possible to establish a probability that adding a profile actually improves the fit (e.g. ). This is sometimes combined with a visual inspection <cit.>. This approach still has the problem that a better fit does not necessarily mean a physically meaningful result and that the light is actually properly associated to bulges and disks (). 

Here we introduce a novel alternative technique based on unsupervised feature learning (deep learning). The main novelty is that the best model to fit a galaxy is set a-priori, instead of by looking at the residuals maps a-posteriori. The objective is then to measure, given a galaxy image, which analytic model among a finite set of possibilities is preferred to describe the surface brightness distribution. Recall that this is different from a morphological classification. We are not aiming at obtaining the true morphology but to assess if a given analytical model is appropriate to describe the galaxy.

We proceed in two main steps described in the following. 



 §.§ Training on simulated analytic galaxies


We first simulate a set of 100.000 synthetic galaxies reasonably spanning all the range of structural parameters expected (see table <ref>) using the GalSim code[<http://galsim-developers.github.io/GalSim/index.html>]. Images are convolved with a real PSF and realistic noise from CANDELS images is added as explained in section <ref>. For this particular application, we only simulate one filter (F160W) that will be used to define the model to be fitted. The H-band filter is chosen as a reference since it is the detection band and also the deepest. 



   


We then define 4 types of profiles among the simulated galaxies:



  * Pure : B/T>0.8 and n_bulge>2.5. These are galaxies for which the surface brightness profile should be  well described with a Single model.

  * Pure Exponential: B/T<0.2 or B/T>0.8 and 0.5<n_bulge<1.5. Objects that are disk dominated thus for the surface brightness profile is well captured with a single exponential profile or a one component profile with a index <2. 



  * Bulge + Exponential: 0.2<B/T<0.8 and n_bulge>2.5. Systems that clearly require two components, one with an exponential profile and another with a large index.

  * Pseudo-bulge + Exponential: 0.2<B/T<0.8 and n_bulge<2. Systems that still require two components, but both with low values of the index. 


Independent Convolutional Neural Networks (CNN) are trained in a binary classification mode to isolate the given type of profile from the others in the simulated galaxies. An introduction to CNNs is out of the scope of this work. For more information, we refer the reader to  and  where more details are given. In this work we train 4 different machines with the same architecture. The input of the network is  a simulated 2D image (with noise and PSF) centered on the galaxy (64×64 pixels) and the output is a probability that the image is described by the model it was trained to identify. The model has 4 convolutional layers of increasing depth (from 16 to 64) and 2 fully connected layers. A 3×3 max pooling is performed after each convolutional layer to reduce the number of parameters and a 10% dropout is applied during training to avoid over-fitting. Additionally, a 1% gaussian noise is added in the first layer to avoid that the network learns features on the noise pattern. The model configuration was established after testing different architectures. Slight modifications do not change the main results. The model is trained until convergence and evaluated on the validation dataset. 

At the end of the training process, each simulated galaxy has 4 associated probabilities. Recall that the probabilities do not add to 1 since they were estimated with four independent CNNs. Since for the simulated galaxies, we know the model that was generated, we can quantify the ability of the CNN to distinguish between different profiles on an independent test dataset which was not used during the training phase. Following a standard procedure, we use the area under the ROC curve. The ROC curve is defined as the relation between the True Positive Rate (TPR) and False Positive Rate (FPR) for different probability thresholds. Ideally one would like to have large TPR and small FPR values. So the largest is the area under the ROC curve, the better. The classification is always a trade off between both. In order to define the optimal threshold needed for the classification of our sample, we computed the following two parameters: Purity (P) and Completeness (C):

    C=TPR=TP/TP+FN


    P=1-FPR=TN/TN+FP

TP and FP stand for true and false positives respectively, while TN and FN are true and false negatives. Specificity is a measurement of how contaminated a selection of a given class is by galaxies not belonging to that class. Sensitivity is a measurement of how well the machine recovers all galaxies belonging to a given class. In figure <ref> we show how these two quantities change depending on the applied probability threshold. As expected, higher probabilities indicate more pure but less complete samples.
The plots confirm that the CNN models are able to distinguish between the 4 different types of profiles. We notice that a probability threshold of p=0.4 results in a reasonable trade-off between purity and completeness, around 80%-90%. In order to test the robustness of this choice, figure <ref> shows how C and P depend on galaxy properties (total magnitude and half light size) for each class. Both quantities are stable, confirming that the chosen threshold can be used over the full parameter space covered by the catalog.











 §.§ Knowledge transfer to real galaxies


The above results are based on simulations. The critical step is to use the 4 machines to classify our real galaxies. Our aim is to provide, for each object, a probability that a given model (i.e. pure bulge, pure disk, bulge+exponential, pseudo-bulge+exponential) is preferred to describe its surface brightness profile. We insist that this is not a morphological classification and it is fundamentally different than the visual morphological catalog presented in <cit.>. It measures, given a set of available analytic profiles, which one is better suited to fit the surface brightness distribution of a given galaxy. This information is intended to be used as input for fitting the light profiles as explained in the following sections.

In transferring the trained machine to real data, there is a risk that the CNN does not provide reliable results since simulations and real galaxies are obviously different. For example, even though our simulations include realistic instrumental effects and noise, we did not include close companions in the training set nor irregular galaxies or bars. The validity of the knowledge transfer between simulations and real data is difficult to evaluate since, as opposed to the simulations, the true profile is not known in the real sample 

We perform some posterior sanity checks In order to verify that the trained model properly represents reality:




  * First, we visually inspect a significant number of color images of different classes provided by the CNN. Although, this is not a quantitative measurement, it allows to identify obvious errors. Examples of the galaxies classified according to the different types of profiles are shown in figure <ref>. Generally speaking, when the algorithm indicates that 2 components are required, a bulge and a disk are observed in the color image.

 
  * In terms of types of models, we find that for the majority of the galaxies (∼ 55%) a bulge+exponential model is preferred.  ∼20% of the remaining galaxies are well fitted with an exponential model, ∼15% are preferentially fitted with 2 low index components and ∼10% are well described with one profile with n>2. The fraction of classical double component galaxies strongly depend on the stellar mass and  on the morphology. Since we are studying massive galaxies (M _*>10.3 M_⊙), we expect that most of them are characterized by two components. Our results is in rather good agreement with the predictions (). The majority of galaxies are described by a classical two component model with a central bulge with large index and an exponential disk.  However, it is worth noticing that if no model pre-selection is done and a 2 component fit is blindly applied to the entire sample, the systematic error reaches ∼30% and thus most probably dominates the error budget. 
 
  
  * We measure that <10% of the galaxies in the real sample have all 4 probabilities below P=0.4, i.e. their surface brightness profile are not properly described by any of the 4 considered models. All other galaxies have at least one of the probabilities above this threshold. This is a good indication that the network did not find drastic differences between the training/test simulated samples and the real data. Otherwise the probability distributions would have been very low for all real galaxies.
 
 
  * For galaxies for which no model is preferred (P<0.4), the visual classification taken from <cit.> indicates that the vast majority of them (>80%) are classified as irregulars with high probability (see also figure <ref>). It is expected that their surface brightness profile is not properly described by any of the considered models. We notice additionally that for the irregulars for which an optimal profile was actually found a pure disk like profile is the best solution.
 

We include in the catalog these probabilities, which indicate how accurately each of the four light profiles fits the light distribution of each galaxy.

This allows to select a model to fit a-priori and that way reduce systematic uncertainties as described in section <ref>. 





§ MULTI-Λ FITS WITH GALFITM


The main tools we used to perform the fits are GalfitM and Galapagos-2 from the MEGAMORPH project <cit.>. They are based on Galapagos and Galfit <cit.>. The main difference is that they can simultaneously fit all images at different wavelengths (as opposed to an independent fit for each band). As shown in the aforementioned works, the advantage of such an approach is that, by combining data from all filters, we effectively increase the S/N and naturally use the color information. Therefore the fit is better constrained down to fainter magnitudes than when considering all bands independently. To accomplish this, the wavelength dependence of the structural parameters of galaxies is parametrized with a family of Chebyshev polynomials. The order of the polynomial for each quantity is a user-configurable parameter which sets the degree of freedom. The fitting algorithm then minimizes the coefficients of the function for each structural parameter. If the degree of freedom is equal to the number of filters, then the parameter is effectively independent in each band, as it is obviously the case for the fluxes. For the other parameters, the choice (of the degrees of freedom allowed) is a trade-off between allowing total independence or setting no variation with wavelength (thus reducing the number of free parameters). More details can be found in <cit.> and <cit.>.

There is no obvious way of selecting the optimal configuration. The wavelength dependence of the structural parameters will certainly vary from galaxy to galaxy. Our approach has been to empirically test different configurations and use them to estimate random uncertainties as discussed in section <ref>. For each galaxy we fit 2 types of models: a 1 component model and a 2 component +Exponential model. Then, for each of the models we adopt three different setups for GalfitM as shown in table <ref>. In all setups, the fluxes of both components are left free, the centroids of galaxies are held constant over wavelength (we assume that the images were properly aligned). The position angles of the galaxy and the axis ratios are also kept constant since these quantities are not expected to present strong wavelength dependence. The most critical parameters are the index and the effective radius. We explore the effect of the wavelength dependence of the size by allowing a quadratic variation in setups 1 and 4 and restricting to constant in setups 2 and 5. Additionally, the maximum degree of freedom is reduced to the number of bands used in each field. For the index of the bulge (+Exponential model) we only allow a linear variation (given that the bulge is normally dominated by old stellar populations, we do not expect a strong wavelength dependence of the index). However, we changed the range from 0-8 in setups 4 and 5 and  to 2.5-8 in setup 6. This is used, as explained in section <ref>, to evaluate our procedure for model selection based on CNNs (see section <ref>). The properties of all runs performed are summarized in table <ref>.



  
  

§ QUANTIFICATION OF UNCERTAINTIES

 
 
 In the following sections, we evaluate the overall accuracy of our final bulge/disk catalog using different methods.
 
 
 

 §.§ Accuracy of model selection

 In order to test the validity of our methodology to select the best model, we compare the outputs delivered by galfitM with the expectations according to the CNN based classes. If both the best model class and the fitting procedure work as expected, one would expect that the best fit model converges towards the expected best profile. In figure <ref>, we show the H-band index distributions of the bulge component for galaxies classified in the 4 profile classes detailed previously (using a probability threshold of 0.4). For obvious reasons, for objects for which a single model is preferred, we plot the global index as well as for pure exponential profiles. Figure <ref> summarizes the criteria used for the selection. We clearly see that the distributions are different for every type of model and follow the expected trends. Pure disks and pseudo-bulges have almost all indices lower than 2. Pure bulges peak at values of n∼ 3-4. The distribution for objects that require a 2 component fit extends to large values as well. However, there is a fraction of objects for which our CNN based model selection technique would have preferred a model with a high index bulge while the fitting procedure converges to a solution with a lower value (dashed orange line in figure <ref>). Given that we expect a contamination of ∼15% (see fig. <ref>), our results is higher than anticipated. In order to test if this is a problem of galfit converging to a local minimum, we use the results of setup 6 (table <ref>) in which the index of the bulge was forced to be larger than 2.5 at all wavelengths. For approximately 50% of the objects, the fitting procedure converged to a new solution with n_b exactly equal to 2.5, i.e. the boundary condition. We considered therefore that for these objects a low index bulge is the best solution. However, for the remaining 50%, the new setup provided a solution with n_b>2.5 in agreement with the CNN classification. The corrected distribution is shown in figure <ref>. As expected, the number of bulges with n_b<2, has been reduced and is now compatible with the expected contamination level. We notice that, if no model selection is applied, half of the bulges have low index values (<2) as shown in figure <ref>.  







In figure <ref>, we show the distribution of the difference between the sizes of disks and bulges for objects identified as requiring two components. As expected, for the vast majority of the objects, the disk component has a larger effective radius than the bulge. For comparison, we also show the same distribution for objects that were classified as pure disks by the CNN model. In that case, roughly half of the population has a bulge larger than the disk. This confirms that a 2 component model is not suited for these galaxies. Including the bulges of these objects in any scientific analysis would definitely introduce a systematic error that could potentially bias the results.







 

 §.§ Statistic systematic uncertainties with simulations



 
We use simulations to estimate the global accuracy of the fits as is commonly done in the literature. In particular, one key quantity to calibrate is the limiting magnitude (or S/N) below which measurements remain unbiased. From <cit.>, we know that for galaxies fainter than F160W=24.5 (at the wide depth), statistical errors on the structural parameters derived from 1 component fits exceed ∼ 20%. This threshold needs to be calibrated for 2 component fits. 

 In order to do so, we follow a standard procedure based on simulations of analytic profiles as done in several previous works (e.g. ). Namely, we generate mock galaxies with 2 analytic profiles (a profile for the bulge component and an exponential for the disk). The structural parameters of each component as well as the total magnitudes and B/T ratios are randomly distributed to span the range of values found in real data. We then convolve each profile with a real PSF and embed the galaxies in a real background. 

A key difference compared to previous approaches is that our fitting procedure takes into account all wavelengths simultaneously. This needs to be captured in the simulations too, i.e. the profile needs to be simulated in every band. A random distribution of all parameters as a function of wavelength is not a good approximation to reality. As done in <cit.>, we simplified the problem by selecting a real low redshift galaxy (z∼0.5) that clearly has two components (from visual inspection and our CNN classification) from which we obtained a bulge/disk decomposition (see figure <ref>). We then use the best SED models for the bulge and the disk as templates for our simulations (see section <ref> for details on the SED fitting). 




We first assign a random B/T (in the i-band rest-frame) and a total magnitude (in the H band) in the range of [18-25] to every simulated galaxy. This allows us to fix the bulge and disk fluxes in the H band. We also associate a random redshift in the range [0.01,3] and shift the SED template described above accordingly.  
We can then associate a realistic magnitude in all other bands with a typical SED of a bulge and a disk while also accounting for the real redshift distribution. The surface brightness profile of the disk is rendered using an exponential profile. For the bulge  component, we draw a random value of the index between 2 and 4, and we assume no wavelength dependence. In order to have a realistic wavelength dependence for the radii, we used random examples from the real data, requiring that the bulge always be smaller than the disk. The final simulated sample contains  ∼ 4000 galaxies.

We then run Galapagos-2 with exactly the same settings used for the real data and compare the input and recovered parameters to assess statistical errors. 



A complete visualization of all the parameter space is complex given the number of parameters. We decided to highlight the dependence of magnitude and size of bulges and disks on galaxy morphology (quantified through B/T), magnitude and redshift. Indeed figure <ref> compares input and output magnitudes and sizes for two different filters: F850 and F160W. The accuracy of the photometry and stellar masses is discussed in section <ref>. As expected, we observe that the errors in the structural parameters of the disk (bulge) increase towards high (low) B/T values. When one of the two components dominates over the other it becomes more difficult to quantify both components properly.  
 For fainter galaxies (F160>23-23.5) the bulge magnitude start deviating from the zero bias line. This suggests that we have reached the S/N limit below which results become biased (as pointed out by  for 1 component fits). We apply a magnitude cut at F160W=23 to keep a zero bias and a scatter lower than 30%. Only galaxies brighter than this magnitude limit will be included in the final catalog released with this work. 
 Similar behavior is observed for the sizes. However, globally the errors remain within ∼10-20%.
 We emphasize that the simulations performed are not representative of the real evolution of galaxy SEDs since we are using a unique template. The effect of the redshift on the measurements is explored in section <ref> and figure <ref>.


 

 §.§ Individual errors

 
The above procedure has provided a statistical analysis of the systematic uncertainties which are on average close to zero, meaning that the method is intrinsically unbiased. In this section we estimate individual uncertainties (random and systematic) for all the galaxies in the catalog. These are particularly important to derive stellar population parameters through SED fitting.
Although Galfit provides error bars, they are known to be underestimated and not very indicative of the true error (see e.g. ). This is generally explained because Galfit assumes a gaussian noise distribution when computing the errors (which is not a good approximation for HST drizzled images) and also does not consider the effects of companions. In GalfitM, the situation is even more dramatic since the fits in all bands are not independent and the constraints in the different parameters are therefore coupled. The galfit manual states that in such a case, the estimation of errors is not reliable. We thus use an alternative approach to estimate the errors. 

The relation between the measured value of a given parameter and the true value is approximated by:

    Measured = True+ b ±σ


The bias (b) or systematic uncertainties is expected to be close to zero as seen from the simulations and therefore the error budget is dominated by the random uncertainties (σ).
We adopt a similar method as the one described in  <cit.>. The method main assumption is that galaxies with similar structural properties share similar errors.

We first use the analytic simulations of section 5.2 to estimate the systematic uncertainties. Since the ground truth is known and the same model is  used to generate galaxies and to fit them, the differences between the input and output in this idealized case can be interpreted as the intrinsic systematic errors.
We thus define for every galaxy in the catalog and the simulation the following vector 𝐩 :


    𝐩=(m/σ_m,log(Re_b)/σ_log(Re_b),log(Re_d)/σ_log(Re_d),n_b/log(n_b),BT/log(BT))



In the above equation, m designates the apparent magnitude in a given filter, n is the index, Re_d and Re_b are the effective radii of the disk and the bulge component and B/T is the ratio between the flux of the bulge and the total flux. Each value is normalized by the dispersion to have similar variation ranges for all parameters. 
We then compute the euclidian distance of each real galaxy to all simulated objects and select the 30 closest objects. The bias for that galaxy is estimated as the 3-sigma clipped median of the difference between the input values of the simulations and the ones recovered by .

To estimate the random errors, we cannot rely on the simulations since they do not capture all the complexity of real data.  Since we have fitted all galaxies with different settings for each type of model (+ Exponential, Sersic) the differences on the resulting fits between the two settings can be used to estimate random uncertainties.  By doing this we assume that both settings have similar systematic biases. Therefore,  for each galaxy in the catalog and every band we compute the following two vectors 𝐩^1𝐂 𝐩^2𝐂 for 1-component and 2-component respectively:

 
    𝐩^1𝐂=(m/σ_m,log(n)/σ_log(n),log(R_e)/σ_log(R_e),f_sph,disk,irr)



    𝐩^2𝐂=(m/σ_m,log(n_b)/σ_log(n_b),log(Re_d)/σ_log(Re_d),log(Re_b)/σ_log(Re_b), 
    BT/log(BT),f_sph,disk,irr)



The meaning of the  different symbols are the same of the previous equations. f_sph, f_disk and f_irr are the probabilities that the galaxy looks like a spheroid, disk or irregular respectively <cit.>. 

Results for systematic and random uncertainties are shown in figures <ref> and <ref>. The systematic errors are close to zero over most of the parameter space. The measurements are essentially dominated by random uncertainties. The typical random uncertainty for the magnitude of the bulge and disk depends mostly on  the total magnitude and the bulge-to-total ratio as already assessed by the simulations from the previous section. Indeed the random errors increase for galaxies magnitude larger then 23 as it can be see in the bottom panels of figure  <ref>. Moreover, for objects with B/T<0.2, the error on the bulge magnitude increases to 0.6-0.7 magnitudes and the error on the disk magnitude rises to ∼0.5 for B/T>0.8. This confirms that the properties of embedded components are estimated accurately only if the dominant component represents less than 80% of the light. There is little or no dependence of the magnitude errors on other parameters such as the size or the index of the bulge. Regarding the errors on sizes of both components, the average error for the bulge is ∼20% and ∼10% for the disks, with again, a dependence on B/T.
Finally, the two panels in the bottom side of figure <ref> explore the dependence of size and magnitude errors with the redshift. We do not observe strong gradients or correlations. This is due to the fact that the majority of our galaxies are concentrated between 0.5<z<1.6 where the angular scale do not drastically change, thus the uncertainties on the size are not affected by this effect. In addition to that the absence of correlation is related on the magnitude selection. Indeed, by selecting bright galaxies (H<23), we are also selecting the largest ones. These two effects dominate the final trends observed.










 
  
  



















 
 
 
 

 Figure <ref> shows the error dependence on wavelength. Both errors on magnitude and size are larger in shorter wavelengths. This is somehow expected since bluer bands are shallower. It is also easy to explain that bulges are more severely affected since they are expected to be redder and therefore fainter in the bluer wavelengths.
 
  
 
 
  

 §.§ Comparison with the literature


Finally, we also perform a comparison with the literature. Namely,  did a 1-component fit to all galaxies in CANDELS down to H=24.5 in two NIR filters (F125 and F160) independently. As described in the previous sections, the method used in this work differs from  in the sense that all bands are fitted simultaneously. Although this technique is intended to benefit from a better S/N, one needs to make sure that no systematics are introduced in the process. At least for bright objects, similar results should be obtained in both works.
Figure <ref> compares the magnitude, index and half light radii from our catalog with the ones from <cit.>. There is reasonable agreement with a systematic difference compatible with zero and a scatter on the order of ∼10% increasing at fainter magnitudes as expected. The scatter is on the order of the error reported in the measurement as discussed by <cit.>. This result confirms that our procedure works as expected at least for a 1 component fit and that no systematic biases are introduced by using all filters jointly.








§ STELLAR POPULATION PROPERTIES OF BULGES AND DISKS



Using the methods described in the previous sections, we obtain 4-7 point SEDs for bulges and disks in our sample together with a measure of uncertainties in the derived fluxes in each band (see section <ref>). We then perform a standard SED fitting with the FAST code <cit.> and derive stellar masses, SFRs, ages and metallicities for all bulges and disks. Additionally, using the best-fit models we also derive rest-frame U,V,J colors for bulges and disks. The input models are grids of Bruzal & Charlot (2003) models that assume a Chabrier (2003) IMF and Calzetti et al. (2000) extinction law. For the star formation history (SFH) we applied an exponentially declining model, with tau in the range of [8:10]. 

Figure <ref> shows some examples of the best-fit templates for galaxies with one or two components according to the CNN based classification. More examples can be checked online in the public release of the catalog[<https://lerma.obspm.fr/huertas/form_CANDELS>]. 
 
 
 
The mass distributions obtained for bulges and disks are shown in the two panels of figure <ref>.

Random uncertainties in the stellar population parameters are estimated by performing MonteCarlo simulations for each galaxy. Assuming the errors on the fluxes as computed in section <ref>, we generate 500 SED realizations for each galaxy that we fit with FAST. The uncertainty on each parameter is derived as the minimum and maximum value within the 68% of realizations with the smallest reduced χ^2 values. Rest-frame colors are also computed for each realization and the uncertainty is computed in an analogous way.

Using the simulations described in <ref>, it is also possible to estimate a global statistical uncertainty on the stellar populations properties. Since the main quantity that will be used in forthcoming works are B/T ratios, we focus here on the stellar masses of both components. We notice that uncertainties in other parameters might be large given the reduced wavelength coverage and therefore should be used with caution. 

 We first run FAST on the simulated SEDs of both components and then perform a second run on the recovered fluxes from GalfitM on the same simulated galaxies.  Figure <ref> shows the comparison of both estimates. We find that the bias is close to zero and the dispersion is on the order of ∼0.2-0.3 dex, which is the typical error expected for SED based stellar masses. This results in an unbiased estimate of the stellar-mass bulge-to-total ratio with a typical scatter of ∼0.2. Notice that this does not mean that the true stellar mass is recovered. It indicates however that GalfitM properly recovers the fluxes of the disk and bulge components without introducing additional systematics as already demonstrated for the sizes in section <ref>. Another additional check is shown in the right panel of figure <ref>, which compares the stellar mass obtained by fitting the photometry of a 1-component fit to the stellar mass obtained by adding the masses of the bulge and disk components. The two estimates agree within a ∼0.2 dex uncertainty as expected. 


Finally, in order to have an independent estimate of the reliability of the stellar masses derived with only 4-7 filters, we compare our values with the stellar masses estimated in CANDELS which cover a larger spectral range. This is shown in the two first panels of figure <ref>. Our estimates are unbiased, but with a scatter of ∼0.4. This scatter should be a combination of model dependence (we did not use the same stellar populations models) and spectral sampling. As a matter of fact, when the sample is divided between galaxies for which we have a spectral sampling of 4 points or less and the others, the results are still unbiased but the scatter increases from ∼0.3 with 4+ filters to ∼0.4 with less than 4 filters.





§ FINAL CATALOG


The catalog release is made in three different tables. The same unique identifier is used to link all tables. 




  * Gold catalog: Contains a table with our best selected configuration for each galaxy.  For each galaxy we provide the structural parameters in all bands derived from what we think is the best setup (see Table <ref>) according to the CNN best class. The selection essentially follows the flow chart of figure <ref>. For each parameter we provide an errorbar computed with the method described in section <ref>. Additionally, we provide an estimate of the stellar mass bulge-to-total ratio as well as rest-frame colors for bulges and disks in the two GOODS fields. For standard use, we recommend this catalog.


  * An additional table with all structural measurements in all bands and setups described in this work is also provided along with the CNN outputs. This table should be used to explore different selections.


  * A third table with the derived stellar population properties for bulges and disks resulting from the SED fitting as well as the uncertainties. This table allows one to compute bulge to total stellar mass ratios. 



 






§ SUMMARY AND CONCLUSIONS


This work presents a catalog of multi-wavelength bulge-disk decompositions of ∼17.600 galaxies in the five CANDELS fields. The dataset is mass complete down to ∼5×10^10 M_⊙  at z∼2. 

Each galaxy is fitted with a 1-component model and a 2 component (bulge) + exponential (disk) model and with three different setups each in which we modify the wavelength dependence of the size and index of the bulge component. We used the galfitM/galapagos2 code from the megamorph project, allowing US to simultaneously fit images at different wavelengths.  One key new ingredient of this work is that we introduced a new method to address the systematic uncertainties arising from the use of a wrong model to describe the surface brightness profile of galaxies. The technique, based on convolutional neural networks, provides a quantitative measurement  of how well a given profile (1 or 2 components) describes a given galaxy. We show that our proposed method can distinguish between different profiles with a ∼80-90% accuracy. It also allows us to reduce the contamination from unphysical components at the ∼10% level and study a clean sample of bulges and disks. 

Through extensive simulations, we show that the structural parameters derived for bulges and disks are globally unbiased. We develop a method based on the comparison of results from different runs to provide individual error bars for each structural parameter. We show that the typical error for the bulge (disk) magnitudes are ∼0.2 (∼0.1) mags. For galaxies with B/T<0.2 (B/T>0.8) the error on the bulge (disk) magnitudes increases to ∼0.5 mags. Sizes are estimated within ∼10-20% uncertainty.

The derived spectral energy distributions of both components with realistic errors are then fitted with stellar population models to estimate stellar masses and bulge-to-total mass ratios. We show that the statistical uncertainties are on the order of ∼20%. We also provide rest-frame colors and SFRs which will be analyzed in forthcoming woks.

The catalog including all derived quantities is made public with the present work. 



§ ACKNOWLEDGMENTS


The present work has been funded by a French Doctoral Contract. The authors are also grateful to Google for the unrestricted gift given to the University of Santa Cruz to carry out the project: "Deep-Learning for Galaxies". 
  
   
mnras

  

