Strong size evolution of disc galaxies since z=1 Readdressing galaxy growth using a physically-motivated size indicator

,


Introduction
Many alternatives for measuring the size of galaxies have been proposed over the years (see a recent review in Chamba 2020).Constrained by the technological limitations of their time, many works delimited galaxy sizes by using isophotal contours, for example 25 mag arcsec −2 in the B band (Redman 1936) or 26.5 mag arcsec −2 (Holmberg 1958).Other choices for estimating the size were based on the radius containing a certain amount of light from the object.The most favoured proposal by far has been that of de Vaucouleurs (1948) to use the effective radius (r e , the radius containing half the light) as a size indicator for spheroidal galaxies.With the advent of wide-area surveys such as the Sloan Digital Sky Survey (SDSS; York et al. 2000), the need to have an automatic and robust way to characterise the sizes of millions of galaxies resulted in the effective radius becoming the most popular (and almost exclusive) galaxy size proxy since then.
Although using the effective radius has some distinct advantages, such as its robustness against image depth (as the galaxy's light profiles decrease exponentially or even more steeply), its Full Table 1 is available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr(130.79.128.5) or via https:// cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/682/A110identification with the galaxy size can lead to confusion (see e.g.Papaderos et al. 1996;Breda et al. 2019;Suess et al. 2019;Chamba et al. 2020;dos Reis et al. 2020).The main problem with the effective radius is that it does not provide a measure of the global extent of the galaxy, but only characterises the size at which 50% of the innermost light is distributed.This can therefore lead to situations where galaxies with very different effective radii (depending on whether or not there is a bright source, e.g. a bulge, bar, ring, etc., at their centre) can have the same global extent (as characterised by the position of their disc edge, for example).To overcome the problems of the previous size characterisations, Trujillo et al. (2020) proposed a physically motivated definition of the size of a galaxy.According to the new definition, the size of a galaxy would be given by the farthest radial location where the gas has been efficiently able to collapse and transform into stars.In this way, this new size proposal could potentially separate the regions of the galaxy that are mostly made up of the stellar component formed in situ from the region at the outskirts, which is likely to be dominated by the accreted material (see e.g.Font et al. 2020).Trujillo et al. (2020) show that using the new size definition, the observed dispersion in the size-mass relation decreases by a factor of about three with respect to that resulting from using the effective radius (for a theoretical discussion about the meaning of a decrease in the dispersion, see e.g.Sánchez Almeida 2020).
An expected outcome of the gas density threshold for star formation is a sharp drop in the surface brightness profiles of the galaxies in the outer parts.These sudden declines on the surface brightness profiles were found in the past in edge-on disc galaxies and were dubbed as truncations (van der Kruit 1979;van der Kruit & Searle 1981a,b).Motivated by both theoretical expectations and observational findings which agree on a value of ∼1 M pc −2 for the stellar mass density at the truncation position in present-day MW-like galaxies galaxies (Martínez-Lombilla et al. 2019;Díaz-García et al. 2022), Trujillo et al. (2020) proposed using the radial location where the stellar mass density profile reaches such a value as a proxy for measuring the size of the galaxies.There is no guarantee, however, that such a proxy based on MW-like galaxies should be able to characterise the radius to which in situ star formation takes place for all galaxies regardless of their mass, morphology, environment, or evolutionary epoch.For this reason, Chamba et al. (2022) explored the edges of galaxies in the local Universe for a sample of about a thousand objects covering a large range in mass and morphologies.Chamba et al. (2022) found that the stellar mass density at the truncation in these objects is moderately mass-dependent, being smaller for dwarf galaxies (∼0.6 M pc −2 ) than for objects of Milky Way mass (∼1 M pc −2 ).We have followed the same methodology here as in Chamba et al. (2022) to identify the size of galaxies but this time with the aim of exploring their evolution with cosmic time.
The purpose of this work is to revisit the size evolution of the disc galaxies with stellar masses similar to the MW using the new physically motivated definition by Trujillo et al. (2020).Since most previous works on galaxy size evolution have used the effective radius as a size indicator (see e.g.Shen et al. 2003;van der Wel et al. 2014;Lange et al. 2015, to name just a few), the use of a radius that is related to the in situ star formation can potentially provide a radically different size evolution of these objects.Indeed, for disc galaxies with the mass of the MW, the observed size evolution using the effective radius as an indicator (see e.g.van der Wel et al. 2014;Mowla et al. 2019;Nedkova et al. 2021;Kawinwanichakij et al. 2021) is found to be relatively mild (with an increase in size since z ∼ 1 of only 20−30%).This mild size evolution of disc galaxies using the effective radius is found no matter how the selection of disc galaxies was made (either by morphology or by their starforming nature), the wavelength where r e was measured (usually taking the V-band restframe, approximately λ ∼ 5000 Å), the facility used to image the galaxies, or the different codes used to derive the galaxy structural parameters (see e.g.Häußler et al. 2022;Euclid Collaboration 2023).Interestingly, this reported mild evolution provided by the use of the effective radius contrasts sharply with simple theoretical predictions based on first principles' arguments suggesting the evolution of nearly a factor of ∼2 of the size of disc galaxies since z ∼ 1 (see e.g.Mo et al. 1998).
To test whether the use of the new size definition produces results more on the line of theoretical expectations, it is necessary to have very deep, high spatial resolution images.For this reason, we take advantage of the Hubble Space Telescope (HST) observations of the CANDELS fields.Deep HST images in the past have shown their power to detect and characterise low surface brightness features (see e.g.Trujillo & Pohlen 2005;Azzollini et al. 2008;Trujillo & Bakos 2013;Buitrago et al. 2017;Borlaff et al. 2019).In order to maximise the possibility of exploring the outer parts of the galaxies with high signal-tonoise ratio and therefore identify their edges reliably, we restrict our study to only massive (M stellar 10 10 M ) disc systems.In addition, we only study galaxies with z < 1 because of the growing surface brightness cosmological dimming -proportional to (1+z)3 when working with flux density (as it is the case in the AB magnitude system) for resolved objects (Giavalisco et al. 1996;Law et al. 2007;Ribeiro et al. 2016).
This work is structured as follows: Section 2 describes the datasets we used and the selection criteria for our galaxy sample.Section 3 explains the methodology we adopted in order to derive the size of our galaxies.Section 4 shows the newlyderived mass-size and mass-stellar density at the edge position relations.Section 5 reports our conclusions.The appendices outline the complementary tests that we have carried out in order to prove the reliability of our results.Hereafter, our assumed cosmology is Ω m = 0.3, Ω Λ = 0.7 and H 0 = 70 km s −1 Mpc −1 .We use a Chabrier (2003) Initial Mass Function (IMF), unless otherwise stated.Magnitudes are provided in the AB system (Oke & Gunn 1983).

Data and sample selection criteria
As our reference dataset, we use the latest available images (v1.0) in the CANDELS survey 1 (Grogin et al. 2011;Koekemoer et al. 2011).In the case of the GOODS-South field, however, we use the data from the Hubble Legacy Field 2 (HLF v2.0; Illingworth et al. 2016) which contains a larger area coverage.Additionally, for the Hubble Ultra Deep Field, we make use of its associated optical imaging (Beckwith et al. 2006) while in the near-infrared we take the ABYSS data (Borlaff et al. 2019).This last survey is a new rendition of the HUDF WFC3 imaging taking special care in preserving the LSB structures 3 .The filters we use throughout this paper are the F606W (V-band) and the F814W (I-band) from the ACS camera, and the F125W (Jband) and the F160W (H-band) from the WFC3 camera, except in the GOODS-South field where we supplement them with the ACS F435W (B-band) and the WFC3 F105W (Y-band).The set of filters used in this work are displayed in Fig. 1.
Our target selection is based on the CANDELS public catalogues4 (Santini et al. 2015;Stefanon et al. 2017;Nayyeri et al. 2017;Barro et al. 2019).We took the galaxy IDs, photometric masses and spectroscopic redshifts (only those flagged as the ones with top quality), discarding those objects with SExtractor stellarity >0.8.We added other high-quality spectroscopic redshifts from other sources, namely the LEGA-C DR2 redshifts with f use = 1 (van der Wel et al. 2016;Straatman et al. 2018) and later from the DR3 data (van der Wel et al. 2021), the ZCOSMOS Final Data Release (Lilly et al. 2009) with Confidence Class 3 and 4 and hCOSMOS (Damjanov et al. 2019) secure redshifts (r-value >5).With this information at hand, we kept in our sample only those galaxies with M stellar > 10 10 M and up to a distance of z spec < 1.1, in order to optimise the signal-to-noise ratio for every object in our sample.This produces a total number of 2192 objects.We supplement these data with single Sérsic structural parameters coming from the CANDELS public catalogs (van der Wel et al. 2012).However, for well-resolved galaxies (in practice those at z < 0.2), these structural parameters (namely effective radii, Sérsic indices, axis ratios and position angles) are sometimes only representative of the luminous central galaxy parts that drive the χ 2 galaxy fits (Papaderos et al. 1996;dos Reis et al. 2020) angle and axis ratio of the galaxies are key to extract the surface brightness profiles, for those objects with z < 0.2, we re-estimate their position angles and axis ratios by fitting an ellipse to the region with surface brightness (in the H band) between 23.5 and 24 mag arcsec −2 .By doing so, we obtain a more representative axis ratio and position angle of the outer parts of the galaxies.Besides, for all galaxies, the axis ratio and position angle values were visually scrutinised and changed if necessary.
To select the disc galaxies in our sample, we gathered machine-learning derived morphologies from Huertas-Company et al. (2015).Specifically, we select those objects that fulfil the conditions for being disc-dominated, namely those catalogued either as pure discs, bulge+disc, or irregular discs (see Sect. 6 of the aforementioned paper).From our initial 2192 objects, our sample was then reduced to 1442 classified as disc galaxies.We also removed: 158 interacting galaxies (as their edges are ill-defined), 101 galaxies lacking imaging in either F606W or F814W bands, 76 galaxies with either noisy imaging or located close to the survey borders, 28 extremely compact galaxies (with a visual extension less than 5 kpc) whose surface brightness profiles does not allow the identification of an obvious edge, 17 visually identified by us as elliptical galaxies, five galaxies with odd pixel values, four visually identified by us as stars, four galaxies with non-realistic colours and one stellar spike, for a total final sample comprising 1048 galaxies.
Summarising, our galaxy sample selection criteria for our targets were: -M stellar > 10 10 M .
-Galaxies not being affected by any observational artefact.
-Non-interacting galaxies.For galaxies fulfilling these conditions we created postage stamps of 150 × 150 kpc, to a maximum of 1000 × 1000 pixels for ACS stamps and 500 × 500 pixels for WFC3 stamps.Some examples for MW-like galaxies in our sample at different cosmic times can be found in Fig. 2. From low to high redshift, their IDs are: J205432.42+000231.09(to show a nearby Universe counterpart from Chamba et al. 2022), 12783_COSMOS, 8902_EGS, 16245_COSMOS, 14358_EGS and 8221_UDS.These IDs are constructed joining the CANDELS IDs adding their respective galaxy fields.
Since our sample of galaxies is selected on the basis of spectroscopic redshifts, a potential bias in the analysis is that the sub-sample of disc galaxies with spectroscopic redshift could be biased towards objects with brighter effective surface brightness than the general disc population.What we have done to investigate this issue is to calculate the effective surface brightness in H-band of all the disc galaxies in the CANDELS sample with z < 1 (i.e.where the redshift has been derived both photometric and spectroscopically).We have checked whether the distribution of the effective surface brightness differs depending on whether the galaxy has a spectroscopic redshift or not.We find that the distribution of effective surface brightness in the original parent sample and our spectroscopically selected sub-sample are very similar, making it unlikely that our selection criterion introduces any potential bias.

Methodology
In this section, we describe the steps we have followed to determine the galaxy edges, and therefore, retrieve the size of these objects using the size definition proposed in Trujillo et al. (2020) and Chamba et al. (2022).

Observed surface brightness profiles
As mentioned in Sect. 1, our proxy for the galaxy size corresponds to the radial location of the gas density threshold that enables star formation.As a consequence of this physical process, we expect this to leave its imprint on the galaxy surface brightness and radial colour profiles5 .To retrieve accurate galaxy surface brightness profiles, it is mandatory to perform a careful masking of the neighbouring objects.We created optical and near-infrared masks according to SExtractor (Bertin & Arnouts 1996) detections in the I-and H-bands, the reddest filters in the HST cameras used in this work.To be conservative, as the ACS data has large spatial resolution and therefore the deblending of the objects is more simple, we joined both masks to create the final near-infrared masks.This was done resampling the optical masks to the near-infrared spatial resolution.Our masking algorithm is an evolution of the one successfully implemented in other publications such as Trujillo et al. (2007) and Buitrago et al. (2008Buitrago et al. ( , 2013)).However, the output masks have been improved for the majority of the galaxies studied to ensure that there is no obvious source of spurious light entering our analysis.This is done by manually masking very faint and small sources that could have an effect, especially at the faintest surface brightness levels.
An accurate background removal is key to determine a reliable surface brightness profile.To perform this task, we followed the strategy described in Pohlen & Trujillo (2006) and Buitrago et al. (2017), measuring the background in regions very close to the galaxies in order to be representative of the local sky background.Using as a first crude size for the galaxy the isophote at 27.5 mag arcsec −2 in the V-band, we estimate the local background using pixels located in elliptical apertures with semi-major axes between 2 and 3 times such rough size.That is  (Trujillo & Fliri 2016), while the rest of images are HST data.These objects were selected to have representative sizes (R edge , golden ellipses) of galaxies at those redshifts.Note how the sizes (using our in-situ star formation criterion) decrease with redshift, while the effective radii (r e , blue ellipses) do not change appreciably.done in every filter, and then we subtracted the mean background flux within that annulus.
The determination of the galaxy size according to the definition used in this work relies on the determination of a sharp change (i.e. a truncation) on the surface brightness (stellar surface mass density) profiles of the galaxies.Traditionally, surface brightness profiles have been derived by azimuthally averaging the galaxy flux in annulli at progressive further galactocentric distances.We call this technique the "ellipse method".It has the advantage of increasing the signal-to-noise of the profiles but, on the other hand, this method is also prone to blur the truncations due to its mixing of light from combining potentially different radial locations of the truncations depending on the azimuthal angle.This is the reason behind computing the surface brightness profiles in two manners in the present work, using both this ellipse method and along the semi-major axis of the galaxies (calling this latter one the "slit method").For the first technique, we derived the galaxy fluxes using resistant mean determinations (i.e.removing 3σ outliers) in one-kpc-wide elliptical annulli.The centre, axis ratio and position angle of the elliptical annulli were fixed to the values describing the light distribution of the galaxies retrieved in the H-band images (see Sect. 2).H-band images are used to avoid dust obscuration effects and to better track the galaxy stellar mass.On the other hand, we used the slit method to trace the surface brightness profiles of the galaxies along their semi-major axis.This is the strategy usually performed for edge-on disc galaxies (e.g.Martín-Navarro et al. 2012).The disadvantage of the slit method with respect to the ellipse method is that the signal-to-noise is poorer due to the lower number of pixels used.To minimise as much as possible this issue, when creating the semi-major axis profiles, our slit apertures had a width of 3 kpc and we summed both sides of the galaxy.This method is superior for dealing with highly inclined galaxies and also for finding sharp drops in the disc light, as only pixels in a given direction contribute to the final statistics.We have studied the location of the edge of galaxies using both the ellipse method and the slit method.We find that the location of the edge can vary by about 10% from one method to the other.Larger deviations of this variation can be expected if the outer parts of the galaxies are irregular.

Sloan-restframe equivalent profiles
To derive the stellar mass density profiles, we converted the observed surface brightness profiles into their equivalent in Sloan restframe filters.To that end, we took advantage of the already existing recipes to derive mass-to-light ratios from Sloan colours (see e.g.Bell et al. 2003;Roediger & Courteau 2015).In A110, page 4 of 22 order to perform this task correctly, we corrected the observed surface brightness profiles using the following steps: -Galactic extinction correction.To do this, we used the values provided by the NASA NED Extinction Calculator6 : with µ corr,1 and µ obs being the corrected and observed surface brightness profiles, respectively.-Cosmological surface brightness dimming correction.We use (1+z) 3 as we are working with flux densities as explained before: -Disc inclination correction.For this last step, we used the expression assuming z 0 /h = 0.12 (for the determination of the α j coefficients see Sect.5.2 in Trujillo et al. 2020), where z 0 indicates the disc's model scale height and h its scale length.We do not correct our surface brightness profiles for internal dust extinction.While this effect is potentially moderate for lowinclination galaxies, it could be relevant for the galaxies in our sample with high inclination.
The Sloan restframe profiles are calculated following the procedure in Buitrago et al. (2017).We interpolated the observed surface brightness profiles linearly between contiguous observed bands in order to obtain the restframe Sloan magnitudes.From them, we also derived radial colour profiles.To further improve our results, we also checked whether the derived Sloan colour values were consistent with the predictions of the E-MILES library7 (Vazdekis et al. 2012(Vazdekis et al. , 2016;;Ricciardelli et al. 2012).We explicitly check that the Sloan restframe colours were consistent with stellar population ages smaller than the age of the Universe at each galaxy's redshift.We use the predictions provided by the Padova+00 isochrones (Girardi et al. 2000) and a Chabrier (2003) IMF for any given metallicity.If the Sloan colours are compatible with ages older than the Universe at that cosmic time, we marked these colour points as unreliable.We also mark as unreliable those points in the radial colour profiles with errors greater than 0.2 mag.

Stellar mass density profiles
We created the stellar mass surface density profiles following the scheme described in Bakos et al. (2008) that utilises the following expressions: where m abs, ,λ is the absolute magnitude of the Sun at the filter8 corresponding to the wavelength λ and µ λ is the corrected surface brightness profile using that filter.The final result is given in units of M pc −2 .The mass-to-light ratio at the filter corresponding to λ is calculated as follows: where the coefficients a λ and b λ are tabulated in Roediger & Courteau (2015), assuming a Chabrier ( 2003) IMF.We obtained all possible combinations of base profile (g, r, i, z) and colour (g − r, g − i, g − z, r − i, r − z, i − z) to create as many as possible stellar mass density profiles.This gives us a total of 24 different stellar mass density profiles.The final stellar mass density profiles we use in this work are the resistant mean combination of each individual stellar mass profile generated above.The estimation of the errors is done by taking into account the individual errors for each individual mass profile (which is basically given by the error in the colour used to compute it) and the intrinsic dispersion of the 24 independent mass profiles.The latter arise from the uncertainty of the Roediger & Courteau (2015) method and reflects the error in the M/L determination.Both errors are summed in quadrature.
As for the Galactic extinction and inclination errors, they are not taken into account since both corrections are systematic and identical for all profiles.
As a consistency test, we derived the total stellar masses of our galaxies using the derived stellar mass density profiles.For galaxies displaying low inclinations (axis ratio >0.3), the total stellar masses were derived using the stellar mass density profiles derived from the ellipse method.For galaxies displaying high inclinations (axis ratio <0.3), the total stellar masses were derived using the stellar mass density profiles derived from the slit method.It is important to note that this is done to better reflect the observed symmetry of the object under study and to minimise the potential effect of the thickness of the disc in our total stellar mass estimations.The integration of the radial profiles to get the total stellar mass cover from zero to the last radial point consider reliable (see Sect. 3.2).This last restriction includes situations when the stellar mass density profile started to rise as a consequence of being artificially affected by the light of a neighbour source.
Our integrated stellar masses are compared to the CANDELS-derived masses (see Sect. 2) and the recent ASTRODEEP GS43 catalog (Merlin et al. 2021).This last one is a state-of-the-art catalog containing consistent photometry for 43 medium-and broad-filters available in the GOODS-South field.There are a total of 159 galaxies overlapping with our sample in this case.The outcome of this comparison is shown in Fig. 3.While both CANDELS and GS43 masses display a similar scatter in ∆mass (=log(M ours )-log(M catalog )) with respect to our determinations, 0.13 and 0.10 dex for CANDELS and GS43 respectively, there is an offset with respect to CANDELS masses (0.19 dex).This is not found in the more comprehensive ASTRODEEP GS43 catalog (0.003 dex).This test gives us confidence on our stellar mass density profiles derivation.Therefore, on what follows, we adopt the total stellar masses derived from the integration of our stellar mass density profiles.

Deriving the size of the galaxies
We use the radial position of the truncations in our profiles as a proxy to derive the physically motivated size radius.Following Fig. 1 in Chamba et al. (2022), the truncation positions (R edge ) are identified as sharp changes in the galaxy surface brightness, colour and/or stellar mass density profiles in the outermost galaxy regions.On what follows, the truncation feature is dubbed as edge to reinforce the idea that our proxy traces the end of the in-situ star forming disc.When the edges are not obvious on the stellar surface mass density and/or surface brightness profiles, the radial colour profiles are particularly useful to retrieve the radial position of the edge.In fact, Pfeffer et al. (2022) and Chamba et al. (2022) show that the rising of the characteristic U-shape in their colour profiles (Azzollini et al. 2008;Bakos & Trujillo 2012;Watkins et al. 2016Watkins et al. , 2019;;Martínez-Lombilla et al. 2019) is connected with a significant drop in the star formation of the galaxy.At that radial position, the disc light diminishes while the potential emergence of a stellar halo gains in importance creating a redder colour.Therefore, using that feature can help us to easily identify the location of the edge of the disc.
Quantitatively, the process we have followed to identify galaxy edges is as follows.Based on the work done in the present-day Universe on Milky Way-like edge-on galaxies such as NGC 4565 and NGC 5907 (see Martínez-Lombilla et al. 2019), we expect the change on slope of the stellar surface mass density profiles to be around 2−3 (i.e.h1/h2 ∼ 2−3, where h1 is the scale length of the disc before the truncation and h2 is the scale length after the truncation).For galaxies with lower inclinations, the change in slope at the edge is expected to be milder for a number of reasons.First, the contrast (and therefore the S/N) between the inner disc and the outer region of the galaxy decreases due to the line-of-sight integration.Second, the use of ellipse profiles (instead of slits) averages light from different angles of the galaxy, and if the edge is not at exactly the same radial distance in all directions, this will also blur the edge detection.We have simulated (using IMFIT; Erwin 2015) a disc galaxy with different inclinations to see what change in the h1/h2 ratio we would expect.The result of this analysis is that even for a perfectly symmetric disc, if we start with a value of h1/h2 = 3 for an edge-on configuration, for a lower inclination we measured values with h1/h2 ∼ 1.5.In real galaxies, taking into account the effects of the point spread function (PSF), deviations from symmetry, etc., 1 < h1/h2 < 1.5 is then expected.
From our initial sample of 1048 disc galaxies, about 60% of them (i.e.634) show a change in the slope in the outer parts of the stellar surface mass density profiles that is compatible with the expectation (i.e.h1/h2 > 1).In cases where h1/h2 ∼ 1 (within the errors in their determination), the use of the stellar surface mass density profiles is not sufficient for our purposes to detect an edge with some confidence.For this reason, our next step is to use the surface brightness profiles in the observed HST bands V, I, J and H and to re-estimate h1/h2 around the potential location of the edge.If h1/h2 is greater than 1 then we use such a value as the location of the edge.This gives a further 8 galaxies with an identified edge.
If the value of h1/h2 is still not greater than 1 with sufficient confidence, then we examine the shape of the rest-frame colours (mainly g − r) to obtain an estimate of the edge position.In this case, we take advantage of the well-known U-shape of the colour profiles for disc galaxies (see e.g.Bakos et al. 2008).This includes a further 287 galaxies with edges identified.Only 119 galaxies (11% of our sample) have no clear edges.In these cases, we identify the edge as the visual limit to which we can see the galaxies in our images.Finally, once the location of the edge (R edge ) is found, the stellar mass surface density at that position (Σ edge ) is what we assign as the stellar mass surface density at the edge.In Appendix A, we show some examples of the expected variation in the determination of the R edge using the different methods.In general, we find good agreement (with uncertainties on the R edge position within 1 kpc uncertainty).In addition, we probe the robustness of the size determination using recently released very deep and superior spatial resolution JWST data (Eisenstein et al. 2023;Rieke et al. 2023).
To test the accuracy of our approach, the authors (F.B., I.T.) independently carried out the determination of galaxy sizes (R edge ) by visually comparing all the profiles and images for the galaxies in our sample.To be sure that the determination of the edge in our profiles is not an artefact of the noise of our images, we calculated the limiting surface brightness of each our images.To make this estimate, we randomly selected pixels on the masked images and calculated the scatter of the pixel distributions in apertures corresponding to 1 × 1 arcsec.Our limiting surface brightness corresponds to a 3σ fluctuation of the background noise on these areas.A detailed description of the calculation of these limits is given in Román et al. (2020).While the background pixel noise is a simple approximation to characterise the true surface brightness limit of the image, especially in the presence of large scale light gradients that may have been left over during data reduction, in the particular case of the galaxies we are working with, which are significantly smaller than the size of the HST detectors, we believe that our local characterisation of the background noise is a fair representation of the local surface brightness limits around such galaxies.
The edge feature and, therefore, the size determination is found in all (except for 11% of the cases, see above) to be located at radial distance where the observed surface brightness profiles are brighter than the limiting surface brightness for each photometric band.Once the radial position of the edge is found, its correspondent value in the stellar mass density profile determines the stellar mass density at the edge (Σ edge ) parameter, using the ellipse method for galaxies with axis ratios ar >0.3 and the slit method for ar <0.3.
We display in Fig. 4, some examples of our derived profiles/quantities in order to determine the size of three MW-like galaxies (columns).There, we show the following panels (rows from top to bottom) -with the in situ star formation radius (R edge ) appearing as a dash-dotted vertical line with a different colour for each galaxy-: -The first row shows colour composite images for the galaxies (A) 9217_COSMOS, (B) 20537_EGS and (C) The first row displays their colour RGB images with their sizes (R edge ) highlighted by dot-dashed lines whose positions are also marked in the plots below.The orientation of the "slit" we used to extract the slit profiles plotted in the panels of the third row is shown in yellow.The second row shows the surface brightness profiles in V-, I-, J-and H-band using azimuthally averaged profiles ("ellipse method") while the horizontal lines denote the images' observed 3σ limiting surface brightness magnitudes in 1 × 1 arcsec 2 boxes.The legend gives the equivalent restframe central wavelength for each observed filter.The third row is similar to the second, but this time following the semi-major axis ("slit method") profiles.The observed limiting surface brightness is written explicitly this time in the legend.For the fourth row, the bottom left plot shows the g − r colour profiles and the bottom right the stellar surface mass density profiles.
A110, page 7 of 22 12887_GOODSS.These RGB images were generated, together with the rest of colour images in the paper, using the algorithm rgb-asinh in the gnuastro 9 software package (Akhlaghi & Ichikawa 2015;Akhlaghi 2019).The edge position appears as a dot-dashed line, the same that also marks its position in the panels below.The slit orientation we utilised for deriving the "slit method" profiles (see third item in this list) is shown in yellow on top of the galaxy images.
-The second row shows the observed surface brightness profiles (i.e. using the "ellipse method").The legend indicates the equivalent restframe wavelength that each filter's central wavelength is probing.The dashed horizontal lines in the bottom are the observed 3σ limiting surface brightness magnitude derived using 1 × 1 arcsec 2 apertures for each image.As expected for a radius related with star formation, edges are more clearly detected in the bluer bands.One can also easily identify the bulge and the spiral structures in these profiles.
-The third row shows the observed surface brightness profiles along the semi-major axis (i.e. using the "slit method").The legend indicates again the observed 3σ surface brightness limits.Note that, as expected the profiles are noisier with respect to the azimuthally averaged ones.However, the edges are more obvious as the irregular shape of the light distribution of these galaxies do not blur the edge features resulting from the ellipse method.
-The fourth row, bottom left position, shows the restframe g − r colour profiles.The individual profiles were corrected by galaxy inclination, Galactic extinction and cosmological dimming as indicated in Sect.3.2.While we have created all colours possible with the derived Sloan-equivalent profiles, we chose to show only g − r for the sake of clarity.They display their characteristic U-shape.When the stellar mass and surface brightness profiles do not clearly show the edge, colour profiles are key to elucidate its radial position.
-The fourth row, bottom right position, shows the averaged stellar surface mass density profiles based on the ellipse method.Using all possible combinations of base profiles and colours, we created resistant mean profiles.The stellar surface mass density (Σ edge ) at the edge position is derived from these profiles.The main parameters of our sample -coordinates, redshift, masses, sizes (R edge ) and the stellar surface mass densities at the edge positions (Σ edge )-could be found in Table 1.Its full version is included in the digital version of this article, where we also provide the surface brightness (both the observed and the Sloan rest-frame equivalent) and the stellar surface mass density profiles from our sample of galaxies.The latter (i.e. the Sloan rest-frame equivalent and the stellar surface mass density profiles) are corrected for Galactic extinction and inclination. 9https://www.gnu.org/software/gnuastro/

The size-mass relation
Figure 5 shows the size of the disc galaxies (using R edge as the size proxy due to their connection with the gas density threshold enabling star formation) versus the galaxy stellar mass for five different redshift bins.
Our local reference are the disc galaxies from the sample by Chamba et al. (2022) whose sizes are determined in a similar way than ours.Although the Chamba et al. (2022) sample covers all morphological types, we only show the local objects fulfilling a similar selection criteria as ours.This means taking galaxies with a visually clear disc component, that is, we select those local galaxies with Hubble type >−2 (i.e. from S0 galaxies to later types).This local sample is displayed with grey colour data points in Fig. 5, while the galaxies analysed in this work are coloured according to their stellar surface mass density at the radial position of the edge (i.e.Σ edge ).R edge error bars come from our spatial sampling when obtaining the surface brightness profiles (i.e. 1 kpc).The more massive disc galaxies are larger in general at any redshift.There is also a progressive decrease of the size of the galaxies with redshift.This decrease of the size of the objects is independent of their stellar mass.Similar trends are also found in the evolution of the size-mass relation when using the effective radius as a proxy for the galaxy size.
To quantify the global evolution of the size-mass relation of our galaxies, we have fitted the size-mass relation with a powerlaw and measure the change on the y-intercept with redshift.A visual inspection of the evolution of the relationship suggests (as a first approximation) that the slope of the size-mass relation does not change much with redshift.Therefore, to simplify the analysis, we fix the slope to the one found in the local sample (i.e.b R = 0.34).The fit to the local sample is always shown as a black dashed line, while the fit at each redshift range appears as a green line.The slope and the y-intercept values can be found in Table 2 with the R subscript.The errors for the y-intercept parameters, as well as the observed scatter values for the size-mass relation at the different redshifts, are derived from bootstrapping half of the sample using 10 4 different realisations.One can see that the scatters are fairly constant (∼0.1 dex) at all redshift bins.
The observed scatter for the size-mass relation at all redshifts is the result of both the intrinsic scatter of the relation plus the observational uncertainties at measuring both the size and stellar mass of the galaxies.To quantify which fraction of the observed scatter is due to uncertainties at the stellar mass determination, we follow the strategy depicted in Trujillo et al. (2020).We place our galaxies in the fitted size-mass relation at a given redshift (i.e. at each stellar mass, we associate a size corresponding to the fitting line of that size-mass relation).Then, we randomly generate artificial samples by allowing a stellar mass variation following a Gaussian distribution with a scatter of 0.13 dex (i.e. the error in our stellar mass determination A110, page 8 of 22  2022) sample and the green solid line depicts the linear regression at each redshift (fixing the slope to the local value).The galaxies in our sample are colour coded based on their stellar surface mass density at the edge position (Σ edge ).Individual error bars in our data points come from our spatial sampling when obtaining the surface brightness profiles (i.e. 1 kpc).A progressive departure from the nearby Universe values towards smaller sizes at higher redshift is noticeable.At a fixed stellar mass, more compact galaxies have larger Σ edge .
Table 2. Best fitting values using least squares fits of the parameters of the size-mass and mass-Σ edge relations using the linear form (in the log-log space) a + bx.Notes.The subscript R stands for R edge and Σ for Σ edge respectively.We also add observed scatters of the relationships.
according to the comparison with CANDELS 10 ).Then we measure the scatter of these artificial samples.The scatter associated with the uncertainty in the stellar mass (for each redshift) appears in Table 3 column σ R,mass .Of course, this uncertainty in the stellar mass has also an impact in the Σ edge -mass relation and it appears in Table 3 column σ Σ,mass .We add the values for the local sample for completeness.The other source of observational scatter in the size-mass relation is the uncertainty at measuring the position of the edge of the galaxies.This effect is not straightforward to model as it depends on the increasingly lower spatial resolution at higher redshift and the spatial binning we have used (1 kpc) to increase the signal-to-noise of the surface brightness profiles.Again, to have a rough estimation of this error, as we did for measuring the observed scatter due to the stellar mass uncertainty, we place our galaxies in the fitted size-mass relation at a 10 Although the comparison with the GS43 sample suggests a lower error on the stellar mass determination (σ = 0.10 dex), we have decided to use the most conservative value of 0.13 dex found at comparing with the CANDELS sample.
given redshift.Then, we randomly generate artificial samples by allowing a size (R edge ) variation following a Gaussian distribution with a scatter of 1 kpc because our inherent spatial discretisation.The contribution to the observed scatter in the R edge -mass and Σ edge -mass relation can also be found in Table 3. There, the intrinsic scatter after correcting quadratically by both the size and mass uncertainty is given by columns σ R,int and σ Σ,int .
Taking into account that we have not modelled all the possible effects that can affect the location of R edge , as the uncertainty in the background removal of the images, the reported intrinsic scatter of our relationships should be considered as upper limits (see a comprehensive discussion in Stone et al. 2021).
in the stellar surface mass density profiles at the R edge position.
As it was shown previously for the size-mass relation, the Σ edgemass relation of the disc galaxies is well described by a power law at all redshift bins.The least-squares fit for the Chamba et al. (2022) sample is the black dashed line, while the green one denotes the fits for the higher-z samples with the slope fixed to the local value (b Σ = 0.32).The fits can also be found in Table 2 -those with the Σ subscript-, and are computed in the same fashion as those for Fig. 5.The evolution of Σ edge with cosmic time is as follows: first, more massive galaxies present larger Σ edge values and, second, there is a progressive departure from the nearby Universe values towards higher densities as the redshift increases.This last trend is even more noticeable than for Fig. 5 since the shift is larger than the one for galaxy sizes.In addition, the colour gradient in our plots shows that at fixed stellar mass, the larger the stellar mass density at the radial position of the edge, the smaller the galaxy size.The observed and intrinsic scatters (σ Σ,mass and σ Σ,int ) are a factor of ∼3 larger than the ones measured in the size-mass relation.As it is the case with the size-mass relation, the observed scatters are also fairly constant with redshift.
Because of the cosmological dimming one would expect that there would be an observational bias against observing the galaxy edges as the redshift of our galaxies grows.That should produce an artificial trend favouring only the detection of galaxy edges that are located at higher stellar densities (i.e.shorter radial distances) as the redshift increases.To explore whether this is affecting our results, we have estimated what is the minimum stellar mass density that we can reliably explore at the different redshifts.By reliable here we mean (based on our experience with galaxies in the local Universe, see e.g.Chamba et al. 2022) that the surface brightness at the edge of the galaxy is about 2 magnitudes brighter than the limiting surface brightness of the data used to identify the boundaries of the galaxy.We have used two different approaches, one where we assume the age of the stellar population at the galaxy edge is fixed and the same at all redshift, and another where the age of the stellar population is becoming older as the cosmic time progresses.We start both approaches by assuming an age of 2 Gyr and a metallicity of [M/H] = −0.4 at z = 1.The reason for taking a subsolar metallicity is because we are representing the stellar populations at the edge of the galaxies.Taking into account our observed surface brightness limits at each redshift, if the stellar populations at all redshift were 2 Gyr (because there is a continuous star formation rejuvenating the stellar population), then we should be able to detect the edge of the galaxies up to stellar mass density of ∼1.3 M pc −2 at z = 1 and ∼0.6 M pc −2 at z = 0.5.This is below the values we observe which at those redshifts are larger than 2 M pc −2 (see Table 4).In the case where the galaxy does not form more stars at redshift lower than 1, then our detection limits would be very similar with the limit for detecting stellar mass density at z = 0.5 being ∼0.9 M pc −2 .
In addition to the previous analysis, we have conducted a large number of tests to explore the robustness of the evolution of the size-mass and stellar surface mass density-mass relations we have found.These tests are all illustrated with figures in the appendices.Here we describe what we have done.First, we have explored whether using edge-on discs we find a similar evolution of the size-mass relation with redshift.While using edge-on disc has some disadvantages as the effect of the internal dust is larger, it has also some important advantages as the location of the edge of the galaxy is easier to spot due to the larger contrast (higher brightness) of the disc due to the line-of-sight integration.We show in Appendix B that the sizes of the edge-on discs are compatible with the sizes derived for disc galaxies with lower inclinations reinforcing our findings.
In Appendix C, we explore what is the effect of having galaxies with different types of inclinations in our sample.In particular, we cut out all the high inclined discs by removing all the galaxies with axis ratio lower than 0.3.This value is motivated by the fact that the inclination correction we have applied is based on a model with a fixed z 0 /h ratio (see Fig. 1 of Trujillo et al. 2020).Only for axis ratio lower than 0.3, the assumed z 0 /h ratio can play a significant role on the inclination correction.Our relationships do not get affected by these removal of highly inclined galaxies.This is reported on Table 4.There we show the R edge and Σ edge values for a MW-like galaxy (i.e.M stellar ∼ 5×10 10 M ) according to the fits to these new filtered relationships.Those values (with the subscript ar) are all within the 1σ uncertainty with respect to the total sample.
Another test that we have conducted is to explore the effect of a potential contamination of spheroid-like galaxies in our sample (i.e.objects with Sérsic index n > 2.5) than the visual selection of our sample can be misidentified as discs.The removal of these objects do not affect either our main results (see Appendix D and Table 4 parameters with the subscript n).
Finally, and connected with our first test, we have probed whether there is any correlation between image depth and galaxy size in Appendix E. To conduct this test, we have taken advantage of the fact that the galaxies in our sample come from five different cosmological fields which are observed with different exposure times (i.e.depth).In Fig. E.1, galaxies are colour coded according to the depth of the images.Redder colours indicate galaxies from the deeper fields (GOODS South and North) while bluer colours correspond to the shallower fields (COS-MOS, UDS and EGS).In the common stellar mass interval (3 to 7 × 10 10 M ) for all redshift bins, we have probed whether the null hypothesis (i.e. that the samples are identical) can be rejected.By running a Kolmogorov-Smirnov test, we do not find any significant statistical evidence for different distributions A110, page 10 of 22 Table 4. Values of R edge and Σ edge parameters derived from our fits to galaxies with the stellar mass of the MW (i.e.∼5 × 10 10 M ).Notes.Columns: (1) Redshift bin, (2) Number of galaxies in each bin, (3) R edge , (4) Σ edge , ( 5)-( 7) same as before but for the axis ratio filtered sample (i.e.excluding objects with axis ratio <0.3), ( 8)-( 10) same as before but for the Sérsic index filtered sample (i.e.excluding objects with n > 2.5).
among the deep and shallow samples.All the above tests give us confidence to conclude that the evolution for the R edge -mass and Σ edge -mass relationships are not affected by any obvious systematic.
For the sake of comparison, we also included in Appendix F the effective radius-mass relation for the galaxies of our sample.This relation is the one that has been used in the last years as the standard for exploring the size evolution of the galaxies.We have coloured the galaxies according to the size as provided by the radial position of the edge (i.e.R edge ).As we have done with R edge as a proxy for galaxy size, we have quantified the evolution of r e by fixing the slope to the local value and fitting the observed relations with a power law (green lines).As previously found in the literature, the size evolution of the disc galaxies, using as proxy for the size the effective radius, is significantly smaller than using R edge .Note also how the scatter significantly grows at large redshift: from 0.15 dex in the local sample to 0.27 dex in the highest redshift bin.There is a weak tendency of galaxies with larger effective radius to be also larger using our size definition.

The cosmic size evolution of disc galaxies
In order to illustrate the size evolution of the disc galaxies since z ∼ 1, we have quantified the size (R edge ) evolution at a fixed stellar mass using two independent ways.In the first case, we use the global fitting to the size-mass relations as a reference.We derive the size of a galaxy with a stellar mass of the MW (i.e.M stellar ∼ 5 × 10 10 M ) according to the best fitting result.The cosmic size evolution of a disc galaxy with such stellar mass is given in the left panel of Fig. 7 (solid blue line).The other way to illustrate the cosmic size evolution is by calculating the average size of all the galaxies in a common stellar mass bin (3×10 10 < M/M < 7×10 10 ) for the different redshift bins.Note that we do not use the entire mass range as the different redshift bins (which covers a different cosmic volume) have galaxies with a different mass distribution.To avoid a bias towards low mass (at low redshift) or high mass (at high redshift) we stick to the previous (and common) stellar mass range at all redshifts.The size evolution of the galaxies in that mass bin is shown in the left panel of Fig. 7 (dashed blue line).We repeat the same exercise for the value of the stellar surface mass density (Σ edge ) at the radial position of the edge.The evolution of the above two quantities is similar (within the error bars) using the two different approaches.
The typical size evolution of a disc galaxy using R edge as a size proxy is a factor of ∼2 since z = 1 (i.e. from 12 kpc some 8 Gyr ago to 24 kpc for a present-day MW-like galaxy).Conversely, the stellar surface mass density at the edge position increases an order of magnitude with redshift.To be more precise, Σ edge increases from a value of ∼1 M pc −2 at z = 0 to ∼13 M pc −2 at z = 1.Disc galaxies with the stellar mass of the MW were fundamentally different 8 Gyr ago than nowadays: the extension of their star forming discs were a factor of 2 smaller and their stellar mass surface densities at the end of the star forming disc were a factor of 10 higher.
In the right panels of Fig. 7, we quantify the galaxy size evolution using both R edge and r e .We parameterized this evolution as a function of the Hubble parameter H(z) (top plot, solid lines) and the inverse of the cosmic scale factor, i.e. (1 + z) (bottom plot, dotted lines).The evolution of the galaxy size is calculated using the results of the fitting to the entire population (i.e. the values shown in the solid blue line in the left plot).The red data points correspond to the evolution of the effective radius published by van der Wel et al. (2014; their Table 1 and Fig. 6) for their star-forming sample.The purple points show the evolution of the effective radius but this time only using the galaxies analysed in this work (OS = Our Sample).
The relative evolution of the effective radius of the galaxy sample selected by their star formation activity (i.e.van der Wel et al. 2014) -red lines -and the effective radii for our morphologically selected sample -purple lines -is similar.However, the effective radius of the former sample is around a factor of 2 larger in r e than the visually classified as discs.There might be a number of reasons for this discrepancy.First, the stellar mass range and the photometric masses assumed in both samples are not the same.Second, the colour-selected (i.e.star-forming) galaxies probably have prominent discs and therefore, the effective radius is larger than for objects that could host more passively evolving populations.
The main message of Fig. 7 is that the size evolution of the disc galaxies using R edge as the proxy for galaxy sizes is stronger than the observed evolution of their effective radius.Another important outcome of the size evolution we have found is that it can be very well represented as a simple function of the evolution of the scale factor (i.e.(1 + z) −1 ).In fact, we find a best fitting function of (1 + z) −1.04±0.03 .This increase by a factor of 2 in the size of the MW-like galaxies since z = 1 is more in line with the expected size evolution of disc galaxies based on simple theoretical predictions (see e.g.Mo et al. 1998).

The disc growth of MW-like galaxies
We can provide a rough estimate of the speed of the radial evolution of the star forming radius (i.e. the in-situ growth of the A110, page 11 of 22    7. Left column: Size evolution of disc galaxies using as a proxy the radial location of the edge of the star forming discs (R edge ; see left blue axis).We also show the evolution of the stellar mass surface density at the edge (Σ edge ; see right green axis).The solid lines show the evolution corresponding to a disc galaxy with the stellar mass of the MW (i.e.M stellar ∼ 5×10 10 M ) according to the best linear fitting to the entire population (see text for details).The dashed lines correspond to the average properties of galaxies within the stellar mass range 3 × 10 10 < M/M < 7 × 10 10 .The redshift of each data point corresponds to the median value of the galaxy redshifts for each redshift bin.MW-like galaxies increase on average their sizes by a factor of ∼2 since z = 1.The stellar mass surface density at the end of their star forming discs (Σ edge ) decreased by an order of magnitude since that epoch.Right column: Evolution of the size of the galaxies in our sample using as a size indicator either the edge of the star forming discs (R edge ; blue points) or its effective radius (r e,OS ; purple points with the subscript OS indicating the galaxies of Our Sample).In the upper panel this evolution is fitted using a power-law function of the Hubble parameter H(z) while in the lower panel we use a power-law function of the inverse of the scale factor (i.e.(1 + z)).We add for comparison the evolution of the effective radius for the general population of star forming galaxies given in van der Wel et al. (2014).A110, page 12 of 22 discs) with cosmic time under the assumption that disc galaxies do not evolve significantly in stellar mass since z = 1.Based on the difference between R edge nowadays and at z = 1, the average speed growth we derive for a disc galaxy with the stellar mass of the MW (∼5 × 10 10 M ) is 1.60 ± 0.45 kpc Gyr −1 .For individual galaxies, which certainly are evolving in stellar mass since z = 1, the above value is a lower limit to the average speed growth of the star forming disc.
Interestingly, the above number can be compared with the estimation of the growth speed for two present-day MW-like galaxies (NGC 4565 and NGC 5907).Martínez-Lombilla et al. (2019) found a present-day growth for these two galaxies of 0.6−1 kpc Gyr −1 .These values are within a factor 2 to 3 the crude estimation we have done in the present work.
Finally, it is worth comparing the size evolution we found in this work with previous attempts of measuring this quantity using break signatures on the surface brightness profiles.Trujillo & Pohlen (2005) and Azzollini et al. (2008) did that exercise using the break radius.It should be noted that at the time these works were carried out, the break radius and the truncation of the profiles were interchangeably called the same.Currently, we know the break radius reflects a signature well inside the star-forming discs (normally associated to the end of the main spiral arms) while the edge of the star forming discs corresponds to a decline on the number of stars outside the star forming radius (see e.g.Martín-Navarro et al. 2012).According to Azzollini et al. (2008) the location of the break radius evolved by 0.4 kpc Gyr −1 .Most of the galaxies in Azzollini et al. (2008) sample are less massive than the ones studied here.Therefore, it is not straightforward to conclude whether the evolution of break and the edge radius is similar or different with cosmic time.The analysis of that issue is beyond the scope of this paper.

Conclusions
We present the size evolution since z = 1 for MW-like galaxies using a physically motivated size proxy based on the radial extension of the in-situ star formation (Trujillo et al. 2020).Following the methodology prescribed by Chamba et al. (2022), to estimate such extension, we look for truncations on the stellar surface mass density and surface brightness profiles of the galaxies which are indicative of a sudden drop on the star formation activity at that radial distance.To carry out this work, we have used the deepest available HST cosmic field observations.We have explored a total of 1048 massive (M stellar > 10 10 M ) disc galaxies with high-quality spectroscopic redshifts.
To characterise the galaxy sizes (R edge ) and the stellar mass surface density at the position where the edge is located (Σ edge ), we have conducted a careful detailed analysis of the galaxy surface brightness profiles (both azimuthally averaged and along the semi-major axis).In addition to the surface brightness profiles in all the bands available for these galaxies, we have extracted colour and stellar mass density profiles.The quantities we provide are corrected by Galactic dust extinction, inclination and cosmological dimming.
By using R edge , we find that disc galaxies with stellar masses similar to MW (i.e.∼5 × 10 10 M ) have increased their size by a factor of 2 since z = 1.The stellar mass surface density at the position of the edge, on the opposite, has decreased by a factor of 10 since that epoch.This decrement by an order of magnitude is potentially linked to the same decrease in the star formation rate since that time (and references therein Madau & Dickinson 2014).The lowering in the star formation activity with cosmic time would decrease the released energy on the available gas, therefore allowing lower gas densities to collapse and form stars progressively farther away in the outer parts of the discs.The strong evolution of the stellar mass surface density at the edge of the disc galaxies since z = 1 contrasts sharply with the evolution of the central regions, which remain essentially unchanged since that time (see Appendix G).
The picture that emerges of the size evolution of disc galaxies by using R edge is significantly different from the one using the effective radius.For instance, there have been recent claims about a negligible size (i.e.effective radius) evolution for MWlike galaxies in the last 10 Gyr (Hasheminia et al. 2022).These results are not in contradiction with our findings.This is because the effective radius is driven by the galaxy light concentration (i.e. the presence or not of a prominent bulge) and does not reflect the intuitive extension of the galaxy (cf.its disc) as R edge does.
There are, of course, a number of caveats in the present study.First, we have not corrected the effect of internal dust obscuration in the surface brightness profiles of our galaxy sample.This could affect our determination of the total stellar mass and the stellar mass surface density at the edge, especially for highly inclined galaxies.We stress, however, that at least for the local MW-like galaxies seen edge-on such as NGC 4565 and NGC 5907, the amount of dust at the edge position is not significant (compared to the internal regions, see e.g.Martínez-Lombilla et al. 2019).Second, our sample is composed of (fairly massive) galaxy discs that, by the way they are selected, avoiding interactions and contamination by luminous neighbouring objects, might be biased towards low-density environments.Finally, even though we make use of HST observations, the effect of the PSF complicates the accuracy at estimating the size and axis ratio (0.2 PSF FWHM in the reddest bands translates into 1.4 kpc spatial resolution at z = 1) for the smaller and higher redshift galaxies in our sample.
We would like also to remark the finding of a large number of massive (M stellar > 5 × 10 10 M ) and compact (R edge < 10 kpc) discs found at z = 1 that are not found in the local Universe (see Fig. 5).These very dense mini-disc galaxies gradually disappear from our sample as the cosmic time progresses.The disappearance of these objects at low redshift remind us that, in order to characterise the size evolution of individual objects, we need to understand in detail the so-called progenitor bias (new galaxies entering the selection criteria at later epochs, i.e.Carollo et al. 2013;Shankar et al. 2015;Zanisi et al. 2021).
Looking ahead to the immediate future, the very deep, highly spatially resolved and near-infrared images provided by the James Webb Space Telescope (JWST) will be instrumental for the purpose of better characterising the size of the galaxies (see some examples in Appendix A) even at redshifts beyond z = 1.Additionally, next generation synoptic facilities, such as the Euclid, Roman and Rubin telescopes, will also deliver the required depths ∼29−30 mag arcsec 2 (3σ in 10 × 10 arcsec boxes, see e.g.Euclid Collaboration 2022a,b; Martin et al. 2022) to conduct these studies with much larger statistics (see a preview of the expected results at such a depth in Trujillo et al. 2021).

Appendix A: Comparison with JWST
Deep JWST data from the GOODS-South region have recently been published (Eisenstein et al. 2023;Rieke et al. 2023).These data (of greater depth and better spatial resolution than those used in the present paper) allow us to make a quick comparison with the HST data and thus to explore the robustness of our estimate of galaxy sizes.They also give us an indication of what to expect in future work to measure the size of distant galaxies using the techniques described in this paper.
In Fig.
A.1, we show the detailed analysis of three galaxies located in the GOODS-South.All these galaxies have a similar stellar mass of M stellar ∼ 5 × 10 10 M .The colour images of the galaxies are produced using the JWST images, which have an exquisite spatial resolution, allowing us to study the internal structure of the discs of these galaxies in great detail.In addition to the images, we show the surface brightness profiles using both the ellipses and the slit along the major axis.To facilitate comparison with the HST data, we show only the surface brightness profile in the F775W (HST) and F090W (JWST) bands due to their wavelength and depth proximity in this particular field.We also show the rest-frame g-r colour using the set of bands from each telescope separately and the surface stellar mass profile.
In each panel, we show the location of R edge using only the information contained in each panel.We also do this separately for HST and JWST.The superior data quality of JWST allows us to identify the galaxy edges more clearly.This is particularly evident in the stellar surface mass density and colour profiles.The increase in spatial resolution (and therefore narrower PSF) of JWST relative to HST at a given wavelength (in this case around the observed 8500Å) reduces the amount of light in the outer part of the galaxy, making the surface brightness and stellar mass density profiles steeper beyond R edge .As a result, the h1/h2 ratio (i.e. the ratio between the inner and outer disc slopes) is larger in JWST than in HST, and the R edge is easier to identify.In general, we find good agreement in the estimation of R edge using both the different methods and facilities.The typical uncertainty in these examples is around ±1 kpc.Comparison between HST and JWST datasets for a sample of three MW-like galaxies with similar stellar masses (M stellar ∼ 5 × 10 10 M ) at different redshifts (from left to right): A) 17058_GOODSS at z = 0.6218 in blue, B) 12634_GOODSS at z = 0.6684 in magenta, and C) 12887_GOODSS at z = 0.997 in dark yellow.The first row shows their colour RGB images (based on JWST data; 70×70 kpc), with their sizes (R edge ) highlighted by dotted dashed lines, whose positions are also marked in the plots below.The orientation of the "slit" we used to extract the profiles plotted in the panels of the third row is shown in yellow.The second row shows the surface brightness profiles in F775W (HST) and F090W (JWST) using azimuthally averaged profiles ("ellipse method").The third row is similar to the second, but this time following the semi-major axis profiles ("slit method").The fourth row shows the g − r colour profiles derived with each telescope, and the fifth row shows the stellar surface mass density profiles based on HST and JWST.In each panel, the vertical dashed lines indicate the position of R edge according to the methodology proposed in this paper.The coloured vertical regions indicate the uncertainty at measuring the position of the edge of the galaxy using independently the information given in each panel.

Fig. 1 .
Fig. 1.Set of HST filters used in this study.The y-axis refers to the filters' transmission (left) and also the position where the restframe SDSS central wavelengths move with redshift (right).Dashed filters are only available for the GOODS-South field.

Fig. 2 .
Fig. 2. MW-like (M stellar ∼ 5 × 10 10 M ) galaxies at six different cosmic epochs in our study.Each individual stamp has a 70 × 70 kpc size.Note that the background for each image appears in black and white in order to highlight the level at which the noise starts to dominate.The lowest redshift galaxy belongs to our Local Universe reference sample (Chamba et al. 2022) based on the IAC Stripe 82 Legacy Project(Trujillo & Fliri 2016), while the rest of images are HST data.These objects were selected to have representative sizes (R edge , golden ellipses) of galaxies at those redshifts.Note how the sizes (using our in-situ star formation criterion) decrease with redshift, while the effective radii (r e , blue ellipses) do not change appreciably.

Fig. 3 .
Fig. 3. Comparison between our photometrically derived stellar masses and those in the CANDELS and ASTRODEEP GS43 catalogues (∆mass = log(M ours )-log(M catalog )).The dashed lines represent the best Gaussian fits (their scatter is also displayed).The arrows with the sigma values are displaced vertically for the blue distribution for the sake of clarity.

Fig. 4 .
Fig. 4. Summary of our analysis for three MW-like galaxies at different redshifts (from left to right): (A) 9217_COSMOS at z = 0.44 in blue, (B) 20537_EGS at z = 0.69 in magenta, and (C) 12887_GOODSS at z = 0.997 in dark yellow.The first row displays their colour RGB images with their sizes (R edge ) highlighted by dot-dashed lines whose positions are also marked in the plots below.The orientation of the "slit" we used to extract the slit profiles plotted in the panels of the third row is shown in yellow.The second row shows the surface brightness profiles in V-, I-, J-and H-band using azimuthally averaged profiles ("ellipse method") while the horizontal lines denote the images' observed 3σ limiting surface brightness magnitudes in 1 × 1 arcsec 2 boxes.The legend gives the equivalent restframe central wavelength for each observed filter.The third row is similar to the second, but this time following the semi-major axis ("slit method") profiles.The observed limiting surface brightness is written explicitly this time in the legend.For the fourth row, the bottom left plot shows the g − r colour profiles and the bottom right the stellar surface mass density profiles.

Fig. 5 .
Fig. 5. Galaxy size (using as a proxy the expected radial location of the end of the in situ star formation, i.e.R edge ) versus galaxy stellar mass for our galaxy sample.The black dashed line represents the linear regression (in log-log plane) to the local values (in grey, visible at all redshifts) using the Chamba et al. (2022) sample and the green solid line depicts the linear regression at each redshift (fixing the slope to the local value).The galaxies in our sample are colour coded based on their stellar surface mass density at the edge position (Σ edge ).Individual error bars in our data points come from our spatial sampling when obtaining the surface brightness profiles (i.e. 1 kpc).A progressive departure from the nearby Universe values towards smaller sizes at higher redshift is noticeable.At a fixed stellar mass, more compact galaxies have larger Σ edge .

Fig. 6 .
Fig. 6.Stellar surface mass densities (Σ edge ) at the radial position of the edge versus galaxy stellar mass for our galaxy sample.The black dashed line represents the linear regression to the local values from Chamba et al. (2022; in grey, visible at all redshifts) and the green solid line depicts the linear regression fitting at each redshift.The data points are coloured according to their size (R edge ) and the Σ edge error bars come from the error bars in the stellar mass density profiles at the R edge position.A progressive departure from the nearby Universe values towards higher densities is appreciable.The data points are colour coded according their galaxy size.At a fixed stellar mass, those galaxies with the larger stellar surface mass densities are those smaller in size.

Fig.
Fig. 7. Left column: Size evolution of disc galaxies using as a proxy the radial location of the edge of the star forming discs (R edge ; see left blue axis).We also show the evolution of the stellar mass surface density at the edge (Σ edge ; see right green axis).The solid lines show the evolution corresponding to a disc galaxy with the stellar mass of the MW (i.e.M stellar ∼ 5×10 10 M ) according to the best linear fitting to the entire population (see text for details).The dashed lines correspond to the average properties of galaxies within the stellar mass range 3 × 10 10 < M/M < 7 × 10 10 .The redshift of each data point corresponds to the median value of the galaxy redshifts for each redshift bin.MW-like galaxies increase on average their sizes by a factor of ∼2 since z = 1.The stellar mass surface density at the end of their star forming discs (Σ edge ) decreased by an order of magnitude since that epoch.Right column: Evolution of the size of the galaxies in our sample using as a size indicator either the edge of the star forming discs (R edge ; blue points) or its effective radius (r e,OS ; purple points with the subscript OS indicating the galaxies of Our Sample).In the upper panel this evolution is fitted using a power-law function of the Hubble parameter H(z) while in the lower panel we use a power-law function of the inverse of the scale factor (i.e.(1 + z)).We add for comparison the evolution of the effective radius for the general population of star forming galaxies given in van derWel et al. (2014).
Fig. A.1.Comparison between HST and JWST datasets for a sample of three MW-like galaxies with similar stellar masses (M stellar ∼ 5 × 10 10 M ) at different redshifts (from left to right): A) 17058_GOODSS at z = 0.6218 in blue, B) 12634_GOODSS at z = 0.6684 in magenta, and C) 12887_GOODSS at z = 0.997 in dark yellow.The first row shows their colour RGB images (based on JWST data; 70×70 kpc), with their sizes (R edge ) highlighted by dotted dashed lines, whose positions are also marked in the plots below.The orientation of the "slit" we used to extract the profiles plotted in the panels of the third row is shown in yellow.The second row shows the surface brightness profiles in F775W (HST) and F090W (JWST) using azimuthally averaged profiles ("ellipse method").The third row is similar to the second, but this time following the semi-major axis profiles ("slit method").The fourth row shows the g − r colour profiles derived with each telescope, and the fifth row shows the stellar surface mass density profiles based on HST and JWST.In each panel, the vertical dashed lines indicate the position of R edge according to the methodology proposed in this paper.The coloured vertical regions indicate the uncertainty at measuring the position of the edge of the galaxy using independently the information given in each panel.

Table 1 .
Galaxy sample main parameters.
Notes.The full table is available at the CDS.

Table 3 .
Expected contributions to the observed scatters for the R edge -mass and Σ edge -mass relations due to different observational errors.