The Dark Energy Survey Supernova Program: Investigating Beyond-Λ CDM

We report constraints on a variety of non-standard cosmological models using the full 5-year photometrically-classified type Ia supernova sample from the Dark Energy Survey (DES-SN5YR). Both Akaike Information Criterion (AIC) and Suspiciousness calculations find no strong evidence for or against any of the non-standard models we explore. When combined with external probes, the AIC and Suspiciousness agree that 11 of the 15 models are moderately preferred over Flat-Λ CDM suggesting additional flexibility in our cosmological models may be required beyond the cosmological constant. We also provide a detailed discussion of all cosmological assumptions that appear in the DES supernova cosmology analyses, evaluate their impact, and provide guidance on using the DES Hubble diagram to test non-standard models. An approximate cosmological model, used to perform bias corrections to the data holds the biggest potential for harbouring cosmological assumptions. We show that even if the approximate cosmological model is constructed with a matter density shifted by ΔΩ m ∼ 0 . 2 from the true matter density of a simulated data set the bias that arises is sub-dominant to statistical uncertainties. Nevertheless, we present and validate a methodology to reduce this bias.


INTRODUCTION
Our understanding of the Universe fundamentally changed in the late 1990s with the remarkable discovery that the expansion of our Universe is accelerating (Riess et al. 1998;Perlmutter et al. 1999).This discovery established ΛCDM as the standard model of cosmology, which asserts that the Universe at late times is dominated by dark energy in the form of a cosmological constant Λ and cold (nonrelativistic), pressureless dark matter (CDM).However, the nature of dark energy remains a mystery.
★ E-mail: uqrcamil@uq.edu.au(RC) In this paper we use the complete photometrically-classified type Ia supernova (SN Ia) data set from the Dark Energy Survey (DES)which represents the largest, most homogeneous SN data set to date -to assess whether the latest SN Ia data prefers any non-standard cosmological models over ΛCDM.
While ΛCDM fits most data well, it lacks a physical motivation and is currently unable to alleviate tensions between early time and late time measurements of parameters such as the current expansion rate of the Universe,  0 (Aghanim et al. 2020;Riess et al. 2022).These two limitations have led to a wealth of exotic cosmological models being proposed (see Di Valentino et al. 2021, for a review).
Non-standard cosmological models attempt to explain observa-tions in a variety of ways, ideally with some physical justification.Models that mimic the late time acceleration include dynamical vacuum energy, cosmic fluids, scalar fields as well as modifications to the theory of general relativity.Other models challenge our assumption of large-scale homogeneity and isotropy, and attribute the dimming of distant supernovae to local spatial gradients in the expansion rate and matter density, rather than due to late time acceleration (Alonso et al. 2010).
Previous analyses have shown that many non-standard models are able to explain the current data (e.g.Davis et al. 2007;Sollerman et al. 2009;Li et al. 2011;Hu et al. 2016;Dam et al. 2017;Zhai et al. 2017;Lovick et al. 2023), although none have shown strong improvement over ΛCDM.In general non-standard models have only been a good fit to the data if they are able to mimic the expansion history of ΛCDM for some choice of parameters.These analyses conclude that new, more statistically powerful data, across a wide range of cosmological observations are required to discriminate between models.
The Dark Energy Survey was designed to provide such data and to reveal in detail both the expansion history and large-scale structure of the Universe.Type Ia supernovae are one of the four pillars of DES science, the others being baryon acoustic oscillations (BAO; DES Collaboration et al. 2024a), galaxy clustering (Rodríguez-Monroy et al. 2022;Porredon et al. 2022;Pandey et al. 2022), and gravitational lensing (Gatti et al. 2021;Amon et al. 2022;Secco et al. 2022).
In this paper we focus on the DES-SN5YR sample (DES Collaboration et al. 2024b) containing 1829 SNe.The DES-SN5YR sample consists of 1635 SNe from the the full five years of the DES survey, of which 1499 have a machine learning probability of being a type Ia larger than 50 per cent and range in redshift from 0.10 to 1.13.This is combined with an external sample of 194 spectroscopically confirmed low- SNe Ia (see Section 6).
Our work builds on previous analyses of non-standard models in two ways.(1) we carefully analyse any cosmological assumptions and approximations that have gone in to the derivation of the information that appears in the Hubble diagram, and estimate their impact.We also provide a prescription for others who would like to use DES SN data to test their own non-standard models, and to provide confidence that there are no hidden assumptions that could bias their result.(2) We constrain a set of non-standard models using the DES-SN5YR sample, with the aim of providing the tightest constraints using SNe Ia measurements alone.
This paper is organised as follows.In Section 2 we describe the cosmology pipeline used to produce a Hubble diagram, focusing on aspects of the pipeline that contain any cosmological model dependence.In Section 3, we introduce a new parameter,   , that can be used as a single-number to summarize supernova cosmology constraints in the -Ω m plane.This new parameter is useful for testing the impact of the reference cosmology used in our simulated bias corrections in Section 4 and the fiducial cosmology used while determining the standardised magnitudes of SN Ia in Section 5. Section 6 describes the data (SN and other external data sets) that we use in this analysis.In Section 7 we describe the theory behind the cosmological models we test and present our results.We discuss our results in Section 8 and conclude in Section 9.

COSMOLOGY PIPELINE
Here, we focus on some areas of the DES-SN5YR baseline analysis described in Vincenzi et al. (2024) -all the way from light curves to cosmology -that are, or may appear to be, subject to cosmological dependencies (highlighted in red in Fig. 1).We aim to provide clarity for others who want to use the DES-SN5YR sample to fit their own models.
The pipeline, illustrated in Fig. 1, is run within the pippin framework (Hinton & Brout 2020), built around several key components including the SALT3 light curve fitting algorithm (Kenworthy et al. 2021), the superNNova photometric classifier (Möller & de Boissière 2020), snana light curve fitting and simulation for bias corrections (Kessler et al. 2009b) producing a bias-corrected Hubble diagram with the "Beams with Bias Corrections" (BBC) formalism (Kessler & Scolnic 2017).We now describe each in turn.

Light curve fitting
To convert sparse light curve observations to SN standardization parameters we use the SALT2 model framework (Guy et al. 2007;Guy et al. 2010) as implemented by the SALT3 software (Kenworthy et al. 2021).SALT3 fits for the time of B-band peak and encapsulates the SN behaviour using three parameters:  0 is the overall amplitude of the light curve;  is related to the  −  colour of the SN at peak brightness; and  1 describes the width of the light curve (stretch).For further details on the light curve fitting used on the DES-SN5YR sample see Taylor et al. (2023).
The SALT3 framework is cosmology independent, except for the assumption that light curves are time-dilated (White et al. 2024) by a factor of (1+ obs ).Note that the observed redshift is used to calculate time dilation, therefore there is no peculiar velocity correction at this stage.

SN Ia distances
The distance moduli,  obs, of SNe Ia are calculated using the modified Tripp equation (Tripp & Branch 1999),  obs, =  , +  1, −   −  host, − M − Δ bias, (1) where   = −2.5log(0 ); M, is a combination of the SN Ia absolute magnitude, , and the Hubble constant  0 , which is marginalised over (see Section 6.3); and  &  are nuisance parameters that represent the slopes of the stretch-luminosity and colour-luminosity relations respectively. is an additional nuisance parameter that accounts for a correlation between standardised SN luminosities and host-galaxy stellar mass,  * .This dependency is modelled as a massstep correction (Conley et al. 2010;Brout et al. 2019).The final term in equation ( 1), Δ bias is applied to each SN to correct for selection effects.

BEAMS with Bias Corrections
The BEAMS with Bias Corrections (BBC; Kessler & Scolnic 2017) framework returns a Hubble diagram from a photometrically 1 identified sample of SNe Ia.It does this by maximising the BEAMS likelihood (Section 2.3.1) that accounts for the probability of the SN event being a core-collapse (CC) contaminant while also incorporating bias corrections (Section 2.3.2) and determining global nuisance parameters, ,  and  from equation (1) (Section 2.3.3).Therefore, along with the Hubble diagram, BBC also outputs the fitted global nuisance parameters, the uncertainty on the estimated distance moduli,  , , and a classifier scaling factor that is introduced in Section 2.3. " ,  populations Figure 1.An overview of the light curve to cosmology pipeline.Here, an emphasis is placed on potential cosmological dependencies (red) and the flow of parameters at each stage.Note that we have also omitted parameter uncertainties and the associated covariances for clarity.However, we have included the final uncertainty term,  , which includes the intrinsic scatter and a contribution based on the probability of the SN event being a CC contaminant (see Section 2.3.1).A subscript  has been added to SN-dependent parameters.Dashed boxes represent external pippin inputs.

The BEAMS likelihood
Photometric SN samples rely on a classifier to provide a probability of each SN being type Ia or else a contaminant such as core collapse SN or peculiar SN Ia.The DES-SN5YR baseline analysis uses machine learning techniques to classify SN via the open-source algorithm superNNova (Möller & de Boissière 2020).2This classification has no cosmological dependence beyond the assumption that the light curves are time dilated by (1 +  obs, ).The predictions of these classifiers,  Ia, are incorporated into the cosmology pipeline by using the 'Bayesian Estimation Applied to Multiple Species' (BEAMS) approach (Kunz et al. 2007;Hlozek et al. 2012;Kunz et al. 2012).The BEAMS approach, involves maximising the BEAMS likelihood, which includes two terms, one that models the SN Ia population and another that models a population of contaminants.Compared to the traditional likelihood used in spectroscopic samples, the BEAMS likelihood adds one fit parameter, the  CC scaling factor  CC .The distance uncertainties are then renormalized to ensure that likely contaminants have inflated distance uncertainties and are down-weighted when fitting cosmology.For detailed descriptions of the BEAMS likelihood see Kunz et al. (2012), Kessler & Scolnic (2017) and Vincenzi et al. (2022).

Bias corrections
The BBC approach uses the BEAMS formalism, and estimates the final term in equation ( 1), Δ bias , using simulations that model the survey detection efficiency, Malmquist bias as well as other biases introduced in the analysis.Simulations of the DES-SN5YR sample are generated using snana3 (Kessler et al. 2019) where light curves are modelled using the SALT3 framework and the 'Dust2Dust' fitting code (Popovic et al. 2023) measures the underlying population of stretch and colour, including their correlations with host properties.
The simulations used for bias corrections within the baseline analysis are performed using a reference cosmology of Flat-ΛCDM with parameters (Ω  , ) ref = (0.315, −1.0).There is an underlying assumption in the BBC framework that the bias correction simulations accurately describe the intrinsic properties of the SNe Ia and selection effects.
The bias correction step thus holds the biggest potential for harbouring cosmological assumptions that could influence the cosmological results.However, the dependence on the reference cosmology has been shown to be weak for models that have similar4 evolution of magnitude versus redshift (Kessler & Scolnic 2017;Brout et al. 2019).Nevertheless, in the analysis of non-standard cosmologies that have the flexibility to deviate significantly from the standard cosmological models, this may no longer be true.In Section 4, we extend on previous work and quantify the cosmological bias resulting from more extreme reference cosmologies in the context of the DES-SN5YR baseline analysis, and provide a prescription for how to fit models that deviate from the reference cosmology significantly in their evolution of magnitude versus redshift.

BBC fit
The global nuisance parameters, ,  and  are used to standardise SN magnitudes and are determined using the BBC fitting algorithm (which has previously been referred to as SALT2mu), following the redshift binning procedure in Marriner et al. (2011) and equation 3 of Kessler & Scolnic (2017).BBC employs a fiducial cosmology5 that provides an arbitrary smooth Hubble diagram in each redshift bin.BBC fits for ,  and  by minimizing the Hubble residuals to the fiducial cosmology among   = 20 logarithmically-spaced redshift bins as well as fitting for a magnitude offset in each bin.
The default fiducial cosmology used in the BBC fit, for the DES-SN5YR analysis, is the Flat-ΛCDM model with parameters ( 0 , Ω m ) = (70, 0.3).This choice may cause confusion within the community regarding a potential cosmology dependence.Fig. 2 provides an exaggerated visualization of the BBC fit to show i) fitting for magnitude offsets in redshift bins allows the data to better resemble its naturally standardized state (with  fit ,  fit consistent with the true values); ii) the magnitude offsets (approximately) map the fiducial cosmology on to the true one by quantifying how much the observations deviate from the reference cosmology in each redshift bin; and iii) that this procedure removes the dependence on cosmological parameters.Marriner et al. (2011) show that the fit for  and  is decoupled from the choice of fiducial cosmology if the number of redshift bins is sufficiently large.Furthermore, Kessler et al. (2023) performs a limited study that looks at the standard deviation of the Hubble residuals of the BBC fit (see Table 1 of Kessler et al. 2023).In Section 5, we re-test this result and extend on the work of Marriner et al. (2011) and Kessler et al. (2023) by explicitly testing extreme cosmologies as well as showing that the impact on cosmology-fitted parameters is negligible.Finally, we present an alternate approach that does not require a fiducial cosmology and achieves consistent fits for  and .

SN Ia distance uncertainties
Following the Pantheon+ analysis (Brout et al. 2022) where  S3fit, includes the uncertainties on the light curve parameters and the associated covariances; while  lens, and  , are uncertainties associated with lensing effects and spectroscopic redshifts, respectively. (  ,   ,  * , ) and  floor (  ,   ,  * , ) are survey-specific scaling and additive factors that are estimated from the BBC simulations.Finally,  vpec, accounts for uncertainties due to peculiar velocities, including both uncertainties in linear-theory modelling and non-linear unmodelled peculiar velocities, as discussed in Sec 2.4.

Modelling peculiar velocities
The redshift that is compared to SN distances should be entirely due to the expansion of the universe.However, in practice the redshift that we measure contains contributions due to peculiar velocities of the SN and its host galaxy.The DES-SN5YR baseline analysis uses peculiar velocities presented by Peterson et al. (2022), which are determined from the 2M++ density fields (Carrick et al. 2015) with global parameters and group velocities used from Said et al.
differ from the reference cosmology used to simulate SNe used for bias corrections.
A) As  offset,i does not in general equal  offset,ii this procedure allows the data to better resemble the true cosmology (black dashed line) approximately mapping the fiducial cosmology on to the true one by quantifying how much the observations deviate from the fiducial cosmology in each redshift bin.In B), the data has been shifted to the fiducial cosmology for illustrative purposes and in C) we shift the data back.Therefore, for this example,  avg would be positive (the data actually sits above our fiducial cosmology), however  offset, and  offset, would be positive and negative respectively (because the data sits above and below the average offset respectively).While this example is exaggerated it is useful to provide insight into BBC and highlight that the method has minimal cosmological dependence.
(2020) and Tully (2015) respectively and a 240 km s −1 uncertainty on these estimates.While the determination of the peculiar velocity corrections includes a fiducial cosmology, the corrections have the largest impact at low redshifts where the cosmology dependence is negligible.Although Peterson et al. (2022) show that the impact of peculiar velocity corrections on  0 and  fits are at the 1 per cent level, the impact of the fiducial cosmology in the derivation of those corrections is negligible compared to the uncertainty in the peculiar velocity map, and therefore we do not consider it further in this work.

THE Ω m − 𝑤 DEGENERACY
There is a degeneracy between the equation of state of dark energy and the matter content of the universe for distance indicators within generalised dark energy models.It has long been known that this degeneracy makes it more difficult to assess systematics on Ω m and  separately.Large shifts in the best fit parameters may not be significant if they occur along the degeneracy direction, but the same size shifts could be very significant if they occur perpendicular to the degeneracy direction.In the DES cosmology analysis we use two methods to account for that degeneracy.The first is setting a prior on matter density 6 and only considering changes in , the other is testing a new parameter   () that allows us to present a single non-degenerate number summarising a SN Ia constraint in the -Ω m plane.
In Fig. 3 we show lines of constant   () overlaid on to the 1 and 2 contours for the DES-SN5YR sample.Since the   () parameter is redshift dependent, it is not as universal as a parameter such as  8 =  8 √︁ Ω m /0.3, which defines a quantity that is relatively independent of the  8 and Ω m degeneracy in lensing studies.Instead, we can select a redshift that matches the degeneracy direction of the sample.In the top right subplot of Fig. 3 we show that   (0.2) makes a good approximation for the -Ω m degeneracy line for the DES-SN5YR sample.Using the   (0.2) parameter, we can therefore use a single number to approximate the DES-SN5YR constraints on the Flat-CDM model and find   (0.2) = −0.340± 0.032 (which includes statistical and systematic uncertainties).
Changes to the analysis that only cause shifts along the degeneracy direction have a very small effect on   even though they can have a misleadingly large effect on Ω m and  (misleading since those shifts are strongly correlated).  is thus an excellent measure by which to evaluate the impact of analysis choices on the supernova cosmology results (see Fig. 4).

REFERENCE COSMOLOGY IN THE BIAS
CORRECTION SIMULATIONS Kessler & Scolnic (2017) show that any dependence on the reference cosmology is weak when the reference cosmology is similar to the true evolution of magnitude versus redshift (see Sec 6.  systematic and also show that using a reference cosmology even 10 away from the true cosmology has less than a 1 shift in the results.We also present an iterative method that can be used to reduce even that small systematic offset.

Testing the impact of the reference cosmology
To examine the impact that the reference cosmology used for the bias correction simulations has on our cosmology fits, we generate and analyze 25 realisations of simulated data.These are created with a Flat-CDM cosmology with parameters ( 0 , Ω m , ) = (70, 0.315, −1.0).We also generate six different BBC simulations, each with a unique reference cosmology.For comparison, in Fig. 4e we plot each reference cosmology (dashed lines) relative to the cosmology used to generate our simulated data (orange).
The average shifts in Ω m and  from the perfect scenario in which the reference cosmology is equal to the true cosmology of our simulated data are shown in Fig. 4a and Fig. 4b respectively. 7n Fig. 4f we plot the results in the  − Ω m plane for a single realisation.The contours and solid orange square are for the ideal case in which the reference cosmology matches the true cosmology.The other symbols show the results when using different reference cosmologies, where the open symbols show the input reference cosmology and the solid symbols show the resulting best fit parameters.
This shows that while the shifts in Ω m and  seem large, when viewed in 2D parameter space they all fall along the Ω m −  degeneracy direction and are thus all well within 1.
The dot-dashed line in Fig. 4f shows the   (0.2) parameter, representing the degeneracy line.Note that the ideal redshift for     to match the degeneracy direction will change depending on the data set.In Fig. 4c we plot the average shift in   (0.2) and in Fig. 4d we plot the shift in  after applying a strong prior on the matter density Ω m = Ω m,true ± 0.001.The fact that the shifts in   (0.2) and | Ω m,true ±0.001 are negligible shows that the impact of the reference cosmology is small and limited to the degeneracy direction, in agreement with the results from Kessler & Scolnic (2017).
We also performed two additional tests that are the inverse of those performed above.Instead of varying the reference cosmology, we fixed the reference cosmology to the baseline cosmology used in the DES-SN5YR analysis and generated 25 realisations of simulated data using both (a) Flat-CDM cosmology with parameters ( 0 , Ω m , ) = (70, 0.350, −0.8) and (b) Flat- 0   CDM cosmology with parameters ( 0 , Ω m ,  0 ,   ) = (70, 0.495, −0.36, −8.8).These cosmologies were chosen to match the ∼ 10 offset brown point in Fig. 4 and the best fit Flat- 0   CDM result in the DES-SN5YR analysis respectively.The results are given in Table 1.For test (a), we again find that the impact of using the incorrect reference cosmology is negligible.For test (b), we see larger shifts in cos- Table 1.Shifts in the best fit parameters using the DES-SN5YR baseline reference cosmology, from the perfect scenario where the reference cosmology is equal to the cosmology used to generate the simulated data.Here, the uncertainties are on the shift in the mean -not the uncertainty on the parameters, which is larger.

𝑎
(0.350, −0.80, 0) 0.02 ± 0.05 0.000 ± 0.008 -(0.495, −0.36, −8.8) -0.18 ± 0.06 −1.6 ± 0.4 * Model used to generate the 25 realisations of simulated data † Determined used a prior on the matter density of Ω m = Ω m,true ± 0.001.mological parameters.However, in this case, their is an additional degeneracy between  0 −   that is not accounted for when applying the prior on Ω m .To visualise this, we plot the 25 realisations in Fig. 5, which, shows that the best fit points are aligned along the degeneracy line and consistent with the truth.We also note that the uncertainties given in Table 1 are on the shift in the mean.The shifts are Δ 0 = 0.18±0.28,Δ  = −1.6±2.2 when using the uncertainty on the parameters.
In summary, this result validates that the BBC baseline approach used in DES Collaboration et al. ( 2024b) is able to return a Hubble diagram that represents the true distance versus redshift relation to within 1 even given a reference cosmology that is ∼ 10 from the truth (brown point in Fig. 4) or varies by ∼ Δ = 0.15 (brown dashed line in Fig. 4e).The apparent bias observed in Fig. 4a and Fig. 4b is due to showing shifts in degenerate parameters separately, without considering the combined influence on the distance versus redshift relation.Importantly, we can be confident in our bias corrections if the expansion history of a non-standard cosmological model falls within the region bounded by the blue and brown dashed lines in Fig. 4e.

The iterative method
Section 4.1 validates the procedure used in the DES-SN5YR baseline analysis, showing that the reference cosmology has a small impact on the cosmological results relative to the statistical uncertainties.However, the BBC reference cosmology may become a dominating systematic for future surveys such as the Rubin Observatory's LSST, which will include hundreds of thousands of well measured SNe Ia (LSST Science Collaboration et al. 2009).Furthermore, Fig. 4 shows that in the case where one finds a tension with other data sets at the extreme ends of the degeneracy direction (e.g. if the CMB contours were at the top left or bottom right in Fig. 4f), it would be beneficial to ensure a close match to the reference cosmology.Since we performed a blind analysis, we did not know whether there would be a discrepancy between the BBC reference cosmology and the final fitted cosmology results.We therefore prepared the following method to correct the reference cosmology if the discrepancy was significant.
It was suggested by Kessler & Scolnic (2017) that an iterative procedure can be applied where  ref is updated with the previous  fit value, to reduce this bias.This procedure is summarised in Fig. 6.In this work, we test the iterative method by applying it to 10 realisations of simulated data created with a Flat-CDM cosmology with parameters ( 0 , Ω m , ) = (70, 0.350, −0.8).This cosmology was selected due to its location in parameter space, which is approximately perpendicular to the Ω m −  degeneracy line in the direction of a general CMB prior and lies outside a 2 region (based on DES-SN5YR simulations) of the default BBC reference cosmology.8Table 2 shows the weighted average shift in cosmological parameters from the truth after 10 realisations.Note that the Ω m prior was only applied on our final results and was not used during the iterative process.We report both Δ| Ω m,true ±0.001 and Δ  (0.2) and find that both are closer to the truth after applying the iterative method.In particular, we find that | Ω m,true ±0.001 has shifted by 0.006 and   (0.2) has shifted by 0.008 closer to the truth.
We note a limitation of this work that we have not explicitly shown † With a prior on the matter density of Ω m = 0.350 ± 0.001 the iterative method converges (because repeatedly redoing the simulations is computationally intensive).However, we performed a third iteration on two random realisations and found that the iterative method remained stable.
The iterative method was not implemented in the current DES results, because after unblinding we found the best fit cosmology to be sufficiently close to the reference cosmology so as to make any bias insignificant (in Sec.4.1 we found Δ ∼ 0.01 given a reference cosmology 10 from the truth).Nevertheless, we conclude that iterating the reference cosmology is a viable method to reduce this bias for future analyses where the reference cosmology may become a dominating systematic.

TESTS OF COSMOLOGY DEPENDENCE WITHIN THE BBC FIT
In this section we validate the baseline analysis assumption that the fit for nuisance parameters is decoupled from the choice of fiducial cosmology using 20 logarithmically space redshift bins (for these tests we restrict ourselves to  and ).
In total, we generated 100 statistically independent realisations that resemble the DES-SN5YR sample in a spatially Flat-ΛCDM universe with parameters ( 0 ,   , Ω m ) = (70, −19.253, 0.3).We ran all 100 realisations through the entire PIPPIN pipeline six times with each run distinguished uniquely by the choice of fiducial cosmology within the BBC fitting procedure.The choice of fiducial cosmologies was chosen such that they vary significantly in the evolution of magnitude versus redshift and are shown in the bottom panel of Fig. 7.
The left panel of Fig. 8 compares the maximum likelihood  and  values for each of the 100 realisations.The top left sub-plot represents the ideal case where the fiducial cosmology is equal to the true cosmology used to simulate the data.Here, we show how the averages of the 100 maximum likelihood values (blue crosses) compared to the true values (black dashed lines).We also make the equivalent comparison after fitting for cosmological parameters, shown in the right panel of Fig. 8.In Fig. 7 we present the shifts in the average of the maximum likelihood , , Ω m ,  and   (0.20) values as a result of varying the fiducial cosmologies within the BBC fit.We also show how the shifts in cosmological parameters impacts the evolution of magnitude versus redshift relative to the ideal case.
We find that the determination of the global nuisance parameters,  and , has a weak dependence on the choice of fiducial cosmology; these results are in agreement with those by Marriner et al. (2011).Extending on the work by Marriner et al. (2011), Fig. 7 shows that the BBC fit is able to recover the ideal cosmological parameters with less than a 1 tension of the standard error given 100 realisations even when using extreme fiducial cosmologies.The two fiducial cosmologies that result in the largest shift in cosmological parameters are unsurprisingly also the two cosmologies that deviate the most in the slope of the distance versus redshift relation ( 0 , Ω m , Ω Λ , ) = (70, 0.0, 1.0, −1.0) and (70, 1.0, 0.0, −1.0).However, both the   (0.2) panel and Fig. 8 show that shift is along the degeneracy direction.
Finally, the lower right of Fig. 7 shows the  differences between the fiducial cosmologies (dashed lines) and even shifts of  up to 0.5 across the -range have negligible impact on the best fit expansion history (solid lines).

Is a fiducial cosmology required?
Often, the role of the fiducial cosmology within the BBC fit causes confusion -both because of perceived cosmology dependence (which we have shown is negligible for any reasonable cosmology in Section 5) and because it is mistaken for the reference cosmology used to generate the BBC simulations that estimates the  bias term in equation (1).
Here, we explore replacing the fiducial cosmology (along with the fitted magnitude offsets in each bin) within the BBC fit with a spline interpolation of the SN magnitudes.To accomplish this, we modify the BBC procedure.Recall that within the current BBC procedure the Hubble residuals are minimized to a fiducial cosmology among 20 independent redshift bins, given a set of global nuisance parameters and 20 offsets in magnitude.Here, we instead minimise the Hubble residuals to a spline interpolation of the SN magnitudes, determined at each fitting step, where we used the weighted average redshift,  avg and distance moduli,  avg in 20 redshift bins as knots.
We compare these two procedures by recreating a simplified BBC fitting procedure that attributes all of the intrinsic scatter to coherent variation at all epochs and wavelengths,  int .9Further complexity is not required as the intrinsic scatter is incorporated into the uncertainties in the same way if we use a fiducial cosmology or a spline and we only need to test consistency between the two methods. 10able 3 compares the fitted nuisance parameters using the same light curve sample when using two different fiducial cosmologies (see Table 3 for model parameters) and a spline that is determined at each fitting step.All parameters are consistent demonstrating the following.First, that the results from our simplified BBC fit are again insensitive to the choice of fiducial cosmology.Second, that a spline is viable alternative to a fiducial cosmology and may reduce confusion as to the role of the fiducial cosmology in future pipelines.

DATA
Having established that the derivation of the DES-SN5YR Hubble diagram is robust to the choice of reference and fiducial cosmological models, we turn to using the Hubble diagram to derive constraints on a range of non-standard models which differ in their background expansion and are therefore sensitive to the DES-SN5YR data.To test the non-standard cosmology fitting code, we generated 25 simulations and ensured that fitted parameters of each model were consistent with the input cosmology.The input cosmology for these simulations used Flat-ΛCDM, for models that could reduce to Flat-ΛCDM for some values of their parameters.Otherwise, we used the model being tested as the input cosmology to generate the 25 realisations.

The DES-SN5YR sample
The DES-SN survey covers ∼27 deg 2 over 10 fields across the DES footprint (see Smith et al. 2020).The survey, which ran for five years using the Dark Energy Camera (DECam; Flaugher et al. 2015).DES detected over 30,000 SN candidates, from these 1635 were deemed SNe Ia-like and included in the DES-SN5YR Hubble diagram with 1499 photometrically classified as type Ia SNe using SuperNNova (Möller et al. 2022;Vincenzi et al. 2024).The DES-SN5YR sample includes publicly available low- SNe Ia from the Harvard-Smithsonian Center for Astrophysics, CfA3 (Hicken et al. 2009) and CfA4 (Hicken et al. 2012), the Carnegie Supernova Project, CSP (Krisciunas et al. 2017) and the Foundation Supernova Survey (Foley et al. 2018).These low- samples span a redshift range of 0.01 to 0.1.However, SNe Ia in the low- sample with redshifts < 0.025 are excluded to minimise the impact of peculiar velocities.With this cut applied, the low- sample comprises 194 SNe, for a total of 1829 SNe in the DES-SN5YR sample; for more details see Möller et al. (2022); Vincenzi et al. (2024) and Sánchez et al. (2024).

External probes
Our data must be interpretable in context of the parameters of the cosmological models that we test.In this work, many of these are defined as modifications to the background expansion and do not describe how the CMB or galaxy power spectrum may change.Additionally, we would like to be agnostic about the pre-recombination history, and in particular the size of the sound horizon   or  * .
Fortunately, as we describe below, we may still combine the DES-SN5YR cosmological constraints with measurements based on observations from the Cosmic Microwave Background (CMB) and Baryon Acoustic Oscillations (BAO) by the use of derived parameters with clear physical meaning.We do not use data from weak lensing surveys in this work.

Cosmic Microwave Background
The CMB data may be expressed in terms of the 'shift parameter'  (Bond et al. 1997), defined in the literature as where  * is the redshift at the surface of last-scattering,  () ≡  ()/ 0 is the normalized redshift-dependent expansion rate and The physical meaning of  in the context of non-standard cosmological models may be understood if the baryon density   = Ω  ℎ 2 is fixed (for example by nucleosynthesis constraints).Although  is sometimes interpreted as set by the location of the peaks and troughs of the CMB power spectrum (if the sound speed is fixed by   and  CMB ), this relies on the absence of additional energy components in the pre-combination era (for example, early dark energy models as reviewed in Poulin et al. 2023).Alternatively,  may also be understood as localised around the surface of last scattering in the following way.During recombination, photons stream out of overdensities and suppress power on small scales in a process known as Silk damping (Silk 1968).Again at fixed   , successive spectral peaks are lower than their predecessors as the multipole  increases, and the rate of suppression   ∝ exp −2(/ Silk ) 2 (see for example Mukhanov 2004) is proportional to the Hubble expansion rate at the time of last scattering.We may therefore define where   () is the transverse comoving distance defined as We see that  ′ ≃  provided the universe is matter-dominated at the time of last scattering.It may be calculated that  ′ ≃ 1.8 × 10 −3  Silk where the prefactor is only sensitive to cosmological parameters by a factor of (1 +  * ) 1/2 and in turn  * does not depend much on the cosmology.Hence  ′ , which is explicitly proportional to  ( * ), connects  to the Silk damping scale which we take as a safe assumption for the range of models we test.Chen et al. (2019) converted the Planck 2018 (Aghanim et al. 2020) TT,TE,EE + lowE measurements to a prior on , finding  = 1.7502 ± 0.0046 for models assuming spatial flatness and  = 1.7429 ± 0.0051 for models that allow curvature.We use these priors in this work.We also note that Lemos & Lewis (2023) remove latetime cosmology dependence from the CMB likelihoods by using flexible templates for late-ISW and CMB-lensing.We convert their baseline results (Early-ΛCDM, see Table 1 of

Baryon Acoustic Oscillations
Baryon acoustic oscillations represent a sharply-defined acoustic angular scale on the sky given by where   (  ) is the transverse comoving distance to the drag epoch, and   is the comoving sound horizon given by and   is the baryon sound speed, while  * and  * are defined in the same way using  * .BAO measurements are given as the ratio of   to either the Hubble distance,   () = / (), transverse comoving distance,   () or a combination of the two termed the dilation scale,   () ≡  2  ()  () . To interpret these in terms of distances,   is needed.However, in this work, we cancel the dependence on the sound horizon scale by using the ratio of the BAO distance with the distance to CMB as, where    = {  ,   ,   }, and we remind the reader that   ( * ) = (/ 0 )/ √ Ω  .In this way, the data represents the ratio of the angular scales of the sound horizon on the surface of last scattering and at the effective redshift of the BAO.The cosmological dependence of  * /  may be neglected.
We use BAO data from the extended Baryon Oscillation Spectroscopic Survey (eBOSS; Dawson et al. 2016;Alam et al. 2021), which is the cosmological survey within SDSS-IV (Blanton et al. 2017).Specifically, we use the BAO-only measurements from SDSS MGS (Ross et al. 2015), SDSS BOSS (Alam et al. 2017), SDSS eBOSS LRG (Bautista et al. 2021) The covariance matrices provided by eBOSS11 have been incorporated into this study with the use of the uncertainties (Lebigot 2009) python package and the final measurements shown in Table 4.Note that although these measurements contain information from the CMB we will refer to these measurements as BAO- * from here on.

Constraining cosmological models
In general, the parameters of an individual cosmological model are constrained by minimizing a  2 likelihood given by χ2 and for DES-SN5YR,   =  model, −   for the   ℎ SN.However, the absolute magnitudes of SNe Ia are degenerate with  0 .For this analysis, no assumption on  0 is presumed and instead  0 is treated as a nuisance parameter that is analytically marginalised over by modifying equation ( 12).The modified  2 likelihood is given by (Goliath et al. 2001), where and and where we sum over all matrix elements, , .For the combined constraints we sum the  2 likelihoods from all data sets as cobaya12 (Torrado & Lewis 2019; Torrado & Lewis 2021), a robust code for Bayesian analysis, was used to minimize equations ( 13) and ( 16).The convergence of MCMC chains was assessed in terms of a generalized version of the  − 1 Gelman-Rubin statistic (Gelman & Rubin 1992), which measures the variance between the means of the different chains in units of the covariance of the chains.For our work, we adopted a more stringent tolerance than cobaya's default value, namely  − 1 = 0.001.

DES Collaboration et al. (2024b) presents cosmological results
for the standard cosmological model and simple variations such as allowing the dark energy equation of state to be other than  = −1 and/or vary with scalefactor.In this work, we extend on that analysis and present constraints on more exotic non-standard cosmological models.
For each of the models we investigate, the same basic theory applies and the theoretical distance moduli can be calculated as, L () is the luminosity distance and follows the relation, where  is the cosmological redshift and  obs is the observed redshift.However, the Friedmann equation (describing how the Hubble parameter changes with scalefactor or redshift) differs.
In the following subsections, we briefly introduce each model and present the associated normalized Friedmann equation  (), used to determine   () (equation 8).We also present the associated parameter constraints using the DES-SN5YR sample alone and after combining the DES-SN5YR sample with the CMB- and BAO- * (summarised in Table 5).For all fits, we report the median of the marginalised posterior and cumulative 68.27 per cent confidence interval.The best fit Hubble diagrams are shown in Fig. 9.

Cosmography
The cosmographic approach is a smooth Taylor expansion of the scale factor,  that makes minimal assumptions about the underlying cosmological model, however retains the assumptions of homogeneity and isotropy (Visser 2004;Zhang et al. 2017;Macaulay et al. 2019).In cosmography, its useful to define the deceleration parameter, the jerk parameter, and the snap parameter, where  is directly related to the accelerated expansion of the universe and  = 1 at all times for a spatially flat universe with a cosmological constant.Here, we Taylor expand the scale factor for a flat universe and take the series expansion to four terms, where 3 2 0 + 3 3 0 − 4 0  0 − 3  0 −  0 and  0 ,  0 and  0 are the current epoch deceleration, jerk and snap parameters respectively.
We fit cosmographic expansion to third (equation 22 excluding the  3 term) and fourth order (equation 22) with our constraints shown in Fig. 10.For the third order fit we find  0 = −0.362+0.067 −0.069 and evidence for an accelerating universe at > 5.When we fit to fourth order, we find  0 = −0.06+0.11  −0.13 , which is consistent with zero however we note that the snap parameter is poorly constrained by the DES-SN5YR alone and find  0 = 1.4 +4.6  −3.3 .This result is analogous to the DES-SN5YR key paper results on the Flat- 0   CDM model who find a  0 consistent with zero when   is included in the fit.We also ensured that our fits were not over influenced by a particular redshift range and found consistent results after (a) removing low- data using the DES SNe alone and (b) removing high- SNe at  > 0.80.

Parametric models for the equation of state
The parametric models we consider here consider time varying dark energy with different functional forms of the dark energy equation of state, .When all components have a constant equation of state, Friedmann's equation is simply where the sum is over matter (  = 0), curvature (  = −1/3), and dark energy with a constant equation of state ( de = constant), which could be a cosmological constant ( Λ = −1).Radiation (  = 1/3) could also be included but is negligible for our redshift range.When testing dark energy with a time-varying equation of state one needs to make the substitution, The simplest parametric model is where  is generalised to an arbitrary constant while retaining spatial flatness (Flat-CDM).This is the baseline cosmological model used within the DES-SN5YR analysis (DES Collaboration et al. 2024b), who also test a flat model with a time varying dark energy in the form of () =  0 +  (1−) (Flat- 0   CDM).While we do not refit these models here, we convert the constraints on Flat- 0   CDM using a linear variation of (), which is anchored to a pivot redshift   instead of  = 0 (Flat- p 0   CDM), such that  p () = and in the case   = 0,  p = 1 the Flat- 0   CDM parameterisation is recovered.
We also test two other parameterisations.Firstly, the DES-SN5YR baseline model with spatial curvature as an additional free parameter (CDM).Secondly, a model where () varies linearly in redshift instead of scalefactor (Flat- 0   CDM), such that () =  0 +    (Weller & Albrecht 2002) and results in a Friedmann equation given by Results for the parametric forms of the equation of state that we test within this work are summarised in Table 5 and the associated contours are plotted in Fig. 11.
Using the DES-SN5YR alone, the CDM model is statistically consistent with a cosmological constant value of  = −1; however both Flat-  CDM and Flat- p  CDM favour a time-varying component to  that increases with time.We note that the Flat- 0   CDM model was constrained in DES Collaboration et al. ( 2024b) finding ( 0 ,   ) = (−0.36+0.36  −0.30 , −8.8 +3.7 −4.5 ) using DES-SN5YR alone.Here, we refit and convert these results to the equation of state at the pivot redshift and find (  0 ,   ,   ) = (−1.00+0.13  −0.14 , −8.6 +3.8 −4.5 , 0.078).When we combine DES-SN5YR with the CMB- and BAO- * our results are still consistent with a time-varying component to  that increases with time with the best fit   and   (we find   = 0.274 for the pivot redshift) both remaining > 1 from a static .
Interestingly, with the combined data sets, all parametric forms of the dark energy equation of state result in a best fit  > 1 from a cosmological constant and all favour a  > −1.

"Thawing" scalar field models
Light scalar fields provide a dynamical model for evolving dark energy inspired by scalar field models for primordial inflation.In the simplest incarnation of these models, the true vacuum energy density (or cosmological constant) of the universe is assumed to be zero, and dark energy is a transient phenomenon arising from the fact that a classically evolving scalar field  with effective mass   < ∼  0 has not yet have reached its ground state.In most particle physics models, light scalars are not technically natural, so it is conventional to consider models in which the small scalar mass is protected by a  weakly broken shift symmetry, as is the case for the pseudo-Nambu-Goldstone boson (PNGB) model introduced by Frieman et al. (1995).
Assuming the canonical Lagrangian for a scalar field, L = (1/2)       −  (), neglecting spatial perturbations the equation of motion of the field in an expanding universe is given by where the expansion rate is given by =  −1/2 is the Planck mass, and the energy density of the field is The time-evolution of   is determined by  and by the equation of - * Cannot reduce to the cosmological constant for any set of parameters.† Best fits are > 2 from the subset of parameters that reduce to the cosmological constant.§ We convert the constraint on the void fraction to the dressed matter density, which is related by state parameter,   =   /  , where the scalar field pressure is For a given form of the potential  () and initial value of the scalar field, (  ) ≡   at some early time   ≪  0 , this dynamical system can be solved to obtain () and thus the expansion history (assuming spatial flatness) where   = 3 2   2 0 /8.For "thawing" scalar field models (the thawing/freezing nomenclature is from Caldwell & Linder (2005)), which include standard potentials of the form  = (1/2) 2   2 +  4 (with  > 0), the PNGB model  () =  2  2 (1 − cos(/  )), and polynomials  () =  =1     with   ≥ 0, at early times the driving term / in equation ( 27) is subdominant compared to the Hubble-damping term 3 .In this limit, the field is effectively frozen at its initial value   , hence (  ) = 0,   (  ) =  (  ), and   (  ) = −1.
Once the expansion rate drops below the curvature of the potential,  < ∼ √︁ | 2 / 2 |, the field begins to roll down the potential, develops non-negligible kinetic energy, and   grows from −1.The parameters of  () and the value of   jointly determine   () and the current scalar energy density, Ω  =   ( 0 )/  .
While there have been a variety of approximate solutions and fits to late-time scalar field evolution (e.g., Dutta & Scherrer (2008); de Putter & Linder (2008); Chiba ( 2009)), numerical experiments show that the redshift-evolution of   for thawing models is very well approximated by where the value of  is only very weakly dependent on  0 and on the form of  () and is generally in the narrow range  = 1.35 − 1.55.As a consequence, these models are characterized by a quasi-onedimensional parameter space that can be taken to be  0 (with  = 1.45 ± 0.1).This approximation holds if the effective scalar mass   is not large compared to  0 (otherwise, the field will begin oscillating around the minimum of its potential by the present epoch.) In Fig. 11, we show constraints on  0 and Ω  marginalized over the narrow thawing-model prior on .For DES-SN5YR alone, we find Ω  = 0.306 +0.041 −0.042 and  0 = −0.83+0.12 −0.14 ; including CMB and BAO measurements, the resulting constraints are Ω  = 0.323±0.007and  0 = −0.867+0.041 −0.040 , i.e., a 3 deviation from  0 = −1.As shown in Table 6, for the combined data sets the thawing model is moderately preferred over ΛCDM based on the AIC.
The current data provide no meaningful constraint on the parameter  that determines the speed with which   grows from its asymptotic value of −1.That is, if we widen the theory prior on  to allow values  ≫ 1, the best-fitting values are very large, with very large uncertainties.Note that for  ≫ 1, () = −1 down to very low redshift  ≪ 1, so cosmic distances vs. redshift should be indistinguishable from those in ΛCDM.

Chaplygin gas models
Chaplygin gas models deviate from ΛCDM by invoking an exotic background fluid with an equation of state  = − − (Kamenshchik et al. 2001;Bento et al. 2002;Fabris et al. 2004) where  is a positive constant.Chaplygin gas models represent pressureless dark matter in the early universe and dark energy in recent times and therefore may also be able to unify dark matter and dark energy (Bilić et al. 2002).
The simplest form of Chaplygin gas, which was introduced by Kamenshchik et al. (2001), has an equation of state  ∝  −1 ( = 1).This model is referred to as the Standard Chaplygin Gas (SCG) model with a Friedmann equation given by SCG has been shown to be inconsistent with other data sets (Bean & Doré 2003;Sandvik et al. 2004;Davis et al. 2007) however will be re-tested within this work.Generalised Chaplygin Gas (GCG), which maintains  as a free parameter, results in a Friedmann equation given by (34) and reduces to ΛCDM for  = 0 and Ω We note that as ΛCDM is recovered for  = 0, the SCG model (which has  = 1) cannot reduce to ΛCDM for any parameter choice.As a result it may not be surprising that, in contrast to the SCG model, the GCG model has been shown to be consistent with the previous data combinations (Davis et al. 2007;Barreiro et al. 2008;Sollerman et al. 2009;Xu & Lu 2010;Zhai et al. 2017) consisting of the ESSENCE, SDSS-II, Constitution and Pantheon SN data sets (Wood-Vasey et al. 2007;Sako et al. 2007;Hicken et al. 2009;Scolnic et al. 2018).Barreiro et al. (2008) suggest that GCG can be thought of as an interacting form of ΛCDM.The analogous interacting form of CDM was proposed by Zhang et al. (2006) termed New Generalised Chaplygin Gas (NGCG).The Friedmann equation for the spatially flat NGCG model is given by and can be reduced to CDM for  = 0.
In Fig. 12 we present the contours for the Chaplygin gas models we investigate in this work.Constrained by the DES-SN5YR alone, the SCG model provides the lowest central value for the matter density of all models tested within this work at Ω m = 0.121 ± 0.035.We note that this is due to the model favouring a high curvature, equivalent to Ω k = 0.43 ± 0.12.When combined with external priors, the SCG model is unable to simultaneously fit the different data sets (see Fig. 12a), which show extremely strong disagreement in the best fit parameters and highlighted by the poor Akaike Information Criterion (AIC) result of ΔAIC= 276.9 relative to Flat-ΛCDM (see Section 8.2 and Table 6).
Using the DES-SN5YR alone, the remaining Chaplygin Gas models FGCG, GCG and NGCG are consistent within 1 ( = 0 and  = −1 for NGCG) of a cosmological constant.When combined with the CMB- and BAO- * both the FGCG and GCG models find  > 1 from a cosmological constant.For the NGCG model, the best fit  is consistent with a cosmological constant, however favouring  > −1.

Cardassian models
Cardassian models, first proposed by Freese & Lewis (2002) deviate from ΛCDM with the following modification to the Friedmann-Lemaítre-Robertson-Walker metric (FLRW) equation, where the usual FLRW equation is recovered for  = 0. Cardassian models invoke no vacuum energy (Λ = 0), instead the additional term in equation (36) (  ) is initially negligible and only begins to dominate in recent times.Once the second term dominates, it causes the universe to accelerate.Therefore, with this modification, pure matter (or radiation) alone can drive an accelerated expansion.Some motivations for the addition of this term have been suggested and include self-interaction of dark matter (Gondolo & Freese 2002), as well as the embedding of our observable three-dimensional brane in a higher-dimensional universe (Chung & Freese 2000).The original power-law Cardassian model results in a Friedmann equation of the same functional form as that of CDM where  = −1 and therefore does not need to be tested separately.Wang et al. (2003) generalises this model by introducing an additional free parameter  > 0. This model is called Modified Polytropic Cardassian (MPC) expansion which follows, and collapses to Flat-CDM for  = 1 where  =  − 1.
Our constraints in the  −  plane for MPC expansion are shown in Fig. 13a.We find  = 13.3 +4.7 −6.5 using DES-SN5YR alone, inconsistent with  = 1 by ∼ 2.This result is inconsistent with previous analyses by Zhai et al. (2017) and Magañ a et al. (2018) however these analyses both include constraints from probes other than SN.Our results are consistent with these previous analyses and  = 1 when we supplement the DES-SN5YR data with external probes, we find  = 1.38 +0.49−0.42 .

Interacting dark energy & dark matter
In typical cosmological models, dark matter and dark energy are assumed to evolve independently.However, dark energy and dark matter provide the largest contribution to the energy budget of the universe so it is worth investigating if these components can interact.Interacting dark energy & dark matter (IDE) models are therefore those which allow for this interaction (Freese et al. 1987) and are desirable as they allow solutions with a constant dark energy to matter ratio, solving the coincidence problem.
In this paper, we consider a popular subset (Barnes et al. 2005;Guo et al. 2007;Li et al. 2009;He et al. 2010;Li & Zhang 2014;Hu et al. 2016;Wang et al. 2016;von Marttens et al. 2019) of IDE models where the total energy density of dark energy and dark matter is conserved however the particular densities evolve as, where  m and   represent the density of matter and dark energy respectively,  is the dark energy equation of state and  is the interaction kernel which indicates the rate of energy transfer between the two components.We investigate three spatially flat IDE models where  has the general form  =   (  ,  m ), the function  (  ,  m ) specifies a particular IDE model and  is the coupling parameter between the dark components.The sign of  describes the energy flow between the interacting components where  < 0 corresponds to a flow of energy from dark matter to dark energy.The paramaterizations of  and the respective Friedmann equations are: where The IDE models in equations ( 39), ( 40) and (41) will be referred to respectively as IDE1, IDE2 and IDE3 throughout this work.The results for the three IDE models we test are summarised in Table 5 and the contours are shown in Fig. 13.
Using DES-SN5YR alone and after combining the DES-SN5YR with priors from the CMB- and BAO- * all of the IDE models tested are consistent within 1 of no interaction between the dark components,  = 0 and  = −1.
We also note that the CMB- puts a stringent constraint on the interaction for the IDE2 model, where find  = 0.000 ± 0.001.The tightness of the constraint on  is expected and in agreement with previous works (Guo et al. 2007;Wang et al. 2016).This is due to the CMB- data not allowing a large deviation from the standard matter-dominated epoch along with the second term in equation ( 40).

Modified gravity
Dvali-Gabadadze-Porrati (DGP) brane world models first introduced by Dvali et al. (2000) arise from a mechanism where the observed 4D gravity is embedded on a brane in 5D Minkowski space.As a result, locally the gravitational potential propagates in 4 dimensions reducing to General Relativity.However, at large distances the gravitational potential propagates in 5D or 'leaks out into the bulk' deviating from General Relativity and causing accelerated expansion.Two branches of cosmological solutions in the DGP model have distinct properties.The solution examined in this work is the so-called self-accelerating branch where the late-time acceleration of the universe occurs without the need of a cosmological constant (Deffayet 2001) and is described by where − Ω  and the length scale for which the 'leaking' takes place is   and Ω   = 1/4 2   2 0 .Therefore, the Flat-DGP and DGP models have the same number of free parameters as Flat-ΛCDM and ΛCDM respectively.
Inspired by the DGP model, Nicolis et al. (2009); Deffayet et al. (2009) introduced Galileon cosmology, which is a scalar field class of models that are invariant under a shift symmetry in field space.Importantly, the Galileon scalar has no effect on the expansion rate during early times due to a natural screening mechanism, the Vainshtein effect in which non-linear effects can effectively decouple the scalar field from gravity (De Felice & Tsujikawa 2011).In late times, there exists a tracker solution (GAL) that is stable and selfaccelerating with a very negative equation of state  < −1.The Friedmann equation for the GAL model has the same number of free parameters as ΛCDM and is given by where Ω  = 1 − Ω  − Ω  .Both the DGP and GAL models provide a good fit to DES-SN5YR alone.However, when we include external probes, our results (summarised in Table 5) are in agreement with previous works (Lombriser et al. 2009;Li et al. 2011;Xu & Zhang 2016;Zhai et al. 2017;Peirone et al. 2018) that show the DGP and GAL models to be inconsistent with multiple data sets, as seen in Figs.13e & 13f.

Timescape cosmology
So far, the models examined all seek to explain the observed acceleration of the universe, assuming a FLRW geometry which is exactly homogeneous and isotropic.However, the local Universe is far from homogeneous and possesses a cosmic web of structures dominated in volume by voids.Timescape cosmology (Wiltshire 2007a(Wiltshire ,b, 2009) ) discards the approximation of a FLRW universe and instead considers a Buchert average (Buchert 2000) over spatially flat wall regions and negatively curved voids.While the Buchert formalism has been investigated in other works, Timescape cosmology also accounts for a geometry difference between the Buchert average and an observer in a gravitationally bound system within the wall regions, for a universe dominated by voids.Wiltshire (2008) shows that this two-scale model results in a difference in clock rates that accumulates over cosmic time.In this work we use the Timescape tracker solution where the luminosity distance is calculated as, where F () ≡ 2 1/3 +  1/3 6 ln is defined implicitly in terms of the redshift by and Note that  v0 is the current epoch void fraction and the only free parameter of the Timescape model (as we treat  0 as a nuisance parameter in this work), which is related to the dressed 13 matter density parameter by The time,  and Hubble parameter, H0 in equations ( 44), ( 45), ( 46) and ( 47) are the volume averaged values, which are related to values we observe in a wall region by and We note that an average expansion law only holds on scales greater than the statistical homogeneity scale, which corresponds to a CMB rest frame redshift of order  ∼ 0.021 − 0.040 (Scrimgeour et al. 2012;Ntelis et al. 2017).In this work we adopt the value used to quote the key results in Dam et al. (2017) of  min = 0.033.We re-run the entire pipeline with this cut, which reduces our low- sample by 68 SNe (see Section 6; from here on we will refer to this modified sample as DES-SN5YR cut ).We also use CMB rest-frame redshifts excluding peculiar velocity corrections of the host galaxy, which are calculated assuming a standard FLRW model to remain consistent with previous work by Dam et al. (2017).Finally, we retest the Flat-ΛCDM model with these same changes to make a consistent comparison between the two models.In addition to the above changes to the DES-SN5YR data, we also note that the conversion of redshift increments to a radial comoving distance involves different assumptions about spatial curvature in the FLRW and Timescape models (see Appendix D2 from Dam et al. (2017) for more details).Therefore, we do not include the CMB- summary statistic as outlined in Section 6.2.1 when constraining the Timescape and Flat-ΛCDM models and include only angular measurements on the BAO scale (BAO- * ⊥ from here on) from the SDSS data (  ( * )/  () constraints from Table 4).
Using DES-SN5YR alone, we find  0 = 0.791 +0.039 −0.034 , equivalent to a dressed matter density of Ω m = 0.292 +0.043  −0.051 and for Flat-ΛCDM find Ω m = 0.362 +0.019 −0.018 .These results are consistent 13 The dressed parameters are defined such that they take numerical values similar to those of cosmological parameters within FLRW models.

Figure 14.
Constraints on the matter density from the DES-SN5YR cut data set only (blue) and BAO- * ⊥ (yellow) as well as the DES-SN5YR cut combined with BAO- * ⊥ prior (black).We show both the constraints from the Timescape and Flat-ΛCDM models (Section 7.8) with the same modifications to the data.In particular, we apply a redshift cut of  min = 0.033 and excluding peculiar velocity corrections.Note that for Timescape cosmology, the void fraction is related to the dressed matter density by with constraints found by Dam et al. (2017) using the JLA catalogue (Betoule et al. 2014).Fig. 16 shows consistent matter density predictions between Flat-ΛCDM in the baseline analysis and Flat-ΛCDM cut after including a redshift cut at  min = 0.033 and excluding peculiar velocity corrections.However both are just outside the 68 per cent confidence interval of the Planck TTTEEE-lowE prediction (Ω m = 0.3166 ± 0.0084; Aghanim et al. 2020).In contrast, the Timescape model has a lower central value for the matter density in agreement with Planck.When combining the DES-SN5YR with BAO- * ⊥ , we find Ω m = 0.446 +0.010 −0.009 for the Timescape model and for Flat-ΛCDM find Ω m = 0.332 +0.011 −0.010 .These results are shown in Fig. 14.It is apparent from the upper panel that the datasets BAO- * ⊥ and DES-SN5YR are in tension in the Timescape model, and this model is therefore disfavoured relative to Flat-ΛCDM by the AIC statistic.

Goodness of fit
To investigate the goodness of fit for each of the models we present the  2 for various data combinations, see Table 6, where  2 = −2 ln L max and L max is the maximum likelihood of the entire parameter space.
The number of degrees of freedom ( dof ) is equal to the number of data points minus the number of cosmological parameters constrained for each model.For DES-SN5YR and DES-SN5YR cut , we approximate the number of data points by summing the BEAMS probability of each SN being Type Ia and find  B(Ia) = 1735 and 1666 respectively.The additional number of data points when including the CMB- , BAO- * or BAO- * ⊥ are 1, 14 and 7 respectively.Using DES-SN5YR alone, we find that all models tested within this work result in good fits to the data.However, the SCG, DGP and GAL models have a poor  2 when combining DES-SN5YR with the CMB- and BAO- * as they are unable to reconcile the additional data sets.To a lesser extent this also afflicts the Timescape model.This can be seen visually in Figs. 12a,13e,13f,& 14 where the parameter space of the combined contours do not share a common region with all probes.

Model comparisons
To assess whether additional parameters invoked in the more complex models are justified given the data, we use the Akaike Information Criterion AIC ≡ 2 − 2 ln L max (Akaike 1974), where  is the number of parameters in the model.We also use the Suspiciousness (Handley & Lemos 2019), which is defined as ln  = ln −ln  where  is the Bayes Ratio and  is the Bayesian information.Handley & Lemos (2019) note that the Bayes ratio is prior-dependent and show that Suspiciousness is prior-independent due to the combination with the Bayesian information.
In Table 6 we quote the 1 2 ΔAIC14 and the difference in the logarithm of the Suspiciousness, Δln  relative to Flat-ΛCDM.To asses the strength of this preference, Trotta (2008) suggests that Δ > 1, Δ > 2.5 and Δ > 5 indicates weak, moderate and strong evidence respectively, against the model with the higher Δ value.In both cases, more negative values indicate that the data prefers the extended model over Flat-ΛCDM.We determine ln  using anesthetic software (Handley 2019) with the nested sampling outputs from polychord (Handley et al. 2015b,a) with 25 ×  live points, 5 ×  repeats and an evidence tolerance requirement of 0.1.
Using the DES-SN5YR alone, both the AIC and Suspiciousness find no strong evidence for or against any of the non-standard models.Both find weak evidence for the third order cosmographic model and moderate evidence for the fourth order cosmographic model.Furthermore, the AIC and Suspiciousness weakly and moderately prefer the Flat- 0   CDM, Flat-   CDM and MPC models over Flat-ΛCDM respectively.The Timescape model, which was fit using the DES-SN5YR cut sample is weakly preferred by both the AIC and Suspiciousness.
When combined with the CMB- and BAO- * both the AIC and Suspiciousness agree that 11 of the 15 non-standard models we investigate are moderately preferred over Flat-ΛCDM.We note that this is not a result of curvature alone with no preference for or against the ΛCDM model.The top performing models include Flat-   CDM with ( 1 2 ΔAIC, Δln ) = (−3.2,−4.17) indicating an evolution of  that increases with time, the thawing model with ( 1 2 ΔAIC, Δln ) = (−3.2,−4.60) and the FGCG model with ( 1 2 ΔAIC, Δln ) = (−3.4,−3.94), which invokes an exotic background fluid.These results suggest that additional flexibility in our cosmological models may be required beyond the cosmological constant. -

Tension metrics
We also use the Suspiciousness to assess whether different datasets are consistent (in contrast to Section 8.2 and Table 6 where the Suspiciousness was used as a model comparison statistic), which is ideal for cases such as ours where we have chosen deliberately wide and uninformative priors (See Lemos et al. 2021, section 4.2).We use the anesthetic software (Handley & Lemos 2019) to determine ln  and produce and ensemble of realisations to estimate sample variance.Using the scale from Trotta (2008), ln  < −5 is considered strong tension, −5 < ln  < −2.5 is considered moderate tension and ln  > −2.5 indicates that the datasets are in agreement.In Fig. 15, we plot the Δln  between the relevant datasets.Note, models already been shown in Section 8.1 to be poor fits to the combined data sets (SCG, DGP and GAL) have been excluded from the plot and all had Δln  << −5.We find a strong tension between the DES-SN5YR with BAO- * ⊥ datasets when fitting the Timescape model.For all other models, we find no indication of tension.

CONCLUSIONS
The DES Supernova survey is the largest, most homogeneous SN data set to date containing 1635 supernovae combined with 194 existing Low- SNe Ia.The statistical power of the DES-SN5YR sample allows us to obtain robust and precise constraints on cosmological models beyond ΛCDM.
We first investigated two important areas of the main DES supernova cosmology analysis that are, or may appear to be subject to cosmological dependencies.
(i) We demonstrated that the assumption of a reference cosmology used to generate simulated light curves and perform selection bias corrections to the data results in a bias that is sub-dominant  to statistical uncertainties.For non-standard models, we also show a region of expansion histories where we are confident in our bias corrections.For the next era of SN experiments, the reference cosmology may become a dominating systematic and as a result, we show that an iterative method (where the reference cosmology is updated in a second iteration based on the best fit cosmology from the first) is viable and can be employed to reduce this bias.
(ii) We demonstrated that the BBC fitting procedure, which uses a fiducial cosmology, is insensitive to that choice of cosmology.We also show that a spline is viable alternative to a fiducial cosmology as it obtains consistent results and may reduce confusion as to the role of the fiducial cosmology in future analyses.
Secondly, we presented constraints on 15 exotic cosmological models using the DES-SN5YR sample alone and after combining the DES-SN5YR with external probes.Using DES-SN5YR alone, we find that all models tested within this work are good fits to the data.This trend continues when we combine the DES-SN5YR with priors from the CMB- and BAO- * except for models that had been previously ruled out.We assessed whether additional parameters invoked in the more complex models are justified given the data by using the Akaike Information Criteria and Suspiciousness.Of the 15 models that we test, we find no strong evidence for or against any of the non-standard models for any of our data combinations.Using the DES-SN5YR alone, the Suspiscousness moderately prefers 3 of the non-standard models along with the fourth order cosmographic model.When combined with the CMB- and BAO- * both the AIC and Suspiciousness agree that 11 models are moderately preferred over Flat-ΛCDM.We show that this is not a result of curvature alone.Our work suggests that additional flexibility in our cosmological models may be required beyond the cosmological constant.

ACKNOWLEDGEMENTS
Author Contributions: RCa developed the code, performed the analysis, drafted the manuscript.TMD conceived the idea, guided the project, ran test code and helped writing.MV contributed to developing part of the analysis, guided the project and commented on drafts.PS commented on drafts, ran test code and helped writing.All authors contributed to this paper and/or carried out infrastructure work that made this analysis possible.Highlighted The Collaborating Institutions are Argonne National Laboratory, the University of California at Santa Cruz, the University of Cambridge, Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas-Madrid, the University of Chicago, University College London, the DES-Brazil Consortium, the University of Edinburgh, the Eidgenössische Technische Hochschule (ETH) Zürich, Fermi National Accelerator Laboratory, the University of Illinois at Urbana-Champaign, the Institut de Ciències de l'Espai (IEEC/CSIC), the Institut de Física d'Altes Energies, Lawrence Berkeley National Laboratory, the Ludwig-Maximilians Universität München and the associated Excellence Cluster Universe, the University of Michigan, NSF's NOIRLab, the University of Nottingham, The Ohio State University, the University of Pennsylvania, the University of Portsmouth, SLAC National Accelerator Laboratory, Stanford University, the University of Sussex, Texas A&M University, and the OzDES Membership Consortium.
Based in part on observations at Cerro Tololo Inter-American Observatory at NSF's NOIRLab (NOIRLab Prop.ID 2012B-0001; PI: J. Frieman), which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation.

Figure 4 .
Figure 4. a) and b): Shifts in Ω m and  (solid points) when using different BBC simulations that are distinguished by a unique reference cosmology (shown by open symbols; listed in the figurelegend).The shifts are measured relative to the perfect scenario (orange square) where the reference cosmology is equal to the true cosmology of our simulated data.c) and d): The associated mean shifts in   (0.20) (with no prior) as well as  determined with a strong prior on the matter density of Ω m = Ω m,true ± 0.001, which minimizes the impact that the Ω m −  degeneracy has on investigating the BBC reference cosmology.For panels a) -d) we have averaged over 25 DES-SN5YR simulations.Note also that the error bars show the uncertainty on the shift in the mean -not the uncertainty on the parameters, which is larger.e): Calculated residual distance moduli of the reference cosmologies (dashed lines) relative to the baseline cosmology (Ω m , ) = (0.315, −1.0) in orange.The solid lines represent the variation in the expansion history from the perfect scenario using the mean of the best fit parameters.f): Best fit parameters (solid points) for 1 realisation of simulated data determined using a unique BBC reference cosmology (shown by open symbols).The 1 and 2 contours shown are for the ideal case (orange square).The grey dotted dashed line represents the   (0.2) parameter.g): Equivalent information to that contained in plot f) but converted to Ω m −   (0.2) space.

Figure 5 .
Figure 5.Comparison of the best fit  0 −   points (with a prior on the matter density, Ω m = Ω m,true ± 0.001) determined using the DES-SN5YR baseline reference cosmology (purple) and when the reference cosmology is set to the input cosmology of the simulations (blue).The points show the maximum likelihood values for each realisation and the crosses represent the averages of the those maximum likelihood values.The ellipses are the 1-and 2 contours representing the dispersion of best fit points.

Figure 6 .
Figure 6.Iterative procedure methodology.During the first iteration, bias corrections are modelled using simulations created using the default reference cosmology with a fixed set of Flat-CDM parameters Ω m,ref = 0.3 and  ref = −1.0.In the second iteration, the simulations are instead created using the maximum likelihood estimates from the first iteration.

Figure 7 .
Figure7.Top panels: Shifts in the average maximum likelihood , , Ω m ,  and   (0.2) values after varying the fiducial cosmology within the BBC fit (Section 5).The error bars used are the standard error of the mean and are therefore much larger for the individual case.The values are shown relative to the ideal case (black dashed line) where the fiducial cosmology is equal to the true cosmology used to simulate the data.Only the model with zero matter density, and pure cosmological constant (plum) shows a more than 1 shift from the fidicual, and comparison with both the   panel and Fig.8shows that shift is along the degeneracy direction.Bottom right: Variation in the evolution of magnitude versus redshift from the ideal case for (a) the input fiducial cosmological parameters (given in the legend) shown as dashed lines and (b) using the mean of the best fit parameter values shown in the zoomed inset axes as solid lines.

Figure 8 .
Figure 8. Left: The best fit  and  for 100 mock realisations for each of six different reference cosmologies as per the legend (see Section 5).The black points show the maximum likelihood values for each realisation and the blue crosses represents the averages of the those maximum likelihood values.The blue ellipses are the 1-and 2 contours representing the dispersion of best fit points.The upper left sub figure represents the perfect scenario where the fiducial cosmology is equal to the true cosmology used to simulate the data.The black dashed lines are used to compare each figure to this ideal case.Right: The equivalent figure after fitting for cosmological parameters, Ω m and .
Lemos & Lewis 2023) into a constraint on the shift parameter and find  = 1.7442±0.0044.Reassuringly, the central value falls between the constraints fromChen et al. (2019).
, SDSS eBOSS ELG (de Mattia et al. 2021), SDSS eBOSS QSO (Hou et al. 2021) and SDSS eBOSS Ly  (du Mas des Bourboux et al. 2020).We note that new BAO measurements from both DES (DES Collaboration et al. 2024a) and the DESI collaboration (DESI Collaboration et al. 2024) were released in the advanced stages of this work and motivates a follow-up analysis with with the inclusion of these data sets.

Figure 9 .
Figure 9. Upper Panel: Hubble diagram of DES-SN5YR with the overlaid best fit Flat-CDM model.We also show the inflated distance uncertainties from likely contaminants.Four lower panels: The difference between the data and the best fit Flat-CDM model from the DES-SN5YR alone.We also overplot the best fit for each model (we exclude the Timescape model as it was fit against a modified Hubble diagram).Spatially flat models are shown as solid lines and models that allow curvature are represented by dashed lines.

Figure 10 .Figure 11 .
Figure 10.Constraints for the 3 rd and 4 th order cosmographic models (Section 7.1) from the DES-SN5YR data set only.The contours represent the 68.3 and 95.5 per cent confidence intervals.

Figure 12 .
Figure 12.Chaplygin Gas models, a) SCG, b) FGCG, c) GCG and d) NGCG (Section 7.4): Constraints from the DES-SN5YR data set only (blue), a prior from the CMB- (green), BAO- * (orange), CMB- + BAO- * (purple) as well as the DES-SN5YR combined with both the CMB- and BAO- * priors (overlaid black contours).The contours represent the 68.3 and 95.5 per cent confidence intervals.The red dashed lines mark the parameters that recover a cosmological constant.

Figure 13 .
Figure 13.Same as Fig. 12 but for the a) MPC model (Section 7.5), the three IDE models (Section 7.6): b) IDE1, c) IDE2 and d) IDE3, as well as the e) DGP model and f) GAL model (Section 7.7).

Figure 16 .
Figure 16.A summary of the best fit matter density for the models constrained by the DES-SN5YR sample.In black are the constraints from DES Collaboration et al. (2024b), in blue and orange are constraints from this work, where the orange points highlight that the Hubble diagram used to constrain the Timescape and Flat-ΛCDM cut models included a redshift cut at  min = 0.033 and excluded peculiar velocity corrections.The purple shaded region represents the TTTEEE-lowE 68 per cent confidence limits for the Flat-ΛCDM model determined by the Planck collaboration (Aghanim et al. 2020).
contributions include: Contributed to project development: DB, DS, JF, RK.Contributed to advising and validation: CL, DS, JF, MSa, PA, PW, RK.Construction and validation of the DES-SN5YR Hubble diagram: AC, AM, BP, BS, CL, DB, DS, GT, HQ, JL, LG, MSa, MSm, MSu, MT, PA, PW, RCh, RK SH.Contributed to the internal review process: CL, GFL, JL, LG, MSu, PW, RK.The remaining authors have made contributions to this paper that include, but are not limited to, the construction of DECam and other aspects of collecting the data; data processing and calibration; developing broadly used methods, codes, and simulations; running the pipelines and validation tests; and promoting the science analysis.TMD, AC, RCa, SH, acknowledge the support of an Australian Research Council Australian Laureate Fellowship (FL180100168) funded by the Australian Government.JF would like to acknowledge Shrihan Agarwal and Yuuki Omori for helpful conversations.MV was partly supported by NASA through the NASA Hubble Fellowship grant HST-HF2-51546.001-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract NAS5-26555.LG acknowledges financial support from the Spanish Ministerio de Ciencia e Innovación (MCIN) and the Agencia Estatal de Investigación (AEI) 10.13039/501100011033 under the PID2020-115253GA-I00 HOSTFLOWS project, from Centro Superior de Investigaciones Científicas (CSIC) under the PIE project 20215AT016 and the program Unidad de Excelencia María de Maeztu CEX2020-001058-M, and from the Departament de Recerca i Universitats de la Generalitat de Catalunya through the 2021-SGR-01270 grant.AM is supported by the ARC Discovery Early Career Researcher Award (DECRA) project number DE230100055.Funding for the DES Projects has been provided by the U.S. Department of Energy, the U.S. National Science Foundation, the Ministry of Science and Education of Spain, the Science and Technology Facilities Council of the United Kingdom, the Higher Education Funding Council for England, the National Center for Supercomputing Applications at the University of Illinois at Urbana-Champaign, the Kavli Institute of Cosmological Physics at the University of Chicago, the Center for Cosmology and Astro-Particle Physics at the Ohio State University, the Mitchell Institute for Fundamental Physics and Astronomy at Texas A&M University, Financiadora de Estudos e Projetos, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro, Conselho Nacional de Desenvolvimento Científico e Tecnológico and the Ministério da Ciência, Tecnologia e Inovação, the Deutsche Forschungsgemeinschaft and the Collaborating Institutions in the Dark Energy Survey. 1.

Table 2 .
Testing the iterative method (Section 4.2): Weighted average (over 10 realisations * ) difference in  and   from the truth for the first and second iterations.
* Ω m = 0.350 and  = −0.8 was used as the true cosmology.

Table 3 .
BBC fitted nuisance parameters for three different fiducial cosmologies, showing the results are stable to the choice of fiducial cosmology or use of a spline (see Section 5.1).

Table 4 .
Summary of the external constraints determined using measurements from eBOSS and Planck. eff   ( * )/  ()   ( * )/  ()   ( * )/  () In this work we use the 'shift parameter'  that is related to the heights of the CMB acoustic peaks and depend on the line of sight distance to the sound horizon.
BAO- * measurements * * The product of the BAO measurements with the CMB acoustic scale.†

Table 5 .
Results for the cosmological models investigated in this work.These are the medians of the marginalised posterior with 68.27 per cent integrated uncertainties ('cumulative' option in ChainConsumer).Cannot reduce to the cosmological constant for any set of parameters.† Best fits are > 2 from the subset of parameters that reduce to the cosmological constant.

Table 5 -
continued Results for the cosmological models investigated in this work.These are the medians of the marginalised posterior with 68.27 per cent integrated uncertainties ('cumulative' option in ChainConsumer).

Table 6 .
Goodness of fit & Model Comparison statistics.A more negative 1 2 ΔAIC and Δln  value indicates a stronger preference over Flat-ΛCDM.
Measurements of Δln  between the DES-SN5YR and the combined CMB- + BAO- * datasets (blue).The modified data sets for the Timescape and Flat-ΛCDM cut are shown in orange.The shaded yellow and red regions represent moderate and strong tension respectively.Note, models already been shown in Section 8.1 to be poor fits to the combined data sets (SCG, DGP and GAL) have been excluded from the plot for clarity and all had Δln  << −5.
of which include ERDF funds from the European Union.IFAE is partially funded by the CERCA program of the Generalitat de Catalunya.Research leading to these results has received funding from the European Research Council under the European Union's Seventh Framework Program (FP7/2007-2013) including ERC grant agreements 240672, 291329, and 306478.We acknowledge support from the Brazilian Instituto Nacional de Ciência e Tecnologia (INCT) do e-Universo (CNPq grant 465376/2014-2).This manuscript has been authored by Fermi Research Alliance, LLC under Contract No. DE-AC02-07CH11359 with the U.S. Department of Energy, Office of Science, Office of High Energy Physics.