**Background:** Climate-resilient grassland systems are critical for stable livestock fodder production, but increasing temperatures and drought frequency threaten yield stability. While short-term studies suggest that higher nutrient inputs increase both yield and interannual yield variance, long-term data (>20 years) on temporal trends in yield variance and their dependence on agronomic management and climatic drivers are lacking. The Park Grass Experiment (PGE) at Rothamsted, initiated in 1856, provides a unique opportunity to examine these relationships over decadal timescales.
**Methods:** The study analyzed total aboveground biomass (sum of 1st and 2nd cuts, 100% dry matter) from 1965 to 2018 across seven contrasting fertilizer treatments (Nil/unfertilized; PKNaMg; N1 PKNaMg with 48 kg N/ha as ammonium sulfate; N2 PKNaMg with 96 kg N/ha; N3 PKNaMg with 144 kg N/ha; farmyard/poultry manure in a 4-year cycle; N*1 with 48 kg N/ha as sodium nitrate), each divided into four liming sub-plots (target pH 7, 6, 5, or no chalk). A mixed model with REML was used, fitting a smoothing spline (k=10 knots, ~5-year intervals) for the mean yield trend and separate residual variances for each sub-period to assess changes in temporal yield variance independently of yield level. Treatment-specific multiple regression analyses tested environmental covariates (air temperature, rainfall, water stress days, N deposition, SO2 emissions, CO2). A novel criss-cross regression approach (extended Finlay–Wilkinson regression) modeled the environmental mean using three selected climatic drivers—accumulated water stress days (March–October), mean air temperature May–June, and mean air temperature July–August—to assess treatment-specific yield sensitivity (slope). Correlations between plant species diversity (Shannon's index, species number, proportion of legumes) and yield metrics were assessed using Kendall's tau. Partial canonical correspondence analyses (pCCA) examined plant community composition dynamics.
**Key Results:** (1) Temporal trends in yield variance showed a consistent unimodal pattern across plots with low to moderate or no N fertilizer: relatively low variance in the 1960s–1970s, a peak around 1995, and declining variance in the past two decades. This pattern was less pronounced at high soil pH and more pronounced at low pH. Exceptions included high N supply or farmyard manure with reduced/no liming, where high variability persisted throughout. (2) Accumulated water stress days (March–October) and mean July–August air temperature were the two main climatic drivers, explaining approximately one-third of yield variance across treatments, and up to 49% (PKNaMg) and 48% (N1 PKNaMg) in specific treatments. Water stress had a stronger impact than temperature (standardized Beta coefficients). (3) The criss-cross regression showed a highly significant treatment × environmental mean interaction (F=3.17; P<0.0001). Slopes ranged from 0.62 (N*1, no chalk) to 1.99 (PKNaMg, pH 6). Liming to pH 6–7 with moderate N supply (N2 PKNaMg, pH 6/7) gave slopes near 1 (0.95–0.99), indicating average yield sensitivity. The unlimed high-N treatment (N3 PKNaMg, no chalk) had a slope of 1.36, indicating high sensitivity to abiotic stress. (4) Mean yields declined over time in high-input treatments: N3 PKNaMg (pH 7/6) and FYM/PM (pH 7/6) decreased from >9 t/ha to nearly 5 t/ha from 1965 to 2018. (5) Higher plant species diversity was significantly negatively correlated with both yield variance (Shannon's index: τ=−0.30; species number: τ=−0.30; legume proportion: τ=−0.20; all P<0.05) and mean yield (species number: τ=−0.50; P<0.05). Treatments with high legume proportions (PKNaMg, N1 PKNaMg) achieved mean yields close to those with higher inorganic N while maintaining higher species richness. (6) pCCA showed that year explained 16.4% of variance in community composition (P<0.001), and treatment explained 64.2% (P<0.001), discriminating grass-dominated high-N unlimed plots from forb-rich unfertilized plots.
**Clinical Implications:** This study provides robust evidence that liming to maintain soil pH 6–7, combined with moderate nitrogen supply (e.g., 96 kg N/ha as ammonium sulfate with PKNaMg), is the most effective agronomic strategy for reducing yield sensitivity to climatic stress and promoting stable grassland productivity over decades. Higher plant species diversity, particularly including legumes, enhances yield stability but is associated with lower absolute yield levels. The novel criss-cross regression approach offers a parsimonious tool for assessing treatment-specific climate resilience that could be applied to other agro-ecological trials. These findings are directly relevant for managing temperate grasslands under climate change, where increasing frequency of combined heat and drought stress is projected to double the area of exposed grassland by 2050 compared to 2021.