© The Authors 2024

1. Introduction

The luminosity and stellar mass functions of galaxies (LMFs) establish the comoving number density of galaxies per luminosity and stellar mass interval, respectively. Luminosity functions are open to be defined according to any photometric band or spectral range, mostly rest-frame luminosities and absolute magnitudes; whereas stellar mass functions are only attached to this stellar population parameter, which is subject to data analysis and scaled according to a stellar population synthesis model set. Owing to the large amount of information contained in the LMFs, these are largely employed for a wide variety of purposes in multiple researches.

Amongst their many applications, the LMFs have been extensively used for setting constraints on cosmological models of galaxy evolution and formation (see e.g. Benson et al. 2003; Lacey et al. 2016) and supporting the development and creation of mock catalogues with realistic galaxy clustering (e.g. Somerville & Davé 2015). In brief, these functions help to test, fine-tune, and calibrate schemes for populating halos from dark-matter-only simulations with galaxies, such as semi-analytic models, the subhalo abundance matching, the halo occupation distribution schemes, and others (Yang et al. 2003; Vale & Ostriker 2004; Zheng et al. 2005; Baugh 2006; Conroy et al. 2006; Benson 2010; Guo et al. 2011; Smith et al. 2017). In the same sense, results from hydrodynamical simulations are often cross-checked with the local stellar mass function of galaxies (Genel et al. 2014; Schaye et al. 2015; Cui et al. 2018). Furthermore, the LMFs can provide a preliminary estimation of the expected number of galaxies included in an imaging survey, as well as for defining their strategy plans. This was also used for example to improve the precision of photometric redshifts (hereafter photo-z) and reduce the number of photo-z outliers at fainter magnitudes by the inclusion of priors (see e.g. Benítez 2000; Ilbert et al. 2006; Brammer et al. 2008; Arnouts & Ilbert 2011; Molino et al. 2014; Hernán-Caballero et al. 2021). This is thanks to the ability of luminosity functions to accounting for the frequency of templates embedded in photo-z codes (which are related to the different spectral-types of galaxies) and that the volume sampled at lower redshifts is smaller (see e.g. Benítez 2000; Brammer et al. 2008).

As the LMFs and their evolution are tightly related to the star formation processes taking place in galaxies (see e.g. Wilkins et al. 2008; Behroozi et al. 2013; Ilbert et al. 2013; Bernardi et al. 2016), these functions have been mainly discussed and included in galaxy evolution studies. In particular, the cosmic evolution of the LMFs or the evolution of the comoving number densities of galaxies with redshift have been a powerful tool to provide hints on the evolutive paths of certain type of galaxies (see e.g. Ferreras et al. 2009a,b; Ilbert et al. 2010, 2013; Brammer et al. 2011; Moustakas et al. 2013; Díaz-García et al. 2019b). In this regard, processes such as quenching (i.e. cessation of star formation Faber et al. 2007; Peng et al. 2015) must be also imprinted in the LMFs. Therefore, these functions have been employed to understand the physical processes leading to the bimodal distribution of galaxies in colour-magnitude, stellar mass-colour, and UVJ-like diagrams and their evolution with redshift (see e.g. Faber et al. 2007; Ilbert et al. 2010, 2013; Moustakas et al. 2013; Schawinski et al. 2014; Díaz-García et al. 2019a). In this regard (and based on the LMFs), some authors figured out that massive quiescent galaxies (log₁₀M_⋆ ≳ 10.7) show a rapid and active increase in number since z ∼ 3 and efficient up to z ∼ 1, that is referred as ‘mass quenching’ and independent of environment (see e.g. Peng et al. 2010; Ilbert et al. 2013). On the other hand, multiple works found that the variation in number density of the LMFs of this type of galaxy at z < 1 is mainly driven by less massive quiescent galaxies (Pozzetti et al. 2010; Ilbert et al. 2010, 2013; Maraston et al. 2013; Moresco et al. 2013; Moustakas et al. 2013; Díaz-García et al. 2019b). Also based on these diagrams, the number density of the so-called ‘green valley’ galaxies can be also used for setting constraints on the typical time scales of quenching processes, where some works apparently point to a fast transition (see e.g. Bell et al. 2004; Muzzin et al. 2013; Schawinski et al. 2014; Angthopo et al. 2019; Díaz-García et al. 2019b). In a series of complementary follow-ups, a number of authors have investigated whether LMFs depend on the environment or the density of galaxies within a projected radius. Overall, it has been widely accepted that the LMFs of galaxies vary according to the local density or environment in which they reside; namely, that more luminous and massive galaxies, especially red luminous galaxies, typically reside in denser environments (Eke et al. 2004; Baldry et al. 2006; Peng et al. 2010; Zandivarez & Martínez 2011; Guo et al. 2014; Nantais et al. 2016; Vázquez-Mata et al. 2020). Interestingly, there is evidence of non-negligible excess of low mass quiescent galaxies and spheroids in high density environments (Peng et al. 2010; Scoville et al. 2013; Moutard et al. 2018) with respect to the field. In fact, this result led to propose a second quenching mechanism dubbed ‘environment quenching’ acting in low mass systems at z < 1, which would result from ram pressure stripping (Gunn et al. 1972; Abadi et al. 1999; Quilis et al. 2000) or strangulation (Larson et al. 1980; Balogh et al. 2000) and may explain the low-mass upturn in the stellar mass function of quiescent galaxies (Drory et al. 2009; Taylor et al. 2015). These facts are also reflected in the fraction of quiescent galaxies in galaxy groups and clusters, where this fraction is higher for denser environments (Woo et al. 2013; Nantais et al. 2016; Liu et al. 2021; González Delgado et al. 2022; Rodríguez-Martín et al. 2022; Sobral et al. 2022).

For tackling the determination of LMFs, a wide variety of methodologies have been proposed and developed during the last decades, which can be used for constraining both functions due to their similarities. Some of them do not assume any functional form for the LMFs and are referred as non-parametric, hence the non-parametric LMF are discretised in stellar mass and/or absolute magnitude bins. One of the most extended non-parametric methods for the determination of the LMFs is the so-called 1/V_max method (Schmidt 1968; Marshall 1985). Alternatively, more non-parametric methods have been proposed, such as the C⁻ (Lynden-Bell 1971), C⁺ (Zucca et al. 1997), and SWML (Efstathiou et al. 1988) methods for the determination of discretised LMFs. On the other hand, when observations are deep enough and the sample is complete in a wide range of stellar mass or luminosity, many authors prefers to develop methods that assume a parametric LMF based on the empirical Schechter function (Schechter 1976). In brief, a Schechter function adopts a power-law function for describing the number density of faint galaxies, whereas the number density of bright galaxies decays exponentially. Thus, a Schechter function is defined according to three parameters: the power or slope of the faint end, the characteristic stellar mass or luminosity, and the normalisation in number density units; parameters that, in general, are correlated or degenerated (Beare et al. 2015; López-Sanjuan et al. 2017). Over the last decades, many works have shown that the LMFs of galaxies are properly described by parametric Schechter-like functions from a very high redshift (see e.g. Pérez-González et al. 2008; Ilbert et al. 2013; Muzzin et al. 2013; Moustakas et al. 2013; Wright et al. 2018; Ishigaki et al. 2018; Bhatawdekar et al. 2019; Bowler et al. 2020; Donnan et al. 2023; Pérez-González et al. 2023; Weaver et al. 2023). Furthermore, recent works have also shown that for a general case, sometimes the use of a double-Schechter function (i.e. a combination of two Schechter functions) is more adequate than the single case, which is commonly interpreted as some galaxy populations are assembled according to different galaxy formation scenarios (see e.g. Drory et al. 2009; Peng et al. 2010; Ilbert et al. 2013). In this sense, the LMFs of quiescent galaxies are often described or fitted by a double Schechter function, while the star-forming ones are parametrised through making use of a single Schechter function (e.g. Li & White 2009; Ilbert et al. 2010; Baldry et al. 2012; Muzzin et al. 2013; Weigel et al. 2016; López-Sanjuan et al. 2017; Ishigaki et al. 2018).

Recent large scale and/or deep multi-filter surveys including or combining broad, intermediate, and narrow-band filters (e.g. COMBO-17¹, COSMOS², ALHAMBRA³, PAU⁴, SHARDS⁵, J-PAS⁶, J-PLUS⁷, and S-PLUS⁸; Wolf et al. 2003; Scoville et al. 2007; Moles et al. 2008; Benítez et al. 2009a, 2014; Pérez-González et al. 2013; Cenarro et al. 2019; Mendes de Oliveira et al. 2019) have proven to be very adequate for the determination of precise galaxy photo-z (σ_z ∼ 0.01; Wolf et al. 2004; Ilbert et al. 2009; Martí et al. 2014; Molino et al. 2014, 2020; Hernán-Caballero et al. 2021). This fact, together with the advantages that this kind of survey offers, makes them very attractive for the study of the LMFs of galaxies since intermediate redshift. Some of these advantages include: (i) there is no pre-selection of the sample since all the galaxies imaged down to the magnitude limit of the detection band are included in the analysis; (ii) images allow for the use of non-fixed aperture photometry for estimating the total flux of each surveyed galaxy, meaning that there is overall no aperture bias in the determination of the total luminosity or stellar mass; and (iii) these surveys can easily map large areas of the sky in a wide redshift range that yield unbiased and statically robust samples of galaxies, meaning that systematic effects as those coming from Poisson errors and cosmic variance may be largely reduced or even neglected in the best cases.

In this paper, we aim at developing a robust methodology for determining the so-called luminosity and stellar mass functions of galaxies, which must be properly adapted to the particularities of large scale multi-filter surveys such as J-PAS, putting an emphasis on the proper determination of these functions up to z ∼ 0.7. Whenever possible, our method must include the outputs or value-added catalogues available from the J-PAS Collaboration. In particular, we specially require for the analysis the inclusion of the typical outputs from our trusted fitting codes for spectral energy distribution (SED) in the J-PAS Collaboration, such as the ones described in González Delgado et al. (2021), in order to perform a solid and statistical analysis accounting for the uncertainties, degeneracies, and correlations amongst all the involved parameters. The ultimate goal consists in checking that all our techniques, J-PAS inputs, and results are self-consistent as to provide proper LMF results that are in agreement with previous studies for a flux-limited sample up to r = 22.0 AB-magnitudes and z = 0.7. As the development of our methods for determining the LMFs is based on a preliminary J-PAS-like data set, dubbed miniJPAS (Bonoli et al. 2021, with further details in our Sect. 2), our conclusions and results are focused on introducing the basics of our techniques, as well as showing the potential of J-PAS data for LMF studies. The latter is motivated by the fact that this survey will scan thousands of square degrees of the northern hemisphere with a unique photometric system including 54 narrow band filters in the optical range. As a consequence of our analysis, we also provide a statistical spectral-type classification of galaxies (quiescent versus star-forming), as well as the stellar mass and luminosity constraints for each of the galaxies that will be made publicly available in a miniJPAS value-added catalogue.

This paper is organised as follows. In Sect. 2, we introduce the miniJPAS survey and the value-added catalogues employed throughout this research. The SED-fitting techniques and methodologies developed to constrain the stellar population properties that are needed to determine the miniJPAS LMFs are summarised in Sect. 3. The process to perform the classification of galaxies according to their spectral type, as well as the luminosity and stellar mass completeness of our sample are detailed in Sects. 4 and 5, respectively. The novel methodology to statistically determine the characteristic parameters defining the LMFs is detailed in Sect. 6. In Sect. 7, we present our main results, which include the miniJPAS LMFs of galaxies up to z = 0.7 and the cosmic evolution of the stellar mass and luminosity densities in the same redshift range. Finally, we discuss our results and conclusions in Sect. 8, while a brief summary of this work is included in Sect. 9.

Throughout this paper we adopt a spatially flat lambda cold dark matter (ΛCDM) cosmology with parameters based on the recent Planck results (Planck Collaboration VI 2020), namely: Hubble constant of H₀ = 67.4 km s⁻¹ Mpc⁻¹ (i.e. h = 0.674), Ω_M = 0.315 (total matter density), and Ω_Λ = 0.685 (cosmological constant density). Stellar masses are quoted in solar mass units [M_⊙] and are scaled according to a universal Chabrier (2003) initial stellar mass function (IMF). Magnitudes are expressed in the AB system (Oke & Gunn 1983).

2. The miniJPAS survey

The miniJPAS survey (Bonoli et al. 2021) was conceived to primarily test the performance of the Javalambre Survey Telescope (JST250) optical system and determine the potential and applications of the J-PAS spectro-photometric data. Thus, the miniJPAS strategy and configuration (Sect. 2.1) was dictated according to the baseline survey strategy of J-PAS (further details in Benitez et al. 2014; Bonoli et al. 2021). For the aims of this work, we build our sample of sources from the miniJPAS photometric catalogue (Sect. 2.2) and we take advantage of the miniJPAS value-added catalogues containing constraints on the galaxy photo-z (Sect. 2.3), and the star/galaxy classification (Sect. 2.4) of each miniJPAS source.

2.1. Observations and instrumentation

The miniJPAS survey was conducted at the Javalambre Astronomical Observatory⁹ (OAJ), a site with an excellent median seeing and darkness that is dedicated to large scale multi-filter surveys to image the northern hemisphere (Moles et al. 2010). In particular, all the miniJPAS observations were done with the JPAS-Pathfinder camera, namely, the first scientific instrument installed at the JST250 telescope. This camera includes a single CCD of 9.2k × 9.2k pixels, with an effective field of view (FoV) of 0.27 deg², and a pixel scale of 0.23″ pix⁻¹. In fact, the greatest difference between miniJPAS and J-PAS resides in the camera employed for the observations: JPCam, a 1.2 Gpixel detector composed of 14 CCDs with a FoV of 4.2 deg² and the same pixel scale than JPAS-Pathfinder. The JST250 telescope is a wide-field telescope of 2.55 m based at the OAJ with an effective collecting area of 3.75 m² and an optimised éttendue of 26.5 m² deg² as to perform large FoV surveys (Cenarro et al. 2018).

One of the main characteristics of the miniJPAS survey is that it was performed using the J-PAS filter system. This unprecedented photometric system comprises 54 narrow-band filters with a full width at half maximum of FWHM ∼ 145 Å (equally spaced every 100 Å), 1 broad-band (u_JAVA, FWHM ∼ 495 Å), and 1 high-pass filter (J1007) extending to the UV and near-infrared ends of the optical range, which yields an effective wavelength range of 3500–9300 Å (further details given in Bonoli et al. 2021). The narrow band filters are properly characterised by transmission curves with very steep side slopes and flat tops, which are usually referred as top-hat filters. As a result, the photometric system delivers spectral energy distributions (SEDs) or J-spectra with an equivalent wavelength resolution of R ∼ 60 for each of the imaged sources. For the miniJPAS case, the filter set was complemented with the ‘standard’ broad-band filters u_JPAS, g_SDSS, r_SDSS, and i_SDSS. Only u_JPAS differs a little from the SDSS-like u band, since it has a redder cut-off.

The miniJPAS field lies on the well-known AEGIS¹⁰ field, centred at (RA, Dec) = (215°, +53°), and is composed of a stripe of four tiles or pointings in such a way that covers the EGS¹¹ field. There is a little overlapping by 3.6′ between contiguous pointings (0.09 deg² of total overlapping area) and the total area observed amounts to ∼1 deg² of the sky. However, after masking bright stars, the window frame, and artefacts (typically from internal light reflections in the system), the effective area of the miniJPAS surveys amounts to 0.895 deg². The average image depths, defined for a 5σ level and a circular aperture of 3″, of the narrow-band filter set range from 23.5 for the bluer bands to 22 AB-magnitudes for the redder parts of the photometric system. Regarding the broad-band filters, the image depths range from 23 to 24 AB-magnitudes. Owing to the observational campaign, the FWHM of the point spread function (PSF) ranges from 0.6″ to 2″, with most of the bands showing FWHMs below 1.5″ (see Fig. 5 in Bonoli et al. 2021).

2.2. Photometry and definition of our miniJPAS sample

The source detection and built of photometric catalogues is based on the well-known Source-Extractor programme (SExtractor; Bertin & Arnouts 1996). As miniJPAS is a survey involving several astronomical topics, catalogues include magnitudes and fluxes for different apertures and magnitude definitions, such as isophotal magnitudes (ISO), PSF-corrected photometry (PSFCOR, based on the method by Molino et al. 2019), as well as Petrosian (PETRO; Petrosian 1976), automatic (AUTO, inspired by Kron’s ‘first moment’ algorithm; Kron 1980), and fixed circular apertures. Except for PSFCOR, all of them were computed in both single- and dual-mode. For the latter, the detection and aperture definition of sources were done with the r_SDSS-band images to subsequently perform forced photometry in the rest of band images.

For our aims, the choice of aperture photometry to be used primarily relies on the compromise of including or approaching the total flux of galaxies with the highest possible signal-to-noise ratio (S/N). With this in mind, the use of dual-mode photometry within the AUTO apertures from SExtractor is specially motivated. In fact, these are intended to estimate the total flux of galaxies where the mean fraction of flux lost is around 6% (Bertin & Arnouts 1996; Graham & Driver 2005) at the same time that circumvents PSF-homogenisation. The ISO, PSFCOR, and PETRO magnitudes were mainly discarded for the LMF analysis addressed in this work because (i) ISO photometry typically loses flux from the ‘wings’ of galaxy light profiles (Bertin & Arnouts 1996; ii) PSFCOR commonly yield lower fluxes by 0.5 magnitudes owing to smaller apertures than AUTO (González Delgado et al. 2021; Hernán-Caballero et al. 2021), meaning that the total luminosity and galaxy stellar mass is underestimated or biased; (iii) PSFCOR colours are typically redder than for AUTO photometry due to the aforementioned smaller aperture, which (in turn) affects the stellar mass determination (see Sects. 3 and 4); and (iv) PETRO uncertainties are on average ∼35% larger (r_SDSS-band magnitude) as a consequence of larger apertures than AUTO, which becomes particularly troublesome when determining the luminosity and stellar mass of fainter miniJPAS galaxies and/or at high redshift. From now on, unless otherwise stated, all the magnitudes in this work refer to AUTO magnitudes.

The parent catalogue for this work is the miniJPAS Public Data Release (PDR201912)¹². As stated in Bonoli et al. (2021), the primary miniJPAS photometric catalogue (table minijpas.MagABDualObj) is complete in terms of detection at r_SDSS = 23.6 and r_SDSS = 22.7 for point-like and extended sources, respectively. However, we define our flux-limited initial sample at r_SDSS ≤ 22.5 for avoiding noisy sources. Moreover, we add extra constraints to guarantee aperture photometry quality: (i) sources must be detected within the image window frame and outside of masked regions (bright stars or artefacts); and (ii) photometry in the detection band (r_SDSS) can only include SExtractor flags equal to FLAGS = 0 and 2, which removes saturated and truncated sources among others. We note that at this point, we have not imposed constraints on point-like sources or redshift, so these issues are statistically addressed during the subsequent LMF methodology and posterior analysis (Sect. 6). As a result, the flux-limited parent catalogue is composed of N_S ∼ 12 600 sources. As we detail in Sect. 3, we imposed extra conditions to perform a proper SED-fitting analysis of this parent catalogue, which ends up slightly reduce the number of sources. Even though this is the parent sample that we use for the present statistical research, it is worth remarking at this point that we based the determination of the miniJPAS LMFs on the cumulative probability of the sample at 0.05 ≤ z ≤ 0.7 and r_SDSS ≤ 22 (largely detailed in Sect. 6). This also guarantees that almost half of the miniJPAS bands (including the detection band) trace the rest-frame bulk of the stellar emission (λ > 4000 Å), thus avoiding any underestimations of the stellar masses (see Davidzon et al. 2017; Weaver et al. 2022, 2023).

2.3. Photometric redshifts

The J-PAS photometric system was specifically designed to measure baryon acoustic oscillations (BAOs; Benítez et al. 2009b; Benitez et al. 2014) along the line of sight, which theoretically requires a redshift precision of σ_z ∼ 0.003 × (1 + z) (Benítez et al. 2009a). To achieve the photo-z precision for this aim, the photo-z group within the J-PAS Collaboration worked on different photo-z codes. For this research, we preferentially focus on the JPHOTOZ package (Hernán-Caballero et al. 2021), which was developed at Centro de Estudios de Física del Cosmos de Aragón¹³ (CEFCA) as part of the data reduction pipeline for J-PAS, also known as JYPE. It is worth mentioning that there are alternative photo-z codes available, which were tested in the collaboration such as Tartu Observatory photo-z (TOPz; details in Laur et al. 2022), which also satisfies the photo-z precision requirements via SED-fitting techniques.

In brief, JPHOTOZ is a set of Python scripts that interface between the database and a customised version of the SED-fitting code LEPHARE (Arnouts & Ilbert 2011), which is preferentially adapted to work with J-PAS-like data. The set of templates used for estimating the J-PAS photo-z is composed of 50 synthetic templates, which were generated by using CIGALE¹⁴ (Boquien et al. 2019). This set of templates was carefully defined through an iterative method based on scores evaluating the combination of templates that optimised the photo-z of miniJPAS galaxies after comparison with spectroscopic redshifts (for further details, see Hernán-Caballero et al. 2021). It is also of note that for a better photo-z determination, a recalibration offset was applied to the input PSFCOR magnitudes, which was recomputed, pointing per pointing, by an iterative process. In addition, the miniJPAS photo-z were determined using the LEPHARE redshift priors from spectroscopic redshifts from the VIMOS VLT survey (VVDS, Le Fèvre et al. 2005), meaning that a prior was applied according to the i-magnitude of the galaxy and the rest-frame colour (g − i) of each template.

To perform a robust statistical analysis, we must have the photo-z probability distribution functions (PDZ in the following), along with its associated ‘odds’ value, for each of the miniJPAS sources. This calls for making use of the PDZ of miniJPAS galaxies, rather than the photo-z with the highest probability (maximum of the PDZ or z_best) for each miniJPAS source (see also Sect. 3). The odds value is the probability enclosed by the PDZ in an interval around its maximum value (Benítez 2000; Molino et al. 2014), which is defined according to z_best ± 0.03 × (1 + z_best) for the miniJPAS case. Certainly, this parameter has proven to be correlated with the quality of each photo-z estimation (Hernán-Caballero et al. 2021, 2023), as well as to be sensitive to catastrophic outliers. Consequently, odds is included in our analysis for a better determination of the LMFs (see Sect. 6). After a comparison with a miniJPAS representative sample of spectroscopic redshifts at r_SDSS < 23, the typical error for the photo-z of miniJPAS galaxies is set at σ_NMAD = 0.013 with an outlier rate of η = 0.39 (further details in Hernán-Caballero et al. 2021). For odds > 0.82 and r_SDSS < 23, the typical error decreases up to σ_NMAD = 0.003 and η = 0.05, which results in ∼5200 galaxies per deg⁻² with this photo-z accuracy. We note that the miniJPAS PDZs range from z = 0 to z = 1.5, that is, a much larger redshift range than the redshift upper-limit (z = 0.7) of the LMFs in this work. All the data related to the photo-z constraints can be downloaded from the CEFCA web portal¹⁵ via Astronomical Data Query Language (ADQL) queries, which are available in the table minijpas.PhotoZLephare_updated with keywords SPARSE_PDF and ODDS.

2.4. Star and galaxy classification

To separate the stars from galaxies, we took advantage of the Bayesian classification by probability distribution function (PDF) analysis initially developed for the J-PLUS survey by López-Sanjuan et al. (2019). This Bayesian analysis accounts for morphological information from the miniJPAS g_SDSS, r_SDSS, and i_SDSS broad-band filters that combines with prior probabilities on the r_SDSS magnitude and parallax information from Gaia catalogues. The star and galaxy classification method assigns a Bayesian probability between unity and zero to determine that a source is classified as a star or galaxy (P_⋆ and P_G, respectively). As this classification only comprises two kinds of sources, it is trivially deduced that P_G = 1 − P_⋆. We find that the distribution of probabilities is strongly bimodal, where the vast majority of the miniJPAS sources present P_G values around unity or null values. For instance, at magnitudes brighter than r_SDSS = 22, only a 3% fraction of miniJPAs sources show probabilities in the range 0.1 < P_G < 0.99. Intermediate values of P_G mainly comprises faint sources. As the majority of the miniJPAS value-added catalogues, the star and galaxy classification assumed in this work can be found at the CEFCA web portal, more precisely, in the table minijpas.StarGalClass and keyword total_prob_star.

It is worth mentioning that additional galaxy and star classifications, such as those of Baqui et al. (2021), were recently added. Machine learning (ML) methods, as those explored in Baqui et al. (2021), have demonstrated their efficacy in performing star and galaxy classifications with a high success rate, specially at faint magnitudes. However, we find that there are little discrepancies between performances of the best ML algorithms from Baqui et al. (2021) and the ‘default’ miniJPAS star and galaxy classification at magnitudes r_SDSS < 22.5 (see Figs. 11 and 12 in Baqui et al. 2021). Therefore, we did not expect that the method chosen to perform the star and galaxy classification (default versus ML) to play a role on the determination of the LMFs up to r_SDSS ≤ 22 at all. However, we also did not investigate the impact of utilising alternative classifications for determining the miniJPAS LMFs since it falls outside the scope of this study.

3. SED-fitting analysis of the miniJPAS sources

Along with the parameters obtained from the value-added catalogues, we obtained the rest of required parameters from SED-fitting techniques (Sect. 3.1) and, more precisely, we based our methods on the PDFs (Sect. 3.2) of the parameters involved to perform a statistical analysis. The SED fitting results allow us to determine not only the rest-frame luminosities and stellar mass, but also the rest-frame colours corrected for extinction for the spectral-type classification (Sect. 4) for each of the sources in our parent sample. Thanks to all these parameters and PDFs, we were also able to determine the stellar mass and luminosity completeness of our flux-limited sample (Sect. 5).

3.1. Updated version of MUFFIT and input ingredients for SED-fitting analysis

The stellar population properties of galaxies needed to perform our analysis were constrained using an updated version of the MUlti-Filter FITting for stellar population diagnostics code (MUFFIT; Díaz-García et al. 2015). In brief, MUFFIT is a SED-fitting code particularly developed and optimised to deal with multi-band photometric data. MUFFIT is based on an error-weighted χ²-test and includes composite models of stellar populations (CSP; a non-parametric star formation history based on mixtures of two simple stellar population models or SSPs) in order to constrain the stellar population properties of galaxies. MUFFIT has proven to be a powerful tool for this aim and easily adaptable to tackle the typical peculiarities that this kind of survey entails (see e.g. Díaz-García et al. 2015, 2019a,b,c; González Delgado et al. 2021).

All the properties of the stellar content of galaxies were established solely on the basis of the stellar continuum of galaxies or colours. For miniJPAS and J-PAS, this is specially relevant since their filter sets mainly include narrow bands that are sensitive to strong nebular or AGN emission lines, which are (in turn) not accounted by SSP models. In this regard, MUFFIT is able to remove bands affected by strong emission lines (e.g. [O II], [O III], Hβ, Hα+[N II], and [S II]) from the SED-fitting analysis; hence, getting a proper determination of galaxy properties via colours. Errors and degeneracies in the parameters due to photon-noise uncertainties are approached via Monte Carlo simulations through the flux uncertainties of each of the bands involved in the SED-fitting analysis. The number of Monte Carlo realisations is set to N_MC = 100, although this number can be increased in a general case at the cost of a longer computational time. From our experience with SED-fitting techniques and multi-filter data, this number is large enough as to perform a statistical analysis of the errors, given the number of sources included in our sample.

Throughout this work, we chose a recent version of the Bruzual & Charlot (2003) SSP models (hereafter, CB17) in order to build our CSP set. In particular, we selected CB17 models for 18 ages ranging from 0.001 to 13.5 Gyr and PARSEC evolutionary stellar tracks with metallicities log₁₀(Z/Z_⊙) = −1.93, −1.53, −0.93, −0.63, −0.33, 0.00, 0.25, and 0.55. In addition, we assumed a universal Chabrier (2003) IMF. The CB17 spectral coverage, λλ 15 Å–36 000 μm, is sufficiently wide as to carry out the SED-fitting analysis of galaxies at the highest redshift of the parent sample (z = 1.5). For this reason, the redshift of our CSPs ranges from z = 0 to z = 1.5, which are equally spaced by Δz = 0.01. Moreover, we add cosmological constraints on the age of the SSP components of our CSPs, in the sense that the age of any component cannot be much older than the age of the Universe at the CSP redshift. Extinctions were added to the CSPs as a foreground screen assuming the same dust attenuation in each of the CSP components. For this purpose, we follow the attenuation law by Calzetti et al. (2000) with values in the range A_V = 0.0–2.5 and assuming a constant of R_V = 3.1. We note that we assume that there is no distinction between extinction and attenuation law. It is worth mentioning that we do not assess potential systematics from the use of a given population synthesis model (see e.g. Díaz-García et al. 2019b), hence we adopt a unique CSP set and/or star formation history.

In the updated version of MUFFIT, photo-z data are treated in a more proper way than in previous versions. Even though it is expected that a large fraction of miniJPAS galaxies have a photo-z precision of σ_NMAD ∼ 0.003 (see Sect. 2.3 and Hernán-Caballero et al. 2021), which would have little impact in the stellar mass determination (details in Díaz-García et al. 2015), some galaxies in the parent sample exhibit very complex PDZs that we must manage properly (see Fig. 1). In brief, the SED-fitting analysis of each source in the sample is restricted to the photo-z values defining the confidence level of 70% of probability in its PDZ (see shaded area in Fig. 1), which is used as a prior probability. Afterwards, each of the Monte Carlo realisations are weighted according the PDZ probability (see red dots in Fig. 1) in order to determine the average stellar population properties and uncertainties (see Eqs. (18)–(20) in Díaz-García et al. 2015) of each miniJPAS galaxy.

Fig. 1. Photo-z probability distribution function (solid black line) of one noisy galaxy in miniJPAS (ID 2406-10239 and odds = 0.52). The shaded area illustrates the photo-z values defining the confidence level of 70% of probability. The red dots show the photo-z values used by MUFFIT to perform the SED-fitting analysis of this source.

Owing to the low S/N seen for some of the galaxies in the parent sample, the χ² minimisation is performed using fluxes instead of magnitudes. Furthermore, to ensure a proper fit of the data, we exclude from the analysis bands flagged with values different to FLAGS = 0, 2 and MASK_FLAGS = 0 and we require a minimum of ten bands with a S/N value higher than 1.5. These constraints result in rejecting a 4% fraction of the sources from the SED-fitting analysis, so our final sample is composed of 12 100 sources. It is worth mentioning that the miniJPAS database includes PDZ constraints for all the sources in the catalogue and no source of the sample was excluded from the analysis because of its photo-z estimate and/or odds value.

As a result of the SED fitting analysis of the final sample, we obtain rest-frame luminosities, stellar mass, rest-frame colours corrected for extinction (usually referred as intrinsic colours, see Díaz-García et al. 2019a), and photo-z (treated as another free parameter during the SED-fitting analysis with values in the PDZ 70% confidence level). Indeed, other stellar population properties such as extinction, age, and metallicity (both luminosity- and mass-weighted) are constrained during the MUFFIT analysis and can be used in future works. We note that stars or point-like sources are included in the SED-fitting analysis, we statistically tackle this issue by making use of the star and galaxy classification (P_G, see Sect. 2.4) in subsequent sections. In Fig. 2, we illustrate a SED-fitting case of one of the brightest galaxies in the miniJPAS survey. This galaxy exhibits an elliptical morphology with an evolved stellar content compatible with a mass-weighted age of 10 Gyr and a sub-solar metallicity of −0.2 dex. It presents red colours that are barely produced by a high dust content (A_V = 0.2) with a stellar mass of 10¹¹ M_⊙. There is no evidence of any presence of emission lines.

Fig. 2. SED-fitting analysis (top panel) and residuals (bottom panel) of a miniJPAS galaxy at z = 0.07 (ID 2470-10239). The black dots and vertical bars in the top panel are the galaxy fluxes and errors as observed by miniJPAS, respectively. The red squares are the best-fitting CSP model. The shaded area shows 2.5 times the photon-noise uncertainty of each band.

3.2. Discretised PDF analysis of galaxies and B-band luminosity

In the following (and due to the nature of the results yielded by our SED-fitting code MUFFIT), we distinguish two sets of results. Firstly, the set of best solutions obtained from the Monte Carlo realisations, which comprises N_MC stellar masses, luminosities, or absolute magnitudes for a set of bands, photo-z, and rest-frame colours corrected for extinction for each of the sources in the final sample. At this point, we assume that the PDF of these stellar population properties is properly described by the Monte Carlo realisations, and therefore, this set of results is equivalent to a discretised version of the PDFs. Secondly, the average or most likely stellar population properties and their errors for each of the sources in the final sample. This second set results from the weighted mean and weighted standard deviation of the Monte Carlo realisations. These values are weighted according to the PDZ and 1/χ² values of each realisation.

Additionally, we computed an extra set of luminosities or absolute magnitudes. In particular, we select the so-called B-broad band from the Classifying Objects by Medium-Band Observations survey¹⁶ (COMBO-17, Wolf et al. 2003). The choice of the B-band for this purpose is based on the fact that it has been used in several works for studying the luminosity functions of galaxies in a wide redshift range. Furthermore, given the observational wavelength range of the miniJPAS survey, rest-frame B-band luminosities (effective or pivot wavelength at λ_pivot ∼ 4570 Å) are directly observed by the miniJPAS photometric system up to z ∼ 1. However, we are able to explore luminosity functions with redder bands thanks to our SED-fitting results, but depending on the redshift and band, the luminosities would actually be an extrapolation of the SED-fitting analysis. For instance, luminosities for the miniJPAS r_SDSS band (λ_pivot ∼ 6250 Å) are not directly observed at redshifts higher than z ∼ 0.5, which set limits in the bands that can be used for determining luminosity functions up to z ∼ 0.7. The process for estimating the B-band luminosity of miniJPAS galaxies lies in the SED-fitting results. Briefly, for each of the sources in our final sample and thanks to the MUFFIT analysis, we have a set of SSP models reproducing the photometric SED or colours at the galaxy redshift. Therefore, from exactly the same combination of SSP models, we are able to reconstruct the rest-frame spectrum of each miniJPAS galaxy, which is subsequently convolved with the B-band and instrument transmission curves¹⁷. We repeat this process for each of the Monte Carlo realisations, getting a discretised PDF of the B-band luminosity for each miniJPAS source. At this point, it is noteworthy that the SED-fitting analysis is performed without including emission lines and, therefore, all the luminosities reported in this work only contain predictions based on the stellar continuum and not about the nebular contribution.

4. Spectral-type classification of galaxies

Over recent decades, many diagrams have been proposed to select or classify galaxies by star-formation activity or spectral type. These classifications mostly involve two main kinds of galaxies, which are unveiled through bimodal distributions of properties. Usually, these two types are referred to as star-forming or quiescent since the discrepancies between them are the result of high versus low levels of star formation activity, respectively (see Díaz-García et al. 2019a, and references therein). For this aim, the most extended diagrams and modern methods confront sets of galaxy features that are sensitive to the spectral-type classification of galaxies such as rest-frame colours versus absolute magnitudes (CMDs, e.g. Wyder et al. 2007; Faber et al. 2007), empirically observed rest-frame colour–colour diagrams (e.g. the so-called UVJ diagrams Williams et al. 2009; Arnouts et al. 2013) and others, including stellar population properties (star formation rate, stellar mass, etc.; Ilbert et al. 2010; Whitaker et al. 2012; Moustakas et al. 2013). However, recent studies have pointed out that samples of quiescent galaxies defined by these diagrams may include a remarkable fraction of dusty star-forming galaxies (5–40% of galaxies, depending on the stellar mass range and redshift; Díaz-García et al. 2019a) and this selection can be significantly improved after correcting colours for extinction or dust reddening (see also Moresco et al. 2013; Schawinski et al. 2014; Díaz-García et al. 2019a,b; Antwi-Danso et al. 2023).

In the present work, we are particularly interested in basing our spectral-type classification of galaxies on the stellar mass-colour diagram corrected for extinction (MCDE, Díaz-García et al. 2019a). The MCDE was proposed as an efficient method in order to build non-biased samples of quiescent galaxies (Díaz-García et al. 2019a) with respect to some previous methods, where the dust correction of the involved rest-frame colours is actually a key step in the process. In the following, we remark some of the advantages that motivated the use of this diagram for the spectral-type classification of galaxies. Firstly, the MCDE diagram includes multiple observables in order to perform the galaxy classification, such as stellar mass, redshift, and rest-frame colours. These aspects will account for the evolution of colours owing to galaxy ageing, as well as low-mass quiescent galaxies exhibit bluer colours on average than the more massive ones. Secondly, these diagrams do not present a strong dependency on the SSP model set employed for the SED-fitting analysis (see Díaz-García et al. 2019a) and/or the star formation history assumed (e.g. parametric or non-parametric) since it is primarily based on rest-frame colours rather than model-dependent stellar population properties such as star formation rates. In fact, the MCDE diagram was shown to yield results that are 98% compatible with the stellar mass versus star formation rate diagram up to z = 1 (see Sect. 5 in Díaz-García et al. 2019a). Thirdly, colours corrected for extinction mildly depend on the choice of the dust extinction law. As showed by Díaz-García et al. (2019a), the sample of quiescent galaxies differs by around 3% in number when the attenuation law by Calzetti et al. (2000) is used instead of the Fitzpatrick (1999) extinction law. Finally, and partly due to the above-mentioned reasons, the limiting values obtained through the MCDE diagram can be easily extended to current and future results from other SED-fitting codes in the J-PAS Collaboration (see González Delgado et al. 2021). It is worth mentioning that we also explored the possibility of using UVJ-like diagrams for the spectral classification of galaxies. However, we found out that rest-frame magnitudes of the near-infrared J-PAS bands (e.g. J1007, which is close to the Y band) were actually extrapolations of the SED-fitting analysis at redshift z > 0.1 as a consequence of the J-PAS photometric system; these, in turn, are tightly related to the star formation history assumption leading to disparate results. For this reason, in this work we discarded UVJ diagrams for the spectral classification of galaxies.

Nevertheless, we found that the sample of miniJPAS sources presented some difficulties in setting limiting values in the MCDE diagram for the spectral classification of galaxies. This is mainly due to two factors. On the one hand, the volume or area imaged in the miniJPAS survey (Ω = 0.895 deg²) limits the number of detected sources, especially in the nearby Universe. On the other hand, the depth of this survey does not allow us to build samples including low-mass quiescent galaxies that are complete in stellar mass at intermediate redshift. For instance, the sample of quiescent galaxies is 95% complete for stellar mass values of log₁₀M_⋆ ≳ 10 dex at z = 0.3 (see Sect. 5). This makes difficult a robust determination of the limiting values as a function of stellar mass and redshift. In other words, the slope of the stellar mass versus the rest-frame colour corrected for extinction show large uncertainties that complicates the galaxy classification. To overcome this problem, we made use of the sample that we also analysed with our SED-fitting code MUFFIT in Díaz-García et al. (2019a) using galaxies from the ALHAMBRA survey, which comprises deeper observations than miniJPAS with F814W < 24.5 across 2.8 deg² of the northern hemisphere (further details in Molino et al. 2014). For this aim, we took all the galaxies from Díaz-García et al. (2019a) and we reconstructed the MCDE diagram using the fittings obtained from ALHAMBRA, but for a rest-frame colour compatible with the miniJPAS photometric system. In particular, we chose the rest-frame colour (u_JAVA − r_SDSS) corrected for extinction, which we refer to as intrinsic colour or (u_JAVA − r_SDSS)_int in the following. For recomputing this colour for the ALHAMBRA galaxies, we follow a similar process than the detailed in Díaz-García et al. (2019a), which (in turn) is the process internally done by MUFFIT for performing the k-corrections of the involved photometric bands (Díaz-García et al. 2015). In brief, we take all the best-fitting models or the same combinations of SSP spectra (mixtures of two SSPs in our case) obtained from the Monte Carlo realisations that reproduces the SED of each galaxy in the sample, but at z = 0 and null extinction (A_V = 0). These sets of rest-frame spectra with null extinction are subsequently convolved with the miniJPAS photometric system (i.e. synthetic photometry) in order to determine the (u_JAVA − r_SDSS)_int colours for each of the ALHAMBRA galaxies and, hence, a MCDE diagram that can be easily extended to miniJPAS galaxies. The next step is to set limiting values for the (u_JAVA − r_SDSS)_int colours according to the stellar mass and redshift of the sample. For this goal, we repeat again the same process detailed in Díaz-García et al. (2019a): (i) after close inspection of the distribution of points in the MCDE, we set straight lines at different redshift bins that split the distribution of quiescent and star-forming galaxies; (ii) for each of the redshift bins, we computed the representative intrinsic colour of quiescent galaxies as a function of stellar mass and redshift; (iii) after removing uncertainty effects of the distributions of colours of quiescent galaxies via a maximum likelihood method (MLE, further details in López-Sanjuan et al. 2014; Díaz-García et al. 2019a), we set the limiting value of quiescent galaxies, ${(u_{\text{JAVA}}-r_{\text{SDSS}})}_{\text{int}}^{\text{lim}}$, as the 3σ limit of the distribution of colours fitted in the previous step.

As a result, the limiting values for the spectral-type classification of galaxies is expressed by an equation of the form:

$$\begin{array}{c}\hfill {(u_{\mathrm{JAVA}}-r_{\mathrm{SDSS}})}_{\mathrm{int}}^{\lim }=a\times ({log}_{10}{M}_{\star }-10)+b\times (z-0.1)+c,\end{array}$$(1)

where a, b, and c are the constants determined through the ALHAMBRA galaxies for the miniJPAS galaxy classification (see values and uncertainties in Table 1). Nevertheless, these coefficients depend on the photometric apertures or magnitude system assumed. In this regard, the miniJPAS PSFCOR magnitudes were computed following the approach detailed in Molino et al. (2014, 2019, namely, based on the same process carried out to perform the ALHAMBRA photometry), meaning that the a, b, and c values obtained for the spectral-type classification of galaxies via ALHAMBRA data are ‘scaled’ according to the rest-frame colours and stellar masses resulted from the miniJPAS PSFCOR photometry. This is important since we are interested in using AUTO magnitudes to constraint the miniJPAS LMFs, which are computed in larger apertures than PSFCOR. As a consequence, stellar masses obtained from SED-fitting analysis by using AUTO magnitudes are systematically larger than the PSFCOR ones and this fact can be also extended to the (u_JAVA − r_SDSS)_int colours. With a view to future results of the J-PAS Collaboration, we applied a second-order correction to adapt the a, b, and c values to the AUTO magnitudes. Therefore, we provided two sets of a, b, and c values (see Table 1) that can be used to perform the classification of miniJPAS galaxies according to results based on PSFCOR or AUTO magnitudes. After a detailed analysis of the results obtained from MUFFIT for AUTO apertures (this work) and those from PSFCOR magnitudes (briefly introduced in González Delgado et al. 2021), we conclude that SED-fitting results of miniJPAS sources with AUTO photometry are in average shifted towards higher stellar mass values and bluer colours by 0.21 dex and 0.07 magnitudes, respectively, with respect to PSFCOR (also in agreement with Figs. 8 and 9 in González Delgado et al. 2021). In addition, we find that these systematic differences are independent of the stellar mass range and galaxy colour and, hence, of the spectral-type. Consequently, the second-order correction of the relation obtained from ALHAMBRA galaxies for the classification of galaxies only comprises the intercept term, c, in Eq. (1). After accounting for the previous offsets, the intercept for galaxy classification for AUTO photometry is ∼0.1 smaller than for PSFCOR (see Table 1).

For the galaxy classification, we assigned probabilities to each of the sources in our sample. For this aim, we make the most of the Monte Carlo realisations obtained from the SED-fitting analysis performed by the MUFFIT code. Essentially, we count the fraction of realisations presenting redder colours than the limiting values obtained by Eq. (1) for AUTO magnitudes (Table 1), denoted as P_Q, for each of the sources in the sample. Under this definition, sources with P_Q ≥ 0.5 are labelled as quiescent galaxies, otherwise these are classified as star-forming galaxies. Actually, P_Q × P_G can be assumed as the probability of a source in our sample of being a quiescent galaxy and this parameter accounts for parameter uncertainties, degeneracies and correlations amongst SED-fitting parameters (e.g. extinction and rest-frame colours, stellar mass and redshift, etc.), the miniJPAS star and galaxy classification, and the redshift and stellar mass dependency of the limiting values for galaxy classification. Likewise, the probability of a source of being a star-forming galaxy can be established as P_SF × P_G, where P_SF = (1 − P_Q). We note that for this work we only classify galaxies as either quiescent or star-forming, that is, extra or alternative classifications such as ‘green valley’ galaxies are not explored in the present paper.

As expected, the distribution of P_Q values strongly correlates with the positions that sources occupy in the MCDE (see Fig. 3). Galaxies with the intrinsic reddest colours (upper parts of the MCDE) present the highest P_Q values (quiescent galaxies, see red dots in Fig. 3), whereas the lower P_Q values lie on the lower parts of the MCDE (star-forming galaxies, see blue dots in Fig. 3). Sources with an uncertain spectral-type classification (i.e. P_Q ∼ 0.5) commonly present intermediate rest-frame colours in the MCDE, which overlaps with the colour range of the so-called green valley galaxies. Overall, the distribution of P_Q values for miniJPAS galaxies is highly bimodal for a magnitude limit of r_SDSS ≤ 22.5 up to z < 0.7. In this regard, we find a low fraction of galaxies with an uncertain spectral-type classification (0.4 ≤ P_Q ≤ 0.6) that ranges from 3% at z ≤ 0.3 to 12% at the higher redshifts in our sample (see insets in Fig. 3) owing to photon-noise uncertainties. Finally, we find that there are correlations between the P_Q values and their positions across the classical stellar mass versus colour diagram, that is, without correcting for dust reddening (see bottom panels in Fig. 3). In fact, the higher (lower) P_Q also lie on the upper (lower) parts of this diagram, but as also found in Díaz-García et al. (2019a), there is a non-negligible fraction of star-forming galaxies (P_Q < 0.5) that present red colours owing to a large dust content or extinction. Moreover, galaxies with 0.4 ≤ P_Q ≤ 0.6 are scattered in the colour range of quiescent galaxies in this diagram.

Fig. 3. Distribution of stellar mass versus rest-frame colour (uJAVA − rSDSS) obtained from miniJPAS sources with rSDSS ≤ 22.5 and PG ≥ 0.5 at different redshift bins. Top: diagram of the (uJAVA − rSDSS) colour after correcting for dust effects or intrinsic colour. Bottom: diagrams without correcting the rest-frame (uJAVA − rSDSS) colour for extinction and normalised histograms of PQ values for each redshift bin (see insets). The dashed black line illustrates the limiting relation for selecting quiescent galaxies and AUTO photometry at the central redshift of each bin (i.e. z = 0.15, 0.4, and 0.6; see Eq. (1)). Dots are colour coded according to the quiescent and star-forming classification (redder and bluer colours, respectively) or quiescent probability (PQ).

5. Luminosity and stellar mass completeness

Owing to the nature of multi-filter surveys like J-PAS, there is no selection criteria other than the photometric depth of the detection band, and therefore, all the sources brighter than that limit in each of the survey pointings are imaged and included in the photometric catalogues. Nevertheless, we distinguish three types of completeness limits that can affect the determination of the miniJPAS LMFs. Firstly, the image depth sets a first unavoidable limitation for defining our galaxy sample. In miniJPAS, it was established that the sample of extended sources is complete up to r_SDSS = 22.7, that is, all the galaxies in the field that are equal or brighter than this magnitude limit are detected in the miniJPAS images and included in the photometric catalogue. As our sample is brighter than this limit, our galaxy sample is not affected by this selection effect. On the other hand, our sample of miniJPAS sources is defined according to a magnitude limit that in turns triggers limitations in stellar mass and B-band luminosity completeness. This would imply that, at a given redshift and close to the magnitude limit of the sample, part of the galaxies below a certain stellar mass and luminosity limit are not observed or included in the photometric catalogue.

It is well known that stellar mass and luminosity completeness of flux-limited samples are properly described by Fermi-Dirac distribution functions (see e.g. Sandage et al. 1979; Ilbert et al. 2010; Díaz-García et al. 2019a). For our aims, we therefore determine two Fermi-Dirac distribution functions for the stellar mass and B-band luminosity completeness (𝒞 hereafter), where each of them is redshift dependent and formally parametrised by two parameters, as follows:

$$\begin{array}{c}\hfill \mathcal{C}(z,M_{\star })=\frac{1}{exp[(M_{\mathrm{F},M_{\star }}(z)-{log}_{10}{M}_{\star })/{\mathrm{\Delta }}_{\mathrm{F},M_{\star }}(z)]+1},\end{array}$$(2)

$$\begin{array}{c}\hfill \mathcal{C}(z,M_B)=\frac{1}{exp[(M_B-M_{\mathrm{F},M_B}(z))/{\mathrm{\Delta }}_{\mathrm{F},M_B}(z)]+1},\end{array}$$(3)

where 0 ≤ 𝒞 ≤ 1 (i.e. the completeness of the sample at certain stellar mass or luminosity at a given redshift), M_F is the stellar mass in dex units or absolute magnitude in the B-band for which the completeness reaches 50% (𝒞 = 0.5), and Δ_F the decrease rate on the fraction of galaxies. It is also of note that both M_F and Δ_F are functions that depend on redshift in a general case, meaning that we do not assume a constant value for Δ_F. From Eqs. (2) and (3), it is easy to obtain the limiting values of stellar mass and absolute magnitude for a given completeness level (${log}_{10}{M}_{\star }^{\mathcal{C}}$ and $M_B^{\mathcal{C}}$, respectively) as:

$$\begin{array}{c}\hfill {log}_{10}{M}_{\star }^{\mathcal{C}}(z)=M_{\mathrm{F},M_{\star }}(z)+{\mathrm{\Delta }}_{\mathrm{F},M_{\star }}(z)\times ln\left(\frac{\mathcal{C}}{1-\mathcal{C}}\right),\end{array}$$(4)

and

$$\begin{array}{c}\hfill M_B^{\mathcal{C}}(z)=M_{\mathrm{F},M_B}(z)-{\mathrm{\Delta }}_{\mathrm{F},M_B}(z)\times ln\left(\frac{\mathcal{C}}{1-\mathcal{C}}\right).\end{array}$$(5)

In addition, we can conclude from Eqs. (4) and (5) that M_F and Δ_F can be determined or constrained as long as we knew two stellar mass and B-band luminosity values at two certain completeness levels.

For the determination of the stellar mass completeness of the miniJPAS flux-limited sample, we adopted the method proposed by Pozzetti et al. (2010). This method relies on the assumption that, at a given redshift, the distribution of the mass-to-light ratio of the fainter galaxies in the sample (M_⋆/L_r) should be similar to the one at the magnitude limit of the sample ($M_{\star }^{\text{lim}}/L_r^{\text{lim}}$). Keeping this in mind, the distribution of stellar mass at a fixed magnitude limit can be reconstructed from the above assumption as:

$$\begin{array}{c}\hfill {log}_{10}{M}_{\star }^{\lim }={log}_{10}{M}_{\star }-0.4\times (r_{\mathrm{SDSS}}^{\lim }-r_{\mathrm{SDSS}}),\end{array}$$(6)

where r_SDSS is the distribution of observed magnitudes for the galaxies with the fainter magnitudes at certain redshift, log₁₀M_⋆ are the stellar mass values for each of these faint galaxies, and $r_{\text{SDSS}}^{\text{lim}}$ the limiting magnitude imposed for the definition of the sample. If our assumption is correct, the distribution of ${log}_{10}{M}_{\star }^{\text{lim}}$ values illustrates the range and frequency of stellar mass values of galaxies close to the magnitude limit. Consequently, all the galaxies with a stellar mass higher than the 95th percentile of the ${log}_{10}{M}_{\star }^{\text{lim}}$ distribution are roughly included in the sample (i.e. 𝒞 = 0.95), whereas galaxies with stellar mass under the 5th percentile of the distribution are not going to be included in the sample (i.e. 𝒞 = 0.05).

Regarding the B-band luminosity completeness, we follow a similar approach, but we assume that the distribution of rest-frame colours of the fainter galaxies in the flux-limited sample at certain redshift is similar to those at the magnitude limit instead. As for the stellar mass completeness case, this results in a distribution of B-band absolute magnitudes at the magnitude limit of the sample formally expressed as

$$\begin{array}{c}\hfill M_B^{\lim }=M_B+(r_{\mathrm{SDSS}}^{\lim }-r_{\mathrm{SDSS}}).\end{array}$$(7)

For this case, the sample would be complete in terms of the B-band absolute magnitude for galaxies brighter than the 5th percentile of the $M_B^{\text{lim}}$ distribution and incomplete for galaxies fainter than the 95th percentile, which we assume that roughly correspond to 𝒞 = 0.95 and 0.05, respectively.

For the determination of the 𝒞 = 0.05 and 0.95 limits of the ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ distributions, we make use of the average results and uncertainties obtained by MUFFIT. We split the sample in redshift bins of width Δz = 0.1 from z = 0.05 to 0.75 for all the miniJPAS sources with P_G ≥ 0.5. Moreover, we split this subsample in quiescent and star-forming galaxies using our galaxy classification probability (i.e. P_Q ≥ 0.5 and P_Q < 0.5, respectively; see Sect. 4). For each of the redshift bins in our subsamples, we select galaxies with observed magnitudes in the range $r_{\text{SDSS}}^{\text{lim}}$ − 0.5 ≤ r_SDSS ≤ $r_{\text{SDSS}}^{\text{lim}}$ for building our distributions of limiting parameters (see Eqs. (6) and (7)). However, for the cases for which the number of galaxies is lower than 30, we extend the bright magnitude limit until this minimum number of galaxies is fulfilled. This particularly happens for quiescent galaxies in low redshift bins since the observed volume and number density of these kinds of galaxies is lower. Finally, we account for the fact that the ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ distributions are built from samples that were determined via SED-fitting techniques and these are affected by uncertainties. If this issue were not addressed, the distribution of limiting magnitudes would be broader than it actually is. Consequently, the upper completeness limit would be overestimated, which would greatly affect and restrict the number of sources during the definition of complete samples, as well as the lower limit would be underestimated. In fact, as detailed in Sect. 7.1, the errors in stellar mass and luminosity are large enough at the fainter magnitudes of the sample as to significantly alter the ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ distributions. To overcome this issue, we applied the MLE method previously mentioned for deconvolving uncertainty effects in the ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ distributions. Interestingly, the ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ distributions of values are properly fitted by Gaussian and/or log-normal distributions, and consequently, the MLE method described in Díaz-García et al. (2019a,b) can be easily and analytically applied here. Then, the 𝒞 = 0.05 and 0.95 limits are obtained from the analytic functions obtained through the MLE (i.e. the intrinsic distributions of ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ after deconvolving the errors in stellar mass and B-band absolute magnitude provided by MUFFIT), which are used to subsequently determine M_F and Δ_F by Eqs. (4) and (5). In this regard, M_F matches the intermediate value between those for 𝒞 = 0.05 and 0.95 (see Eqs. (4) and (5)). In fact, as the ${log}_{10}{M}_{\star }^{\text{lim}}$ and $M_B^{\text{lim}}$ distributions are properly fitted by Gaussian-like functions, it turns out to be a good approach that M_F corresponds to the medians since resulted from the MLE method. We note that as our plan is to study the stellar mass and B-band luminosity completeness for a set of magnitude limit values, we repeat the method detailed above for each of these sample constraints, more precisely, we perform this method for $r_{\text{SDSS}}^{\text{lim}}$ = 21.5, 21.75, 22, 22.25, and 22.5.

With the ultimate goal of obtaining an analytic expression of the M_F and Δ_F parameters as a function of redshift, we fit the ${log}_{10}{M}_{\star }^{\mathcal{C}}$ and $M_B^{\text{C}}$ values obtained in the previous step to a power-law function of the form μ × (z + γ)^ν and δ × (z + ζ)^ϵ, respectively. This process is repeated for the completeness values of 𝒞 = 0.05, 0.5, and 0.95. In fact, these directly determine M_F(z) since ${\log }_{10}{M}_{\star }^{\mathcal{C}=0.5}(z)=M_{\text{F},M_{\star }}(z)$ and $M_B^{\mathcal{C}=0.5}(z)=M_{\text{F},M_B}(z)$. On the other hand, Δ_F, M⋆(z) and Δ_F, MB(z) can be obtained by Eqs. (4) and (5) and the fittings just mentioned. The μ, ν, γ, δ, ϵ, and ζ values for the completeness levels 𝒞 = 0.05, 0.5, and 0.95 and a flux-limited sample of r_SDSS ≤ 22 (i.e. $r_{\text{SDSS}}^{\text{lim}}$ = 22) are shown in Table 2. In Appendix A, we extend our results for alternative magnitude-limited samples (see Tables A.1–A.4), which can be used in future J-PAS studies. Overall, all the values that we obtain for each of the redshift bins are properly fitted by power-law functions independently of the spectral-type of galaxies (see Fig. 4) and the sample is complete (𝒞 = 0.95) for higher luminosity and stellar mass limits at increasing redshift. As previously observed in similar works, samples of quiescent galaxies are complete at higher luminosity and stellar mass limits than star-forming galaxies at same redshift, which is reflected in larger M_F values. On the other hand, we find that Δ_F, M⋆ is smaller for quiescent galaxies than for star-forming galaxies as a consequence of a larger diversity of mass-to-light ratios in star-forming galaxies. In fact, Δ_F, M⋆ is roughly constant for our star-forming galaxy sample with r_SDSS ≤ 22 with a value of 0.16, whereas for quiescent galaxies this value is always lower and mildly decreases from 0.11 at z ∼ 0 to 0.06 at z ∼ 0.7. Regarding Δ_F, MB, our results point out that this value present larger changes with redshift than for the stellar mass case, although this value is approximately constant at z ≳ 0.3 (Δ_F, MB ∼ 0.10) for both star-forming and quiescent galaxies. At lower redshift, Δ_F, MB increases up to ∼0.3 at z ∼ 0.1, independently of the galaxy spectral-type. Consequently, we find that amongst the four parameter defining the completeness limits of the sample, Δ_F, MB and Δ_F, M⋆ are almost independent parameters on the galaxy spectral-type and redshift, respectively. Finally, it is worth mentioning that given the volume observed in miniJPAS, the completeness limits constrained in this work are valid at 0.1 ≲ z ≲ 0.7, and out of this redshift range these should be managed as an extrapolation of the fittings.

Fig. 4. B-band luminosity and stellar mass completeness as a function of redshift (top and bottom panels, respectively) for star-forming and quiescent miniJPAS galaxies (left and right panels, respectively) from the flux-limited sample at rSDSS ≤ 22 magnitudes (i.e. $r_{\text{SDSS}}^{\text{lim}}$ = 22). The dotted, dashed, and solid black lines illustrate the fit of the 5, 50, and 95% completeness levels (i.e. 𝒞 = 0.05, 0.5, and 0.95) obtained from the completeness methodology detailed in Sect. 5 (triangle, square, and circle shape markers), respectively.

Even though this is not an ideal method for determining the completeness curves at the magnitude limit, we can set rough estimations of them that can be used to set limits for the definition of galaxy samples in the miniJPAS survey. An ideal methodology for the determination of the completeness limits and analytical functions would be a direct comparison with respect to the LMFs functions of a similar and deeper survey in an overlapping area (see e.g. Díaz-García et al. 2019a). Other methods involve the injection of synthetic sources or mock catalogues built from simulations in order to check the effect of the introduction of flux limits in the sample. Nevertheless, the use of mock catalogues involve the inclusion of population synthesis models (as well as diverse star formation histories) and a large set of assumptions or distributions (e.g. involving number densities and morphology) typically based on deeper surveys, which may not be trivially accounted for. In any case, the Pozzetti et al. (2010) method has proven to be a good approach for setting constraints on the distributions of stellar mass and luminosity of galaxies at the magnitude limit of the sample (as also concluded by Meneux et al. 2009; Pozzetti et al. 2010; Díaz-García et al. 2019a, after comparison with simulations and results from deeper surveys), which somehow reflects the kind of bias introduced in the selection of the sample in terms of stellar mass and absolute magnitude. With the incoming of J-PAS in a close future, which will involve large overlapping areas with deeper surveys such as ALHAMBRA, we expect to confirm how much precise this method is for the determination of the completeness limits explored in this work.

6. Methodology for determining the LMF parameters

In this research, we depart from the assumption in which the population of both star-forming and quiescent galaxies are properly described by a single Schechter function. This decision is mainly motivated by the fact that our sample at the low-mass galaxy regime is highly reduced in number and only present at the nearby Universe. Specifically, there are few quiescent galaxies of log₁₀M_⋆ < 9.5 at z < 0.2, which makes largely difficult the proper determination of a second Schechter component at the low-mass or faint regime. This inconvenience is a direct consequence of the imaging depth and the small area imaged in miniJPAS, although in a close future with the arrival of J-PAS and its huge observing area, we expect to explore the LMFs via double Schechter functions. Secondly, in this research the LMFs of the whole galaxy population is the sum of the two galaxy spectral-types, meaning that the LMFs of all the galaxies in miniJPAS are actually double Schechter functions. Nonetheless, the LMFs of a single Schechter function at a given redshift are formally expressed as:

$$\begin{array}{cc}\hfill \mathrm{\Phi }(L_B)\mathrm{d}{L}_B& ={\mathrm{\Phi }}_{L_B}^{\star }{\left(\frac{L_B}{{\mathcal{L}}_B}\right)}^{\beta }exp(-\frac{L_B}{{\mathcal{L}}_B})\frac{\mathrm{d}{L}_B}{L_B},\hfill \end{array}$$(8)

$$\begin{array}{cc}\hfill \mathrm{\Phi }(M_{\star })\mathrm{d}{M}_{\star }& ={\mathrm{\Phi }}_{M_{\star }}^{\star }{\left(\frac{M_{\star }}{{\mathcal{M}}_{\star }}\right)}^{\alpha }exp(-\frac{M_{\star }}{{\mathcal{M}}_{\star }})\frac{\mathrm{d}{M}_{\star }}{M_{\star }},\hfill \end{array}$$(9)

where ℒ_B and ℳ_⋆ are the so-called characteristic B-band luminosity and characteristic stellar mass, respectively, (β + 1) and (α + 1) are referred to as the slopes of the faint-end of the LMFs (see also Eqs. (10) and (11)), while ${\mathrm{\Phi }}_{L_B}^{\star }$ and ${\mathrm{\Phi }}_{M_{\star }}^{\star }$ are their normalizations in number density units. Alternatively, these functions can be also found in their logarithmic form as:

$$\begin{array}{cc}& \mathrm{\Phi }(M_B)\mathrm{d}{M}_B=0.4ln10{\mathrm{\Phi }}_{L_B}^{\star }{10}^{0.4({\mathcal{M}}_B-M_B)(\beta +1)}·\hfill \\ \hfill & ·exp[-{10}^{0.4({\mathcal{M}}_B-M_B)}]\mathrm{d}{M}_B,\hfill \end{array}$$(10)

$$\begin{array}{cc}& \mathrm{\Phi }({\overline{M}}_{\star })\mathrm{d}{\overline{M}}_{\star }=ln10{\mathrm{\Phi }}_{M_{\star }}^{\star }{10}^{({\overline{M}}_{\star }-{\overline{\mathcal{M}}}_{\star })(\alpha +1)}·\hfill \\ \hfill & ·exp[-{10}^{({\overline{M}}_{\star }-{\overline{\mathcal{M}}}_{\star })}]\mathrm{d}{\overline{M}}_{\star },\hfill \end{array}$$(11)

where M_B is the B-band absolute magnitude, ℳ_B is the absolute magnitude of the characteristic B-band luminosity, ${\overline{M}}_{\star }$ the stellar mass logarithm (i.e. ${\overline{M}}_{\star }={log}_{10}{M}_{\star }$), and ${\overline{\mathcal{M}}}_{\star }$ the logarithm of the characteristic stellar mass. Accordingly, the LMF shapes are defined according to two sets of redshift-dependent parameters, {ℳ_⋆(z), α(z), ${\mathrm{\Phi }}_{M_{\star }}^{\star }(z)$} and {ℒ_B(z), β(z), ${\mathrm{\Phi }}_{L_B}^{\star }(z)$}, that depend on the galaxy spectral-type.

6.1. Characteristic parameters and slopes of the faint and low-mass ends of LMFs

The characterisation of the miniJPAS LMFs rely on a maximum likelihood method in which we include the probability distribution functions of the involved parameters (see Sect. 3.2). In a sense, our methodology follows the guidelines and includes some aspects of the works by Sandage et al. (1979), Efstathiou et al. (1988), Zucca et al. (1994), Ilbert et al. (2005), López-Sanjuan et al. (2017) and is adapted for the inclusion of all the peculiarities of the miniJPAS survey, the results from MUFFIT and similar codes in the J-PAS Collaboration, the statistical spectral-type classification of galaxies, and the completeness of our flux-limited sample (details in Sects. 3–5). Based on the idea that the probability of observing a galaxy of stellar mass $M_{\star }^{\prime }$ at a given redshift z′ from a sample complete in stellar mass is $\mathrm{\Phi }(M_{\star }^{\prime },z^{\prime })/{\displaystyle {\int }_{M_{\star }^{\text{f}}(z^{\prime })}^{M_{\star }^{\text{b}}(z^{\prime })}}\mathrm{\Phi }(M_{\star },z^{\prime })\text{d}{M}_{\star }$, where Φ(M_⋆, z′) is the stellar mass function (as in Eq. (9)) and $M_{\star }^{\text{f}}(z^{\prime })$ and $M_{\star }^{\text{b}}(z^{\prime })$ are the low and high stellar mass limits of the sample at the source redshift, respectively; the likelihood (ℒ) of our method can be defined as the joint probability of observing our sample that is expressed as:

$$\begin{array}{c}\hfill \mathcal{L}={\displaystyle \prod _{i=1}^{N_{\mathrm{S}}}}{\left(\frac{\mathrm{\Phi }(M_{\star ,i},z_i)}{{\displaystyle {\int }_{M_{\star }^{\mathrm{f}}(z_i)}^{M_{\star }^{\mathrm{b}}(z_i)}}\mathrm{\Phi }(M_{\star },z_i)\mathrm{d}{M}_{\star }}\right)}^{{\omega }_i},\end{array}$$(12)

where N_S is the number of sources in the sample and ω_i is a weight accounting for different factors such as the probability of a source if being a star, the completeness of the sample at the ith source stellar mass and redshift, priors, and so on (further details given in Sect. 6.2). The maximisation of this likelihood makes possible the determination of the {ℳ_⋆(z), α(z), ${\mathrm{\Phi }}_{M_{\star }}^{\star }(z)$} parameters that characterise the most likely LMFs describing the observations of our sample. As in López-Sanjuan et al. (2017), we assume redshift dependent functions for these parameters of the form: α(z) = α₁ + α₂ × z, β(z) = β₁ + β₂ × z, ${log}_{10}{\mathcal{M}}_{\star }(z)={\overline{\mathcal{M}}}_{\star ,1}+{\overline{\mathcal{M}}}_{\star ,2}\times z$, ℳ_B = ℳ_B, 1 + ℳ_B, 2 × z, ${log}_{10}{\mathrm{\Phi }}_{L_B}^{\star }={\overline{\mathrm{\Phi }}}_{L_B,1}^{\star }+{\overline{\mathrm{\Phi }}}_{L_B,2}^{\star }\times z$, and ${log}_{10}{\mathrm{\Phi }}_{M_{\star }}^{\star }={\overline{\mathrm{\Phi }}}_{M_{\star },1}^{\star }+{\overline{\mathrm{\Phi }}}_{M_{\star },2}^{\star }\times z$. Thus, we only need to constrain the sets of {α₁, α₂, ${\overline{\mathcal{M}}}_{\star ,1}$, ${\overline{\mathcal{M}}}_{\star ,2}$, ${\overline{\mathrm{\Phi }}}_{M_{\star },1}^{\star }$, ${\overline{\mathrm{\Phi }}}_{M_{\star },2}^{\star }$} and {β₁, β₂, ℳ_B, 1, ℳ_B, 2, ${\overline{\mathrm{\Phi }}}_{L_B,1}^{\star }$, ${\overline{\mathrm{\Phi }}}_{L_B,2}^{\star }$} values for both star-forming and quiescent galaxies. In this regard, we note that this method requires no photo-z binning of the data other than the upper and lower redshift limits of the sample (see Sect. 6.2). The inclusion of the PDFs of the parameters involved is performed by the Monte Carlo realisations obtained during the SED-fitting analysis and the weights included in Eq. (12) (further details in Sect. 6.2). In practise, this is equivalent to have a sample of N_S × N_MC sources, where the weight in number of each Monte Carlo realisation is 1/N_MC. Thereby, the likelihood to be maximised is expressed as:

$$\begin{array}{cc}& ln\mathcal{L}={\displaystyle \sum _{i=1}^{N_{\mathrm{S}}\times N_{\mathrm{MC}}}}{\alpha }_1{\omega }_{i}ln{M}_{\star ,i}+{\alpha }_2{\omega }_i{z}_{i}ln{M}_{\star ,i}-\frac{{\omega }_i{M}_{\star ,i}}{{10}^{{\overline{\mathcal{M}}}_{\star ,1}+{\overline{\mathcal{M}}}_{\star ,2}{z}_i}}\hfill \\ \hfill & -\frac{1}{{log}_{10}e}[({\alpha }_1{\overline{\mathcal{M}}}_{\star ,1}+{\overline{\mathcal{M}}}_{\star ,1}){\omega }_i\hfill \\ \hfill & +({\alpha }_1{\overline{\mathcal{M}}}_{\star ,2}+{\alpha }_2{\overline{\mathcal{M}}}_{\star ,1}+{\overline{\mathcal{M}}}_{\star ,2}){\omega }_i{z}_i+{\alpha }_2{\overline{\mathcal{M}}}_{\star ,2}{\omega }_i{z}_i^2]\hfill \\ \hfill & -{\omega }_{i}ln[\mathrm{\Gamma }({\alpha }_1+{\alpha }_2{z}_i+1,\frac{M_{\star }^{\mathrm{f}}(z_i)}{{10}^{{\overline{\mathcal{M}}}_{\star ,1}+{\overline{\mathcal{M}}}_{\star ,2}{z}_i}})\hfill \\ \hfill & -\mathrm{\Gamma }({\alpha }_1+{\alpha }_2{z}_i+1,\frac{M_{\star }^{\mathrm{b}}(z_i)}{{10}^{{\overline{\mathcal{M}}}_{\star ,1}+{\overline{\mathcal{M}}}_{\star ,2}{z}_i}})],\hfill \end{array}$$(13)

where Γ is the incomplete Euler gamma function (i.e. $\mathrm{\Gamma }(s,x)={\displaystyle {\int }_x^{\infty }}{t}^{s-1}{e}^{-t}\mathrm{d}t$), and $M_{\star }^{\text{f}}(z_i)$ and $M_{\star }^{\text{b}}(z_i)$ are the lower and upper stellar mass limits (i.e. the faint and bright limits) of the sample at the redshift of the ith realisation (z_i). A similar expression is found for the likelihood of the B-band luminosity function. As described in Sect. 6.2, the faint limits $M_{\star }^{\text{f}}(z)$ and $M_B^{\text{f}}(z)$ in this research are defined according to the completeness level 𝒞 of the sample employed for the determination of the LMFs, whose expressions and values are studied in detail above in Sect. 5. We note that our estimator loses the normalisation of the LMFs during the maximisation of Eq. (13). Therefore, once the {α₁, α₂, ${\overline{\mathcal{M}}}_{\star ,1}$, ${\overline{\mathcal{M}}}_{\star ,2}$} and {β₁, β₂, ℳ_B, 1, ℳ_B, 2} are constrained, the LMF normalisations can be recovered from the number density of sources in the sample (details in Sect. 6.3).

For the maximisation and sampling of the posterior distribution of the likelihood detailed in Eq. (13), we use the Python implementation of the affine-invariant ensemble sampler for Markov chain Monte Carlo code known as emcee¹⁸ (Foreman-Mackey et al. 2013). This emcee analysis is performed with 100 walkers (i.e. 25 walkers per degree of freedom) and 2000 steps per walker for both the ‘burn-in’ phase and sampling the posterior distribution. We check that there is negligible or small impact in the posterior distribution retrieved during the emcee analysis by the use of a larger number of steps. In fact, we find that after a few dozen steps, emcee promptly converges to a stable solution that maximises Eq. (13) and that the posterior distribution can be properly sampled with only ∼20 walkers and after ∼500 steps. It is worth mentioning that a too large number of steps rapidly increases the computational time of this emcee analysis, which may be especially long owing to the calculation of the Γ terms for each of the realisations of each of the sources. In view of the future J-PAS project, we may therefore reduce the number of walkers and steps, where this decision would be based on a consensus to determine the LMFs in a short enough period of time without renouncing to the precision of the LMF parameters.

6.2. Weights and priors

In this research, the weights ω_i in Eq. (13) are one of the key parameters for the determination of the miniJPAS LMFs. These weights can account for, among others, the probability of each of the sources for being a galaxy, the spectral classification of galaxies, the completeness of the sample at certain redshift and stellar mass/absolute magnitude, how relevant is a source owing to the quality of the determination of the parameters involved, and even the statistical number of galaxies in the sample, which is also needed for the determination of the LMF normalisation (see Sect. 6.3). Moreover, as we base the determination of the miniJPAS LMFs on the cumulative probability of the sample at 0.05 ≤ z ≤ 0.7 and r_SDSS ≤ 22, these constraints on the parent sample can be easily included in a statistically rigorous manner through the mentioned weights.

In the following, we distinguish two kinds of weights with different roles. Hence, the weight ω_i included in Eqs. (12) and (13) can be defined as:

$$\begin{array}{c}\hfill {\omega }_i={\kappa }_i\times \frac{{\tau }_i}{⟨\tau ⟩},\end{array}$$(14)

where κ_i accounts for all the terms that contribute to the number of galaxies in the sample and τ_i for any kind of statistical significance meaning that some sources or realisations should have a higher consideration during the maximisation of Eq. (13). As in Zucca et al. (1994), we introduce the median of the distribution of the τ_i weights, ⟨τ⟩, to avoid artificially increasing the error contours of the LMF parameters during the emcee analysis. In fact, τ_i/⟨τ⟩ can act as a Bayesian prior during the sampling of the posterior distribution.

Accordingly, κ_i includes the probability of a source to be a galaxy (P_G) and those sources with null probability to be a galaxy are removed from the analysis. Secondly, κ_i must also reflect the number of galaxies for the different spectral-types. For this aim, we include the P_Q, i and P_SF, i estimated for each of the Monte Carlo realisations (1 or 0 depending whether the realisation is above or below Eq. (1), respectively; see Sect. 4). Hence, galaxies with P_Q values for the spectral classification different than unity or zero partly contribute to both star-forming and quiescent LMFs with those realisations that are compatible with these galaxy spectral-types. Thirdly, we also include the completeness level of each of the realisations. This allows us to extend the analysis below the 𝒞 = 0.95 completeness level, thus making the most of our miniJPAS sample. For instance, if one of the realisations correspond to an observation for a completeness level of 𝒞 = 0.5, this would mean that only half of the sources were observed or included in the sample, and therefore, the LMFs must account for the double of sources at that redshift and stellar mass/absolute magnitude. Fourthly, as we are using the N_MC Monte Carlo realisations in our analysis, we have to take this into consideration to avoid overestimating the number of sources in the sample. Finally, as a consequence of using a flux-limited sample down to $r_{\text{SDSS}}^{\text{lim}}$ = 22 (i.e. r_SDSS ≤ 22), we account for border effects owing to uncertainties in our r_SDSS detection band. In this regard, we also include fainter sources (down to r_SDSS = 22.5) than this magnitude limit with uncertainties that make them compatible with r_SDSS ≤ 22 and we weight its probability (P_r) according to the following equation:

$$\begin{array}{c}\hfill P_{r,i}=\frac{1}{2}-\frac{1}{2}\mathrm{erf}\left\{\frac{2.5}{\sqrt{2}ln10{\sigma }_{r,i}}[{10}^{0.4(r_{\mathrm{SDSS},i}-r_{\mathrm{SDSS}}^{\lim })}-1]\right\},\end{array}$$(15)

where σ_r, i are the magnitude uncertainty included in the miniJPAS photometric catalogue, r_SDSS its magnitude, and erf is the error function, which is defined as $\text{erf}(x)=2 {\pi }^{-1/2} {\displaystyle {\int }_0^x}{e}^{-t^2} \text{d}t$. Considering all these factors, the κ_i weight used for deriving the stellar mass functions of miniJPAS quiescent galaxies is the following:

$$\begin{array}{c}\hfill {\kappa }_i=\frac{P_{\mathrm{Q},i}\times P_{\mathrm{G},i}\times P_{r,i}}{N_{\mathrm{MC},i}\times \mathcal{C}(z_i,M_{\star ,i})}.\end{array}$$(16)

For the stellar mass functions of star-forming galaxies, κ_i is the same than in Eq. (16), except for changing the P_Q, i for P_SF, i. Furthermore, we note that for the luminosity function case, the completeness of the B-band absolute magnitude must be used instead of the stellar mass completeness one, meaning that 𝒞(z_i, M_B, i) (see Eq. (5)). As N_MC, i is the same number for all the sources in the sample, we note that the number of realisations is not relevant for the maximisation of Eq. (13), but it is needed for the LMF normalisation (see Sect. 6.3).

As we mention in Sect. 2.3, odds is a parameter that correlates with the quality of the PDZs of each of the miniJPAS sources and is sensitive to outliers. In fact, Hernán-Caballero et al. (2021, 2023) showed that J-PAS-like sources with low odds values show larger discrepancies with respect to spectroscopic redshifts and more frequently present larger values than these ones. However, we cannot restrict our sample to only sources with high odds values because this would bias our sample and the LMFs are sensitive to the total number density of galaxies. In order to circumvent this issue, we also include the odds parameter as an extra weight during the maximisation of Eq. (13). We find that including odds values in the weight budget improves the constraints on the faint-end slopes of miniJPAS LMFs. This is a consequence of the fact that faint sources in the sample exhibit lower odds in average, and hence overestimated photo-z, that induces higher α and β values. It is also of note that we do not adopt the photo-z probability associated to PDZs, since the JPHOTOZ package used to constrain these probability distributions includes priors based on distributions of the i band from deeper surveys such as VVDS, meaning that its inclusion may require our LMF results to be similar to the VVDS luminosity functions. In addition and owing to the depth of miniJPAS, the faint-end slopes of LMFs are mainly constrained by galaxies at lower redshifts (i.e. z ≲ 0.25), whereas the bulk of the galaxy sample is concentrated around z ∼ 0.5. This means that the maximisation of Eq. (13), and hence α and β, is dominated in number by sources with stellar masses and luminosities above the LMF faint-end limits (i.e. around ℳ_⋆(z) and ℳ_B(z)). In order to mitigate this fact and balance the significance of the nearby galaxies for the determination of the LMF faint-end slopes, we increase the weight of nearby galaxies for the maximisation of Eq. (13) by accounting for the comoving volume (V_C). Therefore, for this work we only account for odds and comoving volume in the τ_i weight, that is, ${\tau }_i={odd{s}_i}^2\times {(\frac{{\mathrm{d}}^2{V}_{\mathrm{C}}(z_i)}{\mathrm{d}z\mathrm{d}\mathrm{\Omega }})}^{-1}$ and ⟨τ⟩ the median of the τ_i distribution of the sources with non-null ω_i weights.

Regarding other priors and constraints, we set lower limits for the stellar mass and luminosity completeness levels that are going to be included in the maximisation of Eq. (13), which directly affects to the definition of $M_{\star ,i}^{\text{f}}$ and $M_{B,i}^{\text{f}}$. This implies that the weight of a realisation with a stellar mass or absolute magnitude under this completeness limit is null (i.e. ω_i = 0) and can be removed from the analysis process. Regardless of the galaxy spectral-type, we assume that the minimum completeness limit imposed for the determination of the B-band luminosity function is 𝒞 = 0.25, whereas for the stellar mass case is 𝒞 = 0.05. Consequently, $M_{B,i}^{\text{f}}=M_B^{\mathcal{C}=0.25}(z_i)$ and ${\log }_{10}{M}_{\star ,i}^{\text{f}}={\log }_{10}{M}_{\star }^{\mathcal{C}=0.05}(z_i)$ (see Eqs. (4) and (5)). The former limit is set to 𝒞 = 0.25 because this helps to prevent overestimation of the α and β values since the miniJPAS sources with the lowest B-band luminosities also exhibit the lowest odds values. This is not the case for stellar mass functions since low stellar mass galaxies do not preferentially show low odds values, but they are more homogeneously distributed and weighted accordingly (see τ_i above). Furthermore, we find that the assumption of a lower limit for completeness prevents the inclusion of some outliers without significantly altering the results. In the same sense, as the Fermi-Dirac function used for the parametrisation of the completeness (see Eqs. (2) and (3)) is a divergent function, we assume that the sample is complete (i.e. 𝒞(z_i, M_⋆, i) = 1 in Eq. (16)) for values above 𝒞 = 0.975.

During the emcee analysis, we also impose limits to characteristic parameters and the faint-end slopes of the LMFs. In particular, we impose that within our redshift limits these parameters must satisfy the following: −4 ≤ α(z)≤4, −4 ≤ β(z)≤4, 9.5 ≤ log₁₀ℳ_⋆(z)≤12.5, and −24 ≤ ℳ_B(z)≤ − 17. In addition, the assumption of fixed or constant α and β values is frequently used and based on previous studies in deeper surveys, where these parameters show little variation with redshift at least since z = 0.7 (see e.g. Faber et al. 2007; López-Sanjuan et al. 2017). Thus, for decreasing the uncertainties in the parameters involved, in this work we opted for assuming α₂ = β₂ = 0 for the stellar mass function of star-forming galaxies and both B-band luminosity functions. Hence, there is no constraint on α₂ for the stellar mass function of quiescent galaxies. This is motivated by the fact that we are using single-Schechter functions for the LMF parametrisation and the number density of low-mass quiescent galaxies rapidly increases at decreasing redshift, while the number of massive quiescent galaxies does not change significantly (see e.g. Díaz-García et al. 2019b, and references therein).

As a result of all these sample restrictions and probabilities, the effective number of galaxies down to r_SDSS = 22 (i.e. ${\displaystyle \sum _i}\frac{P_{\mathrm{G},i}\times P_{r,i}}{N_{\mathrm{MC},i}}$) employed in this work for determining the LMFs of miniJPAS galaxies at 0.05 ≤ z ≤ 0.7 amounts to ∼5000 galaxies.

6.3. Normalisation of the LMFs

Even though the normalisation is lost during the maximisation of the likelihood in our method, this one can be retrieved by the characteristic parameters, faint-end slopes, and the number densities of our galaxy sample. By definition, the LMFs describes the number density of galaxies per unit of comoving volume and luminosity or stellar mass. In other words, if we integrate the LMFs in a certain comoving volume, we must obtain the number density of galaxies in it. With this in mind, the normalisation of our stellar mass function at a redshift bin of width Δz′ and centred at z′ is:

$$\begin{array}{c}\hfill {\mathrm{\Phi }}_{M_{\star }}^{\star }(z^{\prime },\mathrm{\Delta }{z}^{\prime })=\frac{{\sum }_{i=1}^{N_{\mathrm{S}}\times N_{\mathrm{MC}}}{\kappa }_{i}W(z^{\prime }-z_i)}{{\int }_{\mathrm{\Omega }}{\displaystyle {\int }_{z^{\prime }-\mathrm{\Delta }{z}^{\prime }/2}^{z^{\prime }+\mathrm{\Delta }{z}^{\prime }/2}}\frac{{\mathrm{d}}^2{V}_{\mathrm{C}}(z)}{\mathrm{d}z\mathrm{d}\mathrm{\Omega }}\mathrm{\Delta }\mathrm{\Gamma }(z)\mathrm{d}z\mathrm{d}\mathrm{\Omega }},\end{array}$$(17)

where Ω is the solid angle or the survey footprint area (Ω = 0.895 deg² for miniJPAS), V_C is the comoving volume, and ΔΓ is the function that results of integrating the stellar mass function in the stellar range of our sample, namely:

$$\begin{array}{c}\hfill \mathrm{\Delta }\mathrm{\Gamma }(z)=\mathrm{\Gamma }(\alpha (z)+1,\frac{M_{\star }^{\mathrm{f}}(z)}{{\mathcal{M}}_{\star }(z)})-\mathrm{\Gamma }(\alpha (z)+1,\frac{M_{\star }^{\mathrm{b}}(z)}{{\mathcal{M}}_{\star }(z)}).\end{array}$$(18)

W is the window function for selecting galaxies or realisations in a redshift bin:

$$\begin{array}{c}\hfill W(z^{\prime }-z)=\{\begin{array}{ccc}\hfill 1& \mathrm{if}& -\frac{\mathrm{\Delta }{z}^{\prime }}{2}\le z^{\prime }-z\le \frac{\mathrm{\Delta }{z}^{\prime }}{2}\hfill \\ \hfill 0& \mathrm{otherwise}\end{array}.\end{array}$$(19)

The normalisation of the B-band luminosity functions in this work is as Eqs. (17) and (18), but including β(z) and ℒ_B instead. We note that the normalisation of the LMFs as a function of the galaxy spectral-type is implicit in the κ_i term.

However, if we need the normalisation of the LMFs at a given redshift, this expression can result hard to compute. To make it easier to obtain, we compute Eq. (17) in 20 redshift bins of equal width, ranging from z = 0.05 to 0.7, to subsequently fit this set of values to functions of the form ${log}_{10}{\mathrm{\Phi }}_{M_{\star }}^{\star }={\overline{\mathrm{\Phi }}}_{M_{\star },1}^{\star }+{\overline{\mathrm{\Phi }}}_{M_{\star },2}^{\star }\times z$ and ${log}_{10}{\mathrm{\Phi }}_{L_B}^{\star }={\overline{\mathrm{\Phi }}}_{L_B,1}^{\star }+{\overline{\mathrm{\Phi }}}_{L_B,2}^{\star }\times z$. We explicitly checked that a linear function of this kind is enough to fairly reproduce the evolution of the normalisation with redshift since z = 0.7, although quadratic and cubic equations can be also used for the fitting in a general case.

6.4. Eddington bias

As a consequence of uncertainties and the exponential-like shape of galaxy LMFs for values greater than the characteristic stellar mass or luminosity, these functions can be subject to the so-called Eddington bias (Eddington 1913; Caputi et al. 2011; Grazian et al. 2015). This effect can be especially relevant at the high mass and/or bright-end of the LMFs and mainly results in an overestimation of the number of massive or bright galaxies (see e.g. Ilbert et al. 2013; Obreschkow et al. 2018). If we do not account for this effect, other parameters resulting from the LMFs (e.g. the stellar mass and luminosity densities) would be affected in different ways. In our case, the dependency of the errors on redshift and magnitude (see Sect. 7.1 and Fig. 5) makes that the cosmic evolution of the LMFs may be also affected by Eddington bias, and therefore, the cosmic evolution of the stellar mass and luminosity densities as well. The Eddington bias is mathematically expressed as the convolution of the error distribution with the proper distribution of interest. The LMFs determined following the methodology detailed in this paper are actually bidimensional density distributions in the redshift and stellar mass/luminosity space. Consequently, the Eddington bias effect is included in the LMFs according to

$$\begin{array}{c}\hfill {\mathrm{\Phi }}_{\mathrm{Edd}}({\overline{M}}_{\star },z)={\displaystyle {\int }_{-\infty }^{+\infty }}{\displaystyle {\int }_{-\infty }^{+\infty }}\mathrm{\Phi }(x,y)G({\overline{M}}_{\star }-x,z-y)\mathrm{d}x\mathrm{d}y,\end{array}$$(20)

Fig. 5. Normalised distributions of errors of B-band luminosity (top panels) and stellar mass (bottom panels) at different rSDSS magnitude bins (see insets) for star-forming and quiescent galaxies (left and right panels, respectively) in miniJPAS.

where G is the bivariate normal distribution of errors that includes the photo-z uncertainties, the logarithm stellar mass or B-band absolute magnitude uncertainties, and parameter correlations.

We correct for the Eddington bias effect using the average values of parameters obtained during the SED-fitting analysis. For this aim, we split our galaxy sample in bins of redshift, stellar mass, and absolute magnitude for each of the galaxy spectral-types. For each of these bins, we computed the median values and average uncertainties, as well as the correlations between redshift and the other two stellar population parameters. From this map or distribution of values, we constructed the bivariate normal distribution of errors included in Eq. (20) as a function of redshift, stellar mass, and absolute magnitude. For each of the iterations performed during the sampling of the posterior, we find the LMF functions Φ(x, y) that minimise the difference |log₁₀Φ′(M_⋆, z)−log₁₀Φ_Edd(M_⋆, z)| and |log₁₀Φ′(M_B, z)−log₁₀Φ_Edd(M_B, z)|, where Φ′(M_⋆, z) and Φ′(M_B, z) are the LMFs obtained from the maximisation of Eq. (13) (Sects. 6.1 and 6.2) and its subsequently normalization (Sect. 6.3). Therefore, the Eddington bias correction is performed after the emcee analysis and for each of the realisations sampling the posterior distribution.

7. Results

The combination of all the methodologies and analysis techniques detailed in previous sections, which naturally include a proper treatment of the stellar mass and luminosity uncertainties of miniJPAS galaxies (details in Sect. 7.1), result in the LMFs of the miniJPAS galaxies (Sect. 7.2). Furthermore, we can easily take advantage of these functions for additionally deriving the cosmic evolution of the stellar mass and luminosity densities (Sect. 7.3) and the fraction of quiescent galaxies at different redshift, stellar mass, and luminosity (Sect. 7.4).

7.1. Stellar mass and B-band absolute magnitude uncertainties of miniJPAS galaxies

We find that uncertainties of stellar mass and B-band luminosity mainly correlate with magnitude (see Fig. 5) and, by extension, with redshift. For our final sample, uncertainties in the B-band luminosity (σ_MB, see upper panels in Fig. 5) typically cover values of σ_MB < 0.5, whereas for stellar mass this range increases up to σ_log10M⋆ < 0.7 dex. Sources with the lowest luminosity uncertainties are usually brighter and appear more frequently at lower redshifts, while the more noisy estimations of the B-band luminosity are only present at fainter magnitudes (see top panels in Fig. 5). Interestingly, there is a significant fraction of sources exhibiting typical values of σ_MB ≲ 0.15 ranging from 95% for magnitudes r_SDSS ≤ 20 to 50% for 21.5 ≤ r_SDSS ≤ 22. We find that these fractions are ∼5% higher for quiescent galaxies.

Concerning stellar mass uncertainties (σ_log10M⋆, see bottom panels in Fig. 5), these present greater correlations with magnitude than luminosity uncertainties and also rely on the galaxy spectral-type. In this regard, the stellar mass uncertainties of miniJPAS galaxies are systematically smaller for quiescent than for star-forming galaxies (right- and left-bottom panels in Fig. 5, respectively). More precisely, the ∼85% fraction of all the quiescent galaxies present stellar mass uncertainties of σ_log10M⋆ ≲ 0.2 dex down to magnitude r_SDSS = 21, which can be treated as an upper limit for the stellar mass precision of the luminous red galaxies from miniJPAS. However, only a 70% fraction of star-forming galaxies show uncertainties lower than σ_log10M⋆ = 0.3 dex for magnitudes brighter than r_SDSS = 21. At fainter magnitudes, r_SDSS > 21, the stellar mass uncertainties of galaxies from miniJPAS rapidly increases, especially for the star-forming case. At this faint regime, the stellar mass uncertainties of star-forming galaxies ranges from σ_log10M⋆ ∼ 0.25 to 0.55 dex, whereas for quiescent galaxies is from σ_log10M⋆ ∼ 0.1 to 0.4 dex.

This highlights that stellar masses are in average poorly constrained than luminosity, especially at fainter magnitudes. This is a consequence of stellar mass relying on both the brightness-distance and mass-luminosity relation, which results more uncertain than only distance or photo-z. In fact, these results along with the complex PDZ exhibited by some miniJPAS sources motivate the use of the discretised PDF for constraining the miniJPAS LMFs, all this for making the most possible of our sample at fainter magnitudes.

7.2. Stellar mass and luminosity functions up to z = 0.7

In the redshift range explored in this work, the values of the characteristic parameters, the low-mass and faint-end slopes, and normalisations of the miniJPAS LMFs (see Table 3) present median values and posterior distributions (see Figs. B.1 and B.2) that are properly sampled within the range of values imposed or priors (i.e. −4 ≤ α(z)≤4, −4 ≤ β(z)≤4, 9.5 ≤ log₁₀ℳ_⋆(z)≤12.5, and −24 ≤ ℳ_B(z)≤ − 17). In general, even though the assumption of α₂ = β₂ = 0, the parameters describing our LMFs do not present a significant evolution. However, we find that the evolution of these parameters depend on the galaxy spectral-type, meaning that galaxies may follow different evolutionary paths according to their spectral-type.

Secondly, it is well-known that during the process for obtaining the Schechter parameters, degeneracies and correlations amongst the parameters might arise. For instance, there is a strong degeneracy between the faint-end slope and the characteristic parameters. We find that these degeneracies are extended also to the coefficients describing the evolution of these parameters with redshift as shown by the posterior distributions (see Figs. B.1 and B.2). In this regard, more massive or brighter characteristic parameters lead to lower low-mass or faint-end slopes and viceversa. As a consequence of the direct correlation between the normalisations and the rest of Schechter parameters, see Eq. (17), higher characteristic parameters involve lower normalization values. We also find that these degeneracies are independent of the galaxy spectral-type.

The first aspect that we observe is that the B-band luminosity functions of miniJPAS galaxies present different shapes and evolution with redshift depending on the galaxy spectral-type (see Fig. 6). Thus, the slope of the B-band faint-end (i.e. β₁ + 1) of star-forming galaxies is negative down to z = 0.7, while the one for quiescent galaxies is positive (see left and right panels in Fig. 6, respectively). Regarding the bright-end of the luminosity functions, we find that star-forming and quiescent galaxies present similar values. At decreasing redshift, the number density of both bright star-forming and quiescent galaxies decreases. However, the number of faint quiescent galaxies (M_B ≳ −21) is increasing since z = 0.7, as it is reflected in the evolution of the faint-end of its B-band luminosity function. On the other hand, the evolution in number of faint star-forming galaxies exhibit negligible evolution, or at least is compatible with the uncertainties of this function. As a consequence of the low uncertainties in the determination of the B-band absolute magnitudes of the brightest miniJPAS galaxies (see Fig. 5), the effect of the Eddington bias is negligible for the miniJPAS luminosity functions (see Fig. 6).

Fig. 6. Evolution with redshift of the B-band luminosity functions of star-forming and quiescent galaxies (solid lines in left and right panels, respectively) from miniJPAS. The shaded areas exhibit the uncertainties in the luminosity functions obtained with our methodology. The dash-dot lines show the B-band luminosity functions before correcting of Eddington bias effects. The dashed and dotted lines point the absolute magnitude range in which our sample is affected by incompleteness after and before correcting for Eddington bias effects, respectively.

The stellar mass functions of star-forming and quiescent galaxies also show differences in shape and evolution depending on the spectral-type (see Fig. 7). In addition, these functions show larger evolution than the luminosity ones, especially for quiescent galaxies. The high-mass end of star-forming galaxies show evolution since z = 0.7, which points out that the most massive galaxies, both star-forming and mainly quiescent, are present at lower redshifts in our sample. For quiescent galaxies at z < 0.7, we find that the low-mass end evolves the most, whereas the number density of the high-mass end presents a lesser evolution with redshift. Thus, our results points out that there is a more significant increase in the number of intermediate and low mass quiescent galaxies (i.e. log₁₀M_⋆ < 10.7 dex). The effect of Eddington bias is more remarkable for the stellar mass functions than for the luminosity functions. This is expected since the uncertainties of the stellar masses are greater in our sample (see Fig. 5) and mainly affecting the high-mass end of these functions. In particular, the Eddington bias seems to be more relevant for the star-forming case (see Fig. 7). In fact, without the Eddington bias correction, the stellar masses of star-forming galaxies would be largely overestimated in the high-mass end. For quiescent galaxies, the Eddington bias effect is less prominent because the stellar mass uncertainties of quiescent galaxies are in average lower than for the star-forming ones. However, we still find necessary to perform the Eddington bias correction described in Sect. 6.4. Otherwise, the high-mass end of the stellar mass function of quiescent galaxies would agree with a null evolution since z = 0.7. In general, the evolution of the LMFs show that the number density of quiescent galaxies increases at decreasing redshift, but this evolution in number is mainly due to the evolution in number of less massive galaxies, which confirms the results obtained in previous studies involving the evolution of the number density of quiescent galaxies with redshift (see e.g. Díaz-García et al. 2019a, and references therein).

Fig. 7. Same as Fig. 6, but for stellar mass functions.

7.3. Evolution of the stellar mass and luminosity densities since z = 0.7

The stellar mass and luminosity densities (j_M and j_B, respectively) can be directly derived from the LMF parameters as following:

$$\begin{array}{cc}& j_B(z)=\int L_B\mathrm{\Phi }(L_B,z)\mathrm{d}{L}_B\hfill \\ \hfill & ={\mathrm{\Phi }}_{L_B}^{\star }(z){10}^{0.4(M_{B,\odot }-{\mathcal{M}}_B(z))}\mathrm{\Gamma }(\beta (z)+2),\hfill \end{array}$$(21)

$$\begin{array}{cc}\hfill j_M(z)& =\int M_{\star }\mathrm{\Phi }(M_{\star },z)\mathrm{d}{M}_{\star }\hfill \\ \hfill & ={\mathrm{\Phi }}_{M_{\star }}^{\star }(z){\mathcal{M}}_{\star }(z)\mathrm{\Gamma }(\alpha (z)+2),\hfill \end{array}$$(22)

where Γ is the Euler gamma function and M_B, ⊙ is the B-band absolute magnitude of the Sun, for which we adopt a value of M_B, ⊙ = 5.38 (Binney & Merrifield 1998).

As a result, the cosmic evolution of all the miniJPAS galaxies indicates that the global evolution of the B-band density decreases at decreasing redshift (see black line in left panel of Fig. 8). In this sense, the global evolution of the total sample goes from a value of 8.4 to 8.2 dex. Therefore, the luminosity of galaxies per unit of volume is greater at higher redshift. However, this trend is different when we take the spectral-type of galaxies into account. For quiescent galaxies, the B-band luminosity density remains roughly constant with a value of 7.7 dex (see red line in left panel of Fig. 8). Consequently, the decrease of the luminosity density is mainly lead by star-forming galaxies. For these kinds of galaxies, the luminosity density decreases from 8.2 to 8.0 dex (see blue line in left panel of Fig. 8). Consequently, we find that, except for quiescent galaxies, there is a global decrease of 0.2 dex.

Fig. 8. Cosmic evolution of the B-band luminosity and stellar mass densities (left and right panels, respectively) for quiescent, star-forming, and full miniJPAS samples (red, blue, and black lines, respectively). The shaded areas illustrate the uncertainties of these densities. Left panel: B-band absolute magnitude densities from 2dFGRS (cross markers, Madgwick et al. 2002), 2MASS/SDSS (star markers, Bell et al. 2003), COMBO-17/DEEP2 (squares and diamonds, respectively; Faber et al. 2007), AGES (pentagons, Cool et al. 2012), GAMA (plus marker, Loveday et al. 2012), VIPERS (triangle up, Fritz et al. 2014), Boötes field (triangle left, Beare et al. 2015), and ALHAMBRA (circle, López-Sanjuan et al. 2017) surveys are included for comparison. Right panel: stellar mass densities from 2MASS/SDSS (star markers, Bell et al. 2003), S-COSMOS (circles, Ilbert et al. 2010), GAMA (plus markers, Baldry et al. 2012), UltraVISTA (squares, Ilbert et al. 2013), SDSS/GALEX (crosses, Moustakas et al. 2013), COSMOS2015 (diamonds, Davidzon et al. 2017), GAMA/G10-COSMOS/3D-HST (triangle right, Wright et al. 2018), and DEVILS (hexagons, Thorne et al. 2021).

The cosmic evolution of the stellar mass density exhibits larger differences with respect to the luminosity case. The stellar mass density of the global sample of miniJPAS galaxies is compatible with a mild evolution of ∼0.3 dex from z = 0.7 to z = 0 (see black line in right panel of Fig. 8). The global value at z = 0.7 is ∼8.2 dex and 8.5 dex at the very nearby Universe. Unlike the B-band luminosity density, the evolution of the stellar mass density is mainly driven by quiescent galaxies. There is a strong evolution in stellar mass density of quiescent galaxies since z = 0.7 with a value of 7.8 dex that evolves up to 8.3 dex at z = 0 (see red line in right panel of Fig. 8). On the other hand, the stellar mass density related to star-forming galaxies slightly increases at decreasing redshift (see blue line in right panel of Fig. 8), more precisely, from 7.9 to 8.0 dex. Interestingly, the stellar mass density of both spectral-types present a similar value of 8.0 dex at z = 0.4. Therefore, there is a larger fraction of stellar mass in quiescent galaxies at z = 0, whereas the major part of stars reside in star-forming galaxies at the highest redshift explored in this work. This is mainly related to the fact that the high-mass end of the stellar mass function of quiescent galaxies is shifted to higher values than for star-forming galaxies and, as a consequence, the most massive galaxies at the nearby Universe mostly belong to the quiescent galaxy population.

7.4. Fraction of quiescent galaxies since z = 0.7

The fraction of quiescent galaxies, f_r, as a function of the B-band absolute magnitude and stellar mass at different redshift can be directly inferred from the LMFs by

$$\begin{array}{c}\hfill f_{\mathrm{r}}(M_B,z)=\frac{{\mathrm{\Phi }}_{\mathrm{Q}}(M_B,z)}{{\mathrm{\Phi }}_{\mathrm{Q}}(M_B,z)+{\mathrm{\Phi }}_{\mathrm{SF}}(M_B,z)},\end{array}$$(23)

$$\begin{array}{c}\hfill f_{\mathrm{r}}(M_{\star },z)=\frac{{\mathrm{\Phi }}_{\mathrm{Q}}(M_{\star },z)}{{\mathrm{\Phi }}_{\mathrm{Q}}(M_{\star },z)+{\mathrm{\Phi }}_{\mathrm{SF}}(M_{\star },z)},\end{array}$$(24)

where Φ_Q and Φ_SF are the LMFs of quiescent and star-forming galaxies, respectively.

Interestingly, the brightest galaxies (M_B < −22) are mainly classified as star-forming galaxies at z < 0.7 (see left panel in Fig. 9). However, for B-band absolute magnitudes fainter than M_B ∼ −23.5, the fraction of quiescent galaxies starts to increase until reaching a maximum value of f_r ∼ 0.5, where the position of this maximum fraction depends on redshift. Hence the f_r peak moves towards fainter magnitudes at decreasing redshift, more precisely, from M_B ∼ −23 to −21 in the redshift range explored in this work. For fainter magnitudes than the peak previously mentioned, the fraction of quiescent galaxies abruptly decreases reaching values of f_r < 0.2. It is worth mentioning that for the faint regime (M_B ≳ −19), which (in turn) exhibits f_r values below ∼0.2, the fraction of red galaxies are in fact an extrapolation of the B-band luminosity functions owing to the completeness lower limit imposed for the maximisation of the luminosity function likelihoods (i.e. 𝒞 = 0.25). In fact, this magnitude regime is usually reproduced by the second component of a double Schechter function.

Fig. 9. Fraction of quiescent galaxies, fr, as a function of the B-band absolute magnitude (left panel) and stellar mass (right panel) at different redshift (see inset). The dashed lines illustrate the absolute magnitude and stellar mass range for which our quiescent sample is under a completeness level of 𝒞 = 0.25 and 0.05, respectively, meaning that the fraction of red galaxies are obtained through a extrapolation of our results.

The fraction of red galaxies according to the stellar mass of galaxies also reveals interesting results (see right panel in Fig. 9). On the one hand, the majority of galaxies with stellar mass above log₁₀M_⋆ ≳ 10.7 are actually quiescent galaxies at z < 0.7, where the fraction of red galaxies is above 60% in all cases. Moreover, the fraction of quiescent galaxies with stellar masses in the range log₁₀M_⋆ < 10.7 becomes larger at decreasing redshift, meaning that the number of quiescent galaxies is increasing in importance in more recent epochs. Even though this increase in number, the galaxy sample is largely dominated in number by star-forming galaxies for log₁₀M_⋆ < 10.5. It is also of note that the values of f_r for log₁₀M_⋆ ≲ 9.5 should be treated as f_r lower limits, since we are currently using single Schechter functions for constraining the stellar mass functions of quiescent galaxies owing to the limitations of our sample.

8. Discussion

8.1. Comparisons with previous works

During the last decades, the LMFs of galaxies have been widely studied making use of spectroscopic and/or photometric data including galaxies from very early cosmological epochs (z ∼ 7) to the nearby Universe. This allows for an easy comparison of our results with the LMFs obtained in previous works within our redshift range (0.05 ≤ z ≤ 0.7). For a fair comparison, we compare our results with the fitted or inferred LMFs (i.e. adopting a parametric LMF expressed by single or double Schechter functions) obtained by other authors where the Eddington bias effects are corrected and/or the errors due to photon-noise uncertainties are not big enough (≲0.2 dex) as to have a big impact on the shape of the LMFs.

The B-band luminosity functions of the miniJPAS galaxies show number densities similar to those obtained in previous researches (see Fig. 10). In general, we do not find large discrepancies among luminosity functions and both the faint and bright ends of these functions present similar values than in previous works. This conclusion can be extended for quiescent and star-forming galaxies, and hence for the full sample (see top, middle, and bottom panels in Fig. 10, respectively). Interestingly, the miniJPAS B-band luminosity functions are especially similar to those obtained for the ALHAMBRA survey (see circles in Fig. 10 and López-Sanjuan et al. 2017), which were obtained by a statistical Bayesian method with similarities to ours. Except for the ALHAMBRA case, the number density of quiescent galaxies in miniJPAS with absolute magnitudes fainter than M_B ∼ −21 is ∼0.5 dex lower than for the rest of surveys at z ≳ 0.5 (see middle panel in Fig. 10; Faber et al. 2007; Drory et al. 2009; Cool et al. 2012; Beare et al. 2015). However, this discrepancy can be explained as a consequence of the use of the MCDE diagram for the galaxy classification. This is because this diagram accounts for extinction effects and commonly yields a lower number of low mass quiescent galaxies in comparison to other diagrams since it removes dusty star-forming galaxies from the quiescent sample (for a further analysis of this effect, see Díaz-García et al. 2019a). We note that for the ALHAMBRA luminosity functions, as well as its galaxy classification, the authors accounted for extinction and maybe for this reason these luminosity functions are particularly similar to the miniJPAS ones. Finally, we find that the parametric luminosity functions from the VIPERS survey (Fritz et al. 2014) are the most different functions with respect to the miniJPAS ones, but also to the rest of cases included here. It is also of note that some works pointed to the need for a second quiescent component at M_B ≳ −18. However, as explained in Sect. 6, our miniJPAS LMFs are described according to a single Schechter function since it is largely difficult to fit a second quiescent component with the miniJPAS dataset, and therefore, we cannot confirm this fact.

Fig. 10. Comparison of the B-band luminosity function of miniJPAS galaxies (coloured lines from top to bottom: star-forming, quiescent, and full galaxy sample) with those fitted or inferred (Schechter-like functions) from previous works at different redshift. The shaded area illustrates uncertainties in our results. The dashed lines show the B-band luminosity functions for magnitudes fainter than the completeness lower limit of the sample employed for their determination (i.e. $M_B>M_B^{\mathcal{C}=0.25}(z)$). B-band luminosity functions from COMBO-17/DEEP2 (squares and diamonds, respectively; Faber et al. 2007), COSMOS (triangle right, Drory et al. 2009), AGES (pentagons, Cool et al. 2012), VIPERS (triangle up, Fritz et al. 2014), Boötes field (triangle left, Beare et al. 2015), and ALHAMBRA (circles, López-Sanjuan et al. 2017) surveys are included for comparison.

As revealed by many previous works (see e.g. Bell et al. 2003; Bruzual & Charlot 2003; Vazdekis et al. 2012; Ferré-Mateu et al. 2013; Díaz-García et al. 2019a,b), the SSP model set along with the IMF adopted for building stellar population synthesis models can introduce a bias or shift in the stellar mass measurements obtained by a stellar population analysis (e.g. SED-fitting). For instance, we checked that quiescent and star-forming galaxies from miniJPAS are more massive by 0.04 and 0.11 dex, respectively, when Bruzual & Charlot (2003) SSP models (hereafter, BC03) are used instead of the CB17 ones for the MUFFIT SED-fitting analysis, both with the same Chabrier (2003) IMF. Díaz-García et al. (2019a) also found discrepancies between EMILES SSP models (Extended Medium-resolution Isaac Newton Telescope of Empirical Spectra; Vazdekis et al. 2016) with a universal Kroupa (2001) IMF for both Padova00 and BaSTI isochrones (Girardi et al. 2000; Pietrinferni et al. 2004, respectively) and the BC03 models with a Chabrier (2003) IMF, where EMILES stellar masses are in average higher by about 0.11 and 0.15 dex, respectively. Bell et al. (2003) found that the Salpeter (1955) IMF assumption systematically yields higher stellar masses by about 0.4 and 0.15 dex than Kroupa-like and its ‘diet’ Salpeter (1955) IMFs, respectively. In this regard and independently of the galaxy spectral-type, the stellar masses of miniJPAS galaxies determined using MUFFIT with BC03 SSP models and a Salpeter (1955) IMF are systematically 0.24 dex more massive than using the same model set with a Chabrier (2003) IMF instead. Whenever possible and with the only goal of palliating discrepancies amongst stellar masses owing to model assumptions; thus performing a proper comparison of the stellar mass function resulting from previous studies, we present all the stellar mass functions to a common stellar mass framework that is scaled according to the Bruzual & Charlot (2003) model set and a Chabrier (2003) IMF by adding the aforementioned offsets to the characteristic stellar mass (see Fig. 11). Similarly to the luminosity function case, the stellar mass functions of miniJPAS galaxies exhibit similar trends and values to those obtained in previous studies down to z = 0.7. However, there is evidence of an excess of massive star-forming galaxies in the miniJPAS at low redshift (z ∼ 0.05) with respect to results from 2MASS/SDSS and GAMA suveys (see top-left panel in Fig. 11; Bell et al. 2003; Baldry et al. 2012, respectively). Similarly, the number density of massive quiescent galaxies at z ∼ 0.05 is also higher in miniJPAS than in previous studies, but less remarkable than for the star-forming case (see middle-left panel in Fig. 11). Consequently, this excess in number of nearby massive galaxies is also present for the stellar mass function of the full sample and persists when comparing with results based on COSMOS, 3D-HST, and DEVILS results at z ∼ 0.1 (see bottom-left panel in Fig. 11; Wright et al. 2018; Thorne et al. 2021). This excess might be a consequence of assuming that Schechter parameters (e.g. the logarithm of the characteristic stellar mass) evolve linearly with redshift, which may not be the case close to this redshift regime. However, a further detailed study of the stellar mass function at z ∼ 0.05 in the upcoming J-PAS survey would be necessary in order to confirm this. At z > 0.2, discrepancies between the stellar mass functions of both miniJPAS star-forming and the full sample of galaxies and those from other surveys (see top and bottom panels in Fig. 11 and also Pozzetti et al. 2010; Ilbert et al. 2010, 2013; Davidzon et al. 2017) are small and may partly result of the use of different stellar population synthesis models and techniques. Overall, we find that the stellar mass functions of quiescent galaxies are in good agreement with results from other surveys (see middle panels in Fig. 11). Nevertheless, there are slight discrepancies at the low mass regime (log₁₀M_⋆ ≲ 10 dex at z > 0.3), meaning that our parametrised stellar mass functions at the stellar mass range for which our sample is not complete (𝒞 < 0.05) show lower number densities than in other works. Part of these discrepancies may be explained by the use of the MCDE for the galaxy classification in miniJPAS (see Sect. 4), because other diagrams based on rest-frame colours (e.g. the UVJ diagram) and used for the same aim are commonly contaminated by a 20% fraction of dusty star-forming galaxies (Díaz-García et al. 2019a), especially at low mass. As noted above, we are unable to compare the second component of the stellar mass function of quiescent galaxies, log₁₀M_⋆ ≲ 9.5, at any redshift because our analysis is limited to a single Schechter model.

Fig. 11. Same as Fig. 10, but for stellar mass functions. Stellar mass functions from 2MASS/SDSS (star markers, Bell et al. 2003), zCOSMOS (pentagons, Pozzetti et al. 2010), S-COSMOS (circles, Ilbert et al. 2010), GAMA (plus markers, Baldry et al. 2012), UltraVISTA (squares, Ilbert et al. 2013), COSMOS2015 (diamonds, Davidzon et al. 2017), GAMA/COSMOS/3D-HST (triangle right, Wright et al. 2018), and DEVILS (hexagons, Thorne et al. 2021) surveys are included for comparison.

The cosmic evolution of the B-band luminosity and stellar mass densities go hand in hand with the evolution of the LMFs. Thus, the good agreement between the miniJPAS LMFs and those from previous studies is also reflected in the B-band luminosity and stellar mass densities (see left and right panels in Fig. 8, respectively). For both densities, there is a good quantitative and qualitative agreement with previous studies at 0.05 < z < 0.7, especially when uncertainties and abundance variance are taken into account (see also Sect. 8.2). In this regard, the B-band luminosity densities obtained with galaxies from the ALHAMBRA survey (López-Sanjuan et al. 2017) are especially similar to ours, especially the densities of the full sample. For the rest of B-band luminosity densities included here for comparison (including some of the most extended surveys such as 2dFGRS, COMBO-17, DEEP2, GAMA, and VIPERS; Madgwick et al. 2002; Faber et al. 2007; Loveday et al. 2012; Fritz et al. 2014; Beare et al. 2015), we find little numerical differences that in average do not exceed 0.1 dex. At z < 0.2, results involving data from 2MASS, SDSS, and AGES (see Bell et al. 2003; Cool et al. 2012) show larger differences, > 0.2 dex, which is more remarkable for star-forming galaxies but also in disagreement with the trends of other studies (see left panel in Fig. 8). Regarding the stellar mass densities, as expected, the variance or range of values obtained for this parameter in previous researches is wider than for the B-band luminosity density (see right panel in Fig. 8). In any case, we find that estimations from previous works in the literature are below 0.2 dex in the redshift range explored here. In fact, the stellar mass densities for quiescent and star-forming miniJPAS galaxies properly reproduce the evolutionary trends previously constrained in studies involving a galaxy spectral-type classification (including galaxies and data from 2MASS, SDSS, COSMOS, GAMA, UltraVISTA, and GALEX Bell et al. 2003; Ilbert et al. 2010, 2013; Baldry et al. 2012; Moustakas et al. 2013; Davidzon et al. 2017). For the total stellar mass density, there is also evidence of good agreement with the aforementioned and recent works (e.g. Davidzon et al. 2017; Wright et al. 2018). However, we note that the stellar mass densities obtained using data from the DEVILS survey (Thorne et al. 2021) are systematically higher at z < 0.7 than our values (see right panel in Fig. 8), although we do find a similar evolution of 0.4 dex from z = 0.7 to 0 as in Thorne et al. (2021).

8.2. Cosmic variance of the miniJPAS sample

As a consequence of the non-homogeneity of the Universe at small scales (≲1 Gpc, e.g. Davis et al. 1985), the number density of galaxies in small volumes may present variations amongst different pointings much larger than those predicted from Poisson statistics. These variations in the apparent abundance of galaxies, a.k.a. cosmic variance, may be particularly significant for the so-called ‘pencil beam’ surveys (< 1 deg²) to the point that can be one of the largest sources of uncertainty for intermediate and high-redshift surveys (see e.g. Newman & Davis 2002; Driver & Robotham 2010; Moster et al. 2011). Cosmic variance relies on multiple factors such as the aspect ratio of the survey footprint, the imaged cosmological volume (i.e. area and redshift range), the contiguity of the area surveyed (correlated or independent fields), and the galaxy stellar mass or luminosity (see e.g. Driver & Robotham 2010; Moster et al. 2011).

With this in mind and taking into account that the LMFs themselves are a direct measurement of the abundance of galaxies with different properties, we infer that cosmic variance may play a role in our LMF results and must be quantified according to the miniJPAS characteristics. As a consequence, we expect that the uncertainties of the Schechter parameters describing the miniJPAS LMFs (reported in Sect. 7.2) only account for counting errors and may be underestimated. In order to address this, we perform a rough estimation of the cosmic variance affecting our results, meaning that we roughly approach the uncertainties on the number densities of galaxies at different stellar mass, luminosity, and redshift. This in turn can explain part of the discrepancies and differences found when comparing with previous studies (see Sect. 8.1). For this aim, we adopt the equations proposed by Moster et al. (2011) for setting constraints on the cosmic variance of galaxies by using predictions from ΛCMD theory and the galaxy bias obtained via a halo occupation model to predict galaxy clustering. In particular, we choose all the parameters established for the EGS field owing to the geometry and size of this field is close to the miniJPAS one. Regarding the cosmic variance of galaxies for the B-band luminosity, we take advantage of the well-known stellar mass versus mass-to-light ratio correlation (see e.g. Bell et al. 2003). Making use of the SED-fitting results obtained with MUFFIT and miniJPAS data, we construct the distribution of B-band absolute magnitudes for each of the stellar mass bins included in the cosmic variance recipe proposed by Moster et al. (2011), all this at different redshift and for a constant width of Δz = 0.2. Hereafter, we assume that the cosmic variance at a luminosity equal to the mode of each of the B-band absolute magnitude distributions is roughly the same than the cosmic variance of the stellar mass bin employed for constructing each of these distributions.

Overall, the cosmic variance of the miniJPAS survey present theoretical values ranging from 15 to 40% (see Fig. 12), as expected for an EGS-like survey (details in Moster et al. 2011) of similar area. Despite these fractional values are only indicative, they are very useful for understanding how cosmic variance may affect our results. Firstly, for a constant redshift bin width, Δz, the miniJPAS cosmic variance decreases at increasing redshift as a consequence of a larger surveyed comoving volume (see light and dark green lines in Fig. 12). On the contrary, for a constant volume case, the cosmic variance increases at increasing redshift since Δz is smaller (see Moster et al. 2011). As a consequence, our sample at the nearby Universe is subject to a greater cosmic variance than at high redshift. Secondly, the cosmic variance of both bright and massive galaxies is greater than for their faint and less massive counterparts at a given redshift and volume. This implies that at a fixed redshift, the bright or massive ends of LMFs are going to be worse determined than the less massive parts. This fact also reflects that the brightest and most massive galaxies are less frequent in cosmological surveys, and therefore, these are subject to greater cosmic variance and Poisson variations.

Fig. 12. Fractional cosmic variance of our sample (Y-axis) at different redshift (see inset), stellar mass (X-axis, top panel), and B-band absolute magnitude (X-axis, bottom panel) for redshift and stellar mass bins of width Δz = 0.2 and Δlog10M⋆ = 0.5, respectively.

Cosmic variance is a general uncertainty that, a priori, may also affect the values of the Schechter parameters of the LMFs obtained during our likelihood maximisation. A precise determination of how cosmic variance is affecting our results is beyond the scope of this work and would require a more sophisticated analysis, which should account for our adoption of redshift-dependent Schechter parameters¹⁹. However, previous works based on semi-analytic models and real observations showed that different luminosity bins are highly correlated, especially for values below the characteristic luminosity (Smith 2012; López-Sanjuan et al. 2017). This bin-to-bin covariance is mainly due to sample variance arising from the presence of large-scale structures and valid for both volume- and flux-limited samples (Smith 2012). This means that if there is an upward/downward fluctuation in the average number of galaxies in one LMF bin, the other bins will share an upward/downward fluctuation, meaning that cosmic variance is not a random process along the stellar mass or luminosity. Consequently, we expect that the Schechter parameter most affected by cosmic variance is the LMF normalization (see also Trenti & Stiavelli 2008; Smith 2012), whereas the LMF shape may be slightly affected by it (see also Sect. 8.3). For the incoming J-PAS survey (∼8000 deg²), we expect that these LMFs are going to be hardly affected by cosmic variance, which in turn will be useful for determining how cosmic variance affects our methods depending on the geometry and volume used for constructing galaxy subsamples (see López-Sanjuan et al. 2015, 2017).

8.3. Lessons learned for J-PAS and future prospects

The analysis of the miniJPAS sample does not show striking or new results about the LMFs of galaxies down to z = 0.7, since these functions have been largely studied in the last decades. In general, the values of the characteristic parameters, the low-mass and faint-end slopes, and normalisations of the LMFs present values that are in the range of values obtained in recent and previous studies (see Table 3 and Sect. 8.1).

In fact, the good agreement between the B-band luminosity functions and those from other surveys reflects that current star and galaxy classification, completeness limits, photo-z constraints, and aperture photometry are yielding results that are self-consistent and in good agreement with previous works. Moreover, the great degree of similarity with the ALHAMBRA luminosity functions contributes to think that the spectral classification of galaxies along with the completeness of the miniJPAS sample are properly determined for miniJPAS and J-PAS owing to the concordance with the results from this much deeper survey. This also supports the idea that the statistical treatment of the data via probability distribution functions, as the obtained from MUFFIT and the other codes in the J-PAS Collaboration (see González Delgado et al. 2021), is highly recommended for shedding light on the determination of some of the stellar population properties of galaxies. In addition, the inclusion of the PDZs in the SED-fitting analysis is necessary for a proper estimation of the uncertainties affecting some of the stellar population parameters such as stellar mass and luminosity, especially at fainter magnitudes (see also Fig. 5 in this work and Fig. 12 in González Delgado et al. 2021). In particular, J-PAS will be a survey that can accurately constrain the stellar mass of luminous red galaxies (LRGs) with an apparent magnitude of r_SDSS = 22 with an average precision of ∼0.2 dex.

Based on the results obtained in Sect. 7.4 (see also Fig. 9), the LMFs obtained in this work can be used to provide a more robust spectral-type classification of galaxies in future studies. For instance, our results illustrate that, after the Eddington bias correction, the most massive galaxies at low redshift are almost entirely quiescent galaxies. This means that even when a very massive galaxy in the nearby Universe exhibits a value of P_Q ∼ 0.5, this is actually a quiescent galaxy that owing to errors is scattered close to the colour limits used to separate both types (formally expressed by Eq. (1) and Table 1). Therefore, the LMFs describe the probability of observing a galaxy of each type within the MCDE diagram that in turn can be used as a prior to perform a more robust classification by a Bayesian PDF analysis (see e.g. López-Sanjuan et al. 2019). Nevertheless, we note that this refinement should not be addressed before an LMF analysis, since this might introduce a bias in the results. We plan to tackle this more refined galaxy classification in a future work.

On the other hand, we find necessary a revision of the LMFs at the very nearby Universe with J-PAS. The access to many more nearby galaxies (z < 0.2) will allow us to explore the possibility of using double Schechter functions to parametrise the LMFs of quiescent galaxies, as well as to check whether the Schechter parameters are appropriately represented by functions that linearly depend on redshift down to z ∼ 0.1. In addition, a large volume of galaxies will allow us to discard cosmic variance effects, thus reducing the uncertainties in the determination of these functions at any redshift. These facts will facilitate the proper determination of the LMFs at the low mass and faint regimes, which are underestimated in comparison with previous results, hence discarding whether this results of the use of the MCDE diagram for the spectral-type classification of galaxies.

In light of our results and based on our experience for the determination of the miniJPAS luminosity functions, we recommend the inclusion of the odds parameter in future studies involving J-PAS photo-z and/or PDZs, or at least considering whether there may be an impact and/or bias on the results when accounting for odds. This arises from the fact that faint miniJPAS galaxies (r_SDSS ∼ 22) usually present lower odds values that in turn have photo-z biased to higher values after comparing with spectroscopic redshifts, which is also supported by the luminosity functions constrained in this work. Otherwise, the β parameter in the Schechter function for luminosity (see Eq. (8)) would be overestimated and more rapidly evolved after comparing with previous results from the literature. It is also of note that a definition of the sample based on the odds parameter (e.g. by choosing only sources with odds higher than a certain limiting value) can also result in a selection bias of the data set.

9. Summary and conclusions

This work is mainly focused on the development of a robust methodology for the determination of the luminosity and stellar mass functions (LMFs) of galaxies by solely using photometric measurements from multi-filter imaging surveys. Even though this methodology can be easily extended to other kind of surveys, we are especially interested in the application of our methods in the incoming J-PAS survey, putting emphasis in the proper determination of the LMFs of galaxies up to z ∼ 0.7. In addition, we are especially interested in a method involving the value-added catalogues and the typical outputs from the SED-fitting codes used in the J-PAS Collaboration in order to perform a solid and statistical analysis accounting for the uncertainties, degeneracies, and correlations amongst all the involved parameters.

As J-PAS is still an ongoing survey, we use a previous data set referred as miniJPAS (a stripe of ∼1 deg² centred at the EGS field, see Bonoli et al. 2021) that was dictated according to the J-PAS strategy and specifically conceived for testing the performance and potential of the future J-PAS. A portion of the ingredients needed for the determination of the LMFs were previously calculated by the J-PAS Collaboration and included in miniJPAS value-added catalogues. More precisely, we made use of the ‘default’ photometric redshift constraints (photo-z, Hernán-Caballero et al. 2021) and the star and galaxy classification of sources (López-Sanjuan et al. 2019), along with the forced-aperture photometry included in the general miniJPAS catalogues (Bonoli et al. 2021). However, the determination of LMFs requires additional inputs (e.g. stellar mass and luminosity constraints for each of the miniJPAS galaxies, spectral-type classification of galaxies, and completeness of the flux-limited sample) that are constrained by a set of techniques specifically designed and built from scratch for this study.

Firstly, we determined the probability distribution functions of the B-band luminosity and stellar mass via SED-fitting techniques. For this purpose, we used an updated version of the SED-fitting code MUFFIT (Díaz-García et al. 2015), which is carefully developed and optimised to deal with multi-band photometric data. In brief, the SED-fitting analysis performed by MUFFIT is solely based on the stellar continuum of galaxies or colours (i.e. nebular and AGN emission lines are removed during the analysis), using composite models of stellar populations (mixtures of two SSPs, i.e. non-parametric star formation history) for a recent version of the Bruzual & Charlot (2003) SSP models (CB17), attenuation law of Calzetti et al. (2000), miniJPAS fluxes obtained for AUTO aperture photometry, and photo-z probability distribution functions provided for each of the miniJPAS sources. Secondly, we adopted the stellar mass versus rest-frame colour diagram corrected for extinction (MCDE, Díaz-García et al. 2019a) for performing a spectral-type classification of galaxies (quiescent and star-forming); thus, we were able to prevent the inclusion of dusty star-forming galaxies in the quiescent sample. As a result, we assign probabilities of being quiescent galaxies according to the position of sources within this diagram, uncertainties, and degeneracies in the involved parameters. Thirdly, we build parametric functions for the stellar mass and B-band luminosity completeness of our sample as a function of redshift and for different magnitude limits of the r_SDSS detection band.

The LMFs of the miniJPAS galaxies are formally characterised by single-Schechter functions and according to the spectral-type classification of galaxies, while the LMFs of the full sample result of the sum of the quiescent and star-forming ones. The determination of the parameters defining the LMFs are constrained by a novel maximum-likelihood method, which makes the most of our sample and all the statistical results obtained within the J-PAS Collaboration. Thanks to the inclusion of weights, we were able to include the probability distribution functions and correlations of the parameters involved (photo-z, stellar mass, and B-band luminosity), sample completeness, various probabilities (e.g. the star-galaxy classification), and priors (e.g. the range of values for the Schechter parameters) for each of the galaxy spectral-types explored in this work (i.e. quiescent and star-forming) in a robust statistical way. Moreover, the miniJPAS LMFs are subsequently corrected for the so-called Eddington bias effect.

Broadly, the LMFs obtained by miniJPAS galaxies point to smooth evolution with redshift down to z = 0.7, where their shapes or Schechter parameters mainly rely on the galaxy spectral-type. The LMF variations with redshift for star-forming galaxies mainly affects the bright and high-mass ends, while their faint and low-mass ends remain roughly constant. In this sense, we find that there is a decreasing number of bright star-forming galaxies (M_B < −23) at a decreasing redshift, while our results are compatible with a subtle increase of massive star-forming galaxies (log₁₀M_⋆ > 11 dex). Regarding quiescent galaxies, we reached similar conclusions for the bright and high-mass ends of the LMFs than for the star-forming case. However, there is evidence of a greater evolution in the number density of low/intermediate-mass and faint quiescent galaxies (log₁₀M_⋆ < 10 and M_B > −20, respectively). In addition, the cosmic evolution of the global B-band luminosity density is undergoing a slight decrease of ∼0.1 dex from z = 0.7 to z = 0. There is evidence that this density decrement is mainly driven by star-forming galaxies, since the B-band luminosity density of quiescent galaxies is compatible with a non-evolution. On the other hand, the global stellar mass density increases by ∼0.3 dex in the same redshift range, which is mainly driven by the quiescent galaxy population. At z < 0.7, the B-band luminosity functions point out that the brightest galaxies are mainly star-forming galaxies, whereas more than a 60% fraction of massive galaxies (log₁₀M_⋆ > 10.7 dex) are quiescent galaxies. The fraction of red galaxies for the faint and low mass regimes (M_B ≳ −22 and log₁₀M_⋆ < 10.5 dex, respectively) is below 20%, meaning that star-forming galaxies dominate in number and are much more frequent than quiescent galaxies in these regimes. Nevertheless, the fraction of red galaxies with M_B ≳ −22 and/or log₁₀M_⋆ < 10.5 dex is gaining prominence at decreasing redshift.

After comparing with LMFs from previous works involving photometric and/or spectroscopic data at redshift 0.05 ≤ z ≤ 0.7 (e.g. SDSS, DEEP2, zCOSMOS, UltraVISTA, 3D-HST, etc.), we find a reasonable good agreement with the miniJPAS LMFs for both quiescent and star-forming spectral-types. Likewise, the cosmic evolution of the B-band luminosity and stellar mass densities show similar values and trends than those obtained in similar researches. These results highlight the fact that the miniJPAS photo-z constraints, star and galaxy classification, spectral-type, stellar mass, and rest-frame luminosities computed within the J-PAS Collaboration, along with the methodology introduced in this manuscript, are self-consistent with respect to determining the LMFs of galaxies down to z = 0.7. This outcome is in good agreement with previous works.

Acknowledgments

The authors are grateful to the anonymous referee for the valuable and constructive comments, which helped to improve this paper. L.A.D.G., R.M.G.D., R.G.B., G.M.S., J.E.R.M., I.M., and J.M.V. acknowledge financial support from the State Agency for Research of the Spanish MCIU through ‘Center of Excellence Severo Ochoa’ award to the Instituto de Astrofísica de Andalucía (SEV-2017-0709) and CEX2021-001131-S funded by MCIN/AEI/10.13039/501100011033. L.A.D.G., R.M.G.D., R.G.B., G.M.S., and J.E.R.M. are also grateful for financial support from the project PID-2019-109067-GB100. S.B. acknowledges support from the Spanish Ministerio de Ciencia e Innovación through project PID2021-124243NB-C21. A.H.C. receives financial support from the Spanish Ministry of Science and Innovation (MCIN/AEI/10.13039/501100011033) and ‘ERDF A way of making Europe’ through the grant PID2021-124918NB-C44. In addition, I.M. thanks for financial support from the project PID2019-106027GB-C41 and J.M.V. from project PID2019-107408GB-C44. Based on observations made with the JST/T250 telescope and PathFinder camera for the miniJPAS project at the Observatorio Astrofísico de Javalambre (OAJ), in Teruel, owned, managed, and operated by the Centro de Estudios de Física del Cosmos de Aragón (CEFCA). We acknowledge the OAJ Data Processing and Archiving Unit (UPAD) for reducing and calibrating the OAJ data used in this work. Funding for OAJ, UPAD, and CEFCA has been provided by the Governments of Spain and Aragón through the Fondo de Inversiones de Teruel; the Aragón Government through the Research Groups E96, E103, E16_17R, and E16_20R; the Spanish Ministry of Science, Innovation and Universities (MCIU/AEI/FEDER, UE) with grant PGC2018-097585-B-C21; the Spanish Ministry of Economy and Competitiveness (MINECO/FEDER, UE) under AYA2015-66211-C2-1-P, AYA2015-66211-C2-2, AYA2012-30789, and ICTS-2009-14; and European FEDER funding (FCDD10-4E-867, FCDD13-4E-2685).

References

Appendix A: Stellar mass and B-band luminosity completeness for samples with different magnitude cuts

In this appendix, we include the coefficients defining the stellar mass and B-band luminosity completeness of the miniJPAS galaxies for samples with different magnitude cuts following the methodology detailed in Sect. 5. Even though these coefficients were not explicitly used in this research, they can be very useful to define the stellar mass or luminosity completeness of samples in future works involving miniJPAS data. In particular, we present the coefficients for magnitude cuts of r_SDSS ≤ 21.5, 21.75, 22.25, and 22.5 in the detection band (Tables A.1–A.4, respectively).

Appendix B: Constraints on the parameters defining the miniJPAS LMFs

In this section, we include the sampling of the posterior distributions of the coefficients defining the parametric B-band luminosity and stellar mass functions (Figs. B.1 and B.2, respectively) that resulted from the likelihood maximisation of Eq. 13 (largely detailed in Sects. 6.1–6.3) and the subsequent correction for Eddington bias effects (see Sect. 6.4). These results are obtained according to the galaxy spectral-type and constrain the diverse evolutionary tracks followed by the two kinds of galaxies explored in this work (i.e. quiescent and star-forming galaxies), as well as the degeneracies and/or correlations amongst the parameters involved in the single Schechter functions (see Sect. 7.2 for a further description).

Fig. B.1. Constraints on the parameters defining the B-band luminosity functions of miniJPAS star-forming and quiescent galaxies (top and bottom corner panels) according to the methodology developed in Sect. 6. Each dot correspond to one of the solutions contained in the MCMC chains. The dashed and dotted black lines illustrate the confidence regions containing the 68 and 90 % of probability, respectively. The dashed vertical lines establish the 1σ confidence interval and median values of the parameter histograms.

Fig. B.2. Same as Fig. B.1, but for stellar mass functions.

All Tables

All Figures

	Fig. 1. Photo-z probability distribution function (solid black line) of one noisy galaxy in miniJPAS (ID 2406-10239 and odds = 0.52). The shaded area illustrates the photo-z values defining the confidence level of 70% of probability. The red dots show the photo-z values used by MUFFIT to perform the SED-fitting analysis of this source.
In the text

	Fig. 2. SED-fitting analysis (top panel) and residuals (bottom panel) of a miniJPAS galaxy at z = 0.07 (ID 2470-10239). The black dots and vertical bars in the top panel are the galaxy fluxes and errors as observed by miniJPAS, respectively. The red squares are the best-fitting CSP model. The shaded area shows 2.5 times the photon-noise uncertainty of each band.
In the text

Fig. 3. Distribution of stellar mass versus rest-frame colour (uJAVA − rSDSS) obtained from miniJPAS sources with rSDSS ≤ 22.5 and PG ≥ 0.5 at different redshift bins. Top: diagram of the (uJAVA − rSDSS) colour after correcting for dust effects or intrinsic colour. Bottom: diagrams without correcting the rest-frame (uJAVA − rSDSS) colour for extinction and normalised histograms of PQ values for each redshift bin (see insets). The dashed black line illustrates the limiting relation for selecting quiescent galaxies and AUTO photometry at the central redshift of each bin (i.e. z = 0.15, 0.4, and 0.6; see Eq. (1)). Dots are colour coded according to the quiescent and star-forming classification (redder and bluer colours, respectively) or quiescent probability (PQ).

In the text

Fig. 4. B-band luminosity and stellar mass completeness as a function of redshift (top and bottom panels, respectively) for star-forming and quiescent miniJPAS galaxies (left and right panels, respectively) from the flux-limited sample at rSDSS ≤ 22 magnitudes (i.e. $r_{\text{SDSS}}^{\text{lim}}$ = 22). The dotted, dashed, and solid black lines illustrate the fit of the 5, 50, and 95% completeness levels (i.e. 𝒞 = 0.05, 0.5, and 0.95) obtained from the completeness methodology detailed in Sect. 5 (triangle, square, and circle shape markers), respectively.

In the text

	Fig. 5. Normalised distributions of errors of B-band luminosity (top panels) and stellar mass (bottom panels) at different rSDSS magnitude bins (see insets) for star-forming and quiescent galaxies (left and right panels, respectively) in miniJPAS.
In the text

Fig. 6. Evolution with redshift of the B-band luminosity functions of star-forming and quiescent galaxies (solid lines in left and right panels, respectively) from miniJPAS. The shaded areas exhibit the uncertainties in the luminosity functions obtained with our methodology. The dash-dot lines show the B-band luminosity functions before correcting of Eddington bias effects. The dashed and dotted lines point the absolute magnitude range in which our sample is affected by incompleteness after and before correcting for Eddington bias effects, respectively.

In the text

	Fig. 7. Same as Fig. 6, but for stellar mass functions.
In the text

Fig. 8. Cosmic evolution of the B-band luminosity and stellar mass densities (left and right panels, respectively) for quiescent, star-forming, and full miniJPAS samples (red, blue, and black lines, respectively). The shaded areas illustrate the uncertainties of these densities. Left panel: B-band absolute magnitude densities from 2dFGRS (cross markers, Madgwick et al. 2002), 2MASS/SDSS (star markers, Bell et al. 2003), COMBO-17/DEEP2 (squares and diamonds, respectively; Faber et al. 2007), AGES (pentagons, Cool et al. 2012), GAMA (plus marker, Loveday et al. 2012), VIPERS (triangle up, Fritz et al. 2014), Boötes field (triangle left, Beare et al. 2015), and ALHAMBRA (circle, López-Sanjuan et al. 2017) surveys are included for comparison. Right panel: stellar mass densities from 2MASS/SDSS (star markers, Bell et al. 2003), S-COSMOS (circles, Ilbert et al. 2010), GAMA (plus markers, Baldry et al. 2012), UltraVISTA (squares, Ilbert et al. 2013), SDSS/GALEX (crosses, Moustakas et al. 2013), COSMOS2015 (diamonds, Davidzon et al. 2017), GAMA/G10-COSMOS/3D-HST (triangle right, Wright et al. 2018), and DEVILS (hexagons, Thorne et al. 2021).

In the text

Fig. 9. Fraction of quiescent galaxies, fr, as a function of the B-band absolute magnitude (left panel) and stellar mass (right panel) at different redshift (see inset). The dashed lines illustrate the absolute magnitude and stellar mass range for which our quiescent sample is under a completeness level of 𝒞 = 0.25 and 0.05, respectively, meaning that the fraction of red galaxies are obtained through a extrapolation of our results.

In the text

Fig. 10. Comparison of the B-band luminosity function of miniJPAS galaxies (coloured lines from top to bottom: star-forming, quiescent, and full galaxy sample) with those fitted or inferred (Schechter-like functions) from previous works at different redshift. The shaded area illustrates uncertainties in our results. The dashed lines show the B-band luminosity functions for magnitudes fainter than the completeness lower limit of the sample employed for their determination (i.e. $M_B>M_B^{\mathcal{C}=0.25}(z)$). B-band luminosity functions from COMBO-17/DEEP2 (squares and diamonds, respectively; Faber et al. 2007), COSMOS (triangle right, Drory et al. 2009), AGES (pentagons, Cool et al. 2012), VIPERS (triangle up, Fritz et al. 2014), Boötes field (triangle left, Beare et al. 2015), and ALHAMBRA (circles, López-Sanjuan et al. 2017) surveys are included for comparison.

In the text

Fig. 11. Same as Fig. 10, but for stellar mass functions. Stellar mass functions from 2MASS/SDSS (star markers, Bell et al. 2003), zCOSMOS (pentagons, Pozzetti et al. 2010), S-COSMOS (circles, Ilbert et al. 2010), GAMA (plus markers, Baldry et al. 2012), UltraVISTA (squares, Ilbert et al. 2013), COSMOS2015 (diamonds, Davidzon et al. 2017), GAMA/COSMOS/3D-HST (triangle right, Wright et al. 2018), and DEVILS (hexagons, Thorne et al. 2021) surveys are included for comparison.

In the text

	Fig. 12. Fractional cosmic variance of our sample (Y-axis) at different redshift (see inset), stellar mass (X-axis, top panel), and B-band absolute magnitude (X-axis, bottom panel) for redshift and stellar mass bins of width Δz = 0.2 and Δlog10M⋆ = 0.5, respectively.
In the text

Fig. B.1. Constraints on the parameters defining the B-band luminosity functions of miniJPAS star-forming and quiescent galaxies (top and bottom corner panels) according to the methodology developed in Sect. 6. Each dot correspond to one of the solutions contained in the MCMC chains. The dashed and dotted black lines illustrate the confidence regions containing the 68 and 90 % of probability, respectively. The dashed vertical lines establish the 1σ confidence interval and median values of the parameter histograms.

In the text

	Fig. B.2. Same as Fig. B.1, but for stellar mass functions.
In the text