Big-Bang Nucleosynthesis

Big-Bang nucleosynthesis (BBN) turns the expansion history at temperatures near an MeV into a set of primordial light-element abundances. Helium-4 is mainly a clock for weak freeze-out and neutron decay; deuterium is mainly a measure of the baryon density and the competition among nuclear burning rates. The calculation is therefore unusually powerful, but every numerical constraint is conditional on a specified reaction network, abundance sample, and statistical model. Observational and software statements on this page are current through 2026-08-19.

Required background. States, Processes, and the Laws of Thermodynamics supplies Thermodynamics and Statistical Mechanics for equilibrium ratios and thermal rates. Probability and Statistical Inference supplies normalized likelihoods, covariance, priors, and nuisance-parameter reasoning.

Helpful background. The Einstein Equation supplies the gravitational convention, while The Friedmann Equations supplies the radiation-era expansion rate. Relativistic Matter and Compact Objects supplies relativistic-fluid thermodynamics, and Classical Field Modes and Canonical Variables supplies the Flat-Spacetime Quantum Field Theory mode conventions used to label particle species. Free Fields, Creation Operators, and Fock Space then supplies occupation-state counting.

The weak clock and the neutron fraction

When weak reactions are rapid and lepton chemical potentials are negligible, chemical equilibrium fixes the neutron-to-proton ratio to

$$r_{np}^{\rm eq}\equiv\frac{n_n}{n_p} \simeq\exp\left(-\frac{Q}{T}\right), \qquad Q=m_n-m_p\simeq1.293\ {\rm MeV}.$$

This leading nonrelativistic estimate neglects the small neutron–proton mass-prefactor difference and finite-temperature and radiative corrections.

The interconversion rate scales approximately as $\Gamma_{n\leftrightarrow p}\sim G_F^2T^5$, whereas radiation domination gives $H\propto g_*^{1/2}T^2/M_{\rm Pl}$. Cooling therefore drives $\Gamma/H$ downward. Freeze-out is gradual rather than instantaneous, but the crossing near $T\sim1\ {\rm MeV}$ explains why the surviving neutron fraction is sensitive to the weak rates, the expansion rate, and the neutron lifetime. Free neutrons continue to decay before efficient nuclear assembly. Here $g_*$ counts energy-density degrees of freedom and $M_{\rm Pl}=G^{-1/2}$ is the nonreduced Planck mass.

Once almost every surviving neutron is bound into helium-4, counting two neutrons per helium nucleus gives

$$Y_p\equiv\frac{4n_{^4{\rm He}}}{n_b} \simeq\frac{2r_{np}}{1+r_{np}}.$$

Thus $r_{np}\simeq1/7$ gives $Y_p\simeq0.25$. This is an organizing estimate, not the precision algorithm: a network integrates the continuous weak rates, neutron decay, finite-temperature plasma effects, and nuclear reactions. The physical sequence and this estimate are developed in Fields, Molaro, and Sarkar 2026, §24.2, printed pp. 1–3 (Official PDF).

The deuterium gate and the reaction network

At first, energetic photons destroy deuterium faster than nuclei can retain it, even after $T$ falls below the deuteron binding energy. Expansion lowers the high-energy photon tail until the deuterium bottleneck opens. Reactions then rapidly channel neutrons and protons toward helium-4, with residual D, helium-3, and mass-7 nuclei left behind.

For abundances $Y_i=n_i/n_b$, a network has the schematic form

$$\dot Y_i=\sum_r N_{ir}\,\lambda_r(T,\rho_b,\ldots) \prod_j Y_j^{\nu_{jr}}, \qquad \eta\equiv\frac{n_b}{n_\gamma}\propto\Omega_bh^2.$$

Here $N_{ir}$ is the net number of species $i$ produced in reaction $r$, and $\nu_{jr}$ counts its incoming nuclei. Detailed balance supplies reverse rates in the thermal regime. In the baryon-density relation, $\eta$ means the post-electron–positron-annihilation ratio at a fixed present photon temperature, and $h\equiv H_0/(100\,\mathrm{km\,s^{-1}\,Mpc^{-1}})$. During pair annihilation the photon number per comoving volume changes, so the instantaneous $n_b/n_\gamma$ is not that constant late-time ratio. The network follows the conserved baryon density and the evolving plasma temperature. The main controls separate only approximately:

Control	Strongest diagnostic	Physical reason
weak rates and $\tau_n$	$Y_p$	set the neutron fraction and decay interval
$H(T)$ or $\Delta N_{\rm eff}$	$Y_p$ and D/H	change the time available for weak and nuclear processing
$\Omega_bh^2$	D/H	higher density burns deuterium more efficiently
$d(p,\gamma)^3{\rm He}$ and $d+d$ rates	D/H	control the leading deuterium-destruction channels

This division is useful for intuition, but a joint likelihood must evolve one consistent network. Changing $N_{\rm eff}$, for example, can affect both the background expansion and neutrino/plasma physics; changing a laboratory rate can move several nuclides at once.

What the abundance measurements actually select

The numerical example below uses two deliberately explicit inputs. For helium it adopts the LBT $Y_p$ Project IV result

$$Y_p^{\rm obs}=0.2458\pm0.0013.$$

The project began with 54 uniformly observed metal-poor H II regions; 47 passed its 95% fit-quality criterion, 41 survived additional reliability screening, and the quoted weighted mean uses 15 high-signal, low-metallicity targets with ${\rm O/H}\lt4\times10^{-5}$. The spectral model propagates electron temperature and density, collisional excitation, stellar absorption, radiative transfer, ionization, relative-flux and telluric calibration, and model-quality flags. It uses updated hydrogen and helium emissivity grids but assigns them zero quantified uncertainty; that is an assumption, not evidence that atomic emissivities are exact. The quoted error includes a dispersion scale factor $S=1.15$. These choices are documented in Skillman et al. 2026, §§2–4 (Project I, arXiv:2601.22232v1) and Aver et al. 2026, abstract and §§4–7 (Project IV, arXiv:2601.22238v1).

For deuterium the example uses the PDG 2026 twelve-system compilation,

$$D_5^{\rm obs}\equiv10^5({\rm D/H})=2.508\pm0.029.$$

The PDG first combines the twelve most precise recent measurements as $2.508\pm0.027$, retaining the most recent analysis when one absorber is repeated, then inflates the uncertainty by $S=1.08$ for the quoted result; see Fields, Molaro, and Sarkar 2026, §24.3, printed p. 3, eq. (24.2), and Table 24.1, printed p. 5 (Official PDF). This is not the only defensible compression. Kislitsyn et al. obtain $2.533\pm0.024$ from all compiled systems and $2.501\pm0.028$ from a nine-measurement precision subset; an earlier homogeneous seven-system analysis gives $2.527\pm0.030$. Absorber selection must therefore be part of the result. Relevant choices include spectral resolution and signal-to-noise, simplicity of the velocity structure, H I and D I blending, continuum placement, metallicity and ionization corrections, possible astration, and whether the fit was blinded. See Kislitsyn et al. 2024, §§2–4 and Table 4, printed pp. 2–10 and Cooke, Pettini, and Steidel 2018, §§4–5, printed pp. 12–14. The displayed errors are one-dimensional uncertainty compressions under the respective models; correlated continuum, calibration, absorber, and population choices remain sources of systematic uncertainty and model dependence.

Helium-3 is excluded here because stellar processing complicates a primordial extrapolation. Lithium-7 is also excluded from the likelihood: metal-poor stellar atmospheres give an abundance far below the standard-network prediction, and neither stellar depletion nor new early-universe physics has yet supplied a generally accepted resolution. Exclusion prevents that open problem from being hidden inside a two-observable fit. The observational status and lithium discrepancy are reviewed in Fields, Molaro, and Sarkar 2026, §§24.3 and 24.5, printed pp. 3–7 (Official PDF).

A bounded, versioned likelihood

A production calculation should run a named code and nuclear-rate release. PRIMAT includes precision weak and plasma corrections Pitrou et al. 2018, revised preprint 2025, §§II–IV, pp. 7–27 (Open PDF), while PArthENoPE 3 is an independent public network implementation with explicitly revised deuterium rates Gariazzo et al. 2022, §2, preprint pp. 3–6 (Open PDF). Agreement between codes is a valuable check, not a replacement for rate and input versioning. The local emulator below instead reproduces selected derivatives from PRyMordial runs with the PRIMAT and NACRE-II rate compilations reported by Burns.

Let

$$\delta_b=\ln\frac{\Omega_bh^2}{0.02237}, \qquad \delta_\tau=\ln\frac{\tau_n}{878.4\ {\rm s}},$$

and let $p_\gamma,p_n,p_p$ be lognormal, one-standard-deviation shifts of $d(p,\gamma)^3{\rm He}$, $d(d,n)^3{\rm He}$, and $d(d,p)t$. The PRIMAT linear response is centered on the Standard Model value $N_{\rm eff}^{\rm SM}=3.044$, so $\Delta N_{\rm eff}\equiv N_{\rm eff}-3.044$. It is

$$\begin{aligned} m_Y={}&0.24691\left(1+0.039\delta_b+0.731\delta_\tau\right) +0.014\Delta N_{\rm eff},\\ m_D={}&2.445\left(1-1.639\delta_b+0.420\delta_\tau\right) +0.333\Delta N_{\rm eff}\\ &-0.018p_\gamma-0.015p_n-0.012p_p. \end{aligned}$$

With $\boldsymbol d=(Y_p^{\rm obs},D_5^{\rm obs})$ and, as a declared cross-study approximation,

$$C_{\rm obs}=\operatorname{diag}(0.0013^2,0.029^2),$$

the normalized likelihood is

$$\ln\mathcal L=-\frac12\left[ (\boldsymbol d-\boldsymbol m)^{\mathsf T}C^{-1} (\boldsymbol d-\boldsymbol m) +\ln\det(2\pi C)\right], \qquad C=C_{\rm obs}+C_{\rm res}.$$

For the reproducible page-level fit set $C_{\rm res}=0$ and take $\Omega_bh^2\sim\mathcal N(0.02236,0.00015^2)$, $\tau_n/{\rm s}\sim\mathcal N(878.4,0.5^2)$, $p_j\sim\mathcal N(0,1)$ independently, and $\Delta N_{\rm eff}\sim\operatorname{Uniform}(-1,1)$. The last prior is a page-level domain restriction for the local expansion, not an observational measurement; it restricts $N_{\rm eff}$ to $2.044$–$4.044$. A full result must replace the three-rate truncation by the versioned network and test the linear approximation. Observational calibration and sample uncertainty belongs in $C_{\rm obs}$ or explicit observational nuisance parameters; correlated nuclear-rate uncertainty belongs in a correlated prior for the $p_j$ or the full network. Any remaining $C_{\rm res}$ must exclude directions already represented by $\Omega_bh^2$, $\tau_n$, $\Delta N_{\rm eff}$, and $p_j$ to avoid double counting. The fit must not reuse a published CMB+BAO+BBN posterior as a prior on the same BBN data.

Reproducing the local response and baseline comparison

Burns’s sensitivity atlas gives the following derivatives. The first four are logarithmic; the rate entries are absolute changes in $D_5$ per quoted one-standard-deviation lognormal rate shift.

Variation	PRIMAT response	NACRE-II response
$\partial\ln Y_p/\partial\ln\tau_n$	$+0.731$	$+0.731$
$\partial\ln D_5/\partial\ln\tau_n$	$+0.420$	$+0.420$
$\partial\ln Y_p/\partial\ln(\Omega_bh^2)$	$+0.039$	$+0.039$
$\partial\ln D_5/\partial\ln(\Omega_bh^2)$	$-1.639$	$-1.638$
$\partial D_5/\partial p_\gamma$	$-0.018$	$-0.020$
$\partial D_5/\partial p_n$	$-0.015$	$-0.082$
$\partial D_5/\partial p_p$	$-0.012$	$-0.055$

These values are from Burns 2026, Table 3, printed p. 82, and Table 6, printed p. 85 (arXiv:2603.22414v2). The additive $\Delta N_{\rm eff}$ coefficients are from Burns 2026, Table 5, printed p. 84. They immediately reproduce useful checks. Starting from the PRIMAT baseline, a 1% increase in $\Omega_bh^2$ gives $Y_p=0.247006$ and $D_5=2.40513$ to linear order; a 1% increase in $\tau_n$ gives $Y_p=0.248706$ and $D_5=2.45522$. A shift $\Delta N_{\rm eff}=0.1$ gives $Y_p=0.24831$ and $D_5=2.4783$. Individual $+1\sigma$ PRIMAT shifts in the three displayed destruction rates give $D_5=2.427$, $2.430$, and $2.433$. The comparatively low linear-fit quality $R^2=0.9582$ for the baryon response of D/H warns against extrapolating that derivative far from the fiducial point.

At the same weak-rate normalization, the published baselines are

Rate compilation	$Y_p$	$D_5$
PRIMAT	$0.24691\pm0.00012$	$2.445\pm0.038$
NACRE-II	$0.24683\pm0.00015$	$2.517\pm0.104$

Adding the quoted observational and theory standard deviations in quadrature puts the LBT helium value $0.85\sigma$ below the PRIMAT baseline and the PDG deuterium value $1.32\sigma$ above it. The corresponding NACRE-II offsets are $0.79\sigma$ low and $0.08\sigma$ low. These pulls are checks, not a model ranking: NACRE-II carries a much wider deuterium theory error, the two rows are alternative rate compilations rather than independent data, and a normalized comparison must retain covariance and the determinant term. The baselines and linear-propagation uncertainties are tabulated in Burns 2026, Table 7, printed p. 86; the atlas omits off-diagonal covariances in that budget.

Interpretation and limits

Standard BBN is successful because one baryon density and one expansion history describe helium and deuterium observations that arise from very different astrophysical systems. That agreement is conditional, not assumption-free. Nuclear cross sections, the neutron lifetime, plasma and neutrino corrections, stellar absorption and processing, chemical evolution, ionization models, sample heterogeneity, and network implementation all enter before a statement about new physics can be made.

The clean diagnostic is therefore not “BBN measures $N_{\rm eff}$” in isolation. It is: under a named thermal history, reaction network, rate release, abundance selection, covariance, and prior measure, the jointly predicted abundance vector can constrain changes in the expansion rate or particle content. A new neutron-lifetime average, nuclear-rate evaluation, helium emissivity calculation, deuterium absorber compilation, or BBN-code release requires this comparison to be rerun.

Exercises

Suppose all surviving neutrons end in helium-4. Infer $r_{np}$ from the LBT value $Y_p=0.2458$. Then use the PRIMAT response to estimate the helium abundance after a 1% increase in $\tau_n$, holding every other input fixed.

Solution

Inverting $Y_p=2r/(1+r)$ gives

$$r_{np}=\frac{Y_p}{2-Y_p} =\frac{0.2458}{1.7542} =0.14012.$$

For the response check, start from the network baseline rather than the observational central value. With $\delta_\tau=\ln(1.01)$,

$$Y_p\simeq0.24691\left[1+0.731\ln(1.01)\right] =0.248706.$$

Using the equivalent power-law form $0.24691(1.01)^{0.731}$ gives $0.248712$; their small difference is beyond the local linear order retained here.

References

• Aver, E. et al., “The LBT $Y_p$ Project IV: A New Value of the Primordial Helium Abundance,” arXiv:2601.22238v1 (2026), arXiv.
• Burns, A.-K., “Inside the Black Box of Big Bang Nucleosynthesis: A Comprehensive Sensitivity Atlas for the Precision Era,” arXiv:2603.22414v2 (2026), arXiv.
• Cooke, R. J., M. Pettini, and C. C. Steidel, “One Percent Determination of the Primordial Deuterium Abundance,” Astrophysical Journal 855, 102 (2018), DOI.
• Fields, B. D., P. Molaro, and S. Sarkar, “Big Bang Nucleosynthesis,” in F. Takahashi et al. (Particle Data Group), Review of Particle Physics, International Journal of Modern Physics A 41, 2630011 (2026), revised August 2025, Official PDF.
• Gariazzo, S., P. F. de Salas, O. Pisanti, and R. Consiglio, “PArthENoPE Revolutions,” Computer Physics Communications 271, 108205 (2022), DOI.
• Kislitsyn, P. A. et al., “A New Precise Determination of the Primordial Abundance of Deuterium: Measurement in the Metal-Poor sub-DLA System at $z=3.42$ towards Quasar J1332+0052,” Monthly Notices of the Royal Astronomical Society 528, 4068–4081 (2024), DOI.
• Pitrou, C., A. Coc, J.-P. Uzan, and E. Vangioni, “Precision Big Bang Nucleosynthesis with Improved Helium-4 Predictions,” Physics Reports 754, 1–66 (2018), arXiv:1801.08023v3 revised 2025, DOI, Open PDF.
• Skillman, E. D. et al., “The LBT $Y_p$ Project I: An Improved Determination of the Primordial Helium Abundance—Project Description, Sample Selection, Observations, and Methodology,” arXiv:2601.22232v1 (2026), arXiv.