Interpreting and Stabilizing Machine-learning Parametrizations of Convection

1 Introduction

• Sub-grid-scale parameterization is necessary because global climate models cannot both resolve convective-scale motions and run decadal simulations, and it is a major source of climate-change uncertainty.Neural networks can predict unbounded heating and moistening outside the training-data envelope, making coupled instability catastrophic.Source: §1 · pp. 1–2
• Neural networks often become numerically unstable when coupled to atmospheric fluid dynamics, and the causes of these instabilities are not fully understood.Proposed contributors include ambiguous input relationships and network architecture or hyperparameter choices.Source: §1 · p. 2
• The study adapts interpretability methods to convection parameterizations and asks whether their stability can be predicted before coupling them to global models.The work focuses on tailoring partial-dependence plots and saliency maps to this task.Source: §1 · p. 2
• The analysis compares sensitivities of GCRM- and SPCAM-trained neural networks using lower-tropospheric stability, mid-tropospheric moisture, and coupled gravity-wave diagnostics.These environmental controls are physically relevant because moisture affects convection intensity and stability helps determine its depth.Source: §1 · pp. 2–3

2.1 Global cloud-resolving Model

• The GCRM-trained network uses 80 days of three-hourly SAM fields at 4-km resolution, coarse-grained to 160-km grid boxes.Its heating and moistening targets are budget residuals, with coarse-model tendencies estimated using ten 120-second cSAM steps.Source: §2.1 · p. 3

2.2 Super-parametrized Climate Model

• SPCAM embeds eight SAM columns in each CAM grid column, and its neural networks emulate the embedded models’ redistribution of temperature and water vapor.Training uses higher-frequency 30-minute data and a longer, two-year dataset, while the embedded-model scale separation and idealizations constrain aspects of emergent dynamics.Source: §2.2 · pp. 3–4

2.3 Neural network Parameterization

• Both parameterizations map column thermodynamic profiles and auxiliary inputs to heating and moistening, using feed-forward neural networks trained by minimizing a weighted mean-squared-error cost.The GCRM configuration initially includes all-level humidity and temperature inputs and crashes after six days; SPCAM’s NN-stable runs for years, whereas NN-unstable simulations crash within two weeks.Source: §2.3 · p. 4

3.1 Variation of Inputs

• The input-variation method bins GCRM and SPCAM data by lower-tropospheric stability and moisture, then evaluates network sensitivities to systematic changes in the averaged inputs.Because the networks are nonlinear, this differs from averaging their outputs; bin-averaged observed heating and moistening provide a comparison for prediction fidelity.Source: §3.1 · p. 5

3.2 Linear Response Functions

• Linear response functions provide local linear sensitivities around selected base states, but the GCRM analysis uses a tropical-mean state that is not radiative-convective equilibrium.Consequently, some positive modes may represent decay toward a true equilibrium rather than instability about that equilibrium.Source: §3.2 · p. 5

3.3 Coupling to Two-dimensional Linear Dynamics

• The wave-coupling framework combines neural-network response functions with simplified two-dimensional gravity-wave dynamics and computes mode growth rates, phase speeds, and vertical structures.It assumes zero mean winds, hydrostatic anelastic flow, rigid-lid boundaries, and neglects rotation and water-vapor effects on buoyancy.Source: §3.3 · p. 6

4 Regularization

• Removing upper-atmospheric temperature and humidity inputs consistently stabilizes GCRM networks, although the authors describe this ablation as ad hoc.The proposed physical rationale is that upper-atmospheric variables can introduce spurious causality into the parameterization.Source: §4 · p. 7
• The same upper-atmospheric input removal does not stabilize SPCAM networks, so the authors develop input regularization for them.They speculate that SPCAM instability is linked to imperfect network fit, exacerbated by limited hyperparameter tuning.Source: §4 · p. 7
• Direct automatic-differentiation LRFs are noisy and nonlinear, producing unphysical unstable modes above 300 m s−1 even for the “NN-stable” network.The authors note that a single-basic-state LRF therefore does not provide a clean interpretation of the network.Source: §4 · p. 7
• Averaging predictions across randomly perturbed inputs offers a cleaner view of the SPCAM network’s behavior.The authors motivate this ensemble approach because the “NN-stable” network performs stably when interactively coupled to GCM dynamics.Source: §4 · p. 7
• The five-step SPCAM input-regularization workflow perturbs the instability-linked base input with Gaussian noise, averages ensemble LRFs, then couples them to wave dynamics.The base state is traced to the instability, and the earliest timestep with a maximum in convective moistening Q2 is selected.Source: §4 · p. 7
• The authors do not recommend large input spread as a general SPCAM stabilization strategy because it moves inputs outside the training set, where large biases are expected.They use the approach to test whether their interpretability framework is robust against online prognostic failure modes.Source: §4 · p. 8

5.1 The Onset of Deep Convection in ML parameterizations

• Both training datasets and their networks predict sharply increasing net precipitation with mid-tropospheric humidity; stability matters most at intermediate moisture around 15–20 mm.The networks’ predictions are close to bin-averaged net precipitation, with differences around 10 mm/day.Source: §5.1 · p. 8
• As LTS decreases, the GCRM-trained network shifts from shallow to deep heating and moistening, while the SPCAM-trained network does not show equally clear deepening.In the GCRM case, LTS appears to control convective height more than overall moistening; the SPCAM result may depend on the selected moisture bin or model differences.Source: §5.1 · p. 8

5.2 Stabilizing via Ablation of Inputs

• Ablating both upper-atmospheric humidity and temperature removed many unstable propagating modes; ablating either alone left numerous unstable modes, including one near 100 m/s.The remaining modes after both inputs were removed were stationary or very slow; the authors suggest both ablations are needed for stable nonlinear simulations with the 3-hourly GCRM-trained network.Source: §5.2 · p. 9
• The fully ablated case retains a 628 km, zero-phase-speed standing mode with an approximately 3-day e-folding time, primarily coupling vertical motion to moisture.For a network predicting convection without radiation, the standing mode’s interpretation is uncertain because of noisy GCRM budget residuals and training.Source: §5.2 · pp. 9, 14
• An unablated network supports a 6283 km propagating mode at 44 m/s with a fast growth rate of 1.29 d−1 and a structure unlike known convective organization.The authors suggest this spurious mode could explain the network’s instability when coupled to a GCM.Source: §5.2 · p. 14

5.3 Stabilizing Gravity Wave Modes via Regularization of Inputs

• With 1% regularization, NN-unstable had many fast-growing modes near 1 day−1 at 5–20 m/s, while NN-stable had only a few slowly growing modes.The largest propagating growth rates fell to about 10−3 day−1 at 25% regularization.Source: §5.3 · p. 17
• At 20% input spread, longer-term CAM tests became possible; at 15% or less, none ran beyond 21 days, consistent with the wave-coupler predictions.The offline analysis indicated especially prominent shutdown of unstable modes near 16% regularization.Source: §5.3 · p. 17

6 Conclusions

• Both networks predict stronger heating and moistening as lower-tropospheric moisture increases, while lower-tropospheric stability controls convection depth and net precipitation rises in moist, unstable columns.The authors note that moisture and stability covary strongly, so the demonstrated sensitivities are not entirely independent.Source: §6 · p. 19
• Unstable propagating gravity-wave modes are associated with online numerical instability in both training settings, while stabilized schemes have stable propagating modes and more stable prognostic tests.The stabilization strategies were input ablation for SAM and input regularization for SPCAM.Source: §6 · p. 19
• Input regularization smooths SPCAM linear response functions and eases interpretation, but it is only incrementally helpful for prognostic stability, not a solution to instability.It also offers an alternative to ensemble-averaging Jacobians without training multiple networks.Source: §6 · p. 19
• The wave-coupling spectra are evaluated about base states near instability rather than radiative-convective equilibrium, and estimates near zero phase speed or growth rate are hard to assess robustly.Inferring true physical modes will require a more rigorous statistical framework.Source: §6 · p. 19
• The two techniques provide ways to inspect machine-learning parameterizations and assess their physical plausibility, with the potential to support more accurate data-driven parameterizations.The authors present this as a potential, not a demonstrated improvement in accuracy.Source: §6 · p. 19

• The equations use humidity q, dry static energy s, and horizontal and vertical velocities as perturbations from a large-scale state.The anelastic geopotential term is φ = p′/ρ0, with a reference density profile specified for the full nonlinear model.Source: §A.1 · p. 20
• Hydrostatic balance and mass conservation complete the prognostic equations, with buoyancy defined as B = gT/¯T.Source: §A.1 · p. 20
• Combining diagnostic relations with zonal momentum yields a prognostic equation for vertical velocity w, using L = ∂z and H = 1/ρ0 ∂zρ0.Source: §A.1 · p. 20
• The elliptic vertical operator A requires two boundary conditions, taken as a rigid lid and impenetrable surface with w(0) = w(HT) = 0.Source: §A.1 · p. 20

• With w, s, and q vertically collocated, the interior elliptic operator A is discretized as a tridiagonal matrix.Source: §A.2 · p. 20
• The index k runs from 1 to N vertical grid cells, and z denotes height.Source: §A.2 · p. 20
• Rigid-lid conditions use w0 = −w1 and wn+1 = −wn for cell-centered vertical velocity, modifying the boundary rows of A.The boundary rows are (Aw)1 = −a1w1 + b1w1 + c1w2 and (Aw)n = anwn−1 + bnwn − cnwn.Source: §A.2 · p. 21