Dietary restriction impacts health and lifespan of genetically diverse mice

Mice

We enrolled 960 female DO mice in 12 waves, corresponding to birth cohorts from generations 22 to 24 and 26 to 28 with 80 first-parity and 80 second-parity animals from each generation born around 3 weeks apart. No more than one mouse per litter was enrolled in the study. The first cohort entered the study in March 2016 and the study was fully populated in November 2017. This schedule was designed to make efficient use of our phenotyping capacity and minimize the potential for seasonal confounding. The sample size was determined to detect a 10% change in mean lifespan between intervention groups with allowance for some loss of animals due to non-age-related events. We used female mice due to concerns about male aggression. Mice were assigned to housing groups of eight animals in large-format wean boxes with positive pressure ventilation and incoming air temperature of 24.4 to 25.6 °C. Environmental enrichments were provided including nestlets, biotubes and gnawing blocks. Mice were randomized by housing group to one of five dietary interventions. Blinding was not possible due to the different feeding requirements in each intervention group. Of the 960 mice that were entered into the study, 937 mice were alive when interventions were initiated at 6 months of age and only these mice are included in our analysis. All of the procedures used in the study were reviewed and approved under Jackson Lab IACUC protocol 06005.

DR

DR was implemented by controlling the timing and amount of food provided to mice. Feeding schedules for DR were started at 6 months of age. All mice were fed a standard mouse chow diet (5K0G, LabDiet). The AL feeding group was provided with unlimited access to food and water. The IF mice were provided unlimited access to food and water. On Wednesday of each week at 15:00, IF mice were placed into clean cages and food was withheld for the next 24 or 48 h for the 1D and 2D groups, respectively. CR mice were provided with unlimited access to water and measured amounts of food daily at around 15:00, 2.75 g per mouse per day for 20% CR and 2.06 g per mouse per day for 40% CR. These amounts were based on AL consumption of 3.43 g per mouse per day that we estimated based on historical feeding data from DO mice. For the 40% CR protocol, a gradual reduction in food intake was implemented: the mice were first subjected to 20% CR for 2 weeks, then to 30% CR for an additional 2 weeks, before transitioning to the full 40% CR. In the 2D IF protocol, mice were initially acclimatized to the 1D IF regimen for 2 weeks. Mice were co-housed with up to eight mice per pen. Co-housing is standard practice for CR studies³⁹; competition for food was minimized by placing food directly into the bottom of the cage, allowing individual mice to ‘grab’ a pellet and isolate while they eat. CR mice were provided with a 3 day ration of food on Friday afternoon, which resulted in weekly periods of feasting followed by a period of food deprivation of approximately 1 day for the 20% CR mice and 2 days for the 40% CR mice, comparable to the IF fasting periods. The 15:00 feeding time closely approximates the circadian alignment of feeding, starting just before the beginning of the dark cycle, which is the normal active and feeding time of day for mice. This timing has been shown to maximize lifespan extension in mice subjected to 30% CR⁸. The feast–famine cycle induced by the Friday triple feeding has been used in other studies of CR^57,58,59 but direct assessment of the health impacts is lacking.

Food intake (160 mouse independent cohort)

To obtain an accurate estimate of food intake and changes in body weight in response to weekly fasting cycles, we set up an independent cohort of 160 female DO mice. Mice were placed on the same DR protocols as in the main study. Food was weighed daily for a period of 1 week when mice were 30, 36 and 43 weeks of age. Food intake data were normalized to units of g per mouse per day and presented as daily and weekly averages across timepoints by diet. Body composition was determined at 43 and 45 weeks of age using non-imaging nuclear magnetic resonance (NMR) using the Echo MRI instrument, a NMR device with a 5-gauss magnet that is adapted to small-animal studies. NMR data were used to detect changes in body weight and composition before and after fasting. Values of body weight, lean mass, fat mass and adiposity (100% × fat mass/total mass) pre-fasting were co-plotted with the difference between before and after fasting (Friday to Monday for AL and CR; Tuesday to Thursday for 1D IF; Tuesday to Friday for 2D IF).

Phenotyping

We performed three cycles of health assessments of mice at early, middle and late life. These assessments included a 7-day metabolic cage run at around 16, 62 and 114 weeks of age; blood collection for flow cytometry analysis at 24, 71 and 122 weeks; rotarod, body composition, echocardiogram, acoustic startle, bladder function, free wheel running and a blood collection for CBC analysis at around 44, 96 and 144 weeks of age. Furthermore, body weights were recorded weekly and manual frailty and grip strength assessments were carried out at 6 month intervals. All assays were conducted at the Jackson Laboratory according to standard operating procedures.

Body weight

Mice were weighed weekly throughout their lives, resulting in over 100,000 values longitudinally collected for the 937 mice. Body weights were analysed after local polynomial regression fitting within mouse (that is, loess smoothing).

Frailty, grip strength and body temperature

We applied a modified version of the clinically relevant FI³⁸, which was calculated as the average of 31 traits that are indicators of age-associated deficits and health deterioration. Each trait was scored on a scale of 0, 0.5 or 1, where 0 indicated the absence of the deficit; 0.5 indicated mild deficit; and 1 indicated severe deficit. Measurements were taken at the baseline (5 months) and were repeated approximately every 6 months. Simple averaging yielded a raw FI score of between 0 and 1 for each mouse. Frailty scores were adjusted by estimating batch, coat colour and experimenter effects as random factors that were subtracted from raw frailty score values before statistical analysis.

Body composition

We performed dual X-ray absorptiometry analysis using the LUNAR PIXImus II densitometer to collect bone density and body composition (including fat and non-fat lean tissue). Mice were anaesthetized and individually placed onto a disposable plastic tray that was then placed onto the exposure platform of the PIXImus. The process to acquire a single scan lasts approximately 4 min. Measurements were taken at around 44, 96 and 144 weeks of age.

Immune cell profiling using flow cytometry

Peripheral blood samples were analysed by flow cytometry to determine the frequency of major circulating immune cell subsets. Analysis was performed before the start of dietary interventions at 5 months, then at 16 and 24 months of age. These timepoints corresponded to 11 and 19 months of dietary intervention. Red blood cells in PBL samples were lysed and the samples were washed in FACS buffer (Mitenyi, 130-091-222). Cells were resuspended in 25 μl FACS buffer with 0.5% BSA (Miltenyi, 130-091-222 with 130-091-376). Antibodies including Fc block (2.42, Leinco Technologies) were added and incubated for 30 min at 4 °C. Labelled cells were washed and DAPI was added before analysis on the LSRII (BD Bioscience) system. The antibody cocktail contained CD11c FITC, N418 (35-0114-U100, Tonbo Biosciences, 1:100); NKG2D (CD314) PE, CX5 (558403, BD Biosciences, 1:80); CD3e PE-CF594, 145-2C11 (562286, BD Biosciences, 1:40); CD19 BB700, 1D3 (566411, BD Biosciences, 1:40); CD62L PE-Cy7, MEL-14, (60-0621-U100, Tonbo Biosciences, 1:100); CD25 APC, PC61 (102012, BioLegend, 1:80); CD44 APC-Cy7, IM7 (25-0441-U100, Tonbo Biosciences, 1:40); Ly6G BV421, 1A8, (562737, BD Biosciences, 1:80); CD4 BV570, RM4-5 (100542, BioLegend, 1:40); CD11b BV650, M1/70 (563402, BD Biosciences, 1:160); CD45R/B220 BUV496 (RA3-6B2, 564662, BD Biosciences, 1:20); Fc Block, 2.4G2 (C247, Leinco Technologies, 1:100).

Owing to the outbred nature of these mice, flow cytometry markers were limited, and T cell subsets were generally assigned as naive and non-naive by the presence of CD62L and CD44 (immune cell subtype designations are shown in Supplementary Table 7). NKG2D-positive cells were enumerated and may represent memory T cells that accumulate after immune responses⁶⁰. Owing to limitations in flow cytometry markers that identify NK cells and their subsets in the mouse strains contributing to the outbred DO mouse line, NK cells were defined as non-T non-B lymphocytes expressing NKG2D. Within this population, CD11c and CD11b were used to generally define maturation subsets. CD11b expression marks more mature NK cells and CD11c is reduced on the least-mature NK subset⁶¹.

Glucose

At the flow cytometry blood collections at 16, 62 and 114 weeks, mice were fasted for 4 h and glucose was measured using the OneTouch Ultra2 glucose meter from LifeScan along with OneTouch Ultra test strips. At each of the CBC blood collections at 24, 71 and 122 weeks, non-fasted glucose was measured using the glucose meter.

Echocardiogram

Ultrasonography was performed using the VisualSonics (VSI) Vevo 770/2100 high-frequency ultrasound system with 30 and 40 MHz probes. Echocardiography uses pulsed Doppler sonography, applied through the ultrasound probe, to measure blood flow rates and volumes.

Metabolic monitoring cages

Mice were individually housed for 7 days in a metabolic cage (Promethion Model from Sable Systems International) and the activity, feeding and respiration were tracked. Feeding protocols for dietary intervention were maintained. Metabolic cage data were used to assess metabolism, energy expenditure and activity of mice in Y1, Y2 and Y3. Animal-level data were cleaned to remove outliers and instrument failure and summarized as cumulative (food, wheel distance) or median across 5 min intervals (respiratory quotient, energy expenditure). The mean and s.d. were computed at 4 h and 1 h intervals and were plotted as mean ± 2 s.e.m. across timepoint intervals. Animal-level summaries were computed as the average across 7 days of the daily (24 h), light phase (12 h) or dark phase (12 h) median values. Moreover, we computed ‘change’ traits (for example, delta respiratory quotient and delta energy expenditure) as the difference between the 5th to 95th percentiles of all 1 h summaries across the entire 7 day run.

CBC analysis

Blood samples were run on the Siemens ADVIA 2120 haematology analyzer with mouse-specific software as described previously⁶².

Acoustic startle

Startle response was measured in rodents using automated startle chambers, in which a mouse was placed in a clear, acrylic tube attached to a highly sensitive platform that is calibrated to track their startle reflex while being exposed to a series of stimuli at varying decibels and times. Mice were initially exposed to white noise from an overhead speaker, which transitions to a series of randomized, computer-generated stimuli ranging in volume from 70 to 120 decibels at 40 ms in duration and an interval of 9–22 s. The test runs for approximately 30 min.

Rotarod

We used the Ugo–Basile rotarod, which has five lanes evenly spaced along a motorized horizontal rotating rod, allowing for up to five mice to be tested simultaneously. Below each lane is a platform equipped with a trip plate that records the latency for each mouse to fall. At the beginning of the session, mice were placed onto the rod, which began rotating at 4 rpm, slowly increasing to a maximum of 40 rpm, over 300 s. Mice were given three consecutive trials. We reported the mean latency (time to fall) and the slope of latencies across trials, as well as the number of trials with no falls and number of trials with immediate falls. In case a mouse did not cooperate with the test, trials were recorded as missing.

Voiding assay

Cages were prepared by cutting a piece of cosmos blotting paper, 360 gsm, to standard duplex cage dimensions. Shavings were removed from a clean cage, and the paper was taped to the bottom of the cage. Food was provided during this test; however, water was removed to prevent possible leaking onto the blotting paper. Mice were individually housed in a prepared cage for 4 h. At the end of the trial, the mice were returned to their original housing units, and papers were removed and dried for 2–4 h, before being individually bagged. Papers were shipped to Beth Israel Deaconess Medical Center, where they were scanned with ultraviolet light to image and quantify the void spots.

Home cage wheel running

Free-wheel-running data were collected at around 44, 96 and 144 weeks of age. Mice were individually housed for a minimum of 36 h in a special cage suited to house the Med Associate low profile running wheel with a wireless transmitter. The food hopper was removed to allow for seamless movement of the wheel, and food was placed onto the cage floor. The 15.5-cm-diameter plastic wheel sits at an angle on an electronic base, which tracks the revolutions. The battery-powered base allows for continuous monitoring of data, which is then wirelessly transmitted, in 30 s intervals, to a local computer.

Lifespan

Research staff regularly evaluated mice for prespecified clinical symptomology: palpable hypothermia, responsiveness to stimuli, ability to eat or drink, dermatitis, tumours, abdominal distention, mobility, eye conditions (such as corneal ulcers), malocclusion, trauma and wounds of aggression. If mice met the criteria for observation in any of these categories, veterinary staff were contacted. If the clinical team determined a mouse to be palpably hypothermic and unresponsive, unable to eat or drink, and/or met protocol criteria for severe dermatitis, tumours and/or fight wounds, pre-emptive euthanasia was performed to prevent suffering; otherwise, the veterinary staff provided treatment. Both mice euthanized or found dead were represented as deaths in the survival curves. Mice euthanized due to injuries unrelated to imminent death were treated as censored (we recorded a total of 13 censoring events).

Data preparation and analysis

Cloud-based research management software (Climb by Rockstep) was used to track animals, schedule testing and provide a stable repository for primary data collection. Data were regularly reviewed by a statistical analyst during the study for anomalies. Initial data quality control included identifying and resolving equipment miscalibration, mislabelled animals and technically impossible values. If we could not manually correct these using laboratory records, they were removed. Quantitative assays including body weight and temperature were explored for outliers. Quantitative traits other than body weights were corrected for batch effects. To quantify batch effects, we fit a fully random-effects linear mixed model conditioning on diet, body weight at test and age. We adjusted trait values by subtracting the batch model coefficients. Lifespan data were recorded in days but are presented in months (30.4 days per month) for ease of interpretation. Statistical significance for extremely small P values is reported as P < 2.2 × 10^−16 in the main text; non-truncated P values are provided in the Supplementary Information. All analyses were performed using R v.4.2.2 and RStudio v.2022.12.0+353. Data analysis scripts are available (Data availability).

Survival analysis

We performed survival analysis to compare lifespan outcomes for the five study groups. We plotted Kaplan–Meier survival curves and tested the equality of survival distributions across diet groups using log rank tests using an overall test (4 d.f.) and all pairwise comparisons between diets. P values are reported with no multiple testing adjustment, and we considered a comparison to be significant if P < 0.01. We estimated the median and maximum lifespan (90% survival) by diet group with 95% confidence intervals and percentage change relative to the AL group³⁰. Mortality doubling times were estimated from a Gompertz log-linear hazard model with 95% CI and percentage change relative to the AL group (flexsurv R package v.2.2.2).

Longitudinal trait analysis

For traits collected annually or biannually, we were able to explore hypothesized direct and indirect relationships with diet, body weight and age. We used generalized additive mixed models (GAMMs) with a gaussian/identity link to analyse these effects by fitting a series of nonlinear relationships between trait response and covariates. GAMMs for a combination of fixed and random effects (formulae below) on trait response pre-adjusted for batch were fit using the gam() function of the mgcv package in R, with the Newton optimizer and default control parameters. Age was rescaled to proportion of life lived (PLL = age at test/lifespan). The PLL scale removes artifacts due to survivorship bias across groups with different lifespans. All continuous variables except for PLL were rank normal scores transformed (RZ) prior to model fitting.

1. (1)

   RZ(T) ~ D + s(BW) + s(PLL) + (1|ID)

2. (2)

   RZ(T) ~ D + s(PLL) + (1|ID)

3. (3)

   RZ(T) ~ D + s(BW) + (1|ID)

4. (4)

   RZ(T) ~ s(BW) + s(PLL) + (1|ID)

5. (5)

   RZ(T) ~ D + s(BW) + s(PLL|diet) + (1|ID)

where T is trait, D is dietary assignment, BW is body weight at the date closest to T’s collection date, PLL is the proportion of life lived as of T collection date and s() is the smoothing parameter. Each mouse had multiple datapoints across the T collection date. This clustering was accounted for with a random intercept for ID, specified as (1|ID) above. We performed hypothesis tests related to the GAMM fits to explore trait sensitivity to body weight (model 1 (M1) versus M2), PLL, (M1 versus M3), diet (M1 versus M4) and diet-by-trait interaction (M1 versus M5). Using the models specified above and a conservative false-discovery rate (FDR < 0.01, one-step Benjamini–Hochberg method), we therefore identified traits that responded additively to body weight, traits that responded additively to diet, traits that responded additively to PLL (scaled age) and traits that responded interactively to diet and PLL. Traits were categorized as health, metabolism, haematology or immune. For each trait category, bar plots were generated to show the number of traits with significant (FDR < 0.01) associations with body weight, diet, PLL and diet × PLL interactions. We repeated the same analysis with age (age in months) in place of PLL. We applied FDR adjustment to each test across traits and timepoints (Benjamini–Hochberg FDR method).

The complete results of the longitudinal trait analysis are provided in Supplementary Table 6.

Trait associations with lifespan

To identify traits that are associated with lifespan, we performed regression analysis on lifespan with traits at each timepoint after adjusting for effects of diet and body weight. We fit linear models:

1. (1)

   Lifespan ~ diet + BW6 + BW_test,

2. (2)

   Lifespan ~ diet + BW6 + BW_test + trait

3. (3)

   Lifespan ~ diet + BW6 + BW_test + trait + diet:trait

where BW6 is the last preintervention body weight, and BW_test is the body weight at time of testing. For body composition and change-in-body weight traits, we did not include the body weight terms. All continuous variables were rank normal scores (RZ) transformed before model fitting. We performed likelihood ratio tests for the diet and body weight adjusted association (M2 versus M1) and for the diet × trait interaction (M3 versus M2). We applied a FDR adjustment to each test across traits and timepoints (one-step Benjamini–Hochberg method). Traits were categorized and significant (FDR < 0.01) results were tabulated as above.

We estimated the overall strength of association between trait and lifespan as the regression coefficient of the trait term in M2. As the traits are standardized, the regression coefficients are equivalent to diet- and body-weight-adjusted partial correlations (designated as r) and can be compared across traits. Moreover, we stratified the data to estimate diet-specific body-weight-adjusted partial correlations that are shown in the scatter plot panels (Fig. 2f).

The complete results of lifespan association analysis are provided in Supplementary Table 7.

Network modelling

We performed a multivariate network analysis to decompose the effects of DR on lifespan across measured traits. An empirical covariance matrix was estimated using the nonparanormal SKEPTIC estimator derived from the pairwise Kendall–Tau correlation using pairwise-complete data. The estimated covariance matrix was projected to the nearest positive definite matrix by truncating eigenvectors with negative eigenvalues. A sparse low-rank Gaussian graphical model was fitted to this covariance estimate using the ggLASSO Python package with the parameters lambda = 0.1, mu = 100. The result was normalized to obtain an inferred partial correlation matrix, which we use as our phenotype network for downstream analysis.

To cluster the phenotypes, we constructed a weighted k-nearest neighbours graph using the absolute partial correlation from above as similarity weights (that is, we retained the k-largest weights per node) using k = 10. The resulting k-nearest neighbours graph was clustered using spectral clustering to obtain 20 clusters, which were labelled by hand. Diet and lifespan were broken out into their own univariate clusters and the body weight cluster was split into lean tissue mass and adiposity (fat mass/total mass) to account for body composition effects known from the literature. This resulted in a summarized representation of 23 clusters (diet, lifespan and 21 groups of traits, one per cluster plus FTM).

We obtained a covariance decomposition by recomputing a sparse low-rank network (as above) on the reduced representation of 23 features. To obtain the relative importance of different effects of DR on lifespan, we performed covariance path decomposition⁴⁸ of the covariance between diet and lifespan using this reduced representation graphical model. The absolute values of path scores were used to rank the relative importance of each path and normalized to sum to 1 to estimate the fraction of the covariance between DR and lifespan explained. To obtain network visualizations, the position of nodes and orientation of edges were determined by computing max-flow through the path network. We define the path network as the graph formed by taking the partial correlation network but reweighting edges to be the sum of all (absolute) path scores of paths that contain that edge.

QTL mapping

Genetic mapping analysis of all continuous-valued traits including lifespan was performed using the R package qtl2⁶³. Additive covariates for the genome-scan models were diet and body weight variables as indicated in Supplementary Table 5. Founder haplotype effects were estimated treating 8-state genotypes as random effects. Trait heritability was estimated using the R/qtl2 function est.herit() to fit a linear mixed model including an additive kinship matrix. Confidence intervals were obtained by parametric bootstrap. Significance thresholds for QTL mapping were estimated from 1,000 permutations of the trait and covariate data; significant QTL exceed the genome-wide multiple-test adjusted 0.05 threshold (7.5) and suggestive QTL exceed the unadjusted 0.05 threshold (6.0)⁶⁴.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.