Constraining multi-field inflation using the SPHEREx all-sky survey power spectra

We investigate how well the SPHEREx all-sky survey can constrain local primordial non-Gaussianity beyond the parameter f NL using galaxy power spectra. We forecast joint constraints on the parameters f NL, g NL and τ NL obtained assuming a simple two-field curvaton model of inflation. The parameters f NL and g NL characterise the squeezed limits of the primordial bispectrum and trispectrum respectively, and lead to a characteristic scale-dependence of the galaxy bias that increases out to arbitrarily large scales. Values of the parameter τ NL > (6/5f NL)2 cause the galaxy power spectrum to have a stochastic component which also increases out to arbitrarily large scales. Our MCMC forecasts indicate that SPHEREx can provide joint constraints on any two of the three parameters f NL, g NL and τ NL. Due to strong degeneracies among these parameters, measurements of the galaxy power spectra alone may not be sufficient to jointly constrain all three. Constraints on f NL, g NL and τ NL obtained from galaxy power spectrum observations depend on the modelling of underlying nuisance parameters. We study the robustness of our forecast constraints to modelling choices and note that even with relatively conservative modelling assumptions, SPHEREx galaxy power spectra can provide strong evidence of local non-Gaussianity, even if the particular values of f NL and g NL cannot be measured precisely.


Introduction
A major goal of cosmology is to determine the source of the initial seeds of structure in our Universe [1,2].In the simplest scenario, the primordial curvature perturbations were sourced by quantum fluctuations in a single inflaton field, whose energy density is also responsible for driving the accelerated expansion of the Universe during inflation [3][4][5].In this case, the primordial perturbations are nearly Gaussian, with departures from Gaussianity only arising if the inflaton has strong self-interactions [6][7][8] or unusual scenarios in which the inflaton field is not on the slow-roll attractor solution to the equations of motion [9,10].If additional light (m ≪ H, where H is the Hubble rate during inflation) fields are present during inflation, then these fields will acquire a spectrum of primordial fluctuations in addition to the fluctuations from the inflaton.Such fluctuations can be converted into fluctuations in curvature through processes that are generically nonlinear but local in configuration space.Examples of such scenarios are the curvaton model and modulated reheating [11][12][13][14][15]. Curvature perturbations generated in this manner have distinct non-Gaussian statistical distributions that can be used to infer their presence.
The non-Gaussian primordial curvature perturbation, ζ, arising from fluctuations in additional light fields can typically be written as, where Z is a Gaussian random field and the coefficients f NL and g NL characterise nonlinear terms which are subdominant in the limit f NL ⟨Z2 ⟩1/2 , g NL ⟨Z 2 ⟩ ≪ 1. 1 The expansion in eq.(1.1) is said to describe Local Primordial non-Gaussianity (LPnG) because it is generated through a mapping of a Gaussian field that is local in configuration space.The primordial JCAP05(2024)094 bispectrum of ζ from eq. (1.1) is proportional to f NL , while the primordial trispectrum contains two terms, one proportional to g NL and another proportional to f2 NL .A characteristic of LPnG initial conditions is that the higher-order correlation functions of ζ peak in the so-called squeezed-limit, when one of the wavenumbers of the ζ fields is taken to zero.This is a feature that is distinct from primordial curvature generated by fluctuations from a single inflaton field [16].The detection of LPnG would therefore be of tremendous importance, and imply that the initial conditions of our Universe were sourced by a new light field, over and above the inflaton [17].The best constraints on the parameters f NL and g NL currently come from the Planck cosmic microwave background (CMB) datasets [18], according to which f NL = −0.9± 5.1 and g NL = (−5.8± 6.5) × 10 4 at 68% confidence.Typical multi-field inflationary models (including the aforementioned curvaton models) predict f NL ≈ 1 [1].A robust detection of f NL ≈ 1 is therefore an important theoretical target and would conclusively establish the existence of one or more additional light fields during the epoch of inflation.
A generalisation of the non-Gaussian initial conditions in eq.(1.1) is to allow for an observable contribution to the primordial curvature ζ that is Gaussian and sourced by the inflaton, in addition to that sourced the curvaton (see, e.g., refs.[19,20]). 2In this case, the initial curvature perturbations can be written as a sum of terms coming from inflaton, ( Above, ζ I is a Gaussian random field and ζ C is given by an expansion such as eq.(1.1).As we shall see, these scenarios have suppressed non-Gaussianity, relative to a model where ζ I fluctuations are negligible, yet have an O(f 2 NL ) trispectrum that is enhanced relative to the bispectrum in single-source initial conditions from eq. (1.1).This trispectrum peaks in the collapsed limit and its amplitude is characterised by a parameter τ NL that is bounded from below as τ NL ≥ ( 6 5 f NL ) 2 [21,22].Initial conditions that exceed this bound can be thought of as having a collapsed trispectrum that is not fully correlated with the square of the primordial bispectrum.Current bounds on τ NL from Planck CMB data are τ NL ≤ 2800 at 95% confidence and therefore far above the minimum threshold set by f NL constraints.
Although at present the best constraints on primordial non-Gaussianity come from CMB datasets, major advances are expected to come from large-scale structure [1].For scenarios with LPnG, the squeezed-limit N + 1-point primordial correlation functions modulate N -point functions by the large-scale curvature field ζ(⃗ x).This in turn manifests as a modulation of the galaxy distribution by ζ(⃗ x), or a scale-dependent galaxy bias with respect to the density field [23].The scale-dependent galaxy bias is a key signal (in addition to the galaxy bispectrum) used by large scale structure surveys to constrain f NL [1,[24][25][26][27][28][29][30][31][32].Current constraints on f NL from the BOSS galaxy survey are f NL = 9 +33 −35 from the galaxy power spectrum alone and f NL = −33 ± 28 with the addition of the galaxy bispectrum, both at 68% confidence [30].Yet, the scale-dependent bias is the combined effect of the squeezed limits of arbitrary higher point primordial correlation functions.In other words, all the non-linear JCAP05(2024)094 coefficients f NL , g NL . . . in eq.(1.1) each contribute to the scale-dependent galaxy bias (see, e.g., refs.[33,34]).As these non-linear coefficients all contribute to the galaxy power spectrum, there is a strong degeneracy among them when constrained using the galaxy power spectrum alone.On the other hand, as we will discuss, scale-dependent bias constraints can be interpreted as constraints on all higher-order terms appearing in eq.(1.1).The galaxy power spectrum therefore provides broad constraining power on LPnG.The tightest constraints on g NL from scale-dependent bias come from quasar power spectra, −4 < g NL × 10 −5 < 4.7 at 95% confidence.Joint constraints on f NL , g NL have been performed in [26,27].Relatedly, for initial conditions of the form eq. (1.2), the halo bias also acquires a stochastic component that is scale-dependent on horizon scales, and has been used to constrain the relative contributions of inflaton and curvaton fluctuations [20].Generalising future galaxy power spectrum constraints on primordial non-Gaussianity to include these terms beyond f NL is important, as information about these higher order statistics may be helpful in distinguishing between different multi-field inflation models, thus helping to further constrain inflationary physics.
The SPHEREx All-sky galaxy survey aims to achieve constraints on the parameter f NL at a level σ f NL ≲ 0.5 [28].Reaching this threshold would represent an order-of-magnitude improvement over CMB constraints and the first instance of reaching values of f NL predicted by the simplest two-field inflation scenarios [11,12].SPHEREx will achieve this target through joint analysis of galaxy auto-and cross-power spectra in combination with galaxy bispectra.While the signatures of f NL , g NL , and τ NL on the galaxy bispectrum and trispectrum are distinct, owing to the different shapes of their primordial correlation functions, these parameters have nearly common signatures in the galaxy power spectrum.As we will review, both f NL and g NL contribute to scale-dependent galaxy bias, while τ NL ≥ ( 6 5 f NL ) 2 generates a stochasticity between galaxies and matter that is also probed by galaxy autoand cross-power spectra.Therefore, the galaxy power spectra measured by SPHEREx can be used to put broad constraints on primordial non-Gaussianity beyond just those set by interpreting data in terms of the f NL parameter.Our goal in this paper is to forecast constraints achievable by SPHEREx from galaxy power spectra on LPnG where the primordial trispectrum, characterised by g NL and τ NL , is significant.In this, we will pay special attention to degeneracies among these parameters and possible limitations to interpretation due to uncertainties in modelling their signatures on the galaxy bias.While additional information on these parameters will be available from galaxy bispectra and trispectra, in particular information that could break degeneracies among them, the analysis is considerably more complicated and model-dependent.Galaxy power spectra, therefore, could play an important role in establishing the existence of LPnG, that could subsequently be studied in more detail through higher-order galaxy correlation functions.
The rest of this paper is organised as follows.In section 2, we review the derivation of scale-dependent bias from LPnG, allowing for LPnG that is sourced by arbitrarily high orders in an expansion such as eq.(1.1).In the same section, we review the derivation of scale-dependent stochastic bias sourced by multi-field initial conditions such as eq.(1.2).Here, we also allow the curvaton field to have polynomial terms of arbitrary high order, before specialising to the case where terms after g NL Z3 can be neglected.In section 3, we present our likelihood and forecasting pipeline, as well as models we use for the SPHEREx galaxy power spectra.In section 4, we present our forecast constraints on f NL , g NL , and JCAP05(2024)094 τ NL , along with a study of the robustness of these constraints to bias modelling choices.Our conclusions are presented in section 5.

Galaxy bias from local primordial non-Gaussianity
In this section, we review the derivation of the effect of LPnG on the galaxy power spectrum measured by galaxy surveys [23], with particular focus on the impact of higher-order primordial correlation functions, beyond the bispectrum.Throughout our analysis, we assume a weak primordial non-Gaussianity, such that the dimensionless cumulants ⟨ζ N ⟩ c /⟨ζ 2 ⟩ N/2 ≪ 1, for N > 2 and where c refers to the connected part of the correlation function. 3We will, however, allow dimensionless cumulants of different orders to be of comparable magnitude, i.e. we will consider examples where For this to happen with the initial conditions in eq.(1.1), higher-order coefficients, such as g NL , must be much larger in magnitude than lower-order coefficients, such as f NL .For the bispectrum and trispectum specifically, one would require that g NL ∼ f NL A −1/2 S for them to have comparable magnitude, where A s characterises the scalar power spectrum of ζ.Initial conditions such as this can be thought of as having a primordial kurtosis that is small, but as large as the primordial skewness, in comparison to the variance.Non-Gaussianity of this sort is completely consistent with current datasets (and in fact explains the relative magnitude of the current constraints on f NL and g NL ).
In section 2.1, we review how LPnG leads to a characteristic scale-dependence in the galaxy bias that diverges quadratically at the largest scales.In the following section 2.2, we review how LPnG gives rise to a stochastic contribution to the galaxy power spectrum in the case of a two-field curvaton model of inflation wherein both the inflaton and the massless curvaton contribute comparably to the primordial curvature perturbation [20,35].In both sections 2.1 and 2.2, we first present a more general result for scale-dependent bias and stochasticity that holds for non-trivial squeezed limits of arbitrary higher point primordial correlation functions and then provide concrete expressions for primordial non-Gaussianity with only a primordial bispectrum and trispectrum.We work exclusively in the comoving synchronous gauge.This has the advantage that it is the unique gauge in which the linear bias relation (δ g = b g δ cb ) remains valid (to linear order) at large scales [36]. 4

Scale-dependent galaxy bias from local primordial non-Gaussianity
LPnG is characterised by the primordial curvature perturbation having a non-trivial N -point correlation function (N ≥ 3) that peaks in the squeezed limit, i.e. when one of the external momenta becomes soft.It is well understood that to lowest order, the squeezed-limit of a primordial (N + 1) correlation function encodes the response of the small-scale N -point correlation function in the presence of a long wavelength mode of the field whose momentum becomes soft [37,38].In particular, for the squeezed-limit of the primordial correlation functions, one can show that for N ≥ 2, where the linear CDM+baryon density fields are evaluated at typical redshifts probed by galaxy surveys.Under initial conditions with LPnG, the left-hand side of the above equations approaches a non-trivial value implying that the fluctuations of CDM+baryon densities acquire a dependence on an ambient large-scale curvature perturbation (ζ L ) that does not vanish even at the largest scales.The equivalence principle forbids such a dependence from being generated from nonlinear gravitational evolution [17].
To analyse the effect of LPnG on the halo bias, we adopt the peak-background-split (PBS) approach [39].In this approach, the matter density fluctuations can be separated into a large-scale component that does not vary significantly over scales as large as the size of a typical halo and a small scale component which undergoes gravitational collapse to form dark matter haloes.The comoving number density of dark matter tracers, e.g.galaxies or quasars, n g is therefore a function of the connected moments of the linear matter density field δ cb,R smoothed over some small scale R corresponding to the typical Lagrangian radius of dark matter haloes. 5For Gaussian initial conditions, the only non-trivial such moment is the smoothed variance of the linear density field σ 2 R = ⟨δ 2 cb,R ⟩.For non-Gaussian initial conditions, one can have non-trivial higher order moments S N (R; z) = ⟨δ N cb,R ⟩ given by: where W (kR) is the window function used to calculate the smoothed matter density field δ cb,R and T cb (k, z) is the linear transfer function of matter density fluctuations Due to the coupling between small-scale and large-scale modes induced by LPnG following eq.(2.1), the moments of smoothed linear matter density fluctuations S N (R; z) acquire a dependence on a background large scale curvature perturbation ζ L , through the dependence of S N on the squeezed-limit of the primordial N + 1-point correlation function,6

JCAP05(2024)094
The quantity n g will generically depend on the small-scale statistics of the matter field such as its variance and the higher-order cumulants S N .Within a large-scale matter density perturbation δ L , n g [δ L ; {S N (ζ L )}] acquires an additional dependence on the large scale curvature perturbation ζ L .At linear order in ζ L , we have, where S o i are the moments of the smoothed linear density field in the absence of an ambient large-scale curvature perturbation.For typical models with weak non-Gaussianity, only finitely many terms contribute to the summation in eq.(2.6).Relating ζ L to δ L using the transfer function T cb (k, z), one can write a linear bias relation for fluctuations in the abundance of matter density tracers We define b L = ∂ log n g /∂δ L as the scale-independent linear Lagrangian halo bias which acquires a scale-dependent correction in the presence of LPnG, The CDM+Baryon transfer function T cb (k, z) can be obtained by inverting the Poisson equation and is typically written as, where the function T (k) approaches 1 as k → 0, 7 and D(z) is the linear growth factor normalised to be the scale factor during the matter-dominated era.It is convenient to write the scale-dependent contribution to the bias (eq.(2.8)) in the following form where (2.12) Note that in eq.(2.12), B(z) will depend on properties of the tracer population n g (for instance, the mass or luminosity of galaxies), but we have suppressed that argument here.In this way, we see that LPnG induces a scale-dependent galaxy bias that diverges as 1/k 2 at the largest scales and that its magnitude is parameterised by a weighted sum over the squeezed limits JCAP05(2024)094 of primordial correlation functions (given by the factors ∂S N /∂ζ L , see eqs.(2.4) and (2.5)) weighted by the marginal contribution of small-scale statistics of the linear density field to the abundance of dark matter tracers (given by the factors ∂ log n g /∂S N ).This fact makes scale-dependent bias a broad probe of any type of local primordial non-Gaussianity.This is in stark contrast to probes such as N -point functions of the CMB, which only constrain N -point functions of ζ at the same order.We now specialise to the simplest model of LPnG constrained by the CMB datasetsin this model, the primordial curvature fluctuation can be expressed as a polynomial (at most cubic) in a Gaussian random field as Current constraints on the parameters f NL and g NL satisfy the condition for weak non-Gaussianity, i.e.
ζ .Under these conditions, the only relevant primordial correlation functions in addition to the power spectrum are the bispectrum and trispectrum, proportional to f NL and g NL , In our studies of g NL , we are working in the limit that so that f NL ≪ g NL and the second term in eq.(2.15) can be dropped.In the next section, we will generalise to two-field initial conditions, which can enhance the f 2 NL term.The only relevant non-trivial moments of the smoothed linear matter density field (S N defined in eq.(2.3)) are then proportional to f NL and g NL and are as follows, ) where the factors Σ 3 and Σ 4 can be obtained using eq.( 2.3) and are given by: Eqs. (2.4) and (2.5) show that the response of the moment S N to an ambient long-wavelength curvature perturbation is proportional to the smoothed squeezed limit of the (N + 1) linear matter correlation function.Since the only non-trivial linear matter correlation functions (for the model given by eq.(2.13)) are the bispectrum and the trispectrum respectively, -7 -

JCAP05(2024)094
it follows that the only non-trivial responses are those of S 2 and S 3 (i.e. of the smoothed variance and skewness of the matter density field respectively).From eqs. (2.4), (2.16) and (2.17), we obtain: ) where Σ 3 is given by eq.(2.19).Substituting eqs.(2.21) and (2.22) in eq. ( 2.12), we obtain the following equation for the scale-dependent bias in this simple LPnG model: (2.23) Using eq. ( 2.17), we can express Σ 3 (∂ log n g /∂S 3 ) = (5/3) × ∂ log n g /∂f NL and thus rewrite the non-Gaussian contribution to the galaxy bias in its most often quoted form: where the quantities β f (z) = 2∂ log n g /∂ log σ and β g (z) = 3∂ log n g /∂f NL denote the responses of the abundance of matter density tracers to the variance and skewness of the smoothed matter density field, respectively.They are dependent on the complicated physics of galaxy (or other tracer) formation and may not be expressible analytically.In a galaxy power spectrum forecast, we need to make modelling choices for the parameters β f and β g .A strong modelling choice is to assume that the mean tracer density n g has the form of a universal halo mass function -wherein n g is only a function of the ratio δ c /σ, where δ c is the threshold density for spherical collapse in an Einstein-deSitter universe.In this case, one can show that β f = 2b L δ c .Modelling β g requires further assumptions about how the tracer density n g depends on f NL .For our SPHEREx forecasts (under strong modelling assumptions), we use the following fitting formula derived in ref. [34] (see their eq.( 51)): where κ 3 = 0.000329(1 + 0.09z)b −0.09L and ν = (1 + 0.5β f ) 1/2 .The universality assumption has been widely used to model the scale-dependent galaxy bias induced by LPnG.However, hydrodynamical separate universe simulations (see refs.[43,44]) have shown that the universality assumption may not hold true for realistic galaxy populations.In particular, ref. [43] has shown that for galaxies selected by stellar mass, β f can be modelled more accurately as β f = 2δ c (b E − p), where 0.4 ≤ p ≤ 0.7, and b E is the (scale-independent) Eulerian bias.Inspired by this result, we also obtain MCMC forecasts with a non-universal modelling choice whereby β f = 2δ c (b E − 0.5) and β g is given by eq.(2.25), with ν = (1 + δ c (b E − 0.5)) 1/2 .In section 4.3, we study the robustness of our forecast constraints to these choices by allowing p to be a free parameter.
If f NL and g NL are of comparable magnitude, then our expressions for β g and β f show that the contribution of g NL to the scale-dependent bias signal is suppressed w.r.t.contribution of

JCAP05(2024)094
f NL by four orders of magnitude.In this scenario, the primordial correlation functions scale hierarchically, e.g. S N +1 ∼ S N σ, and scale-dependent bias is primarily probing the amplitude of the largest non-trivial cumulant, f NL .In our forecasts with non-zero g NL we focus on the opposite limit, in which g NL ∼ 10 4 f NL , leading to a primordial kurtosis of comparable magnitude to the primordial skewness, S 3 ∼ S 4 , so that both terms contribute comparably to the scale-dependent bias.This situation is consistent with current constraints on both parameters and is also the only scenario in which g NL can be detected in the foreseeable future.
As b E is typically taken as a free parameter for each tracer population, eq. ( 2.24) shows that the parameters f NL and g NL are degenerate with respect to their impact on the scaledependent bias.In fact, for a single redshift and single tracer, f NL and (β g /β f )g NL cannot be distinguished from each other by galaxy power spectrum measurements alone.However, observations of multiple tracer populations, which typically have different values of β g (z) and β f (z), may break this near-perfect degeneracy.In this paper, we perform MCMC forecasts to quantify the possible degradation in the joint constraints on both f NL and g NL that can be obtained from the data collected by the SPHEREx all-sky survey [28].

Stochastic galaxy bias from local primordial non-Gaussianity
As stated before, any model of cosmic inflation that produces LPnG in the primordial curvature perturbations necessarily needs to include at least one light field in addition to the inflaton.A simple model that can produce non-Gaussian primordial fluctuations compatible with current measurements of f NL and g NL incorporates a single massless field called the curvaton [11][12][13] in addition to the inflaton, such that the inflationary expansion of the universe is driven mainly by the dynamics of the inflaton, whereas the primordial fluctuations are completely dominated by fluctuations in the curvaton.If both the inflaton and the curvaton contribute comparably to generating the primordial curvature perturbations, the primordial statistics are changed.In the simplest scenario, the contribution of the inflaton remains Gaussian, while the curvaton contribution remains non-Gaussian and admits an expansion like eq. (2.13), which is local in real space in terms of a Gaussian random variable, i.e.:

JCAP05(2024)094
a scale-independent constant (for generalisations of this, see, e.g., ref. [19]).Assuming this simple model, the bispectrum and trispectrum can be computed as, (2.29) From above, we identify the parameters f NL and g NL that characterize the bispectrum and trispectrum of ζ as, (2.31) (2.32) We see then that for a given level of non-Gaussianity in the ζ C field, allowing ζ to have a contribution from Gaussian inflaton fluctuations reduces the overall amplitude of non-Gaussianity by factors of ).On the other hand, the collapsed (O(f 2 NL )) part of the trispectrum is enhanced relative to the bispectrum by a factor of (1 + ξ 2 ).Defining we see that τ NL is enhanced by contributions from the inflaton [21].The most recent Planck datasets [45] lead to the constraint τ NL < 2800 at 95% confidence -this constraint is consistent with the assumption of weak primordial non-Gaussianity (⟨ζ N ⟩/⟨ζ 2 ⟩ N/2 ≪ 1).Under these conditions, the only non-trivial higher point functions of the primordial curvature perturbations (eq.(2.26)) are the bispectrum and the trispectrum of ζ C , which have the functional forms given by eqs.(2.14) and (2.15), respectively.
To analyse the impact of a primordial curvature perturbation of the kind given by eq.(2.26) on the galaxy number density fluctuations, we note that small-scale correlation functions of the curvature perturbation ζ couples only to an ambient long-wavelength mode of the non-Gaussian part of primordial curvature perturbation, i.e. ζ C .In particular, eq.(2.1) holds in a slightly modified form where now only the non-Gaussian contribution to the curvature, i.e. ζ C , modulates the small-scale correlation functions, whereas for N ≥ 2 the curvature perturbations sourced by the inflaton do not modulate small-scale correlation functions, -10 -

JCAP05(2024)094
The moments of the smoothed small-scale linear matter density field S N (R, z) (defined in eq. ( 2.3)) in the presence of a long-wavelength curvature perturbation From the arguments of the previous section, it follows therefore that the galaxy number density fluctuations obey a relation similar to eq. (2.7), namely where ) is the large scale (linear) density perturbation.Eq. (2.38) shows that the galaxy number density contrast δ g = δn g /n g is not fully correlated with the ambient large scale matter density field δ L -it is also partially correlated with ζ C,L , which has fluctuations that cannot be completely specified from the realisation of ζ or δ L alone.This introduces stochasticity in the distribution of galaxies, relative to the matter field, and is known as stochastic bias [20].From eq. (2.38), the galaxy bias b g and the galaxy power spectrum P gg are given by where, as before, b L = ∂ log n g /∂δ L is the scale-independent, Lagrangian bias, while the non-Gaussian, scale-dependent correction to the galaxy bias is where α(k, z) is defined in the previous section (eq.(2.9)) and where where b E = 1 + b L is the Eulerian galaxy bias, P cc (k, z) is the linear CDM+Baryon power spectrum, and the reader is reminded that ξ 2 represents the ratio of inflaton-to-curvaton power, eq.(2.28), so that τ NL = (6 5 f NL ) 2 (1 + ξ 2 ).Equation (2.46) is the main result of this section, which shows that the galaxy power spectrum acquires a stochastic contribution proportional to the relative contribution of the inflaton to the primordial curvature perturbation w.r.t. the curvaton.We note here an interesting feature of this scenario is that the different k-dependences of the ∆b NG and ∆b 2  NG terms in eq.(2.46) allow the combination (1 + ξ 2 ) to be extracted from P gg even if the coefficient B C in eqs.(2.41) & (2.42) cannot be determined, for instance due to the challenges associated with modelling n g (S N ).
For the stochastic bias signal to be large enough for a viable measurement of τ NL , we need the stochastic part of the galaxy power spectrum to be comparable to its non-stochastic part.In other words, a robust measurement of τ NL can be made using the stochastic halo bias if (following eq.(2.46)), ξ 2 ∆b NG /(b E + ∆b NG ) 2 ∼ O(1) in some wavelength and redshift bin which has a significant constraining power.Our choices of the fiducial values for τ NL will be motivated by imposing this requirement at z = 0 and at the largest scales accessible to the SPHEREx survey (i.e.k min = 0.001 h/Mpc).
In general, the parameter τ NL , which parameterises the collapsed limit of the primordial 4-point function, can be shown to satisfy the Suyama-Yamaguchi inequality [21,22]: τ NL ≥ ( 6

JCAP05(2024)094
f NL , g NL and τ NL around fiducial values consistent with state-of-the-art constraints from CMB datasets.This is a standard Bayesian Markov Chain Monte Carlo (MCMC) method used to obtain parameter forecasts for cosmological surveys -in particular, we use a multitracer galaxy power spectrum likelihood (developed in ref. [46]) for a SPHEREx-like galaxy survey constructed within the MontePython v3.4 MCMC sampling package [47,48].In this section, we give a broad overview of our MCMC analysis method used for the forecasts obtained in this paper and refer the reader to ref. [46] for more details.
Our MCMC analysis involves the following steps: • Modelling the observables (in this case the galaxy auto-and cross-power spectra at different scales and redshifts) for a given set of cosmological parameters.
• Obtaining the posterior likelihood for a given set of cosmological parameters (using the modelling framework defined above) given a 'mock data' provided by observables modelled in the fiducial cosmology.
• Obtaining confidence intervals for the cosmological parameters by sampling the posterior likelihood function using Markov Chain Monte Carlo (MCMC) methods.Alternatively, one could obtain (potentially) less accurate confidence intervals by computing a Fisher matrix.In this paper, unless otherwise specified, we sample the posterior likelihood function using the Metropolis Hastings MCMC sampling algorithm within the MontePython v3.4 MCMC sampling package.
Obtained in this way, the mean cosmological parameters should be the same (within numerical uncertainties) as the fiducial parameters, whereas the confidence intervals provide the forecast sensitivity to the cosmological parameters.The SPHEREx all-sky survey [28] is designed to minimise systematics in observations of galaxy clustering at large scales (going to k min = 0.001 h/Mpc).It will map a large cosmic volume and measure spectroscopic redshifts of galaxies belonging to the all-sky catalogues of the WISE [49], Pan-STARRS [50] and DES [51] galaxy surveys.It is also a deep survey -observing galaxies in eleven redshift bins z min i , z max i i=11 i=1 with mean redshifts z i going from z min = 0.1 to z max = 4.3.Due to its low spectral resolution (as compared to a true spectroscopic galaxy survey), SPHEREx is expected to measure spectroscopic galaxy redshifts to varying degrees of precision and will divide its observed galaxies into five galaxy populations according to their redshift uncertainties.In this way, the data collected by SPHEREx will lend itself to the use of multitracer techniques which can be used to obtain significantly tighter constraints on primordial non-Gaussianity.The observables to be obtained from SPHEREx are the galaxy power spectra and galaxy bispectra in redshift space -in this paper, we study the constraining power of the SPHEREx all-sky survey in constraining LPnG parameters using the galaxy power spectra alone.As discussed, galaxy power spectra alone could provide a first detection of any LPnG that could later be studied in more detail through higher-order correlation functions.Note that since SPHEREx will divide its observed galaxies into five galaxy populations, our observables include the cross-power spectra as well the auto-power spectra of different galaxy samples.
As the first step in our likelihood analysis, we need to model the galaxy power spectra for a given set of cosmological parameters.To do so, we use the linear bias relation (δ g = b g δ cb ) -13 -

JCAP05(2024)094
that is valid (at the linear level) at late times and in the comoving synchronous gauge (see eq. (2.46)).According to the linear bias model, the theoretical galaxy power spectrum in the presence of LPnG between the i th and j th SPHEREx galaxy samples is given by where b g,i and ∆b NG,i are respectively the total galaxy bias (including non-Gaussian corrections) and the non-Gaussian contribution (due to both f NL and g NL ) to the galaxy bias of the i th SPHEREx galaxy population.As mentioned in the previous section, P cc is the linear CDM+Baryon power spectrum.The parameter ξ 2 is defined such that τ NL = (6/5f NL ) 2 (1+ξ 2 ) and encodes the relative contribution of the inflaton w.r.t. the curvaton in the primordial curvature fluctuations.As shown in ref. [46] the sensitivity of SPHEREx to f NL is entirely driven by low-k Fourier modes where use of the linear power spectrum and linear bias factors are justified.Since the k-dependence of the bias is identical for f NL , g NL and τ NL , this is also justified here.We use the Boltzmann code CLASS v3.2.0 [52][53][54] to compute the linear CDM+Baryon power spectrum P cc (k, z), and model the non-Gaussian correction ∆b NG to the galaxy bias according to the analysis of section 2 -in particular following eq.( 2.24) for a cosmology with non-zero f NL and g NL and eq.(2.46) for τ NL ̸ = ( 6 5 f NL ) 2 .For the galaxy power spectrum forecasts obtained in this paper, we model the nuisance parameter β f as Our first sets of forecasts will fix either p = 1 and p = 0.5 with b E a free parameter for each galaxy sample and each redshift.We subsequently explore consequences of treating p as a free parameter for each SPHEREx galaxy sample.The choice p = 1 holds in the particular case when the mean number density of observed galaxies has the same form as a universal halo mass function.If this universality assumption does not hold true, the value of β f (which in general is tracer-dependent) is difficult to estimate from first principles as it encodes all the intractable details of the non-linear physics of galaxy formation.However, the results of hydrodynamical simulations within the separate universe framework conducted in refs.[43,44] suggest that β f is related to the galaxy bias as β f = (b E − p)δ c with a tracer-dependent number p. Our choice of p = 0.5 is representative of modelling for β f following the observations in refs.[43,44] that for galaxies selected by stellar mass, 0.4 ≤ p ≤ 0.7 is a better fit to the results of hydrodynamical separate universe simulations.SPHEREx galaxy samples are not selected by stellar mass, so we forecast two choices of p to illustrate the p-dependence of constraints.
For modelling β g , we use the fitting formula eq.(2.25) derived in ref. [34] which relates β g to β f through ν = (1+0.5βf ) 1/2 , with β f as in eq.(3.2).Note that our first-round forecasts are obtained using these strong modelling choices that assume a perfect knowledge of the functional dependence of the nuisance parameters β f and β g on b E .We also study forecasts under a more conservative modelling choice where p a free, but redshift-independent, parameter.If p were completely free (i.e.allowed to vary independently at different redshifts), then in the limit of large uncertainty on p, β f and β g would become free-parameters independent of b E .In the case where there is large uncertainty on the form of β f and β g , f NL and g NL cannot be distinguished from scale-dependence bias alone.For our free-p analysis the parameters -14 -

JCAP05(2024)094
become sufficiently degenerate that we restrict to forecasting f NL and τ NL only.In this limit constraints on f NL can be thought of as representing constraints on an effective amplitude that combines LPnG and tracer-dependent bias parameters, such as the B parameter in eq.(2.12), while (1 + ξ 2 ) = τ NL /( 65 f NL ) 2 should, in principle, still be constrained even with completely free β f .Throughout our analysis, we use δ c = 1.42 because it shows better agreement between the above fitting formulas and N-body simulations [34].
The galaxy power spectrum observed by a given galaxy survey differs from the theoretical galaxy power spectrum modelled according to eq. (3.1) in several ways.The peculiar velocities of galaxies along the line-of-sight break the isotropy of the matter power spectrum and cause the observed galaxy power spectrum to only be isotropic around the line of sight.This effect is known as redshift-space distortion and can be modelled at the linear and mildly non-linear scales relevant for the SPHEREx survey by making the replacement b g,i → b g,i + f (z)µ 2 and b g,j → b g,j + f (z)µ 2 in eq.(3.1);where f (z) is the linear growth rate of small-scale matter density perturbations f (z) = d log δ cb /d log a, and µ is the cosine of the angle between the wave vector ⃗ k and the line of sight (this is the well-known Kaiser term [55]).In addition to redshift-space distortion effects, there are additional multiplicative factors needed to map eq.(3.1) to a model of the observed galaxy power spectra.Our final expression for the observed galaxy power spectrum between the i th and j th galaxy samples of SPHEREx is (see ref. [46]) where f AP (z) denotes the contribution of the Alcock-Paczynski effect and f BF (k, z) encodes the effect of non-linear bulk flows.The Alcock-Paczynski term corrects for the error in the galaxy power spectrum measured from data assuming a fiducial cosmology which is different from the true cosmology.The term f AP (z) relates the comoving volume elements in the fiducial and true cosmologies and if given by: (3.4) The term f BF (k, z) models the smearing of (small-scale) BAO features in the galaxy power spectrum due to the effect of non-linear bulk flows [28].Following ref. [28], The Lagrangian displacement fields Σ ⊥ and Σ || are given by [28] Σ where D(z) is the linear growth factor and f (z) is the linear growth rate.We set the parameter c rec = 0.5 [28] and Σ 0 = 11 h/Mpc for a fiducial σ 8 = 0.8.The terms f σ z,i (k, z), f σ z,j (k, z) are exponential damping factors that suppress the observed galaxy power spectrum signal at small scales along the line of sight due to the uncertainty in the determination of galaxy redshifts (namely σ z,i and σ z,j which denote redshift uncertainties in the i th and j th galaxy populations to be observed by SPHEREx).For the i th SPHEREx galaxy sample, , where H(z) is the Hubble rate at redshift z and σ z,i is the redshift uncertainty of the i th galaxy sample.We refer the reader to refs.[28,46] for a detailed discussion of these additional terms and their importance in the context of the SPHEREx all-sky survey.
Proceeding further, we treat the observed linear galaxy number density fluctuations at a given (fiducial) k-mode and redshift as Gaussian random variables [46] with covariance (assuming Poisson shot noise) given by where P ij denotes the observed the cross-power spectrum of the i th and j th galaxy populations of SPHEREx (see eq. ( 3.3)) and n g,i is the mean number density of observed galaxies in the i th galaxy population of SPHEREx at redshift z.This translates into an exponential likelihood for the galaxy power spectra at a given k-mode and redshift given by where the proportionality constant J is independent of the sampled cosmological parameters.In eq.(3.9), both C and D are the covariance matrices of galaxy number density fluctuations (at a given k-mode and redshift) that respectively depend on the sampled, s , and fiducial, fid , cosmological parameters in an MCMC run: (3.11) The entire SPHEREx power spectrum likelihood is the product of L k,z over all redshifts bins and independent k-modes.As stated before, the SPHEREx survey is a large-scale survey [28] with the observed k-modes between k min = 0.001 h/Mpc and k max = 0.2 h/Mpc and eleven redshift bins with mean redshifts going from z min = 0.1 to z max = 4.3.Note that for the purposes of obtaining power-spectrum forecasts, we assume that different redshift bins and independent k-modes are uncorrelated -this is a standard approximation used in forecasting the sensitivity of cosmological surveys (see refs.[28,46,56]).An MCMC run used to obtain cosmological forecasts typically maximises the likelihood by minimising χ 2 defined as where V (z i ) is the comoving volume surveyed in the i th redshift bin.
Having thus obtained the likelihood, we use the boosted MCMC sampler MontePython v3.4 to conduct MCMC runs to fit for the following cosmological parameters: {ω b , ω cdm , 100θ s , A s , n s , z reio , M ν }, as well as combinations of {f NL , g NL , τ NL }.In addition, we marginalise over the scale-independent, Gaussian parts of the galaxy biases (given by b 0,i = b g,i − ∆b NG,i ) as well as the bulk flow parameter Σ 0 which parametrises the effect of non-linear bulk flows - 16  given by f BF (k, z) in eq. ( 3.3) and eqs.(3.5)-(3.7)(see refs.[28,46] for more details).During our MCMC analysis, we combine our power spectrum likelihood (eq.(3.12), (3.9)) with a mock CMB power spectrum likelihood named 'fake_planck_realistic' in MontePythonthis likelihood is designed to mimic the sensitivity of the Planck 2018 dataset in constraining cosmological parameters (see ref. [57] for more details).This mock CMB likelihood does not include the CMB lensing power spectrum -this ensures that it is not significantly correlated with our SPHEREx galaxy power spectrum likelihood.Moreover, the mock CMB likelihood does not include the effect of primordial non-Gaussianity and does not by itself provide any constraints on LPnG parameters -this makes it equivalent to a Planck-based prior imposed in order to resolve degeneracies between νΛCDM cosmological parameters (as is done in ref. [28]).
We use a νΛCDM cosmology with fiducial parameters ω b = 0.02218, ω cdm = 0.1205, 100 × θ s = 1.04113,A s = 2.032692 × 10 −9 , n s = 0.9613, z reio = 7.68, M ν ≡ Σ i m νi = 0.06 eV.For simplicity neutrinos are implemented as three degenerate species (see, e.g., ref. [58]).In this paper, we forecast constraints on non-Gaussian initial conditions with primordial trispectra that contribute comparably in magnitude to galaxy power spectra as primordial bispectra with f NL ∼ O(1).To achieve this scenario in which the primordial kurtosis is comparable to the skewness we choose fiducial values of f NL = 1.0 and g NL = 1.0 × 10 4 .For the two-source initial conditions described in section 2.2 that produce a τ NL not specified by f NL alone, we choose to obtain forecasts around two fiducial values of τ NL , namely τ NL = 1.3 × 10 2 and τ NL = 1.3 × 10 3 , well above our fiducial value of f NL = 1 and consistent with CMB datasets [45].

Results
We start by presenting single-parameter forecasts for f NL and g NL which parametrise the squeezed-limit primordial bispectrum and trispectrum respectively.Table 1 shows the SPHEREx multitracer MCMC forecasts for the individual parameters f NL and g NL around their respective fiducial values f NL = 1.0, and g NL = 1.0 × 10 4 obtained within the framework of the two-field inflationary model.Note that for the two-field inflationary model considered in this paper, it doesn't make sense to obtain individual forecasts for τ NL because τ NL cannot be non-vanishing unless f NL ̸ = 0.The forecasts in table 1 are obtained under strong modelling assumptions (see section 3, eqs.(2.24) and (2.25) and surrounding text), where 1/2 with p = 0.5 or p = 1 for all galaxy samples, yet b E is a free parameter for each galaxy population at each redshift.As expected, the forecast constraints on f NL and g NL depend on the assumed value of p.Our individual forecasts for f NL or g NL JCAP05(2024)094 alone are also consistent with the observation at the end of section 2.1 that the contribution of g NL to the scale-dependent bias signal is suppressed w.r.t. that of f NL by four orders of magnitude, due to the relative magnitudes of the dimensionless cumulants ⟨δ 3 cb ⟩/⟨δ 2 cb ⟩ ∼ f NL σ and ⟨δ 4 cb ⟩ c /⟨δ 2 cb ⟩ ∼ g NL σ 2 for f NL and g NL of similar size.Having established the forecast constraints on either f NL or g NL alone, we now consider combined constraints.We conduct MCMC runs to obtain joint forecasts for the pairs of parameters (f NL , g NL ) and (f NL , τ NL ).Expanding this to include all three parameters is a challenging and computationally intensive sampling problem due to strong degeneracies between the parameters.We expect the sensitivity would be washed out and barely more illuminating than the two-parameter case, since we would only see notable constraining power for the parameter(s) with the greatest contribution to the observable signal.The exception to this is the somewhat fine-tuned case where all three parameters contribute comparably, in which case we would expect decreased sensitivity due to the strong three-way degeneracy.For the subcase with two parameters giving comparable contributions to the scale-dependent bias, we expect the result would largely asymptote to that of the two-parameter case, possibly with slightly decreased sensitivity to the parameters.Given this expectation, we do not consider it worth the significant computational resources required (and associated carbon footprint) and omit the three-parameter case.

Joint forecasts on f NL and g NL
To obtain joint forecasts for f NL and g NL , we conduct MCMC runs with the galaxy bias modelled as eq.((2.24)), with β f and β g as in eq.(3.2) and eq.(2.25), and simultaneously fit for f NL and g NL .Table 2 shows our joint MCMC forecasts for f NL and g NL as well as their covariance around fiducial values f NL,fid = 1.0 and g NL,fid = 1 × 10 4 .Our fiducial values for f NL and g NL are consistent with the most recent CMB constraints [18] and are such that the dimensionless skewness and kurtosis of the primordial curvature (⟨ζ 3 ⟩/⟨ζ 2 ⟩ 3/2 and ⟨ζ 4 ⟩/⟨ζ 2 ⟩ 2 ) are of similar order of magnitude.As mentioned in section 2.1, both f NL and g NL contribute comparably to the observed signal of primordial non-Gaussianity (i.e. the scale-dependent bias) in this regime.We again consider two representative values of p, p = 1 and p = 0.5.Controlling for the difference in the order of magnitude of f NL and g NL , table 2 shows that we obtain nearly equal forecast constraints for f NL and g NL , analogous to the individual forecasts in table 1.However, these forecasts are significantly weaker than their respective individual forecasts.This degradation in forecasts for f NL and g NL can be attributed to the strong degeneracy between f NL and g NL , which is evident from eq. (2.24).Both terms produce the same k-dependence in the bias, the terms differ only in how the amplitude of the k-dependent terms vary with tracer population through β f and β g , which are themselves degenerate through the dependence on b E .This degeneracy is reflected in figure 1, which shows the joint 68% and 95% confidence regions for f NL and g NL obtained from our MCMC analysis.The degeneracy between f NL and g NL is also reflected in their covariance (see table 2), which is ∼ −0.9 for both the p = 1 and p = 0.5, with only slight variations in the covariance between the two p-values (visible in figure 1) due to the small changes in the relationship between β g and b E in each case.This is consistent with the Fisher matrix analysis done by ref. [59], where it was found that while multitracer analysis using JCAP05(2024)094 scale-dependent halo bias can distinguish between f NL and g NL , their covariance is always close to −1 due to their degeneracy w.r.t. the scale-dependent bias.
As shown in table 1, our SPHEREx multitracer forecasts under non-Gaussian initial conditions for only f NL ̸ = 0 or only g NL ̸ = 0 are such that ∆f NL /f NL ∼ ∆g NL /g NL ∼ 1 for fiducial values f NL = 1.0 and g NL = 1.0 × 10 4 .This means that with systematic uncertainties eliminated and in the scenario where either f NL ̸ = 0 or g NL ̸ = 0, SPHEREx power spectrum measurements alone can potentially detect LPnG (at the small amplitudes here) with at most 1σ significance (with the modelling assumptions made here).It is worth noting that this continues to hold even with both f NL and g NL being as large as their fiducial valuesin fact, the prospect of SPHEREx detecting LPnG with 1σ significance is slightly improved when both are non-zero.For forecasts obtained with p = 1, the point f NL = g NL = 0 is now disfavoured with ∆χ 2 = 4.6, which puts it far outside the 68% confidence region contour (as compared to the forecasts made for a single parameter in which case the point f NL = 0 is nearly at the boundary of the 68% confidence interval).On the other hand, the point f NL = g NL = 0 is disfavoured with ∆χ 2 = 6.5 for forecasts obtained using the modelling assumption p = 0.5 -suggesting that SPHEREx power spectrum measurements alone can potentially detect LPnG at the (small) amplitudes considered here with as much as 2σ significance if SPHEREx galaxies are better characterized by p = 0.5. 11

Joint forecasts on f NL and τ NL
We now consider joint constraints on f NL and τ NL , with g NL = 0. To obtain joint multitracer forecasts for f NL and τ NL , we conduct two MCMC runs with fiducial values f NL = 1, τ NL = 1.3 × 10 2 and f NL = 1, τ NL = 1.3 × 10 3 respectively.As mentioned in section 2, these fiducial 11 ∆χ 2 = 2.3(6.17) for the boundary of the 68.3%(95.4%)confidence region.
- 19    values are well within the current constraints, yet produce a large stochastic component of the galaxy bias.Note that while conducting our MCMC runs, we need to impose the condition τ NL ≥ (6/5f NL ) 2 following the Suyama-Yamaguchi inequality [21].Table 3 shows the joint MCMC forecast around fiducial values f NL = 1, τ NL = 1 × 10 2 and f NL = 1, τ NL = 1 × 10 3 , respectively, for p = 1.0 and p = 0.5.Figure 2 shows the joint 68% and 95% contours for f NL and τ NL around their respective fiducial values for examples with p = 1 and p = 0.5.In all cases, the covariance between τ NL and f NL is ≈ −0.6 -this is consistent with the fact that f NL and τ NL are only partly degenerate.This is because the contribution of τ NL to the galaxy bias (unlike that of g NL ) goes as ∆b 2 NG (∝ 1/k 4 on large scales, see eq. (2.46)) and does not have the same scale-dependence as the contribution of f NL (∝ 1/k 2 on large scales).The contour plots in figure 2 as well as the covariances in table 3 indeed show that while f NL and τ NL are partly degenerate, they are not as degenerate as f NL and g NL .We also see that the choice of p does not appear to effect the degeneracy between f NL and τ NL .Comparing the left and right panels of figure 2  we see that constraints on f NL are weakened in the presence of larger τ NL .Our MCMC forecasts lead to the observation that in the absence of systematics, SPHEREx can potentially constrain large values of τ NL more precisely at the cost of degraded constraints on f NL .In particular, in the scenario with τ NL = 1 × 10 3 , SPHEREx could detect non-trivial τ NL at ≈ 5σ significance while not, at the same time, having unambiguously detected f NL ̸ = 0.This means that SPHEREx could infer the presence of an additional light field during inflation (and thus indirectly rule out f NL = 0) through a robust measurement of scale-dependent stochastic bias, if the amplitude is sufficiently large.

Impact of modelling choices
Our single-parameter, as well as joint, MCMC forecasts for f NL , g NL , and τ NL (see tables 1, 2, 3) make it amply clear that the constraining power of any galaxy survey in constraining LPnG parameters using the scale-dependent galaxy bias varies with the modelling of the nuisance parameters β f and β g (see subsection 2.1).Given that β f and β g are strongly dependent on tracer population, they are difficult to model from first principles and present a big limiting factor to fully exploiting galaxy power spectrum data to constrain LPnG.The non-trivial relation between β f and the Eulerian galaxy bias with Gaussian initial conditions has been investigated with hydrodynamical separate universe simulations (see refs.[43,44]), which provide some evidence that supports modelling β f as with p being a tracer-dependent number.However, modelling β g in terms of other observable parameters, such as b E is a much more intractable problem and represents a bigger modelling ambiguity.The is due both to the complicated nature of modelling n g from first principles, and the observation that even in the simplest models, β g depends non-linearly on b E , and the β g for each sample will therefore depend on how tracers are selected and grouped together in data analyses [34].
In this paper, we have so far obtained forecasts for the constraining power of the SPHEREx galaxy surveys under the strongest modelling choices, which assume a complete knowledge of the relation between β f , β g and b E .In particular, we assume (see section 3) and β g is given by eq.(2.25) with ν = (1 + 0.5β f ) 1/2 .Different tracer populations could have different values of p, which can change our MCMC forecasts proportionally (i.e.smaller p gives tighter constraints) given that the contributions of f NL and g NL are proportional to β f and β g respectively.We have so far presented constraints for p = 1.0 and p = 0.5 as representative values that reflect the universal halo mass function assumption and some -21 -JCAP05(2024)094 results of hydrodynamical separate universe simulations, respectively.A more conservative approach is to make β f and β g free parameters and impose priors on them in our analysis.It is important to marginalise over the nuisance parameters β f and β f to obtain a more complete picture of how well a given galaxy survey can constrain LPnG using the galaxy power spectrum.We will now discuss how generalising modelling assumptions on β f and β g affects our forecast constraints.
If β f and β g are free parameters, eq.(2.24) illustrates that they are completely indistinguishable from scale-dependent bias measurements alone.In other words, without predictions for the values of β f and β g , the galaxy power spectrum loses its utility as a signal that can potentially distinguish between a primordial bispectrum and trispectrum.This situation is not all bad, however: it is still true that in the presence of squeezed-limits of higher-point primordial correlation functions, SPHEREx could potentially detect a non-trivial LPnG at high significance, yet with uncertainty in which squeezed-limit correlation functions were contributing.This is because the scale-dependent bias measurement is a measurement of a redshift-dependent weighted combination of the parameters f NL , g NL together with all the higher order coefficients in the expansion (see eq. (2.13)).Any measurement of f NL from observations of the scale-dependent bias of galaxies should therefore be interpreted as a measurement of a weighted combination of f NL and g NL (and potentially other higher-order squeezed-limit non-Gaussianity parameters) that depends on the observed tracer populations via the true values of β f and β g (as well as their counterparts for higher-order squeezed limit non-Gaussianity parameters).In particular, non-trivial squeezed limits of higher-point primordial correlation functions can potentially contribute to increasing the observable signal of LPnG -our results show that this happens around the fiducial values of f NL and g NL consistent with current constraints from CMB datasets chosen in this paper. 12n the scenario in which multiple higher-order statistics contribute to scale-dependent bias, one should then consider constraining an effective amplitude of LPnG determined by B(z) in eq.(2.12).If an effective amplitude B(z) is determined to be non-zero, this would already definitively prove the existence of LPnG and thereby additional light fields during inflationeven if the individual coefficients in eq.(1.1) were undetermined.The challenge is that B(z) depends on both the primordial non-Gaussian parameters and the tracer populations.An individual measurement of B from one redshift and tracer population will have much lower signal-to-noise than a measurement of LPnG that combines data from all redshift bins and tracer populations.Indeed, Fisher matrix computations for the most constraining galaxy sample of SPHEREx performed using our likelihood pipeline indicate that a non-trivial LPnG (with f NL ∼ 1, g NL ∼ 10 4 ) cannot be ruled out with any significance by fitting for the individual B(z i ) when these parameters are free to vary independently at each redshift z i .Consequently, one needs a model for B(z) -for instance, a model relating B(z) to the galaxy bias b E (z) or a parameterization of the redshift dependence of B(z) -to be able to harness the full power of the data to unambiguously detect LPnG at these threshold amplitudes.

Discussion
The simplest inflationary models that generate initial conditions with LPnG have a primordial curvature perturbation which has non-Gaussianity predominantly characterised by a bispectrum (parametrised by f NL ) and which is sourced wholly by the fluctuations of a light field separate from the inflaton.Upcoming galaxy surveys such as SPHEREx show good promise in being able to constrain such models from observations of the galaxy power spectrum (and bispectrum) (see ref. [28]).In this paper, we analyse how well the upcoming SPHEREx all-sky survey can constrain LPnG models beyond this simplest class of models.In particular, we consider a two-field inflation model in which the primordial curvature perturbations have a non-trivial trispectrum whose squeezed limit is parameterised by g NL (in addition to a squeezed-limit bispectrum parametrised by f NL ).We also consider the scenario in which primordial curvature fluctuations are sourced from fluctuations of both the inflaton and the curvaton.In this case, the contribution of inflaton fluctuations to the primordial curvature perturbations is Gaussian whereas the contribution of curvaton fluctuations shows LPnG.A substantial inflaton contribution to the primordial curvature perturbations decorrelates the collapsed limit of the primordial trispectrum from the square of the primordial bispectrum -thus introducing an additional parameter in the model (called τ NL ) which parametrises the collapsed limit of the primordial trispectrum.Thus, the non-Gaussianity in this model is characterised by three parameters namely f NL , g NL and τ NL .In this paper, we present -24 -

JCAP05(2024)094
joint MCMC forecasts for the LPnG parameters f NL , g NL , τ NL from galaxy power spectra obtained in the context of the SPHEREx all-sky survey.
In section 2.1, we show that the squeezed limits of arbitrary higher-point correlation functions of the primordial curvature perturbations all impart the same characteristic scaledependence to the bias of matter density tracers -which is an important signal targeted by large-scale galaxy surveys that intend to constrain LPnG using galaxy power spectra.SPHEREx galaxy power spectra could therefore constrain or detect LPnG due to arbitrarily high-order correlation functions.In the context of the two-field inflation model assumed in this paper, this means that we can obtain joint MCMC forecasts for the pair of parameters f NL and g NL using the scale-dependent galaxy bias.However, it also means that f NL and g NL are strongly degenerate w.r.t.their signature on the galaxy power spectrum, causing their joint MCMC forecasts to be substantially degraded compared to MCMC forecasts for the individual parameters f NL and g NL .Figure 1 and table 2 show our joint MCMC constraints on f NL and g NL around fiducial values compatible with state-of-the-art constraints from CMB datasets.As expected, the joint forecasts on f NL and g NL are degraded compared to their single-parameter forecasts (in table 1).The strong degeneracy between f NL and g NL is due to the fact that their contributions to ∆b NG have the same k-dependence and differ only in their dependence on the Eulerian bias, b E .This degeneracy is reflected in their covariance, ∼ −0.9, which also depends on modelling choices for ∆b NG , the scale-dependent bias (eqs.(2.24), (2.25), and (3.2).
Unlike f NL and g NL , the parameter τ NL parametrises the collapsed limit of a primordial correlation function and leads to a stochasticity in the galaxy power spectrum which increases at the largest scales (see section 2.2). Figure 2 and tables 3 show our joint MCMC forecasts for f NL and τ NL around fiducial values consistent with constraints obtained from CMB datasets.In all cases, the covariance between f NL and τ NL is ∼ −0.6.This is a reflection of the fact that while f NL and τ NL are quite degenerate in their impact on the galaxy power spectrum, they are not as degenerate as f NL and g NL owing to the fact the f NL terms in ∆b NG scale as 1/k 2 while the τ NL terms scale as 1/k 4 .Our results also indicate that the SPHEREx survey can potentially constrain a large, but CMB-consistent, value of τ NL tightly though with degraded constraints on f NL .As noted in section 4.2, a robust measurement of τ NL would also indicate the presence of an additional light field during inflation even if the parameter f NL could not be distinguished from zero.
In section 4.3, we note that the constraints on f NL , g NL and τ NL obtainable from galaxy power spectra depend on how one chooses to model the nuisance parameters β f and β g (see eq. (2.24)).Modelling ambiguity in β f as well as β g reduces the utility of the galaxy power spectrum in jointly constraining f NL and g NL .This is because the scale-dependent bias effectively probes a weighted (redshift-dependent) combination of the parameters f NL and g NL , where the relative weights are set by β f and β g .In a general LPnG scenario (e.g.eq.(1.1) to all orders in Z), scale-dependent bias measurements can conservatively only help us ascertain the presence/absence of an additional light degree of freedom (over and above the inflaton); modelling ambiguities in the bias expansion (eq.(2.12) precludes determining which non-trivial squeezed limits of higher-point functions are sourcing scale-dependent bias).To obtain model-independent constraints on different inflationary models that generate LPnG -25 -

JCAP05(2024)094
we need to exploit higher-order galaxy correlation measurements, which provide information complementary to measurements of galaxy power spectra.
Including the galaxy bispectrum (as well as higher-point galaxy correlation functions) in addition to the galaxy power spectrum may break the perfect degeneracy between the scale-dependent galaxy bias and the parameter β f and can tighten the constraint on f NL .However, theoretically modelling galaxy bispectra introduces additional bias parameters (such as the b ϕδ (z) parameters in [30]) over and above the parameter β f defined in eq.(2.24) (see for example [30]) -this presents additional modelling challenges and may necessitate modelling assumptions which significantly affect constraints on f NL , g NL and τ NL . 14Moreover, constraints on the LPnG parameter f NL obtained from the bispectrum really come from information at non-linear scales unlike constraints obtained from the scale-dependent bias which are dominated by information from linear/quasi-linear scales (see [46]).On the other hand, as discussed previously, our joint forecasts show that galaxy power spectra alone can provide a powerful test of single-field inflation free from such complications and more agnostic to which statistics are sourcing the scale-dependent bias.We leave a detailed analysis of the constraining power of galaxy bispectrum (as well as the galaxy trispectrum and other higher-order galaxy correlation functions) obtained by SPHEREx in constraining LPnG to future work.
As for f NL and τ NL , we note in section 4.3 that our joint forecasts are altered by ≲ 20% under more conservative modelling choices than the ones we have assumed in our primary analysis.This means that the conclusion that SPHEREx can potentially measure τ NL more precisely using galaxy power spectra at the cost of degraded constraints on f NL likely holds true even under more conservative modelling choices.Adding galaxy bispectrum data can only make the constraints better [32].
Since modelling ambiguities represent a big limitation in using galaxy power spectrum data to constrain LPnG, it would be interesting to investigate optimal ways of mitigating their impact.In particular, the fact that they are late-time parameters and encode the non-linear physics of galaxy formation leads to the observation that they should couple to primordial non-Gaussianities as well as non-Gaussianities generated through gravitational non-linearities in the same way.This approach might help suggest novel ways to break their degeneracy with LPnG parameters.We plan to investigate along these lines in the future.
.27) where ζ I is the Gaussian contribution to the primordial curvature perturbation that is sourced from the inflaton fluctuations, whereas the non-Gaussian ζ C is sourced from the fluctuations of the curvaton.The quantities f NL and g NL characterise the amplitudes of the bispectrum and trispectrum of the ζ C field.The correlation between the fields ζ I and ζ C is dependent on the particular features of the inflationary model.For this work, we consider a simple model in which ζ I and ζ C are uncorrelated and have proportional power spectra.We define the contribution of the inflaton to the primordial curvature perturbation w.r.t. the curvaton contribution by ξ 2 = P ζ I /P ζ C , (2.28) -9 -

Figure 1 .
Figure 1.Joint MCMC forecasts around fiducial f NL = 1.0 and fiducial g NL = 1 × 10 4 obtained from the galaxy auto-and cross-power spectra of all five SPHEREx galaxy samples.The two sets of contours illustrate constraints for two modelling choices for the scale-dependent bias, p = 1 and p = 0.5.

Figure 2 .
Figure 2. Joint MCMC forecasts for f NL and τ NL (obtained around the fiducial values assumed in this paper) from the scale-dependent bias of all five SPHEREx galaxy samples for p = 1 and p = 0.5 modelling choices.

Figure 3 .
Figure 3. Joint Fisher forecasts for f NL and τ NL around f NL = 1 with τ NL = 1.3 × 10 2 (left) and τ NL = 1.3 × 10 3 (right) and β f = δ c (b E − p) as a function of priors on p.The red plots show forecasts obtained with fiducial p = 0.5 and the black plots show forecasts obtained with fiducial p = 1.0.
.42) Equation (2.38) can now be written in a concise form to give the modified linear bias relation for (Lagrangian) galaxy number density contrasts, L is the part of the density field sourced by ζ C .The only non-trivial responses of the smoothed moments of the density field (for the model given by eqs.(2.26) and (2.27)) can be obtained using eq.(2.37) and are given by: Σ 3 is given by eq.(2.19).Substituting eqs.(2.44) and (2.45) in eq.(2.42), we see that the non-Gaussian correction to the galaxy bias ∆b NG in the two-field initial conditions given by eq.(2.26) is the same as eq.(2.24).The galaxy power spectrum P gg = ⟨δ g δ g ⟩ on the other hand can be shown to be,8 where

Table 1 .
Single-parameter LPnG forecasts for f NL and g NL obtained from the SPHEREx multitracer likelihood for the p = 1 and p = 0.5 modelling choices (see eq. (2.25) and eq.(3.2)).

Table 2 .
Joint MCMC forecast for f NL and g NL obtained from the SPHEREx multitracer likelihood for p = 1 and p = 0.5 modelling choices.

Table 3 .
or the two sub-tables in table 3, we see that larger fiducial values of τ NL lead to tighter constraints on τ NL .This is, perhaps, unsurprising as a large τ NL should lead to a larger signal.On the other hand, Joint MCMC forecast for f NL and τ NL obtained from the SPHEREx multitracer likelihood.For each fiducial value of τ NL , we consider two example values of p = 1 and p = 0.5, entering the model for β f as in eq.(3.2).