MegaMorph – multiwavelength measurement of galaxy structure: physically meaningful bulge–disc decomposition of galaxies near and far

Abstract

Bulge–disc decomposition is a valuable tool for understanding galaxies. However, achieving robust measurements of component properties is difficult, even with high-quality imaging, and it becomes even more so with the imaging typical of large surveys. In this paper, we consider the advantages of a new, multiband approach to galaxy fitting. We perform automated bulge–disc decompositions for 163 nearby galaxies, by simultaneously fitting multiple images taken in five photometric filters. We show that we are able to recover structural measurements that agree well with various other works, and confirm a number of key results. We additionally use our results to illustrate the link between total Sérsic index and bulge–disc structure, and demonstrate that the visually classification of lenticular galaxies is strongly dependent on the inclination of their disc component. By simulating the same set of galaxies as they would appear if observed at a range of redshifts, we are able to study the behaviour of bulge–disc decompositions as data quality diminishes. We examine how our multiband fits perform, and compare to the results of more conventional, single-band methods. Multiband fitting improves the measurement of all parameters, but particularly the bulge-to-total flux ratio and component colours. We therefore encourage the use of this approach with future surveys.

1 INTRODUCTION

The spatial distribution of light within a galaxy is a key observable, with which we can constrain models of galaxy formation and evolution. The typical sizes, surface-brightness profiles and ellipticities of galaxies have been essential in determining the main physical mechanisms at work in producing the galaxy population (e.g. Silk & Norman 1981). We have grown to understand that these properties are the result of multiple competing processes, including rapid collapse, ongoing gas accretion, disc instabilities and the merging of existing stellar systems (e.g. Benson 2010). As a consequence, galaxies are often separable, at least to a degree, into components with distinct spatial structure, kinematics and stellar populations.

Observations often integrate over these components, e.g. aperture photometry, to estimate overall properties for each galaxy. With such quantities, one can gain a general picture of the merger and star formation history of a given galaxy. Typically, however, a variety of histories can produce similar integrated properties. Considering the properties of a galaxy's components separately enables a much more detailed account of its lifetime to be constructed.

The simplest approach to separating the properties of the main galaxy structures is bulge–disc decomposition. This can be applied to imaging data alone, and hence to the largest samples of galaxies available. Although conceptually simple, bulge–disc decomposition remains a challenging task, due to the variety of structures that galaxies display, not to mention the usual observational limitations of resolution and signal to noise (S/N).

A common method is to study and model the one-dimensional (1D) light profile, along the major or minor axis of the galaxy, or azimuthally averaged. These 1D profiles are usually obtained by fitting a set of ellipses to the isophotes in the (2D) image. However, 1D representations of the radial surface-brightness distribution suffer from strong systematic uncertainties since they neglect the differing intrinsic shapes of the disc and bulge components.

A solution to this problem is two-dimensional (2D) decomposition (Byun & Freeman 1995), which utilizes all the spatial information in the images (for more details about standard 1D and 2D methods, see Peng et al. 2010). On the other hand, fitting in 2D is usually more sensitive to features, such as bars and spiral arms, which are difficult to model. The usual procedure for 2D bulge–disc decomposition is to fit a parametric model to the image, accounting for the point spread function (PSF), pixelization and noise properties of the image. The projected surface-brightness profile of each component is typically modelled using an analytic function, the most common choice being the Sérsic profile (Sersic 1968).

Separating galaxy components is supposedly easier for small samples of nearby galaxies where a more interactive fitting process can be applied. Multiple studies have applied 2D decomposition to examine the correlations between bulge and disc properties at optical to infrared wavelengths (e.g. Noordermeer & van der Hulst 2007; Méndez-Abreu et al. 2008; Barway et al. 2009; Tasca & White 2011), to study the coevolution of supermassive black holes and their host galaxy (e.g. Kim et al. 2008; Vika et al. 2012), to investigate the evolution of structure over cosmic time (e.g. Bruce et al. 2012, 2014) and environment (Hudson et al. 2010; Head et al. 2014), to study the structural properties of isolated late-type galaxies (e.g. Durbala et al. 2008), and to measure quasar host galaxy parameters (e.g. McLure, Dunlop & Kukula 2000).

Some studies go a step further and attempt to decompose a third component, usually a bar (e.g. Laurikainen, Salo & Buta 2005; Gadotti 2009, 2011; Weinzirl et al. 2009). In addition to providing measurements of bar properties for study, including a potential bar in the model helps to avoid any such feature from contaminating measurements of the bulge and disc.

A significant issue lies in identifying which components are present, and hence which model parameters are to be trusted. This amounts to choosing the appropriate complexity of model for a given galaxy. Fitting a more complex model usually results in a significantly improved goodness-of-fit statistic (e.g. chi-squared), irrespective of whether or not the model parameters are physically meaningful. This problem is complicated by the presence of galaxy features that are not included in the model, such as cores, non-elliptical and twisted isophotes, dust lanes, etc. Many studies ultimately resort to selecting the most appropriate model by visual inspection of the original images and their fit residuals.

Elliptical galaxies are usually regarded to be one-component systems, and hence they are usually chosen to be modelled by a single-Sérsic profile. However, it is far from clear whether this is physically the best way to describe these systems. Taking a different approach, Huang et al. (2013) fitted three components to each member of a sample of elliptical galaxies, finding that these galaxies can be well described by the combination of three Sérsic profiles, each with low Sérsic index but different effective radii. Huang et al. argue that these components have physically meaningful interpretations. The intermediate-size component is the original, built from early collapse and major mergers. The largest component is comprised of stars accreted in more recent minor mergers. Finally, the most compact component is attributed to central star formation following the dissipative accretion of gas brought in by some of those recent minor mergers.

For large samples of galaxies, more automated approaches to deciding how many components a galaxy comprises are essential. For example, Allen et al. (2006) employed a logical filter to decide whether the results of fitting a bulge–disc model were physically plausible, or whether their single-Sérsic fit should be preferred. They showed that the routine structural decomposition is an important for understanding the bimodality of galactic properties. Simard et al. (2011) have created the largest catalogue of multicomponent galaxy structure to date. They fit one million galaxies with three different models, and used F-tests with a calibrated probability threshold to choose the best model for each galaxy. Lackner & Gunn (2012) expanded the model options five, selecting between them using a logical filter. These studies have provided the first complete estimates of the bulge and disc properties for the local Universe.

To date, most studies have measured structural properties of a galaxy using only one image, in a single waveband. However, modern surveys provide images of the same galaxies in many different bands. In some cases, models are fitted to each band independently. This does not produce reliable colours, however, so more often an initial model is fit to one preferred band, then the structural parameters are fixed during fits to the other bands. Simard et al. (2011, following Simard et al. 2002) take a more consistent approach by fitting their models to images in two bands simultaneously, while Mendel et al. (2014) use a hybrid procedure to produce bulge and disc colours in five optical bands.

Until recently, no method was available that could fit models to an arbitrary number of images at different wavelengths. Driven by a determination to make more effective use of the multiwavelength imaging available from modern surveys, the MegaMorph project (Häußler et al. 2013; Vika et al. 2013; Bamford, in preparation) developed and tested a new version of 2D photometric analysis which constrains a single, wavelength-dependent model using multiple images simultaneously.

This paper is one of a series that investigates the benefits of this multiwavelength approach to measuring galaxy structural properties. In Bamford et al. (in preparation), we present this new tool in detail, describing the new features and demonstrating its use through some specific examples. In Vika et al. (2013, hereafter V13), we test our new method by fitting single-Sérsic models to original and artificially redshifted image of 163 nearby galaxies. In Häußler et al. (2013), we demonstrate our approach on a large data set from the GAMA (Driver et al. 2009; Hill et al. 2011) survey, automating both the preparation of the data and the fitting process itself. The resulting measurements – in particular the variation of structural parameters with wavelength – are studied further in Vulcani et al. (2014). The objective of this paper, is to investigate the ability of galfitm to perform bulge–disc decomposition on galaxy images with a wide range of resolution and S/N. This is achieved by analysing the same sample as V13: large, nearby galaxies in the Sloan Digital Sky Survey (SDSS; Abazajian et al. 2009), with both original images and versions that have been convolved and resampled in order to simulate the galaxies’ appearance at a range of redshifts. A complementary analysis of multiband bulge–disc decomposition, using the same GAMA sample as Häußler et al. (2013), will be presented in a forthcoming paper (Haeussler et al., in preparation).

This paper is structured as follows: in Section 2, we present our data set, give a brief description of galfitm, and then explain how we fit our sample and identify reliable components. In Section 3, we present the distributions of structural parameters obtained from the original SDSS imaging, and examine the stability of these distributions with respect to the effects of distance. In Section 4, we present correlation between structural parameters and a way of separating elliptical from lenticular galaxies in our sample. We provide a summary in Section 5.

2 DATA

2.1 Sample selection and imaging

In this paper, we use the same set of 4026 galaxy images as in V13. These images comprise a sample of 163 nearby galaxies with imaging from SDSS in the u, g, r, i and z passbands. Our galaxies typically extend over more than one SDSS frame, so to create them we employ montage (Jacob et al. 2010), which performs the transformations, rebinning and background adjustment necessary to combine the individual frames into a single mosaic.

In addition to the original images of our galaxies, we use 3863 further images in which the galaxies have been artificially redshifted. We use ferengi (Barden, Jahnke & Häußler 2008) to create a set of ugriz images mimicking the appearance that each of the 163 nearby galaxies would have if they were observed by SDSS at a range of redshifts. The artificial redshifting algorithm applies cosmological changes in angular size, surface brightness and, optionally, shifting of the rest-frame passband (k-correction), to simulate the observation of a given galaxy at a greater distance. We produce images for redshifts 0.01–0.25, in steps of 0.01. Note that not all our galaxies have images for every one of these redshifts, either because they originally have a redshift higher than 0.01, or because the galaxy effectively becomes a point source. As in V13, to avoid confusion with genuine redshift biases, we disable the k-correction feature of ferengi in this work.

Full details of the redshifting process and the data preparation for both original and redshifted images have been given in V13. In this paper, we use the same masks, PSFs and sky estimates.

The morphological breakdown of our galaxy sample is given in the top row of Table 1. All classes, except elliptical, also include barred types. Note that, while we ensure a broad range of morphologies are included, the distribution of Hubble types in this sample is not representative of the local Universe. As described in V13, our sample is comprised of galaxies which have had their structure carefully measured by previous studies. We can then compare these with our semi-automated, multi and single-wavelength results.

2.2 Structural parameters

2.2.1 Fitting galaxies with galfitm

We use a modified version of galfit3 to fit 2D analytic models to our galaxy images. A detailed description of galfit is given by Peng et al. (2002, 2010). We have adapted galfit (version 3.0.2) for the requirements of this project, as briefly described below. To differentiate it from the standard release, we refer to our modified version as galfitm.¹ All the work in this paper uses version galfitm-0.1.3.1.

Standard galfit3 accepts only a single input image with which to constrain the model fit. To utilize multiband data, it was therefore necessary to make a number of significant modifications. However, we have endeavoured to retain the original code unchanged, wherever possible. galfitm is therefore backward compatible and produces almost identical results to galfit3 when used with single-band data (see section 4.1 in V13).

Our modified code can accept an arbitrary number of (pixel-registered) images of the same region of sky at different wavelengths. To these images, galfitm fits a single, wavelength-dependent, model. As for galfit3, this model may comprise one or more components, each with a number of parameters. For example, for a single-Sérsic function, the parameters are centre position (x_c, y_c), magnitude (m), effective radius (r_e), Sérsic index (n), axial ratio (b/a) and position angle (PA). To enable these model components to vary with wavelength, each of their standard parameters is replaced by functions of wavelength. For convenience, these are chosen to be Chebyshev polynomials (see Häußler et al. (2013) and Bamford et al. (in preparation) for more details). Instead of directly fitting the standard parameters, galfitm optimizes the coefficients of these polynomials to best match all the multiband data.

For each standard parameter, the user may select the order of the polynomial that describes its wavelength dependence, and thereby the freedom that parameter has to vary. Some parameters may be entirely fixed; for example in this paper, we set n = 1 for the second Sérsic function, in order to model an exponential disc component. Other parameters may be allowed to vary as a constant with wavelength; e.g. one might allow the central x and y coordinates to vary during the fit, but require that they are the same in every band. Further parameters may be permitted to vary with wavelength as linear, quadratic, or higher order functions; e.g. one could choose to allow the axis ratio b/a to vary linearly with wavelength, in order to account for changes in ellipticity at different wavelengths. Ultimately, specifying a polynomial with as many coefficients as there are input bands allows that parameter to vary freely with wavelength. The user therefore has great flexibility to achieve a compromise between the freedom of the model, physical expectations, and the number of free parameters that must be constrained by the data.

A key element of our technique is that parameter values at the wavelengths of low-S/N images can be partially interpolated or extrapolated from the higher S/N data. However, any significant signal present in those images should have an appropriate influence on the fit. Systematic biases will be thus be reduced in comparison to an extrapolation based only on the high-S/N bands. The risk of such systematics may be further reduced by giving the model more freedom to vary with wavelength (e.g. linear or quadratic variation of r_{e, b}, r_{e, d} and n_b). The cost is an increase in the statistical uncertainties of the parameters in the low-S/N bands.

In this work, we have held most structural parameters fixed versus wavelength. This corresponds to assuming a simplified picture of galaxy structure, in which galaxies comprise only bulge and exponential disc components, each without colour gradients and not departing from an elliptical projected Sérsic profile. Real galaxies may not obey these assumptions. Therefore, while our structural constraints increase the stability of the fits, there is also a risk of introducing systematic biases in cases where the true wavelength dependence of the profile does not correspond to that assumed. Lower S/N images would be most susceptible to such systematics, as their parameters will be influenced by any higher S/N data.

Our approach makes the assumed variation of galaxy structure with wavelength explicit, and allows one to relax these assumptions in a selective and gradual manner. This flexibility allows the user to balance systematic and statistical uncertainties, using independent observational results, physical insight and knowledge of their data set.

We plan to explore the variation of structural parameters with wavelength in detail in a future paper. However, based on preliminary results, we find that the vast majority of results present in this paper do not change by allowing small wavelength variation of the structural parameters.

2.2.2 Model choices

Most galaxies are considered to be primarily two-component systems, comprising a disc, with an exponential (n = 1) profile, and a bulge, typically well represented by a Sérsic function with n ∼ 0.5–4 (Graham 2001). On the other hand, elliptical galaxies are generally considered to be single-component systems, describable by a single-Sérsic profile.

In V13, we performed single-Sérsic fits to all of our images. In this paper, we supplement these with bulge–disc decompositions, performed using galfitm to fit two superimposed elliptical Sérsic models. For the first component, we fit all standard parameters (effective radius r_e, apparent magnitude m, Sérsic index n, axis ratio b/a and position angle PA), while for the second, we fix the Sérsic index to be equal to 1.

In order to carry out a blind test of our analysis, and allow us to draw conclusions relevant for large surveys, we assume that we do not know the Hubble classification for our galaxies. We therefore fit them all (even known ellipticals) with two functions. In Section 4, we will explore what can be learned from this, including the potential for separating single-component systems, i.e. ellipticals and pure disc galaxies, from two-component systems, using structural parameters diagnostics.

We run two sets of fits, each of which is performed on the u-, g-, r-, i-, z-band images for all of our original and redshifted galaxies. For the first (single-band fitting; or SM, to reuse the nomenclature of V13), we fit each band individually. For the second (multiband fitting; MM), we fit each galaxy using all five bands simultaneously.

We allow the magnitudes (for both bulge and disc) to vary completely freely between bands. For the multiband method this amounts to setting the wavelength dependence of magnitude to be described by a quartic polynomial, with as many free coefficients as the number of bands. We allow full freedom as we wish to avoid any potential biases on the recovered magnitudes, and hence colours, which may result from assuming a lower order polynomial dependence.

For the effective radius and Sérsic index, we choose to not permit any variation with wavelength. This effectively ignores colour gradients within each component, but keeps the overall number of free model parameters down, hopefully improving the reliability of the decomposition process. Any measurements of the wavelength dependence of individual components will be noisy and are unlikely to provide significant evidence to contradict the reasonable default position of a constant value. This is therefore what we assume. Our decision is supported by previous results in the literature. For instance, both MacArthur, Courteau & Holtzman (2003) and McDonald et al. (2011) find that the Sérsic index of the bulge, as well as the effective radius of the bulge and the disc, show no significant variation (or a slightly linear relation in rare cases) across optical and NIR wavelengths.

In Fig. 1, we show our results (the effective radius, Sérsic index and spectral energy distribution), for the bulges of six example galaxies fit in our original and artificially redshifted images. For most of our galaxies, the results of the single-band fits (solid lines) show substantial fluctuations with wavelength, which worsen with increasing redshift. The results of our multiband fits (dashed lines), with n and r_e constant with wavelength, recover reasonable values that are close to the average of the higher S/N bands (gri) for single-band fitting. The multiband results are more resilient for different redshifts (e.g. black, red and blue lines in different columns).

We also assume that the shapes of our galaxy components do not change with wavelength, so we set the axis ratio, position angle and galaxy centre to be constant with wavelength.

In both runs (SM, MM), we use the same initial parameters for galaxy centre (x_c, y_c), magnitude (m), Sérsic index (n), effective radius (r_e), position angle (θ), axis ratio (b/a) and sky background value (although different values are used for each galaxy image, see V13 for more discussion of the sky estimate). We experimented with various different schemes for choosing initial parameters values, before selecting the following approach.

The sample of 163 original images is fit first. The initial magnitudes are determined using the MM single-Sérsic model results found in V13. The initial value for the bulge magnitude was set to m_ss + 0.75, where m_ss is the single-Sérsic magnitude. The initial disc magnitude was set to m_ss + 0.65 in order to start with a slightly fainter bulge than disc. The initial effective radius of the bulge was chosen to be 0.5r_{e, ss}, where r_{e, ss} is the effective radius from the single-Sérsic fit. Similarly, the initial effective radius for the disc was set to be equal to the r_{e, ss} We therefore use the observation that bulges are typically smaller than their host discs. We found that starting with an equal bulge and disc there are more chances the bulge to fit parts of the disc component.

The initial Sérsic index of the bulge was chosen to be n_ss, while the Sérsic index of the disc was fixed to unity, for an exponential profile. The initial values of disc axis ratio, disc position angle, disc and the bulge centre (x, y) were set equal to the equivalent single-Sérsic value. In cases where a parameter was variable with wavelength in the single-Sérsic fit, but constant in the current paper, we took the median of the five values. Finally, the initial value of axis ratio for the bulge was set equal to 0.8, and the initial value of position angle for the bulge was arbitrarily set to 10°. We confirmed that the final results of the fits do not depend on the these values.

All the parameters are allowed to vary during the fitting process, with the exception of the disc Sérsic index and the sky background, which were kept fixed. galfit and galfitm give the option to constrain the range of values for each parameter in order to avoid unphysical results. We make use of this option by applying the following constraints. We require both magnitudes (bulge and disc) to vary within the range of 5–35 mag. Similarly, both the effective radius of the bulge and the disc were allowed to vary within the range of 0.04–600 arcsec. We constrain the freedom of the Sérsic index by allowing it to vary within the range of 0.1–15. However, we exclude any bulge with n_b < 0.3 from the final sample. Finally in the case of the centre (x, y), we applied two constraints, one to fix the bulge and disc to have the same centre, and a second to restrict their variation, with respect to the single-Sérsic fit, to be no more than |$\sqrt{s/8}$|⁠, where s is the size of the image.

For each artificially redshifted image, we repeat the same procedure as above to estimate the initial parameter values, but use the MM single-Sérsic result obtained for the same redshifted image. In cases where a single-Sérsic magnitude was unphysically faint, we calculated the initial parameter values by cosmologically adjusting the values obtained for the lowest-redshift artificial image. We apply the same constraints as for with the original galaxies.

In addition to SM and MM, we perform another set of fits to the artificially redshifted images, which we refer to as ‘aperture fits’. For these, we take the structural parameters from the SMr-band results and keep these fixed while performing single-band fits to the u-, g-, i- and z-band images. Only the magnitudes are allowed to vary freely during the fit. In this way, we apply an identical model in all the bands and ensure we only measure the variations in the flux for a fixed ‘aperture’. This approach is commonly applied to ensure meaningful colours.

2.3 Inspection of individual galaxies

Fig. 2 shows the original ugriz images, residuals from the single-Sérsic and bulge–disc MM fits, and the separate MM bulge and disc model components, for an example spiral galaxy The figure also includes various useful numbers for each fit. Using similar figures, we have visually examined the fitted models – and their residuals – for all of our 163 galaxies, to ensure that their shape and size correspond to the real galaxy.

In addition to checking the images, we also inspect all the recovered parameters for both the original and the artificially redshifted images for each galaxy. In a similar manner to V13, Figs 3 and 4 present a summary of the bulge and disc results for two example galaxies. Equivalent plots are available for all the 163 galaxies. The left-hand panel shows the SM results and the right-hand panel the MM results. At redshift zero, we plot the results from the original images. The first row of panels shows the absolute magnitude (M), the second row shows the effective radius, the third row shows the Sérsic index in the case of the bulge panel and the bulge-to-total flux ratio in the case of the disc panels. The last two rows show the axis ratio and the position angle.

In these figures, we determine the absolute magnitude and the effective radius assuming distances simply derived from the observed redshift and adopted cosmology. Therefore, the values shown for the original images in Figs 3 and 4 could differ slightly from later figures, for which we use more directly determined distances when they exist.

Figs 3 and 4 illustrate some of the behaviours seen for our fits as the galaxies are simulated as they would appear at higher redshifts. In Fig. 3, we present the recovered structural parameters for the galaxy NGC 2776. We can see that the single-band results fit different structures in each band, particularly for the bulge, and that the parameters of these structures vary substantially with small changes in the simulated redshift. In contrast, the multiband results are much more stable as a function of redshift, although some small trends are seen before even these results become noisy and unreliable at z ≳ 0.15. The systematic decline in Sérsic index with simulated redshift appears to be a consequence of diminishing spatial resolution, and was also a fairly common feature of the single-Sérsic component fits in V13. For both methods SM and MM, the disc parameters are much more stable than those of the bulge, presumably due to the disc's larger size and less steep inner profile. Note that, even though the bulge and disc magnitudes are completely free to vary between bands, constraining the wavelength variation of n, r_e, PA and b/a through multiband fitting leads to much more stable measurements of the magnitudes, and hence colours.

Fig. 4 presents another set of recovered structural parameters, this time for the galaxy NGC 6314. The first thing to notice is that, for the single-band fit, the effective radius of the bulge is much larger than the disc for the g- and r-bands. This is an indication that, in these bands, the Sérsic function is fitting the disc and the exponential function fitting the bulge, especially given the n behaviour. The usual solution to this problem is to apply constraints on the fit, e.g. insist that the bulge be smaller than the disc. However, using the same set of constraints for the whole sample may introduce biases in other galaxies. Turning to multiband solves this problem without requiring constraints, now the same structural components are fitted in all the bands. In addition, again we see a reduction in the variations caused by small changes in simulated redshift, and that the fit remains more reliable to higher redshifts.

So far we have shown, via specific examples, that that our multiband approach can measure the fluxes and sizes of galaxy bulge and disc components more reliably than if each band is fit individually, at least when we allow no freedom for the r_{e, b}, r_{e, d}, n_b, b/a_b, b/a_d, PA_b and PA_d parameters to vary with wavelength. Substantial variations in the recovered parameters with relatively small changes data quality (redshift) are dramatically reduced. The improved stability is particularly noticeable at low S/N. As a result, it increases the distance out to which meaningful bulge–disc information can be recovered for a galaxy of a given luminosity. In Section 3, we demonstrate these improvements in a more general manner, by considering the average trends of various parameters versus redshift, for our whole galaxy sample.

2.4 Obtaining reliable structural measurements

Before studying the distribution of the galaxy component parameters, we must select a sample with reliable bulge–disc measurements. In this section, we describe the various controls we apply to determine if our fits are physical meaningful. For those galaxies where the fitted model is a poor match to the original image, given our physical expectations, we repeat the fit again with different initial parameter values or additional model components. We aim for our procedure to be applicable in an automated manner, that could be used for large surveys. However, as we have a small sample of galaxies, we still use the tools described in Section 2.3 to inform our choices, to check if the automated selection agrees with visual inspection, and to aid the interpretation of our results.

We start by identifying unphysical models where re-attempting the fit may produce a better outcome. Examining the results of the multiband fits to the original images, we find that 50 galaxies (out of 163) have bulges with r_{e, b}/r_{e, d} > 0.9. After further investigation we separate these galaxies into two cases. The first group consists of 17 galaxies for which the n = 1 component is fitting the inner structure of the galaxy and the free-n component is fitting the disc. For these the measured Sérsic index varies between 0.5 and 1.4. All these galaxies are late type spirals (Sc, Sd, Sm) and 13 of them have a bar. For these cases we believe that the initial parameter values were far from optimal. We choose to fit these galaxies again, using a different set of initial values, with a brighter flux and larger size for the intended disc and lower flux and smaller size for the bulge. The new fit corrects the problem for the vast majority of the 17 galaxies.

The second group consists of 33 galaxies where the bulge fits an inner structure, but also dominates the outer region of the galaxy. In this group, we find 10 early-type galaxies (E/S0), where 9 have bulge Sérsic values 2.8–5.7 and 1, NGC 4458, has n_b = 11.3. The remaining 23 galaxies are spiral galaxies (Sa, Sb, Sc) with bulge Sérsic values between 4 and 11, usually accompanied by a high bulge-to-total ratio. These n_b and B/T values are unusually high for late type spirals. For this second group of galaxies, we initially re-fit the galaxies in the same way as for the first group. This approach corrects the fits for almost one-third of the cases. For the remaining galaxies, we attempt another fit with an addition of a third component (in the form of a central point source), together with the second set of initial parameters. We choose to add a PSF function to account for any extra flux in the centre of these galaxies that could be responsible for the high values of Sérsic index and effective radius.

These new fits return smaller bulges with lower Sérsic indices for another one-third of the cases, and they reduce the Sérsic index without reducing the effective radius for a further four galaxies. For the remainder, we do not adopt the fit results with additional-PSF component, because either the value of the bulge effective radius or Sérsic index was larger than before, or the PSF magnitude was negligible (<30 mag). One case where the addition of a PSF component failed to improve the fit is the elliptical galaxy NGC 4458 (see Appendix A for more details on this galaxy).

After refitting both groups we then select our final samples of trustworthy bulges and discs. We reject bulges in the multiband, original-image results with fit parameters on any of the constraints, insignificant bulge components (at least 3 mag fainter than the disc or below the SDSS point source detection limit), and bulges with effective radius smaller than 5 pixels. We perform these checks in each individual band. The final numbers of galaxies with acceptable bulge measurements are shown in Table 1. A large fraction of late-type spiral galaxies do not possess a reliable bulge measurement. For all these cases, we trust the disc parameters but exclude the bulge properties from further analysis.

In the case of the disc component, we find a few occasions where the effective radius of the disc has taken unreasonably large values or its total brightness is more than 3 mag fainter than the bulge. For these cases, the bulge component returns very similar results to the single-Sérsic fit (Paper II), the disc only accounts for a minor details in the residual. Four of these galaxies (NGC 4360, NGC 4378, NGC 4486 and NGC 4621) are elliptical galaxies and consequently may indeed lack a disc. In addition, we fail to fit a disc component for NGC 4459, even though it is classified as S0. For these objects we trust the bulge parameters but exclude the disc properties from further analysis.

Two problematic cases that were discovered through the above checks are the galaxies NGC 4378 and NGC 4450. For these two galaxies, both the bulge and disc components have been removed. While both these galaxies are two component galaxies, some image distortions hinder the bulge–disc decomposition process.

Around 110 galaxies have both bulge and disc components that we deem as trustworthy, depending on the band. We should stress that this sample includes 19 elliptical galaxies that were fitted with two functions. For these galaxies, we do not know if these two components correspond to truly distinct structures with different kinematics, but we choose to keep them in the further analysis. So far we have not found any indication that their fits are inappropriate, other than their visual classification,² which would not be available for a large automated sample. The number of galaxies that have both trustworthy bulge and disc measurements are broken down by morphology and band in Table 1.

We attempt to identify elliptical galaxies by comparing the goodness-of-fit of our one- and two-component fits. In Fig. 5, we plot the difference in reduced-χ² between the single-Sérsic and bulge–disc fits as a function of bulge-to-total flux ratio. The improvement in reduced-χ² for the elliptical galaxies is in the same range as the other Hubble categories. We find that 18 per cent of our elliptical galaxies and 33 per cent of our lenticular have a dramatic reduction (⁠|$\chi _{{\rm SS}}^2-\chi _{{\rm BD}}^2>0.1$|⁠) in their reduced-χ² by adding an extra exponential function. Only 13 per cent of our elliptical galaxies and 17 per cent of the S0 show a negligible change in reduced-χ² (⁠|$\chi _{{\rm SS}}^2-\chi _{{\rm BD}}^2<0.01$|⁠). In cases where a PSF function has been included we use the |$\chi _{{\rm BD}}^2$| for that fit. The addition of the PSF function in all the cases improved the |$\chi _{{\rm BD}}^2$| by less than 0.01.

Chi-squared should always decrease when a model is given more freedom. The Bayesian Information Criterion (BIC) is based on χ², but penalises additional parameters, in an attempt to provide a guide to whether the additional freedom is warranted by the data. However, applying the BIC to our data finds only nine cases where the single-Sérsic fit is deemed better than the bulge–disc model. All these cases have already been identified as having an insignificant bulge or disc by our above selection criteria. We conclude that, by using the reduced-χ² or BIC, we cannot select a clean sample of elliptical galaxies.

Quantifying the uncertainties of the bulge and disc structural measurements is a challenging task. In V13, we provided the following uncertainties for our single-Sérsic fits: (u, g, r, i, z) for m (±0.13, ±0.09, ±0.1, ±0.11, ±0.12), r_e (±12 per cent, ±11 per cent, ±12 per cent, ±14 per cent, ±15 per cent) and n (±9 per cent, ±11 per cent, ±14 per cent, ±15 per cent, ±17 per cent). These were based on plausible systematic uncertainties in the sky estimation, which typically dominates the error budget. As bulge–disc decomposition is a more complicated task, we expect the uncertainties on our bulge and disc measurements to be even larger. Overall, as we will also see from the further analysis, the bulge parameters are more dependent on the initial conditions, while the disc parameters are more robust. We attempt to determine indicative uncertainties on our fit parameters by refitting a randomly selected sample of 10 galaxies with different sky values. We alter the sky values by our estimated systematic sky uncertainties, as before. We find that both n_b and r_{e, b} can change by up to ∼25 per cent. The parameters of the disc are less strongly affected and are on the same level as the single-Sérsic uncertainties.

As a simple check of our fits, Fig. 6 shows the difference in the recovered r-band total magnitude between our single-Sérsic and two-component models. We colour-code the galaxies based on their Hubble classification. For the vast majority of the galaxies the difference between the total magnitudes is smaller than 0.2 mag. Bright galaxies m < 11 show a small systematic trend, where the single-Sérsic magnitude (m_ss) is brighter than the sum of the bulge and disc magnitude (m_sum). The outliers in this plot are NGC 5850, NGC 4725, IC 0724, NGC 4636, NGC 4303 and NGC 5806 with M_sum − M_ss = 0.58, 0.42, 0.32, 0.28, 0.25 and 0.25, respectively.

We compare our derived bulge and disc parameters with those measured by previous studies, specifically those that analysed small samples of galaxies and carefully fitted each galaxy individually. Initially, we compare with Fisher & Drory (2008), which presents structural parameters for 18 galaxies common to our sample, using V-band imaging from various sources. For comparison we convert our magnitudes to V-band central surface brightnesses. For the remaining parameters, r_{e, b}, r_{e, d} and n_b we compare our multiband results, which do not vary with wavelength, directly with the V-band measurements from Fisher & Drory (2008). Fig. 7 shows these comparisons. The agreement is satisfactory for most of the structural parameters. Both sets of bulge measurements display large error bars. The discrepancy between the two different studies appears to be a result of the method used to model the galaxy (1D versus 2D) and of different masking theme. Fisher & Drory decomposition is based on the major axis of each galaxy while ours utilizes all the spatial information of the image. Various studies, e.g. Ferrari et al. (2004) and Peng et al. (2010), have found that the parameters derived with 1D fitting methods and 2D are not always in agreement. Even different approaches of 1D fitting e.g. major versus minor axis fitting can change the results. The second possible reason is that we fit the entire galaxy while in Fisher & Drory (2008) they manually exclude the inner and outer part of the galaxy. This practice has the advantage that you can exclude for instance the inner part of the galaxy that may not follow the Sérsic function but has the disadvantage that the fitting process depends on personal choices (or personal experience) of what should be included in the fitting and what should be excluded. The four galaxies that show significant offset are discussed further in Appendix A.

Additionally, we compare our results with Möllenhoff (2004) for our 12 common galaxies, some of which are also in the Fisher & Drory (2008) comparison. In order to perform a sensible comparison, we first recalibrated the absolute magnitudes found in Möllenhoff (2004) to the distances used in this paper, and the AB zero-point system. Finally, we converted UBVRI magnitudes to ugriz using the transformations provided by Blanton & Roweis (2007). Our results are similar to Möllenhoff (2004), with ΔM < 0.5 mag, except for the u-band measurements of NGC 2841, NGC 3521, NGC 4274 and NGC 4303. Further information about some of these cases can be found in Appendix A.

3 STRUCTURAL PROPERTIES VERSUS MORPHOLOGICAL TYPE AND REDSHIFT

Our primary aim in this paper is to demonstrate that our multiband decomposition method is able to determine physically meaningful bulge and disc parameters, both for nearby galaxies with high-quality imaging and more distant galaxies with noisier and less well-resolved images. In this section, we therefore study the behaviour of various bulge and disc parameters for galaxies with different morphological types. For each parameter, we first present the distributions as measured on our original SDSS images, using our multiband (MM) approach. We do not show SM results for the original imaging because at such high resolution both methods return similar results (see V13 and below). We then investigate the stability of our measurements on the artificially redshifted images, by examining how the median parameters for each morphological group vary with redshift, and comparing the multiband (MM) and single-band (SM) methods.

In Figs 8 –10 and 12–16, we divide our galaxies into four morphological groups: E, S0–Sa, Sb–Sc and Sd–Irr. Where we present the results of artificial redshifting, we only plot up to a redshift of 0.15 (in contrast to 0.25 in previous figures), as beyond this redshift neither approach produces useable structural measurements.

3.1 Bulge-to-total

First, we consider the relative fluxes of the components. We remind the reader that magnitude is the only parameter in the MM fits that is entirely free to vary with wavelength. In Fig. 8, we show the distribution of r-band bulge-to-total flux ratio, using trustworthy measurements from MM fits to the original images.

Only a few elliptical galaxies have B/T close to one, contrary to expectations. Instead, most have B/T ∼ 0.7, while four have B/T < 0.6. The four cases with small B/T are NGC 4636 (B/T = 0.36), NGC 4649 (B/T = 0.49), IC 3653 (B/T = 0.55), NGC 4473 (B/T = 0.58). The S0–Sa galaxies display a broad distribution with a peak around B/T = 0.5. In the lower panel, we show intermediate and late disc galaxies. Sb–Sc galaxies show an extended distribution with about three-quarters having B/T < 0.5, and with a noticeable peak at B/T < 0.1, while all Sd-Irr galaxies have B/T < 0.3.

In Table 2, we give the average values of the B/T for each Hubble-type bin, as measured using the original images. The uncertainties on the median are estimated as |$1.253 \sigma /\sqrt{N}$|⁠, where σ is the standard deviation and N is the number of galaxies in each Hubble-type bin. We also measure median values when excluding galaxies with r_{e, b} > r_{e, d}. For these galaxies, we suspect that the bulge component fits part of the disc component leading to high B/T values. We only show this second set of results for Sb–Sc galaxies. The remaining morphological bins either contain no galaxies with r_{e, b} > r_{e, d} or the median values do not change.

In Fig. 9, we investigate our ability to recover the bulge-to-total flux ratio as our sample becomes more distant. For easier readability of the plot, we only show the r-band results for all images created with ferengi. We include only those galaxies used in Fig. 8. To facilitate the comparison with the results for the original images, we add the median value of the original images at redshift zero. Following the same colour coding of different Hubble classifications as in Fig. 8, we plot the median B/T value for each group.

Both SM and MM show similar behaviour; B/T decreases at high redshifts for the elliptical and early-disc samples, whereas the later-discs show the opposing trend. As the image resolution decreases and our ability to distinguish two components diminishes, there appears to be a tendency for the two functions to split the total flux equally. However, the MM median values show smaller fluctuations and are more stable to significantly higher redshifts than SM. For instance, the late-types appear well recovered out to z ∼ 0.09 with MM, while the SM results are only stable to z ∼ 0.05.

We can compare our average B/T values (from the original images) with those found in Weinzirl et al. (2009) and Laurikainen et al. (2010). Using a sample of 143 galaxies observed in the H band, Weinzirl et al. (2009) found that ∼69 per cent of bright spirals have B/T < 0.2 and 76 per cent of the bulges have n < 2. We find that 67 per cent of our spiral sample (Sa–Sd) have B/T < 0.2 and n < 2. Similarly, Laurikainen et al. (2010) used a sample of ∼300 S0–Sm galaxies with images in the K band found to determine B/T for each Hubble category: B/T_S0 = 0.39 ± 0.13, B/T_Sa = 0.26 ± 0.12, B/T_Sb = 0.12 ± 0.09 and B/T_Sd = 0.06 ± 0.13. Despite some differences in our sample selections and the wavelength considered, our B/T distributions for different morphologies are highly consistent with these two independent studies.

3.2 Colours

Galaxies display a range of optical colours, which correlate strongly with morphology and structure (e.g. see Willett et al. 2013). However, the total colour of a galaxy averages over any distinct stellar populations it may contain. The bulge and disc components of a galaxy may be expected to comprise contrasting stellar populations due to their different formation histories. Considering their colours individually is thus a sensible first step towards better understanding the distribution of stellar populations within galaxies. For example, quantifying the differences and correlations between the colours of bulge and disc components can help us differentiate between proposed bulge formation scenarios.

Previous work has suggested that there are substantial variations in the colours of bulges and discs between galaxies, while the colours of the two components within a given galaxy are often similar (Peletier & Balcells 1996), though significantly offset (MacArthur et al. 2004; Cameron et al. 2009). Here we briefly present the results for our initial sample, and show the advantage of the multiband technique method for measuring bulge and disc colours.

In Fig. 10, we plot the distribution of the colour difference between the two components. In the top panel, we see that early-type galaxies contain bulges and discs with similar colours. In contrast, the late-types possess bulges that are significantly redder than their discs. The average values of the component colour difference can be found in Table 2. For all spiral galaxies (Sa–Sm), we find 〈Δ(g − i)〉 ∼ 0.3 mag, irrespective of their more detailed Hubble type. Even S0s, considered alone, typically possess bulges that are slightly redder than their discs, with 〈Δ(g − i)〉 = 0.05 ± 0.05 mag.

These values compare very well with previous measurements in the literature. MacArthur et al. (2004) found an average bulge−disc colour difference of 〈Δ(B − R)〉 = 0.29 ± 0.17 mag for a sample of 172 low-inclination disc galaxies (S0–Irr), while Hudson et al. (2010) find 〈Δ(B − R)〉 = 0.23 ± 0.02 mag for L_* discs in eight low-redshift clusters. Similarly, Cameron et al. (2009) reported a colour difference of 〈Δ(u − r)〉 = 0.27 ± 0.04 mag (without their average dust correction) using ∼1500 two-component galaxies extracted from the Millennium Galaxy Catalogue. The bulge–disc colour difference we find for S0s is also consistent with the 〈Δ(g − i)〉 = 0.09 ± 0.01 mag found by Head et al. (2014) for S0s in Coma.

To examine this behaviour in more detail, Fig. 11 presents the colour–magnitude distribution for each component, colour-coded by Hubble type. Elliptical and S0 galaxies typically have both their components on a red-sequence, with (g − i) ∼ 1.2 mag, resulting in the distributions centred around zero in the top panels of Fig. 10. The disc colours of early-spirals (Sa/Sab) are also typically on this red sequence, while the discs of late-spirals (Sb–Sm) inhabit a blue cloud, with later types being fainter (though our heterogeneous sample selection may be somewhat responsible for this). The bulge colours of spirals show a considerable scatter. Some, particularly those of types Scd–Sm, lie in the blue cloud, whereas the bulges of Sab–Sc galaxies are often above than the red sequence. Dust extinction may be responsible for these very red bulges. However, we do not see any significant trend in bulge colour with disc inclination, as one might expect if this were the case.

We now consider the behaviour of the bulge−disc colour difference versus apparent redshift. In the top panel of Fig. 12, we show the MM results. The early-type galaxies have a median colour difference very close to zero, which remains almost constant out to z ∼ 0.15. For Sb–Sc galaxies, the median value is stable till redshift 0.03, after which it is overestimated with respect to the original measurement, but at least relatively stable and differentiated from the early types. Sd–Irr galaxies show a greater degree of variation beyond z ≳ 0.05, although note that this sample contains intrinsically fainter galaxies than the other sets.

In the middle panel of Fig. 12 we show the results of fitting using the aperture method, for which structural parameters are fixed to the r-band results and colours obtained by fitting for the bulge and disc fluxes in the each other band. The initial behaviour is similar to that in the top panel, but with greater variation, such that the different Hubble types are less clearly differentiated beyond z ≳ 0.06. However, we notice that MM median colours for the Sb–Sc galaxies beyond redshift 0.08 are maintained to higher values, compared to the colour at redshift zero, while the aperture median colour drops again close to the dashed line.

Finally, in the bottom panel, we show the SM results, from independent fits to each band. In this case the variations in structural parameters between bands make the colours of each component very noisy, and sensible values cannot be obtained beyond z ≳ 0.03. This emphasizes that colours for the bulges and discs of individual galaxies cannot be directly obtained via independent fits to multiple bands. Even using such measurements in a statistical fashion (e.g. to estimate the average colours of bulges) would be highly unreliable.

3.3 Sérsic index

We now move on to study how the structural parameters (n_b, r_{e, b}, r_{e, d}) are distributed for different Hubble types, and investigate the performance of the multiband fitting in measuring these values. In Fig. 13, we plot the distribution of the bulge Sérsic index, as measured in the original images using the MM method. Elliptical galaxies present a peak around 4. S0–Sa galaxies display a bimodality, with peaks around 1–2 and 4. Intermediate spiral types (Sb–Sc) have a wide range of n_b values, mostly in the range ∼1–4. The five Sb–Sc galaxies with n_b > 4 are NGC 5430, NGC 2841, NGC 3521, NGC 3642 and NGC 4698, with n_b = 4.2, 5.7, 6.7, 6.8 and 7.1, respectively. See Appendix A for further discussion of some of these galaxies. All our Sd–Irr galaxies have bulges with n_b < 2, except for NGC 4653 and NGC 4108B, with n_b = 2.7 and 4.3. The average values are given in Table 2. Typical spiral bulges with n_b ∼ 2 and a progression to lower bulge Sérsic index for later Hubble types corresponds very well to expectations from the literature (e.g. Graham & Worley 2008; Laurikainen et al. 2010; McDonald et al. 2011; Head et al. 2014).

Fig. 14 uses our artificially redshifted images to examine how well we are able to recover the bulge Sérsic index with increasing redshift. This plot is complicated by the bimodal distribution of the S0–Sa class, which results in the median being unstable. Overall, both MM and SM methods recover similar median n_b values. For the elliptical galaxies, the median n_b is well recovered at all redshifts probed. For the spiral classes, the median n_b is more variable, particularly for z ≳ 0.05. However, in general, the MM fits appear to be rather more stable.

3.4 Effective radius

In Fig. 15, we investigate the relationship between the sizes of the bulge and disc and Hubble type, by plotting histograms of the ratio of bulge effective radius to disc effective radius. Note that we do not constrain our bulges to be smaller than our discs, and neither do we subsequently exclude galaxies based on r_{e, b}/r_{e, d}. Consequently, in Fig. 15 we find seven galaxies³ with r_{e, b}/r_{e, d} > 1. Most of these galaxies have peculiarities that interfere with the fit. They are discussed individually in Appendix A.

Disregarding the few galaxies with r_{e, b}/r_{e, d} > 1, we find very little difference between the distributions for different Hubble types. Average values are listed in Table 2. Other studies also tend to find little dependence of the bulge-to-disc size ratio on morphology (e.g. Graham & Worley 2008).

We find that bulges are typically around one-quarter of the size of their accompanying discs. Rather than the ratio of effective radii, other studies typically quote r_{e, b}/h, where h = r_{e, d}/1.678 is the exponential disc scalelength. Furthermore, h is often corrected for inclination-dependent projection and extinction effects, which complicates comparisons. Finally, given the strong wavelength dependence of galaxy effective radius found by Vulcani et al. (2014), measurements at optical versus near-infrared wavelengths may be expected to differ significantly, even when an average extinction correction is applied. Nevertheless, we attempt an approximate comparison.

Assuming some average corrections, our median optical r_{e, b}/r_{e, d} ≈ 0.25 translates into an extinction-corrected, face-on r_{e, b}/h ≈ 0.35. This agrees well with the values found by Khosroshahi, Wadadekar & Kembhavi (2000), Noordermeer & van der Hulst (2007) and Méndez-Abreu et al. (2008), but is a factor of 1.5 larger than that found by the careful analysis of multiple data sets from the literature (Graham & Worley 2008) and twice that found by Laurikainen et al. (2010). The latter study, and some of the works that were included in Graham & Worley (2008), included additional central components in their models, such as bars or nuclei. This may have led to the smaller bulge sizes they measure. Given the care taken in these studies, we suspect that our bulge r_{e, b} may be somewhat overestimated. However, remember that our aim in this work is to perform fits to our nearby galaxies in a simple, automated manner, suitable for large surveys of relatively distant galaxies, and ascertain the performance of this approach.

With this in mind, Fig. 16 shows the median r_{e, b}/r_{e, d} for several Hubble type bins as a function of simulated redshift. For the multiband (MM) fits, we again see that for low-redshift data there is little difference with morphology. The average size ratios remain fairly constant to z ∼ 0.04. After this, as the data quality becomes substantially poorer, a strong trend to increasing r_{e, b}/r_{e, d} is seen, particularly affecting galaxies with lower B/T. Single-band fits perform reasonably similarly (neglecting the ellipticals, for which the reality of the disc is unclear). However, they show a stronger and noisier bias, which sets in at slightly lower redshifts.

Generally, we observe that the lower the data quality, the harder it is to separate the two components and the more similar their properties become. However, it is usually the bulge fit which is most affected, and hence biased. The effective radius of the disc components tend to remain stable for almost the entire redshift range considered, particularly for our multiband fits.

4 CORRELATION OF STRUCTURE PARAMETERS

4.1 Bulge-to-total

For disc galaxies, the overall Sérsic index (of a single-Sérsic model fit) is often regarded as an indication of the bulge-to-total ratio (B/T). Indeed, we adopt this interpretation in Häußler et al. (2013), V13 and Vulcani et al. (2014). With our bulge–disc decompositions, we are now in a position to test this.

In Fig. 17, we plot the bulge-to-total ratio as a function of the single-Sérsic index (n_SS) measured in V13. We see that, as the overall Sérsic index increases, the bulge is responsible for a greater proportion of the galaxy flux, confirming our expectations. Galaxies with a low overall Sérsic index typically contain two components, a bulge and a disc, and the more dominant the bulge component, the higher the overall Sérsic index.

Galaxies of type Scd and later generally have low B/T, while earlier spirals (Sa–Sc) span a wide range of B/T. Interestingly, earlier types tend to have greater n_SS at a given B/T, suggesting that n_SS is also dependent on other aspects of galaxy structure. Most galaxies with n_SS ∼ 1 have B/T < 0.1, with the exception of four galaxies.⁴

For the original images, the resolution and S/N is sufficiently good that fitting bands individually is comparable to our multiband approach. However, for more distant galaxies our multiband method gives more robust measurements of B/T. We demonstrate this in Fig. 18, where we again plot the bulge-to-total ratio as a function of the single-component Sérsic index (n_SS), but now using galaxies artificially redshifted to z = 0.1. This figure shows the results of fitting each band independently (SM), in addition to our multiband measurements (MM). The clearer correlation for MM clearly illustrates the advantage of our multiband method.

A number of studies have presented a correlation between bulge-to-total ratio and the Sérsic index of the bulge (e.g. Graham 2001), particularly as a diagnostic for distinguishing so-called pseudo- and classical-bulges (Durbala et al. 2008; Weinzirl et al. 2009; Laurikainen et al. 2010). Fig. 19 confirms this relationship for our measurements. The scatter is relatively large, probably as a result of the difficulty in constraining the bulge properties, as discussed in the previous section. Nevertheless, it is clear that the more bulge-dominated a galaxy is, the higher its bulge Sérsic index.

We also see a weak correlation between the ratio of bulge and disc sizes, r_{e, b}/r_{e, d}, and the bulge-to-total flux ratio, in Fig. 20. There are some indications that the relation depends on morphological type, but the scatter and incompleteness of our sample prevent us from making definitive conclusions.

4.2 Component axis ratios and the division of ellipticals and lenticulars

Fig. 21 explores the relationship between the axis ratio of the bulge (top panel) and disc (bottom panel) versus the bulge Sérsic index. As before, we include galaxies with elliptical morphologies in these plots for two reasons. First, recent work has blurred the lines between elliptical and lenticular galaxies, with many ellipticals found to contain faint disc components when studied in detail (Krajnović et al. 2013). Secondly, our aim is to inform work on large surveys, which may not have morphological classifications available. For these ellipticals, and despite our nomenclature, we do not go so far as to assume that the exponential component of our model represents a disc, but consider it to be an indication of additional structure that cannot be well modelled by a single-Sérsic component.

Considering all the points in the top panel of Fig. 21, there is a clear correlation such that bulges with higher Sérsic index tend to be more circular (b/a ∼ 1). The vast majority of galaxies with elliptical morphologies are found in the upper-right region, with n_b > 2 and b/a > 0.5, as might be expected. Lenticulars slightly separate out from ellipticals in this plot, generally being limited to slightly lower b/a and a wider range of n_b > 2.

Moving our attention to the lower panel, we first see little correlation between disc axial ratio and bulge Sérsic index. Note that our sample of spiral galaxies is seriously incomplete for inclined systems, due to selection restrictions applied by the studies from which V13 obtained their sample.

Focusing on ellipticals and S0s, we see a surprisingly strong separation between the two morphologies in disc b/a. The vast majority of galaxies with classified as elliptical have (b/a)_disc > 0.5, while the lenticulars have mostly (b/a)_disc < 0.5. We also see an offset of the inclined lenticulars to higher bulge Sérsic indices. We are not certain whether this reflects reality, or is a bias in the decomposition process. Simulations suggest that only small variations in n_b are expected from decomposition effects (Pastrav et al. 2013). S0s are generally not expected to contain significant amounts of dust, so extinction should not play a significant role. In any case, the effects are typically ≲ 0.1 in n_b.

In order to explore the separation of ellipticals and lenticulars in Fig. 21 further, we highlight early-type galaxies using their kinematic classification from Emsellem et al. (2011). We note that almost all the early-type galaxies with low (b/a)_disc are fast rotators, including both of the elliptical galaxies which fall in this region of parameter space dominated by S0s. The early types with (b/a)_disc > 0.5 are a mixture of fast and slow rotators; however, most (or all with n_b > 2) have been classified as elliptical galaxies (see Cappellari et al. 2011 and Krajnović et al. 2013 for a more thorough study of this topic).

Obviously a plot of (b/a)_disc versus n_b, or even just (b/a)_disc alone, is an effective, automated way of separating galaxies with elliptical and lenticular visual morphologies. However, this raises the question of whether such a separation is a physically meaningful thing to do.

The difficulty of distinguishing face-on S0s from ellipticals is a well-known problem. The result in Fig. 21 illustrates this issue in terms of quantitative structural measurements. When a fast-rotating early-type galaxy appears to have an inclined disc ((b/a)_disc ≲ 0.5, it is usually classified as S0. If the same galaxy were viewed closer to face-on, it would be classified as elliptical. In our (non-representative) sample, this amounts to about half of true S0 galaxies (discy, fast rotators that that would have been visually classified as S0 if viewable from other angles) being misclassified as elliptical. We presently do not have a reliable morphological or structural way of recovering these objects, but instead must resort to kinematic information (e.g. Krajnović et al. 2013). However, we remain hopeful that with additional work we can make further progress on an image-based solution to this long-standing problem.

5 CONCLUSIONS

All previous studies that have utilized 1D or 2D photometric bulge–disc decompositions have performed their fits using, at most, two bands simultaneously. In this work, for the first time, we have performed bulge–disc decompositions simultaneously on five-band imaging. To achieve this, we have used galfitm, a modified version of galfit3 which enables a single, wavelength-dependent model to be fit to multiple images of the same galaxy.

We have evaluated the performance of our multiband method by applying it to SDSS ugriz images of 163 nearby galaxies, as well as to another 3863 artificially redshifted images of the same objects. For our models, we use a combination of two Sérsic profiles: one for the disc, with n fixed to 1, and another for the bulge, with free-n.

Using the original images, we have shown that our fitting results generally agree well with structural parameters obtained from the literature, both when we compare specific galaxies and averages for bins of Hubble type. We confirm the standard picture that disc galaxies of earlier morphology have a larger fraction of their flux in a bulge component. However, the sizes of these bulges, with respect to their accompanying disc, do not vary significantly with Hubble type. We find that bulges in spiral galaxies typically have a Sérsic index n ∼ 2, except in the latest types, where it is more usually exponential n ∼ 1. The dominant component of ellipticals has n ∼ 3.5. Puzzlingly early-type discs appear to display a bimodality with peaks around n ∼ 1 and 4.

The observed colours of discs display a classic colour–magnitude diagram, with a well-defined red sequence, inhabited by E–Sb galaxies, and a blue cloud corresponding to later Hubble types. The colour–magnitude diagram for bulges is more complex. The bulges of E–Sa galaxies lie on a red sequence similar to their discs, with only a small average difference in the colours of their bulges and discs. The bulges of Sb–Sc galaxies are often even redder than the early-type red sequence, indicative of dust reddening. On the other hand, for many late-type disc galaxies we find bulges with fairly blue colours, suggestive of recent star formation. Despite this complexity, the average difference in the colours of bulges and discs within the same galaxy is constant for all spiral galaxies, 〈Δ(g − i)〉 ∼ 0.3 mag.

Our fits permit a quantitative illustration of the notorious difficulty of distinguishing between galaxies with elliptical and lenticular morphology. Most early types are well fitted by a combination of a Sérsic profile and an exponential profile. It is not clear how often this exponential profile represents a genuine disc component, or whether it reveals the presence of an extended halo or some other substructure. Nevertheless, a significant inclined disc component is strong indication that an early-type galaxy will be visually classified as S0 and possess fast-rotator kinematics. Unfortunately, the lack of S0 morphologies among galaxies with face-on disc components, despite kinematic evidence indicating the discs are real, suggests that such systems are usually misclassified as ellipticals.

Using our artificially redshifted images we have investigated the range of redshift over which our fit parameters remain stable, and hence are inferred to be reliable. We have demonstrated that the results of our multiband fits show less variation and are reliable to higher redshifts than the results of fitting each band independently. This is true even in the common single-band approach where one first performs a full fit on a preferred band, then fixes the structural parameters in all subsequent fits to those obtained in that preferred band.

Our approach produces somewhat more reliable measurements of the bulge Sérsic index and effective radius, although these are both difficult quantities to measure with consistent accuracy, especially when it is not feasible to fit each galaxy interactively. More promising is our method's performance in recovering the bulge-to-total flux ratio and in differentiating between the colours of the bulge and disc. We therefore recommend the adoption of multiband bulge–disc decomposition, allowing studies to reliably probe to greater distances and lower luminosity galaxies.

This publication was made possible by NPRP grant # 08-643-1-112 from the Qatar National Research Fund (a member of Qatar Foundation). The statements made herein are solely the responsibility of the authors. BH and MV were supported by the NPRP grant. SPB gratefully acknowledges an STFC Advanced Fellowship. We would like to thank Carnegie Mellon University in Qatar and The University of Nottingham for their hospitality. We would like to thank the referee for the constructive comments that improved the paper.

Funding for the SDSS and SDSS-II has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, the US Department of Energy, the National Aeronautics and Space Administration, the Japanese Monbukagakusho, the Max Planck Society and the Higher Education Funding Council for England. The SDSS website is http://www.sdss.org/. The SDSS is managed by the Astrophysical Research Consortium for the Participating Institutions. The Participating Institutions are the American Museum of Natural History, Astrophysical Institute Potsdam, University of Basel, University of Cambridge, Case Western Reserve University, University of Chicago, Drexel University, Fermilab, the Institute for Advanced Study, the Japan Participation Group, Johns Hopkins University, the Joint Institute for Nuclear Astrophysics, the Kavli Institute for Particle Astrophysics and Cosmology, the Korean Scientist Group, the Chinese Academy of Sciences (LAMOST), Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, Ohio State University, University of Pittsburgh, University of Portsmouth, Princeton University, the United States Naval Observatory and the University of Washington.

REFERENCES

APPENDIX A: NOTES ON INDIVIDUAL GALAXIES

Here, we briefly describe some galaxies with peculiar results. In Fig. A1, we present the images and models for some of the following galaxies in order to show the visual appearance of the galaxies that return peculiar results. In the same figure, we have also added three additional galaxies with what we consider good result in order to support a comparison.

1. NGC 3521[SAB(rs)bc] – a high-inclination galaxy with strong dust features that make it extremely difficult to fit a two-component model. Our two-component model measures a large bulge with r_{e, b}/r_{e, d} = 1.94 and n_b = 6.7. This may not be realistic, but it has been retained in the analysis.

2. NGC 3642[SA(r)bc] – a late-type galaxy for which we measure r_{e, b}/r_{e, d} = 2.5 and n_b = 6.7. However, by examining the equivalent of Fig. 3 for NGC 3642, we noticed that that both r_{e, b}/r_{e, d} and n_b decrease with redshift until z = 0.03. After that, the values remain constant until z = 0.12. The recovered values for this range of redshifts are r_{e, b}/r_{e, d} = 0.2 and n_b = 1.5, which seem much more reasonable. These dramatic changes in parameters could be due to the detailed substructure that is visible in the high-resolution images. As the galaxy becomes more distant these substructure are less pronounced and, as a result, the fit parameters better reflect the properties of the bulge and disc.

3. NGC 3893[SAB(rs)c] – a similar case to NGC 3642.

4. NGC 4123[SB(r)c] – a late-type galaxy for which our two-component model measures n_b = 10.5. However, if we add a PSF into the model, the bulge Sérsic index is reduced to 1.86 and r_b decreases by 50 per cent. The corresponding change in disc effective radius is ∼2 per cent. We therefore use the model including the PSF component.

5. NGC 4452[S0(9)] – an almost perfectly edge-on galaxy, which contains a very thin disc (for which our model measures b/a = 0.09). The bulge component fits a second elongated component (b/a = 0.37) with n_b = 1.08 and effective radius 20 per cent larger than the thin disc. The properties of the second component are more consistent with a thick disc, rather than a bulge. However, in the spirit of avoiding special cases, we retain both components in the analysis.

6. NGC 4458[E] – an elliptical galaxy with an internal structure. As a result, the exponential component fits the inner part of the galaxy and the free-Sérsic component measures n_b = 11. Due to the peculiar fitting results, we attempted to fit this galaxy with various initial parameter values. However, the final results always remained the same, and agree for both SM and MM methods.

7. NGC 4459[S0] – after various attempts with different initial values, we failed to fit a significant second component for this lenticular galaxy. In analysis, we therefore only include the bulge component.

8. NGC 4698[SA(s)ab] – a spiral galaxy with peculiar structure: the bulge is elongated perpendicular to the main disc. We measure structural parameters of r_{e, b}/r_{e, d} = 1.5 and n_b = 7.1, consistently using both the SM and MM methods, and the redshifted images. However, these are not in agreement with other studies which focus on the unique structure of this galaxy.

9. UGC08041[SB(s)d] is a similar case to NGC 4123. The two-component model measures n_b = 9.37, while the inclusion of a PSF component reduces this to n_b = 1.17. Both the effective radius of the bulge and disc component change less than 10 per cent. We use the model with the PSF component.