Study Overview

The objective of this study was to determine whether energetic state constrains tail regeneration in adult female brown anoles (Anolis sagrei) and whether supplementation with insulin-like growth factor 1 (IGF1) or insulin-like growth factor 2 (IGF2) modifies this relationship.

  • Energetic state was quantified as percent body mass loss and analyzed as a continuous predictor of regenerative performance.
  • The primary response variable was cumulative percent tail regeneration measured weekly over eight weeks following tail autotomy.
  • Secondary analyses evaluated absolute regenerated tail length and weekly regeneration rate.
  • Linear mixed-effects models included animal identity as a random effect to account for repeated measurements.

The analysis proceeds in the same order as the manuscript:

  1. Verify successful randomization and dietary manipulation.
  2. Establish the relationship between energetic state and regeneration.
  3. Test whether IGF supplementation modifies that relationship.
  4. Evaluate robustness, alternative response variables, dose effects, and late-stage observations.
  5. Generate publication figures and supporting outputs.

1. Reproducibility Setup

1.1 Load packages

# Data manipulation
library(dplyr)
library(tidyr)

# Statistical modeling
library(lme4)
library(nlme)
library(car)
library(lmtest)
library(MuMIn)
library(bestNormalize)
library(multcomp)
library(emmeans)
library(performance)

# Model predictions and summaries
library(ggeffects)
library(Rmisc)

# Visualization
library(ggplot2)
library(patchwork)
library(plotly)
library(viridis)
library(hrbrthemes)

packageVersion("lmerTest")
## [1] '3.2.1'
search()
##  [1] ".GlobalEnv"            "package:hrbrthemes"    "package:viridis"      
##  [4] "package:viridisLite"   "package:plotly"        "package:patchwork"    
##  [7] "package:ggplot2"       "package:Rmisc"         "package:plyr"         
## [10] "package:lattice"       "package:ggeffects"     "package:performance"  
## [13] "package:emmeans"       "package:multcomp"      "package:TH.data"      
## [16] "package:MASS"          "package:survival"      "package:mvtnorm"      
## [19] "package:bestNormalize" "package:MuMIn"         "package:lmtest"       
## [22] "package:zoo"           "package:car"           "package:carData"      
## [25] "package:nlme"          "package:lme4"          "package:Matrix"       
## [28] "package:tidyr"         "package:dplyr"         "package:stats"        
## [31] "package:graphics"      "package:grDevices"     "package:utils"        
## [34] "package:datasets"      "package:methods"       "Autoloads"            
## [37] "package:base"
# Reproducible random-number generation for bootstrap analyses
set.seed(123)

1.2 Shared plotting objects and color palettes

# Position adjustments used in several figures
pd <- position_dodge(0.02)
dodge <- position_dodge(width = 0.6)

# Treatment colors used throughout the analysis
treatment_colors <- c(
  "AdLib"   = "grey40",
  "DR"      = "grey70",
  "DR.IGF1" = "#E87722",
  "DR.IGF2" = "#0C2340"
)

igf_colors <- c(
  "No_IGF" = "grey40",
  "IGF1"   = "#E87722",
  "IGF2"   = "#0C2340"
)

2. Data Import and Curation

2.1 Pre-experimental body-mass dataset

diet <- read.csv(
  "WL.analysis.csv",
  stringsAsFactors = FALSE,
  check.names = FALSE
)

# The source file contains a trailing space in the Week column name.
names(diet)[trimws(names(diet)) == "Week"] <- "Week"

# Standardize factor order for consistent output.
diet <- diet %>%
  mutate(
    Treatment = factor(
      Treatment,
      levels = c("AdLib", "DR", "DR.IGF1", "DR.IGF2")
    )
  )

str(diet)
## 'data.frame':    228 obs. of  14 variables:
##  $ Animal_ID         : chr  "BC039" "BC048" "BC099" "BC158" ...
##  $ ToeClip           : chr  "1455" "1513" "11-Jan" "2330" ...
##  $ BC_Cage           : int  32 28 9 28 8 36 39 39 28 15 ...
##  $ Treatment         : Factor w/ 4 levels "AdLib","DR","DR.IGF1",..: 3 1 1 4 1 3 2 2 3 4 ...
##  $ Diet              : chr  "DR" "Adlib" "Adlib" "DR" ...
##  $ Tail              : num  73 69 82 71 74 79.5 85 60 75.5 79 ...
##  $ SVL               : int  46 48 48 44 42 47 51 46 45 47 ...
##  $ Week              : int  -4 -4 -4 -4 -4 -4 -4 -4 -4 -4 ...
##  $ Mass              : num  2.57 2.71 3.67 2.64 2.46 ...
##  $ Bcond             : num  0.0558 0.0565 0.0766 0.06 0.0586 ...
##  $ Perc_BC.change    : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Delta Mass_Initial: num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Perc_WL           : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ Perc_Orig         : num  100 100 100 100 100 100 100 100 100 100 ...

2.2 Eight-week regeneration dataset

data <- read.csv(
  "R.analysis.currated.csv",
  stringsAsFactors = FALSE,
  check.names = FALSE
)

# Remove observations from four individuals following substantial secondary
# tail loss unrelated to the original experimental autotomy.

excluded_animals <- c("BC121", "BC144", "BC383", "BC561")

#Separate treatments out into single variables (Ex. DR.IGF2 is Restricted and IGF2 at Diet and IGF_Type Variables)
data <- data %>%
  filter(!Animal_ID %in% excluded_animals) %>%
  mutate(
    Treatment = factor(
      Treatment,
      levels = c("AdLib", "DR", "DR.IGF1", "DR.IGF2")
    ),
    Diet = case_when(
      Treatment == "AdLib" ~ "AdLib",
      Treatment %in% c("DR", "DR.IGF1", "DR.IGF2") ~ "Restricted",
      TRUE ~ NA_character_
    ),
    Diet = factor(Diet, levels = c("AdLib", "Restricted")),
    IGF = if_else(
      Treatment %in% c("DR.IGF1", "DR.IGF2"),
      "IGF",
      "No_IGF"
    ),
    IGF = factor(IGF, levels = c("No_IGF", "IGF")),
    IGF_type = case_when(
      Treatment %in% c("AdLib", "DR") ~ "No_IGF",
      Treatment == "DR.IGF1" ~ "IGF1",
      Treatment == "DR.IGF2" ~ "IGF2",
      TRUE ~ NA_character_
    ),
    IGF_type = factor(IGF_type, levels = c("No_IGF", "IGF1", "IGF2"))
  )

str(data)
## 'data.frame':    608 obs. of  18 variables:
##  $ Animal_ID        : chr  "BC039" "BC039" "BC039" "BC039" ...
##  $ ToeClip          : chr  "1455" "1455" "1455" "1455" ...
##  $ BC_Cage          : int  32 32 32 32 32 32 32 32 28 28 ...
##  $ Treatment        : Factor w/ 4 levels "AdLib","DR","DR.IGF1",..: 3 3 3 3 3 3 3 3 1 1 ...
##  $ Tail             : num  73 73 73 73 73 73 73 73 69 69 ...
##  $ SVL              : int  46 46 46 46 46 46 46 46 48 48 ...
##  $ Week             : int  1 2 3 4 5 6 7 8 1 2 ...
##  $ Mass             : num  2.1 2.07 2.1 2.1 2.24 ...
##  $ Dose             : num  2.38 2.38 2.38 2.38 2.38 ...
##  $ Perc.WeightLoss  : num  18 19.3 18.2 18.2 12.7 ...
##  $ Regeneration     : num  0 0.15 6 11 16 20 22 24 0 0.1 ...
##  $ Rate_Reg         : num  NA 0.15 5.85 5 5 4 2 2 NA 0.1 ...
##  $ Rate_Perc        : num  NA 0.217 8.478 7.246 7.246 ...
##  $ Perc.Regeneration: num  0 0.217 8.696 15.942 23.188 ...
##  $ Eggs             : int  0 1 0 0 0 0 0 0 0 0 ...
##  $ Diet             : Factor w/ 2 levels "AdLib","Restricted": 2 2 2 2 2 2 2 2 1 1 ...
##  $ IGF              : Factor w/ 2 levels "No_IGF","IGF": 2 2 2 2 2 2 2 2 1 1 ...
##  $ IGF_type         : Factor w/ 3 levels "No_IGF","IGF1",..: 2 2 2 2 2 2 2 2 1 1 ...
table(data$Treatment, useNA = "ifany")
## 
##   AdLib      DR DR.IGF1 DR.IGF2 
##     160     136     160     152
length(unique(data$Animal_ID))
## [1] 76

3. Verification of Randomization

Animals should not differ systematically in body size or original tail length before treatment implementation.

3.1 Baseline subset

#Week -4 is the beginning of the study, prior to diet implementation. Week 0 is the time of tail autotomy. 
pre_diet <- diet %>%
  filter(Week == -4)

3.2 Initial snout–vent length

model_svl <- lm(SVL ~ Treatment, data = pre_diet)
anova(model_svl)
## Analysis of Variance Table
## 
## Response: SVL
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treatment  3   6.984  2.3279  0.6542 0.5829
## Residuals 72 256.214  3.5585
summary(model_svl)
## 
## Call:
## lm(formula = SVL ~ Treatment, data = pre_diet)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -3.900 -1.200 -0.200  1.100  4.706 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       45.9000     0.4218 108.816   <2e-16 ***
## TreatmentDR        0.3941     0.6223   0.633    0.529    
## TreatmentDR.IGF1   0.3000     0.5965   0.503    0.617    
## TreatmentDR.IGF2   0.8368     0.6043   1.385    0.170    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.886 on 72 degrees of freedom
## Multiple R-squared:  0.02653,    Adjusted R-squared:  -0.01403 
## F-statistic: 0.6542 on 3 and 72 DF,  p-value: 0.5829
plot_svl <- ggplot(
  pre_diet,
  aes(x = Treatment, y = SVL, fill = Treatment)
) +
  geom_violin(trim = FALSE, scale = "width") +
  geom_boxplot(width = 0.10, outlier.shape = NA) +
  geom_jitter(width = 0.10) +
  scale_fill_manual(values = treatment_colors) +
  labs(x = "Treatment", y = "Snout–vent length") +
  theme_classic() +
  theme(legend.position = "none")

plot_svl

ggsave(
  filename = "svl.png",
  plot = plot_svl,
  width = 5,
  height = 4,
  dpi = 600
)

Manuscript result: Initial snout–vent length did not differ among treatment groups, \(F_{3,72}=0.65\), \(p=0.58\).

3.3 Initial body mass

model_mass <- lm(Mass ~ Treatment, data = pre_diet)
anova(model_mass)
## Analysis of Variance Table
## 
## Response: Mass
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treatment  3  0.3761 0.12538  0.5547 0.6467
## Residuals 72 16.2743 0.22603
summary(model_mass)
## 
## Call:
## lm(formula = Mass ~ Treatment, data = pre_diet)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.8368 -0.3412 -0.0487  0.3016  1.2372 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       2.83280    0.10631  26.647   <2e-16 ***
## TreatmentDR       0.19896    0.15684   1.269    0.209    
## TreatmentDR.IGF1  0.12180    0.15034   0.810    0.421    
## TreatmentDR.IGF2  0.09694    0.15231   0.636    0.527    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4754 on 72 degrees of freedom
## Multiple R-squared:  0.02259,    Adjusted R-squared:  -0.01814 
## F-statistic: 0.5547 on 3 and 72 DF,  p-value: 0.6467
plot_mass <- ggplot(
  pre_diet,
  aes(x = Treatment, y = Mass, fill = Treatment)
) +
  geom_violin(trim = FALSE, scale = "width") +
  geom_boxplot(width = 0.10, outlier.shape = NA) +
  geom_jitter(width = 0.10) +
  scale_fill_manual(values = treatment_colors) +
  labs(x = "Treatment", y = "Initial body mass") +
  theme_classic() +
  theme(legend.position = "none")

plot_mass

ggsave(
  filename = "mass.png",
  plot = plot_mass,
  width = 5,
  height = 4,
  dpi = 600
)

Manuscript result: Initial body mass did not differ among treatment groups, \(F_{3,72}=0.55\), \(p=0.65\).

3.4 Original tail length

model_tail_baseline <- lm(Tail ~ Treatment, data = pre_diet)
anova(model_tail_baseline)
## Analysis of Variance Table
## 
## Response: Tail
##           Df Sum Sq Mean Sq F value Pr(>F)  
## Treatment  3  420.0 139.994  2.3623 0.0795 .
## Residuals 64 3792.7  59.262                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model_tail_baseline)
## 
## Call:
## lm(formula = Tail ~ Treatment, data = pre_diet)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -19.000  -4.250   0.875   4.292  20.500 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        73.000      1.814  40.232   <2e-16 ***
## TreatmentDR         1.500      2.743   0.547   0.5864    
## TreatmentDR.IGF1    5.250      2.566   2.046   0.0449 *  
## TreatmentDR.IGF2    5.778      2.566   2.252   0.0278 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7.698 on 64 degrees of freedom
##   (8 observations deleted due to missingness)
## Multiple R-squared:  0.09969,    Adjusted R-squared:  0.05749 
## F-statistic: 2.362 on 3 and 64 DF,  p-value: 0.0795

Manuscript result: Original tail length did not differ significantly among treatment groups, \(F_{3,72}=2.36\), \(p=0.08\).

4. Verification of the Dietary Manipulation

4.1 Percent body-mass loss after acclimation

#Week -1 is the measurement following three weeks of restricted feeding, and one week prior to tail autotomy at "Week 0". 
diet_final_acclimation <- diet %>%
  filter(Week == -1)

model_diet_loss <- lm(Perc_WL ~ Diet, data = diet_final_acclimation)
anova(model_diet_loss)
## Analysis of Variance Table
## 
## Response: Perc_WL
##           Df  Sum Sq Mean Sq F value    Pr(>F)    
## Diet       1  608.22  608.22  13.854 0.0004068 ***
## Residuals 67 2941.39   43.90                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model_diet_loss)
## 
## Call:
## lm(formula = Perc_WL ~ Diet, data = diet_final_acclimation)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -18.1313  -4.6672   0.7973   4.0723  20.1000 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    1.104      1.562   0.707 0.482225    
## DietDR         6.761      1.817   3.722 0.000407 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.626 on 67 degrees of freedom
##   (7 observations deleted due to missingness)
## Multiple R-squared:  0.1713, Adjusted R-squared:  0.159 
## F-statistic: 13.85 on 1 and 67 DF,  p-value: 0.0004068
diet_group_summary <- diet_final_acclimation %>%
  group_by(Diet) %>%
  summarise(
    n = sum(!is.na(Perc_WL)),
    mean_percent_loss = mean(Perc_WL, na.rm = TRUE),
    sd_percent_loss = sd(Perc_WL, na.rm = TRUE),
    mean_percent_original_mass = mean(Perc_Orig, na.rm = TRUE),
    .groups = "drop"
  )

diet_group_summary
## # A tibble: 2 × 5
##   Diet      n mean_percent_loss sd_percent_loss mean_percent_original_mass
##   <chr> <int>             <dbl>           <dbl>                      <dbl>
## 1 Adlib    18              1.10            7.58                       98.9
## 2 DR       51              7.86            6.27                       92.1

Manuscript result: Diet-restricted females lost significantly more body mass than ad libitum-fed females, \(F_{1,67}=13.85\), \(p=0.0004\). The average DR individual lost ~7.8% of their body mass, while those on the Adlib diet lost an average of 1.1%.

4.2 Figure 1A: energetic manipulation

plot_diet <- ggplot(
  diet_final_acclimation,
  aes(x = Diet, y = Perc_Orig)
) +
  geom_violin(
    fill = "grey75",
    color = "black",
    trim = FALSE,
    scale = "width",
    alpha = 0.8
  ) +
  geom_boxplot(
    width = 0.12,
    outlier.shape = NA,
    fill = "white",
    color = "black"
  ) +
  geom_jitter(
    width = 0.08,
    size = 1.8,
    shape = 21,
    fill = "white",
    color = "black",
    alpha = 0.8
  ) +
  labs(
    x = "Diet",
    y = "Percent of original body mass"
  ) +
  theme_classic(base_size = 12)

plot_diet

ggsave(
  filename = "diet.png",
  plot = plot_diet,
  width = 5,
  height = 4,
  dpi = 600
)

Depiction: Violin plots illustrate the distribution of final body mass, expressed as a percentage of initial body mass, for ad libitum and diet-restricted animals. Dietary restriction shifted the overall distribution toward lower body mass, confirming successful induction of an energetic deficit. However, both dietary treatments produced a spectrum of body mass changes, resulting in substantial overlap between groups. This continuous variation in energetic state provided the opportunity to evaluate regeneration as a function of individual body mass loss rather than as a simple categorical comparison between dietary treatments.

5. Establishing the Energetic Trade-off

5.1 Primary cumulative-percent-regeneration model

model_weightloss <- lmerTest::lmer(
  Perc.Regeneration ~ Week * Perc.WeightLoss + Tail +
    (1 | Animal_ID),
  data = data
)

summary(model_weightloss)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Perc.Regeneration ~ Week * Perc.WeightLoss + Tail + (1 | Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 2747.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.1896 -0.5652  0.0844  0.5744  3.1799 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 15.69    3.961   
##  Residual              21.76    4.665   
## Number of obs: 443, groups:  Animal_ID, 65
## 
## Fixed effects:
##                       Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)           13.15695    5.28435  68.47432   2.490  0.01521 *  
## Week                   4.84774    0.19428 399.39842  24.953  < 2e-16 ***
## Perc.WeightLoss        0.10361    0.05420 436.49430   1.912  0.05657 .  
## Tail                  -0.25817    0.06936  68.87769  -3.722  0.00040 ***
## Week:Perc.WeightLoss  -0.03020    0.01021 403.58612  -2.958  0.00327 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Week   Prc.WL Tail  
## Week        -0.144                     
## Prc.WghtLss -0.047  0.632              
## Tail        -0.977 -0.003 -0.123       
## Wk:Prc.WghL  0.135 -0.848 -0.742 -0.012
anova(model_weightloss)
## Type III Analysis of Variance Table with Satterthwaite's method
##                       Sum Sq Mean Sq NumDF  DenDF  F value    Pr(>F)    
## Week                 13550.0 13550.0     1 399.40 622.6418 < 2.2e-16 ***
## Perc.WeightLoss         79.5    79.5     1 436.49   3.6544 0.0565742 .  
## Tail                   301.5   301.5     1  68.88  13.8548 0.0003996 ***
## Week:Perc.WeightLoss   190.5   190.5     1 403.59   8.7518 0.0032750 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Manuscript results:

  • Week: \(\beta=4.85\), \(t=24.95\), \(p<0.0001\)
  • Percent body-mass loss: \(\beta=0.10\), \(t=1.91\), \(p=0.06\)
  • Week × percent body-mass loss: \(\beta=-0.030\), \(t=-2.96\), \(p=0.0031\)

A negative Week × Percent Weight Loss interaction indicates that animals experiencing greater energetic deficits regenerated progressively less tail tissue over time.

5.2 Alternative regeneration metrics

These models retain the same fixed- and random-effect structure while replacing the response variable.

model_weightloss_rate <- lmerTest::lmer(
  Rate_Reg ~ Week * Perc.WeightLoss + Tail +
    (1 | Animal_ID),
  data = data
)

model_weightloss_absolute <- lmerTest::lmer(
  Regeneration ~ Week * Perc.WeightLoss + Tail +
    (1 | Animal_ID),
  data = data
)

summary(model_weightloss_rate)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Rate_Reg ~ Week * Perc.WeightLoss + Tail + (1 | Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 1857.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.4619 -0.7154 -0.2113  0.5768  6.0827 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 0.01256  0.1121  
##  Residual              7.32539  2.7065  
## Number of obs: 380, groups:  Animal_ID, 64
## 
## Fixed effects:
##                        Estimate Std. Error         df t value Pr(>|t|)   
## (Intercept)            4.621180   1.490139  64.189058   3.101  0.00286 **
## Week                  -0.017850   0.123697 359.320996  -0.144  0.88534   
## Perc.WeightLoss        0.061746   0.035109 374.665893   1.759  0.07944 . 
## Tail                  -0.028967   0.018671  48.287823  -1.551  0.12732   
## Week:Perc.WeightLoss  -0.007433   0.006688 359.705927  -1.111  0.26716   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Week   Prc.WL Tail  
## Week        -0.382                     
## Prc.WghtLss -0.248  0.732              
## Tail        -0.902 -0.007 -0.117       
## Wk:Prc.WghL  0.337 -0.817 -0.900 -0.020
anova(model_weightloss_rate)
## Type III Analysis of Variance Table with Satterthwaite's method
##                       Sum Sq Mean Sq NumDF  DenDF F value  Pr(>F)  
## Week                  0.1525  0.1525     1 359.32  0.0208 0.88534  
## Perc.WeightLoss      22.6578 22.6578     1 374.67  3.0931 0.07944 .
## Tail                 17.6323 17.6323     1  48.29  2.4070 0.12732  
## Week:Perc.WeightLoss  9.0474  9.0474     1 359.71  1.2351 0.26716  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model_weightloss_absolute)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Regeneration ~ Week * Perc.WeightLoss + Tail + (1 | Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 2330.5
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.90525 -0.61605  0.07588  0.59909  2.96424 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 6.839    2.615   
##  Residual              8.366    2.892   
## Number of obs: 442, groups:  Animal_ID, 65
## 
## Fixed effects:
##                        Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)            0.005838   3.446978  68.149511   0.002    0.999    
## Week                   2.964134   0.120670 396.631276  24.564   <2e-16 ***
## Perc.WeightLoss        0.017058   0.033812 434.412764   0.505    0.614    
## Tail                  -0.045651   0.045234  68.663813  -1.009    0.316    
## Week:Perc.WeightLoss  -0.005878   0.006342 400.910832  -0.927    0.355    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Week   Prc.WL Tail  
## Week        -0.137                     
## Prc.WghtLss -0.044  0.629              
## Tail        -0.979 -0.003 -0.120       
## Wk:Prc.WghL  0.128 -0.847 -0.738 -0.011
anova(model_weightloss_absolute)
## Type III Analysis of Variance Table with Satterthwaite's method
##                      Sum Sq Mean Sq NumDF  DenDF  F value Pr(>F)    
## Week                 5047.8  5047.8     1 396.63 603.3913 <2e-16 ***
## Perc.WeightLoss         2.1     2.1     1 434.41   0.2545 0.6142    
## Tail                    8.5     8.5     1  68.66   1.0185 0.3164    
## Week:Perc.WeightLoss    7.2     7.2     1 400.91   0.8593 0.3545    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Manuscript results: As a sensitivity analysis, the relationship between energetic state and regeneration was evaluated using two alternative response variables: weekly regeneration rate and cumulative regenerated tail length. Weekly regeneration rate showed no significant effects of week, percent body mass loss, or their interaction, although percent body mass loss exhibited a weak positive trend (p = 0.079). Similarly, when regeneration was expressed as cumulative regenerated tail length, regeneration increased significantly over time (β = 2.96, p < 0.0001), but neither percent body mass loss nor its interaction with week was significant. Together, these analyses indicate that expressing regeneration as cumulative percent regeneration provides the greatest sensitivity for detecting changes in the relationship between energetic state and regenerative performance while accounting for differences in initial tail length among individuals.

5.3 Exploratory visualization by energetic-state category

The primary analyses treat percent body-mass loss continuously. The following three-category plot is retained only as a descriptive visualization.

data <- data %>%
  mutate(
    WL_group = cut(Perc.WeightLoss, breaks = 3)
  )

ggplot(
  data,
  aes(
    x = Week,
    y = Perc.Regeneration,
    color = WL_group,
    group = WL_group
  )
) +
  stat_summary(fun = mean, geom = "line", linewidth = 1.2) +
  stat_summary(fun = mean, geom = "point") +
  labs(
    x = "Week",
    y = "Cumulative percent regeneration",
    color = "Weight-loss category"
  ) +
  theme_classic()

Depiction: For visualization purposes, animals were divided into tertiles based on percent body mass loss, although energetic state was analyzed as a continuous variable in all statistical models. Regeneration trajectories were similar during the early stages of regeneration but began to diverge after approximately Week 5. Animals that maintained the highest percentage of their original body mass exhibited the greatest cumulative regeneration, whereas those experiencing greater body mass loss regenerated more slowly over time. This pattern visually reinforces the primary mixed-model result that the effect of energetic state on regeneration increases as regeneration progresses. This was part of the reasoning for using percent weight loss as a continuous variable rather than the diet treatments.

6. Testing Rescue by IGF Supplementation

6.1 Primary IGF1-versus-IGF2 rescue model

model_igf_split <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss * IGF_type + Tail +
    (1 | Animal_ID),
  data = data
)

summary(model_igf_split)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type + Tail +  
##     (1 | Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 2717.5
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -3.05240 -0.58328  0.07647  0.53288  3.14152 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 15.52    3.940   
##  Residual              19.87    4.457   
## Number of obs: 443, groups:  Animal_ID, 65
## 
## Fixed effects:
##                                    Estimate Std. Error        df t value
## (Intercept)                         8.95632    5.37444  66.82767   1.666
## Week                                5.57153    0.22507 390.98913  24.755
## Perc.WeightLoss                     0.16173    0.06336 418.47176   2.553
## IGF_typeIGF1                        6.23221    2.90106 308.60824   2.148
## IGF_typeIGF2                        4.86116    3.88987 361.49614   1.250
## Tail                               -0.22402    0.07251  66.89851  -3.090
## Week:Perc.WeightLoss               -0.05150    0.01184 388.69919  -4.350
## Week:IGF_typeIGF1                  -1.60695    0.50288 411.98772  -3.196
## Week:IGF_typeIGF2                  -2.61979    0.61694 399.53649  -4.246
## Perc.WeightLoss:IGF_typeIGF1       -0.25717    0.14183 423.73347  -1.813
## Perc.WeightLoss:IGF_typeIGF2       -0.23525    0.17665 422.61347  -1.332
## Week:Perc.WeightLoss:IGF_typeIGF1   0.03935    0.02919 424.77894   1.348
## Week:Perc.WeightLoss:IGF_typeIGF2   0.11330    0.02987 393.97138   3.793
##                                   Pr(>|t|)    
## (Intercept)                       0.100302    
## Week                               < 2e-16 ***
## Perc.WeightLoss                   0.011041 *  
## IGF_typeIGF1                      0.032472 *  
## IGF_typeIGF2                      0.212219    
## Tail                              0.002919 ** 
## Week:Perc.WeightLoss              1.74e-05 ***
## Week:IGF_typeIGF1                 0.001503 ** 
## Week:IGF_typeIGF2                 2.70e-05 ***
## Perc.WeightLoss:IGF_typeIGF1      0.070503 .  
## Perc.WeightLoss:IGF_typeIGF2      0.183676    
## Week:Perc.WeightLoss:IGF_typeIGF1 0.178340    
## Week:Perc.WeightLoss:IGF_typeIGF2 0.000172 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(model_igf_split)
## Type III Analysis of Variance Table with Satterthwaite's method
##                               Sum Sq Mean Sq NumDF  DenDF  F value    Pr(>F)
## Week                          5315.5  5315.5     1 405.93 267.5724 < 2.2e-16
## Perc.WeightLoss                  0.0     0.0     1 419.69   0.0011 0.9735115
## IGF_type                       107.4    53.7     2 336.71   2.7031 0.0684535
## Tail                           189.6   189.6     1  66.90   9.5458 0.0029187
## Week:Perc.WeightLoss             0.0     0.0     1 413.66   0.0021 0.9634575
## Week:IGF_type                  486.0   243.0     2 406.02  12.2327 6.939e-06
## Perc.WeightLoss:IGF_type        87.4    43.7     2 422.06   2.1996 0.1121193
## Week:Perc.WeightLoss:IGF_type  296.9   148.5     2 412.47   7.4732 0.0006484
##                                  
## Week                          ***
## Perc.WeightLoss                  
## IGF_type                      .  
## Tail                          ** 
## Week:Perc.WeightLoss             
## Week:IGF_type                 ***
## Perc.WeightLoss:IGF_type         
## Week:Perc.WeightLoss:IGF_type ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Manuscript interpretation: IGF2 significantly altered the increasingly negative relationship between energetic deficit and regeneration ( P = 0.00017), whereas IGF1 did not (p =0.178) .

6.1.2 Alternative regeneration metrics

These models retain the same fixed- and random-effect structure while replacing the response variable.

model_weightloss_rate <- lmerTest::lmer(
  Rate_Reg ~ Week * Perc.WeightLoss *IGF_type + Tail + 
    (1 | Animal_ID),
  data = data
)

model_weightloss_absolute <- lmerTest::lmer(
  Regeneration ~ Week * Perc.WeightLoss  *IGF_type + Tail + 
    (1 | Animal_ID),
  data = data
)

summary(model_weightloss_rate)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Rate_Reg ~ Week * Perc.WeightLoss * IGF_type + Tail + (1 | Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 1876.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.3778 -0.6999 -0.2014  0.5595  6.1085 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 0.02198  0.1482  
##  Residual              7.35988  2.7129  
## Number of obs: 380, groups:  Animal_ID, 64
## 
## Fixed effects:
##                                     Estimate Std. Error         df t value
## (Intercept)                         3.983893   1.592815  62.654178   2.501
## Week                                0.022785   0.150929 344.637628   0.151
## Perc.WeightLoss                     0.078408   0.043764 366.933538   1.792
## IGF_typeIGF1                        0.174182   1.704939 355.773477   0.102
## IGF_typeIGF2                       -1.023781   2.207736 325.606939  -0.464
## Tail                               -0.018409   0.020299  42.917758  -0.907
## Week:Perc.WeightLoss               -0.013048   0.008160 336.923711  -1.599
## Week:IGF_typeIGF1                  -0.121004   0.326578 366.982012  -0.371
## Week:IGF_typeIGF2                  -0.151883   0.409571 366.830455  -0.371
## Perc.WeightLoss:IGF_typeIGF1       -0.051345   0.097548 357.354390  -0.526
## Perc.WeightLoss:IGF_typeIGF2       -0.006471   0.110176 344.373932  -0.059
## Week:Perc.WeightLoss:IGF_typeIGF1   0.011185   0.019099 360.490004   0.586
## Week:Perc.WeightLoss:IGF_typeIGF2   0.015650   0.020760 364.824718   0.754
##                                   Pr(>|t|)  
## (Intercept)                          0.015 *
## Week                                 0.880  
## Perc.WeightLoss                      0.074 .
## IGF_typeIGF1                         0.919  
## IGF_typeIGF2                         0.643  
## Tail                                 0.370  
## Week:Perc.WeightLoss                 0.111  
## Week:IGF_typeIGF1                    0.711  
## Week:IGF_typeIGF2                    0.711  
## Perc.WeightLoss:IGF_typeIGF1         0.599  
## Perc.WeightLoss:IGF_typeIGF2         0.953  
## Week:Perc.WeightLoss:IGF_typeIGF1    0.558  
## Week:Perc.WeightLoss:IGF_typeIGF2    0.451  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(model_weightloss_rate)
## Type III Analysis of Variance Table with Satterthwaite's method
##                                Sum Sq Mean Sq NumDF  DenDF F value Pr(>F)
## Week                           1.2234  1.2234     1 366.99  0.1662 0.6837
## Perc.WeightLoss               11.8232 11.8232     1 346.80  1.6064 0.2058
## IGF_type                       1.8351  0.9176     2 340.62  0.1247 0.8828
## Tail                           6.0531  6.0531     1  42.92  0.8224 0.3695
## Week:Perc.WeightLoss           1.5280  1.5280     1 366.17  0.2076 0.6489
## Week:IGF_type                  1.7282  0.8641     2 366.60  0.1174 0.8893
## Perc.WeightLoss:IGF_type       2.0494  1.0247     2 350.20  0.1392 0.8701
## Week:Perc.WeightLoss:IGF_type  5.7791  2.8896     2 362.92  0.3926 0.6756
summary(model_weightloss_absolute)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Regeneration ~ Week * Perc.WeightLoss * IGF_type + Tail + (1 |  
##     Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 2318.5
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.63135 -0.56996  0.08751  0.58321  3.08045 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 6.837    2.615   
##  Residual              7.847    2.801   
## Number of obs: 442, groups:  Animal_ID, 65
## 
## Fixed effects:
##                                     Estimate Std. Error         df t value
## (Intercept)                        -1.926997   3.529452  66.015831  -0.546
## Week                                3.286170   0.141627 388.192452  23.203
## Perc.WeightLoss                     0.049363   0.039961 415.191120   1.235
## IGF_typeIGF1                        2.509958   1.855201 301.514022   1.353
## IGF_typeIGF2                        1.628220   2.477492 361.831869   0.657
## Tail                               -0.028452   0.047614  66.202786  -0.598
## Week:Perc.WeightLoss               -0.016560   0.007448 386.035519  -2.223
## Week:IGF_typeIGF1                  -0.529278   0.317018 409.527055  -1.670
## Week:IGF_typeIGF2                  -1.338341   0.388458 396.681548  -3.445
## Perc.WeightLoss:IGF_typeIGF1       -0.112795   0.089886 423.495021  -1.255
## Perc.WeightLoss:IGF_typeIGF2       -0.104107   0.111946 424.229848  -0.930
## Week:Perc.WeightLoss:IGF_typeIGF1   0.005680   0.018437 423.020369   0.308
## Week:Perc.WeightLoss:IGF_typeIGF2   0.062718   0.018801 391.245932   3.336
##                                   Pr(>|t|)    
## (Intercept)                       0.586923    
## Week                               < 2e-16 ***
## Perc.WeightLoss                   0.217422    
## IGF_typeIGF1                      0.177091    
## IGF_typeIGF2                      0.511467    
## Tail                              0.552170    
## Week:Perc.WeightLoss              0.026772 *  
## Week:IGF_typeIGF1                 0.095772 .  
## Week:IGF_typeIGF2                 0.000631 ***
## Perc.WeightLoss:IGF_typeIGF1      0.210217    
## Perc.WeightLoss:IGF_typeIGF2      0.352911    
## Week:Perc.WeightLoss:IGF_typeIGF1 0.758172    
## Week:Perc.WeightLoss:IGF_typeIGF2 0.000932 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(model_weightloss_absolute)
## Type III Analysis of Variance Table with Satterthwaite's method
##                                Sum Sq Mean Sq NumDF  DenDF  F value    Pr(>F)
## Week                          2166.04 2166.04     1 403.24 276.0261 < 2.2e-16
## Perc.WeightLoss                  1.96    1.96     1 421.07   0.2495  0.617681
## IGF_type                        15.79    7.90     2 333.50   1.0061  0.366741
## Tail                             2.80    2.80     1  66.20   0.3571  0.552170
## Week:Perc.WeightLoss             4.31    4.31     1 411.44   0.5492  0.459063
## Week:IGF_type                  103.05   51.53     2 403.37   6.5661  0.001563
## Perc.WeightLoss:IGF_type        16.65    8.33     2 422.70   1.0610  0.347023
## Week:Perc.WeightLoss:IGF_type   87.74   43.87     2 410.23   5.5902  0.004024
##                                  
## Week                          ***
## Perc.WeightLoss                  
## IGF_type                         
## Tail                             
## Week:Perc.WeightLoss             
## Week:IGF_type                 ** 
## Perc.WeightLoss:IGF_type         
## Week:Perc.WeightLoss:IGF_type ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Manuscript results: As a sensitivity analysis, the primary model was repeated using alternative measures of regeneration, including weekly regeneration rate (mm/week) and cumulative regenerated tail length (mm). For weekly regeneration rate, neither week, percent body mass loss, IGF supplementation, nor any interaction terms significantly predicted regeneration rate, although percent body mass loss exhibited a weak positive trend (β = 0.078, p = 0.074). Similarly, when regeneration was expressed as cumulative regenerated tail length, regeneration increased significantly over time (β = 3.29, p < 0.0001), and the interaction between week and percent body mass loss remained significantly negative (β = −0.0166, p = 0.027), indicating that the energetic trade-off became stronger as regeneration progressed. Consistent with the primary analysis, IGF2 significantly attenuated this negative relationship over time (Week × IGF2: β = −1.34, p = 0.0006; Week × Percent Weight Loss × IGF2: β = 0.0627, p = 0.0009), whereas IGF1 showed no significant effects. Together, these analyses demonstrate that the primary conclusions are robust to alternative measures of regeneration, with cumulative regeneration yielding results that closely mirror those obtained using percent regeneration.

6.2 Model-based predictions

predicted_all_weeks <- ggpredict(
  model_igf_split,
  terms = c(
    "Perc.WeightLoss",
    "IGF_type",
    "Week [2,4,6,8]"
  )
)

predicted_all_weeks
## # Predicted values of Perc.Regeneration
## 
## IGF_type: No_IGF
## Week: 2
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     -0.54 | -7.95,  6.87
##             -35 |      0.93 | -4.17,  6.02
##              -5 |      2.69 |  0.13,  5.25
##              50 |      5.92 |  2.18,  9.66
## 
## IGF_type: No_IGF
## Week: 4
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     16.78 | 10.10, 23.46
##             -35 |     15.67 | 11.04, 20.30
##              -5 |     14.35 | 11.96, 16.73
##              50 |     11.91 |  8.62, 15.21
## 
## IGF_type: No_IGF
## Week: 6
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     34.10 | 26.42, 41.78
##             -35 |     30.42 | 25.11, 35.73
##              -5 |     26.00 | 23.33, 28.68
##              50 |     17.91 | 14.19, 21.62
## 
## IGF_type: No_IGF
## Week: 8
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     51.42 | 41.52, 61.33
##             -35 |     45.17 | 38.37, 51.96
##              -5 |     37.66 | 34.36, 40.97
##              50 |     23.90 | 19.12, 28.68
## 
## IGF_type: IGF1
## Week: 2
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     13.18 | -0.91, 27.28
##             -35 |     10.19 |  0.49, 19.90
##              -5 |      6.60 |  1.95, 11.25
##              50 |      0.01 | -6.19,  6.21
## 
## IGF_type: IGF1
## Week: 4
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     22.57 | 10.39, 34.75
##             -35 |     18.97 | 10.59, 27.35
##              -5 |     14.65 | 10.60, 18.70
##              50 |      6.73 |  1.18, 12.28
## 
## IGF_type: IGF1
## Week: 6
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     31.96 | 17.01, 46.91
##             -35 |     27.75 | 17.61, 37.89
##              -5 |     22.70 | 18.08, 27.32
##              50 |     13.44 |  6.29, 20.59
## 
## IGF_type: IGF1
## Week: 8
## 
## Perc.WeightLoss | Predicted |       95% CI
## ------------------------------------------
##             -60 |     41.34 | 20.75, 61.94
##             -35 |     36.53 | 22.67, 50.38
##              -5 |     30.75 | 24.72, 36.78
##              50 |     20.16 | 10.19, 30.13
## 
## IGF_type: IGF2
## Week: 2
## 
## Perc.WeightLoss | Predicted |        95% CI
## -------------------------------------------
##             -60 |     -0.40 | -20.36, 19.55
##             -35 |      0.85 | -13.01, 14.71
##              -5 |      2.35 |  -4.34,  9.04
##              50 |      5.11 |  -2.43, 12.65
## 
## IGF_type: IGF2
## Week: 4
## 
## Perc.WeightLoss | Predicted |        95% CI
## -------------------------------------------
##             -60 |     -1.92 | -17.92, 14.09
##             -35 |      2.43 |  -8.71, 13.56
##              -5 |      7.64 |   2.20, 13.07
##              50 |     17.19 |  11.02, 23.36
## 
## IGF_type: IGF2
## Week: 6
## 
## Perc.WeightLoss | Predicted |        95% CI
## -------------------------------------------
##             -60 |     -3.43 | -19.66, 12.80
##             -35 |      4.00 |  -7.25, 15.26
##              -5 |     12.92 |   7.47, 18.37
##              50 |     29.28 |  22.83, 35.72
## 
## IGF_type: IGF2
## Week: 8
## 
## Perc.WeightLoss | Predicted |        95% CI
## -------------------------------------------
##             -60 |     -4.94 | -25.43, 15.54
##             -35 |      5.58 |  -8.57, 19.73
##              -5 |     18.21 |  11.48, 24.94
##              50 |     41.36 |  33.15, 49.56
## 
## Adjusted for:
## *      Tail = 76.41
## * Animal_ID = 0 (population-level)
plot_predictions <- ggplot(
  predicted_all_weeks,
  aes(
    x = x,
    y = predicted,
    color = group,
    fill = group
  )
) +
  geom_ribbon(
    aes(ymin = conf.low, ymax = conf.high),
    alpha = 0.15,
    color = NA
  ) +
  geom_line(linewidth = 1.2) +
  facet_wrap(~facet) +
  scale_color_manual(values = igf_colors) +
  scale_fill_manual(values = igf_colors) +
  labs(
    x = "Percent body-mass loss",
    y = "Predicted cumulative percent regeneration",
    color = "IGF treatment",
    fill = "IGF treatment"
  ) +
  theme_classic()

plot_predictions

Depiction: Predicted regeneration trajectories demonstrate that the relationship between energetic state and regeneration strengthened over time. By Weeks 6–8, IGF2 reversed the negative association between body mass loss and regeneration observed in untreated and IGF1-treated animals, indicating progressive attenuation of the energetic trade-off.

6.3 Week-specific energetic-state slopes

week_specific_slopes <- emtrends(
  model_igf_split,
  ~ IGF_type | Week,
  var = "Perc.WeightLoss",
  at = list(Week = c(2, 4, 6, 8))
)

summary(week_specific_slopes, infer = TRUE)
## Week = 2:
##  IGF_type Perc.WeightLoss.trend     SE  df lower.CL upper.CL t.ratio p.value
##  No_IGF                  0.0587 0.0491 430  -0.0377   0.1552   1.197  0.2321
##  IGF1                   -0.1197 0.0915 392  -0.2996   0.0602  -1.308  0.1915
##  IGF2                    0.0501 0.1260 366  -0.1976   0.2978   0.398  0.6911
## 
## Week = 4:
##  IGF_type Perc.WeightLoss.trend     SE  df lower.CL upper.CL t.ratio p.value
##  No_IGF                 -0.0442 0.0437 414  -0.1301   0.0416  -1.013  0.3116
##  IGF1                   -0.1440 0.0796 408  -0.3004   0.0124  -1.810  0.0710
##  IGF2                    0.1737 0.1010 307  -0.0255   0.3729   1.716  0.0872
## 
## Week = 6:
##  IGF_type Perc.WeightLoss.trend     SE  df lower.CL upper.CL t.ratio p.value
##  No_IGF                 -0.1472 0.0503 422  -0.2461  -0.0484  -2.928  0.0036
##  IGF1                   -0.1683 0.1000 430  -0.3650   0.0284  -1.682  0.0933
##  IGF2                    0.2973 0.1030 367   0.0943   0.5003   2.880  0.0042
## 
## Week = 8:
##  IGF_type Perc.WeightLoss.trend     SE  df lower.CL upper.CL t.ratio p.value
##  No_IGF                 -0.2502 0.0654 430  -0.3787  -0.1217  -3.828  0.0001
##  IGF1                   -0.1926 0.1390 425  -0.4665   0.0813  -1.382  0.1677
##  IGF2                    0.4209 0.1310 430   0.1639   0.6780   3.218  0.0014
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95
week_specific_slopes_df <- as.data.frame(
  summary(week_specific_slopes, infer = TRUE)
)

Manuscript Result: Week-specific slope estimates demonstrated that energetic state had little influence on regeneration during the early stages of tail regrowth (Weeks 2–4), as none of the treatment-specific slopes differed significantly from zero. By Week 6, untreated animals exhibited a significant negative relationship between body mass loss and regeneration (β = −0.147, p = 0.0036), whereas IGF2-treated animals exhibited a significant positive relationship (β = 0.297, p = 0.0042), indicating attenuation of the energetic trade-off. This divergence became even more pronounced by Week 8, with untreated animals showing a stronger negative relationship (β = −0.250, p = 0.0001) and IGF2-treated animals showing a stronger positive relationship (β = 0.421, p = 0.0014). Although IGF1 exhibited a modest positive trend at Week 4 (p = 0.071), it did not differ significantly from zero at any time point.

6.4 Pairwise comparisons of slopes

The manuscript describes Tukey-adjusted comparisons among supplementation groups. Both the full pairwise comparisons and treatment-versus-control contrasts are retained.

# All pairwise comparisons within week, Tukey adjusted
pairwise_slopes_tukey <- pairs(
  week_specific_slopes,
  adjust = "tukey"
)

pairwise_slopes_tukey
## Week = 2:
##  contrast      estimate     SE  df t.ratio p.value
##  No_IGF - IGF1  0.17847 0.1040 408   1.713  0.2014
##  No_IGF - IGF2  0.00865 0.1350 379   0.064  0.9977
##  IGF1 - IGF2   -0.16982 0.1560 376  -1.090  0.5206
## 
## Week = 4:
##  contrast      estimate     SE  df t.ratio p.value
##  No_IGF - IGF1  0.09976 0.0911 412   1.095  0.5176
##  No_IGF - IGF2 -0.21795 0.1100 325  -1.980  0.1189
##  IGF1 - IGF2   -0.31771 0.1290 350  -2.467  0.0375
## 
## Week = 6:
##  contrast      estimate     SE  df t.ratio p.value
##  No_IGF - IGF1  0.02106 0.1120 430   0.188  0.9808
##  No_IGF - IGF2 -0.44454 0.1150 379  -3.876  0.0004
##  IGF1 - IGF2   -0.46560 0.1440 412  -3.238  0.0037
## 
## Week = 8:
##  contrast      estimate     SE  df t.ratio p.value
##  No_IGF - IGF1 -0.05765 0.1540 426  -0.374  0.9258
##  No_IGF - IGF2 -0.67114 0.1460 430  -4.594 <0.0001
##  IGF1 - IGF2   -0.61349 0.1910 429  -3.210  0.0041
## 
## Degrees-of-freedom method: kenward-roger 
## P value adjustment: tukey method for comparing a family of 3 estimates

Manuscript Result: Tukey-adjusted pairwise comparisons of week-specific slopes showed no significant differences among treatments at Week 2. By Week 4, the slope for IGF2 differed significantly from IGF1 (p = 0.038) but not from untreated controls (p = 0.119). Significant differences between IGF2-treated and untreated animals emerged by Week 6 (p = 0.0004) and became even more pronounced by Week 8 (p < 0.0001), indicating that IGF2 progressively attenuated the energetic trade-off over time. IGF1 did not differ from untreated controls at any time point.

7. Bootstrap Evaluation of the Rescue Effect

The analysis contains two bootstrap procedures. Both are retained here because they address different purposes:

  1. A parametric bootstrap of the full rescue model for visualizing fixed-effect distributions.
  2. A late-stage bootstrap used in the sensitivity analysis of weeks 5–8.

7.1 Parametric bootstrap of the full rescue model

boot_fixed_effects <- function(fit) {
  fixef(fit)
}

bootstrap_full <- bootMer(
  model_igf_split,
  FUN = boot_fixed_effects,
  nsim = 500,
  type = "parametric"
)

bootstrap_full_df <- as.data.frame(bootstrap_full$t)
names(bootstrap_full_df) <- names(fixef(model_igf_split))

bootstrap_full_long <- bootstrap_full_df %>%
  select(
    `Week:Perc.WeightLoss`,
    `Week:Perc.WeightLoss:IGF_typeIGF1`,
    `Week:Perc.WeightLoss:IGF_typeIGF2`
  ) %>%
  pivot_longer(
    cols = everything(),
    names_to = "term",
    values_to = "estimate"
  )

bootstrap_plot <- ggplot(
  bootstrap_full_long,
  aes(x = estimate, fill = term)
) +
  geom_density(alpha = 0.9) +
  geom_vline(xintercept = 0, linetype = "dashed") +
  scale_fill_manual(
    values = c("grey40", "#E87722", "#0C2340")
  ) +
  labs(
    x = "Bootstrap estimate",
    y = "Density",
    fill = "Model term"
  ) +
  theme_classic() +
  theme(legend.position = "none")

bootstrap_plot

ggsave(
  filename = "bootstrapConfid.png",
  plot = bootstrap_plot,
  width = 2.5,
  height = 2,
  dpi = 600
)

8. Additional and Sensitivity Analyses

8.1 Administered hormone dose

igf_data <- data %>%
  filter(Treatment %in% c("DR.IGF1", "DR.IGF2")) %>%
  mutate(
    IGF_type = factor(
      Treatment,
      levels = c("DR.IGF1", "DR.IGF2")
    )
  ) %>%
  group_by(Animal_ID) %>%
  fill(Dose, .direction = "downup") %>%
  ungroup()

model_dose <- lmerTest::lmer(
  Rate_Reg ~ Week * Dose * IGF_type +
    Tail +
    (1 | Animal_ID),
  data = igf_data
)

summary(model_dose)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Rate_Reg ~ Week * Dose * IGF_type + Tail + (1 | Animal_ID)
##    Data: igf_data
## 
## REML criterion at convergence: 964.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.3156 -0.6620 -0.1643  0.4927  6.1166 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 0.000    0.000   
##  Residual              7.389    2.718   
## Number of obs: 200, groups:  Animal_ID, 33
## 
## Fixed effects:
##                            Estimate Std. Error        df t value Pr(>|t|)
## (Intercept)                -0.34343    6.62219 191.00000  -0.052    0.959
## Week                        0.37140    1.17613 191.00000   0.316    0.753
## Dose                        2.27038    2.74356 191.00000   0.828    0.409
## IGF_typeDR.IGF2             4.14492    7.97901 191.00000   0.519    0.604
## Tail                       -0.01755    0.03851 191.00000  -0.456    0.649
## Week:Dose                  -0.23512    0.54290 191.00000  -0.433    0.665
## Week:IGF_typeDR.IGF2       -0.78700    1.57735 191.00000  -0.499    0.618
## Dose:IGF_typeDR.IGF2       -1.96201    3.65967 191.00000  -0.536    0.593
## Week:Dose:IGF_typeDR.IGF2   0.38027    0.72186 191.00000   0.527    0.599
## 
## Correlation of Fixed Effects:
##             (Intr) Week   Dose   IGF_DR Tail   Wek:Ds W:IGF_ D:IGF_
## Week        -0.810                                                 
## Dose        -0.884  0.911                                          
## IGF_DR.IGF2 -0.644  0.684  0.738                                   
## Tail        -0.445 -0.029 -0.010 -0.049                            
## Week:Dose    0.807 -0.993 -0.920 -0.680  0.023                     
## W:IGF_DR.IG  0.606 -0.745 -0.679 -0.923  0.016  0.740              
## D:IGF_DR.IG  0.645 -0.684 -0.750 -0.992  0.047  0.690  0.914       
## W:D:IGF_DR. -0.609  0.746  0.692  0.916 -0.013 -0.752 -0.992 -0.923
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
anova(model_dose)
## Type III Analysis of Variance Table with Satterthwaite's method
##                    Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Week               0.0058  0.0058     1   191  0.0008 0.9777
## Dose               3.6730  3.6730     1   191  0.4971 0.4816
## IGF_type           1.9939  1.9939     1   191  0.2699 0.6040
## Tail               1.5350  1.5350     1   191  0.2077 0.6491
## Week:Dose          0.1147  0.1147     1   191  0.0155 0.9010
## Week:IGF_type      1.8394  1.8394     1   191  0.2489 0.6184
## Dose:IGF_type      2.1237  2.1237     1   191  0.2874 0.5925
## Week:Dose:IGF_type 2.0505  2.0505     1   191  0.2775 0.5989

Interpretation: Modest individual variation in administered dose did not detectably influence regeneration rate or explain the different responses to IGF1 and IGF2.

model.rescue2 <- lmer(
  Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type + Tail +
    (1 | Animal_ID),
  data = data,
  REML = FALSE
)

model.rescue.dose <- lmer(
  Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type +  Tail + Dose +
    (1 | Animal_ID),
  data = data,
  REML = FALSE
)

AIC(model.rescue2, model.rescue.dose)
##                   df      AIC
## model.rescue2     15 2715.726
## model.rescue.dose 16 2717.558

Result: Adding IGF dose as an additional covariate did not improve model fit (ΔAIC = 1.91), indicating that variation in administered dose did not explain additional variation in regeneration beyond treatment identity.

8.2 Endpoint covariance between energetic state and regeneration

final_regeneration_data <- data %>%
  group_by(Animal_ID) %>%
  summarise(
    final_regeneration = max(Regeneration, na.rm = TRUE),
    weight_loss = first(Perc.WeightLoss),
    Treatment = first(Treatment),
    .groups = "drop"
  )

cor.test(
  final_regeneration_data$weight_loss,
  final_regeneration_data$final_regeneration
)
## 
##  Pearson's product-moment correlation
## 
## data:  final_regeneration_data$weight_loss and final_regeneration_data$final_regeneration
## t = 0.56258, df = 61, p-value = 0.5758
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1791079  0.3140206
## sample estimates:
##       cor 
## 0.0718454
model_endpoint_tradeoff <- lm(
  final_regeneration ~ weight_loss * Treatment,
  data = final_regeneration_data
)

summary(model_endpoint_tradeoff)
## 
## Call:
## lm(formula = final_regeneration ~ weight_loss * Treatment, data = final_regeneration_data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -13.0612  -3.8851  -0.1858   3.8499  13.4936 
## 
## Coefficients:
##                              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  17.65763    3.43266   5.144  3.7e-06 ***
## weight_loss                   0.03591    0.21771   0.165    0.870    
## TreatmentDR                  -5.52370    4.36111  -1.267    0.211    
## TreatmentDR.IGF1              1.13593    5.18537   0.219    0.827    
## TreatmentDR.IGF2             -9.61449    6.61079  -1.454    0.152    
## weight_loss:TreatmentDR       0.14270    0.24926   0.572    0.569    
## weight_loss:TreatmentDR.IGF1 -0.36936    0.28830  -1.281    0.206    
## weight_loss:TreatmentDR.IGF2  0.38573    0.33340   1.157    0.252    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.379 on 55 degrees of freedom
##   (13 observations deleted due to missingness)
## Multiple R-squared:  0.2204, Adjusted R-squared:  0.1212 
## F-statistic: 2.222 on 7 and 55 DF,  p-value: 0.04612
anova(model_endpoint_tradeoff)
## Analysis of Variance Table
## 
## Response: final_regeneration
##                       Df  Sum Sq Mean Sq F value  Pr(>F)  
## weight_loss            1   14.82  14.821  0.3642 0.54868  
## Treatment              3  320.43 106.809  2.6245 0.05954 .
## weight_loss:Treatment  3  297.71  99.238  2.4385 0.07418 .
## Residuals             55 2238.30  40.696                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

This endpoint analysis is retained because it demonstrates that the longitudinal mixed-effects model is more sensitive to changes in regenerative dynamics than a single final-time-point comparison.

8.4 Final percent body-mass loss among treatment groups

final_weight_loss <- data %>%
  ungroup() %>%
  group_by(Animal_ID) %>%
  summarise(
    final_weight_loss = max(Perc.WeightLoss, na.rm = TRUE),
    Treatment = first(na.omit(Treatment)),
    .groups = "drop"
  ) %>%
  filter(
    !is.na(final_weight_loss),
    is.finite(final_weight_loss),
    !is.na(Treatment)
  ) %>%
  mutate(
    Treatment = droplevels(as.factor(Treatment))
  )

model_final_weight_loss <- lm(
  final_weight_loss ~ Treatment,
  data = final_weight_loss
)

anova(model_final_weight_loss)
## Analysis of Variance Table
## 
## Response: final_weight_loss
##           Df Sum Sq Mean Sq F value  Pr(>F)  
## Treatment  3  602.9 200.972   2.953 0.03951 *
## Residuals 61 4151.5  68.057                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
final_weight_loss_emmeans <- emmeans(
  model_final_weight_loss,
  pairwise ~ Treatment,
  adjust = "tukey"
)

final_weight_loss_emmeans
## $emmeans
##  Treatment emmean   SE df lower.CL upper.CL
##  AdLib       17.9 1.94 61     14.0     21.7
##  DR          25.6 2.29 61     21.0     30.2
##  DR.IGF1     23.1 2.00 61     19.1     27.1
##  DR.IGF2     24.7 2.00 61     20.7     28.7
## 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast          estimate   SE df t.ratio p.value
##  AdLib - DR          -7.765 3.00 61  -2.586  0.0571
##  AdLib - DR.IGF1     -5.257 2.79 61  -1.884  0.2456
##  AdLib - DR.IGF2     -6.841 2.79 61  -2.452  0.0781
##  DR - DR.IGF1         2.508 3.04 61   0.825  0.8423
##  DR - DR.IGF2         0.924 3.04 61   0.304  0.9901
##  DR.IGF1 - DR.IGF2   -1.584 2.83 61  -0.560  0.9435
## 
## P value adjustment: tukey method for comparing a family of 4 estimates

Manuscript interpretation: Final percent body-mass loss did not differ among the three diet-restricted groups (p > 0.65), supporting the conclusion that the IGF2 rescue effect was not explained by reduced energetic deficit.

ggplot(
  final_weight_loss,
  aes(
    x = Treatment,
    y = final_weight_loss,
    fill = Treatment
  )
) +
  geom_violin(alpha = 0.3, trim = FALSE) +
  geom_jitter(width = 0.1, alpha = 0.7, size = 2) +
  stat_summary(
    fun = mean,
    geom = "point",
    size = 4,
    color = "black"
  ) +
  stat_summary(
    fun.data = mean_cl_normal,
    geom = "errorbar",
    width = 0.15,
    color = "black"
  ) +
  scale_fill_manual(values = treatment_colors) +
  labs(
    x = "Treatment",
    y = "Final percent body-mass loss"
  ) +
  theme_classic() +
  theme(legend.position = "none")

9. Evaluation of Weeks 5–8

9.1 Visual inspection of individual regeneration trajectories

regeneration_raw <- read.csv(
  "R.analysis.currated.csv",
  stringsAsFactors = FALSE
)

regeneration_after_week_2 <- regeneration_raw %>%
  filter(Week > 2)

plot_length_loss_all <- ggplot(
  regeneration_raw,
  aes(
    x = Week,
    y = Regeneration,
    color = Animal_ID,
    group = Animal_ID
  )
) +
  geom_line() +
  geom_text(
    aes(label = Animal_ID),
    hjust = 0,
    vjust = 0,
    size = 2
  ) +
  facet_wrap(
    ~Treatment,
    scales = "free",
    ncol = 2
  ) +
  theme(legend.position = "none")

plot_length_loss_late <- ggplot(
  regeneration_after_week_2,
  aes(
    x = Week,
    y = Regeneration,
    color = Animal_ID,
    group = Animal_ID
  )
) +
  geom_line() +
  geom_text(
    aes(label = Animal_ID),
    hjust = 0,
    vjust = 0,
    size = 2
  ) +
  facet_wrap(
    ~Treatment,
    scales = "free",
    ncol = 2
  ) +
  theme(legend.position = "none")

plot_length_loss_all

plot_length_loss_late

ggsave(
  filename = "length_loss.png",
  plot = plot_length_loss_all,
  width = 8,
  height = 8,
  dpi = 600
)

ggsave(
  filename = "length_loss_sub.png",
  plot = plot_length_loss_late,
  width = 8,
  height = 8,
  dpi = 600
)

9.2 Early-only versus full eight-week model

model_early <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss * IGF_type +
    Tail +
    (1 | Animal_ID),
  data = data %>%
    filter(Week <= 4)
)

model_full <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss * IGF_type +
    Tail +
    (1 | Animal_ID),
  data = data
)

summary(model_early)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type + Tail +  
##     (1 | Animal_ID)
##    Data: data %>% filter(Week <= 4)
## 
## REML criterion at convergence: 1398.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.3603 -0.6766  0.0364  0.5755  3.6846 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept)  2.109   1.452   
##  Residual              13.084   3.617   
## Number of obs: 250, groups:  Animal_ID, 65
## 
## Fixed effects:
##                                    Estimate Std. Error        df t value
## (Intercept)                        -1.03249    3.28353  88.26216  -0.314
## Week                                6.51706    0.61407 214.47589  10.613
## Perc.WeightLoss                     0.13154    0.07686 203.94379   1.711
## IGF_typeIGF1                        5.00021    3.16573 232.96815   1.579
## IGF_typeIGF2                        2.41231    4.08500 229.06414   0.591
## Tail                               -0.12144    0.04048  63.90994  -3.000
## Week:Perc.WeightLoss               -0.04154    0.03320 216.92163  -1.251
## Week:IGF_typeIGF1                  -2.72621    1.30759 221.65810  -2.085
## Week:IGF_typeIGF2                  -2.28603    1.51661 194.62876  -1.507
## Perc.WeightLoss:IGF_typeIGF1       -0.17276    0.16381 230.81261  -1.055
## Perc.WeightLoss:IGF_typeIGF2       -0.04365    0.19071 228.74930  -0.229
## Week:Perc.WeightLoss:IGF_typeIGF1   0.08865    0.07167 232.22843   1.237
## Week:Perc.WeightLoss:IGF_typeIGF2   0.06108    0.07022 197.01077   0.870
##                                   Pr(>|t|)    
## (Intercept)                        0.75393    
## Week                               < 2e-16 ***
## Perc.WeightLoss                    0.08852 .  
## IGF_typeIGF1                       0.11558    
## IGF_typeIGF2                       0.55542    
## Tail                               0.00385 ** 
## Week:Perc.WeightLoss               0.21227    
## Week:IGF_typeIGF1                  0.03822 *  
## Week:IGF_typeIGF2                  0.13335    
## Perc.WeightLoss:IGF_typeIGF1       0.29270    
## Perc.WeightLoss:IGF_typeIGF2       0.81917    
## Week:Perc.WeightLoss:IGF_typeIGF1  0.21739    
## Week:Perc.WeightLoss:IGF_typeIGF2  0.38547    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(model_full)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type + Tail +  
##     (1 | Animal_ID)
##    Data: data
## 
## REML criterion at convergence: 2717.5
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -3.05240 -0.58328  0.07647  0.53288  3.14152 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 15.52    3.940   
##  Residual              19.87    4.457   
## Number of obs: 443, groups:  Animal_ID, 65
## 
## Fixed effects:
##                                    Estimate Std. Error        df t value
## (Intercept)                         8.95632    5.37444  66.82767   1.666
## Week                                5.57153    0.22507 390.98913  24.755
## Perc.WeightLoss                     0.16173    0.06336 418.47176   2.553
## IGF_typeIGF1                        6.23221    2.90106 308.60824   2.148
## IGF_typeIGF2                        4.86116    3.88987 361.49614   1.250
## Tail                               -0.22402    0.07251  66.89851  -3.090
## Week:Perc.WeightLoss               -0.05150    0.01184 388.69919  -4.350
## Week:IGF_typeIGF1                  -1.60695    0.50288 411.98772  -3.196
## Week:IGF_typeIGF2                  -2.61979    0.61694 399.53649  -4.246
## Perc.WeightLoss:IGF_typeIGF1       -0.25717    0.14183 423.73347  -1.813
## Perc.WeightLoss:IGF_typeIGF2       -0.23525    0.17665 422.61347  -1.332
## Week:Perc.WeightLoss:IGF_typeIGF1   0.03935    0.02919 424.77894   1.348
## Week:Perc.WeightLoss:IGF_typeIGF2   0.11330    0.02987 393.97138   3.793
##                                   Pr(>|t|)    
## (Intercept)                       0.100302    
## Week                               < 2e-16 ***
## Perc.WeightLoss                   0.011041 *  
## IGF_typeIGF1                      0.032472 *  
## IGF_typeIGF2                      0.212219    
## Tail                              0.002919 ** 
## Week:Perc.WeightLoss              1.74e-05 ***
## Week:IGF_typeIGF1                 0.001503 ** 
## Week:IGF_typeIGF2                 2.70e-05 ***
## Perc.WeightLoss:IGF_typeIGF1      0.070503 .  
## Perc.WeightLoss:IGF_typeIGF2      0.183676    
## Week:Perc.WeightLoss:IGF_typeIGF1 0.178340    
## Week:Perc.WeightLoss:IGF_typeIGF2 0.000172 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Reported comparison:

  • Weeks 1–4: IGF2 interaction \(\beta=0.088\), \(p=0.213\)
  • Weeks 1–8: IGF2 interaction \(\beta=0.114\), \(p=0.000138\)

9.3 Within-individual repeatability

icc(model_full)
## # Intraclass Correlation Coefficient
## 
##     Adjusted ICC: 0.439
##   Unadjusted ICC: 0.114

Reported result: Adjusted ICC approximately 0.47, indicating moderate within-individual repeatability across the full eight-week dataset.

9.4 Residual variation across weeks

full_model_residual_data <- model.frame(model_full) %>%
  mutate(
    residual = residuals(model_full),
    absolute_residual = abs(residual)
  )

model_absolute_residuals <- lm(
  absolute_residual ~ factor(Week),
  data = full_model_residual_data
)

anova(model_absolute_residuals)
## Analysis of Variance Table
## 
## Response: absolute_residual
##               Df  Sum Sq Mean Sq F value    Pr(>F)    
## factor(Week)   7  279.36  39.909  6.1601 6.937e-07 ***
## Residuals    435 2818.18   6.479                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

9.5 Late-stage bootstrap

late_data <- data %>%
  filter(Week >= 5)

model_late <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss * IGF_type +
    Tail +
    (1 | Animal_ID),
  data = late_data
)

late_boot_function <- function(fit) {
  fixef(fit)["Week:Perc.WeightLoss:IGF_typeIGF2"]
}

bootstrap_late <- bootMer(
  model_late,
  FUN = late_boot_function,
  nsim = 1000
)

late_bootstrap_interval <- quantile(
  bootstrap_late$t,
  c(0.025, 0.5, 0.975),
  na.rm = TRUE
)

late_bootstrap_interval
##       2.5%        50%      97.5% 
## 0.06210499 0.15405441 0.24375908

Reported result: Bootstrap resampling of the late-stage regeneration data confirmed that the IGF2 rescue effect was both positive and highly reproducible. The median estimated rescue coefficient was 0.153 (95% bootstrap CI = 0.069–0.244), with the confidence interval remaining entirely above zero, indicating that the attenuation of the energetic trade-off by IGF2 was robust to repeated resampling of the data.

10. Evaluation of Original Tail Length as a Covariate

Models are fitted with maximum likelihood for valid AIC comparison.

model_rescue_without_tail <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss * IGF_type +
    (1 | Animal_ID),
  data = data,
  REML = FALSE
)

model_rescue_with_tail <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss * IGF_type +
    Tail +
    (1 | Animal_ID),
  data = data,
  REML = FALSE
)

AIC(
  model_rescue_without_tail,
  model_rescue_with_tail
)
##                           df      AIC
## model_rescue_without_tail 14 2723.345
## model_rescue_with_tail    15 2715.726
summary(model_rescue_with_tail)
## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
##   method [lmerModLmerTest]
## Formula: Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type + Tail +  
##     (1 | Animal_ID)
##    Data: data
## 
##       AIC       BIC    logLik -2*log(L)  df.resid 
##    2715.7    2777.1   -1342.9    2685.7       428 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.0681 -0.5929  0.0779  0.5390  3.1805 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 14.17    3.764   
##  Residual              19.44    4.409   
## Number of obs: 443, groups:  Animal_ID, 65
## 
## Fixed effects:
##                                    Estimate Std. Error        df t value
## (Intercept)                         8.95617    5.17142  72.03898   1.732
## Week                                5.57428    0.22248 400.21150  25.056
## Perc.WeightLoss                     0.16235    0.06253 429.96459   2.596
## IGF_typeIGF1                        6.18235    2.84050 331.02945   2.177
## IGF_typeIGF2                        4.88263    3.81574 378.73719   1.280
## Tail                               -0.22414    0.06977  72.01932  -3.213
## Week:Perc.WeightLoss               -0.05169    0.01170 397.72469  -4.417
## Week:IGF_typeIGF1                  -1.60184    0.49654 422.19869  -3.226
## Week:IGF_typeIGF2                  -2.62240    0.60958 409.23308  -4.302
## Perc.WeightLoss:IGF_typeIGF1       -0.25435    0.13957 438.29668  -1.822
## Perc.WeightLoss:IGF_typeIGF2       -0.23644    0.17381 436.11440  -1.360
## Week:Perc.WeightLoss:IGF_typeIGF1   0.03902    0.02879 436.27083   1.355
## Week:Perc.WeightLoss:IGF_typeIGF2   0.11350    0.02952 403.28248   3.844
##                                   Pr(>|t|)    
## (Intercept)                        0.08758 .  
## Week                               < 2e-16 ***
## Perc.WeightLoss                    0.00974 ** 
## IGF_typeIGF1                       0.03022 *  
## IGF_typeIGF2                       0.20147    
## Tail                               0.00197 ** 
## Week:Perc.WeightLoss              1.29e-05 ***
## Week:IGF_typeIGF1                  0.00135 ** 
## Week:IGF_typeIGF2                 2.12e-05 ***
## Perc.WeightLoss:IGF_typeIGF1       0.06908 .  
## Perc.WeightLoss:IGF_typeIGF2       0.17443    
## Week:Perc.WeightLoss:IGF_typeIGF1  0.17609    
## Week:Perc.WeightLoss:IGF_typeIGF2  0.00014 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(model_rescue_with_tail)
## Type III Analysis of Variance Table with Satterthwaite's method
##                               Sum Sq Mean Sq NumDF  DenDF  F value    Pr(>F)
## Week                          5341.3  5341.3     1 415.87 274.7259 < 2.2e-16
## Perc.WeightLoss                  0.0     0.0     1 433.80   0.0003 0.9860797
## IGF_type                       108.4    54.2     2 356.60   2.7887 0.0628429
## Tail                           200.7   200.7     1  72.02  10.3213 0.0019675
## Week:Perc.WeightLoss             0.1     0.1     1 424.05   0.0042 0.9484235
## Week:IGF_type                  487.0   243.5     2 415.99  12.5230 5.233e-06
## Perc.WeightLoss:IGF_type        87.2    43.6     2 436.27   2.2438 0.1072814
## Week:Perc.WeightLoss:IGF_type  298.2   149.1     2 422.81   7.6681 0.0005355
##                                  
## Week                          ***
## Perc.WeightLoss                  
## IGF_type                      .  
## Tail                          ** 
## Week:Perc.WeightLoss             
## Week:IGF_type                 ***
## Perc.WeightLoss:IGF_type         
## Week:Perc.WeightLoss:IGF_type ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Reported result: Inclusion of original tail length improved model fit by \(\Delta AIC=6.9\), so Tail was retained as an additive covariate in the primary rescue model.

11. Publication Figures

11.1 Figure 1: dietary restriction and the energetic trade-off

## 11.1 Figure 1: dietary restriction and energetic constraint

library(dplyr)
library(ggplot2)
library(lme4)
library(lmerTest)
library(emmeans)
library(patchwork)

# ============================================================
# PREPARE DATA
# ============================================================

figure1_weight_data <- diet_final_acclimation %>%
  mutate(
    Diet_display = case_when(
      Treatment == "AdLib" ~ "High",
      Treatment %in% c("DR", "DR.IGF1", "DR.IGF2") ~ "Low",
      TRUE ~ NA_character_
    ),
    Diet_display = factor(
      Diet_display,
      levels = c("High", "Low")
    )
  )

figure1_regeneration_data <- data %>%
  mutate(
    Diet_display = case_when(
      Treatment == "AdLib" ~ "High",
      Treatment %in% c("DR", "DR.IGF1", "DR.IGF2") ~ "Low",
      TRUE ~ NA_character_
    ),
    Diet_display = factor(
      Diet_display,
      levels = c("High", "Low")
    )
  )


# ============================================================
# PANEL A
# DIETARY RESTRICTION INCREASED BODY-MASS LOSS
# ============================================================

panel_1a <- ggplot(
  figure1_weight_data,
  aes(x = Diet_display, y = Perc_WL)
) +
  geom_violin(
    fill = "grey75",
    color = "black",
    trim = FALSE,
    scale = "width",
    alpha = 0.8
  ) +
  geom_boxplot(
    width = 0.12,
    outlier.shape = NA,
    fill = "white",
    color = "black"
  ) +
  geom_jitter(
    width = 0.08,
    size = 1.8,
    shape = 21,
    fill = "white",
    color = "black",
    alpha = 0.8
  ) +
  labs(
    title = "Dietary restriction increased weight loss",
    x = "Diet",
    y = "Percent body-mass loss"
  ) +
  theme_classic(base_size = 12) +
  theme(
    plot.title = element_text(face = "bold"),
    axis.title = element_text(face = "bold"),
    axis.text = element_text(face = "bold")
  )


# ============================================================
# PANEL B
# DIETARY RESTRICTION REDUCED REGENERATION
# ============================================================

panel_1b <- ggplot(
  figure1_regeneration_data,
  aes(
    x = Week,
    y = Perc.Regeneration,
    shape = Diet_display,
    linetype = Diet_display,
    group = Diet_display
  )
) +
  stat_summary(
    fun = mean,
    geom = "line",
    linewidth = 1.2,
    color = "black"
  ) +
  stat_summary(
    fun = mean,
    geom = "point",
    size = 3,
    color = "black"
  ) +
  stat_summary(
    fun.data = mean_se,
    geom = "errorbar",
    width = 0.15,
    color = "black"
  ) +
  scale_shape_manual(
    values = c("High" = 16, "Low" = 17)
  ) +
  scale_linetype_manual(
    values = c("High" = "solid", "Low" = "dashed")
  ) +
  labs(
    title = "Dietary restriction reduced regeneration",
    x = "Week",
    y = "Percent regeneration",
    shape = "Diet",
    linetype = "Diet"
  ) +
  theme_classic(base_size = 12) +
  theme(
    plot.title = element_text(face = "bold"),
    axis.title = element_text(face = "bold"),
    axis.text = element_text(face = "bold"),
    legend.position = "bottom"
  )


# ============================================================
# DATA FOR PANELS C AND D
# NO IGF TREATMENT
# ============================================================

tradeoff_no_igf <- data %>%
  filter(
    IGF_type == "No_IGF",
    !is.na(Perc.WeightLoss),
    !is.na(Perc.Regeneration)
  )


# ============================================================
# PANEL C
# WEEK 6 RELATIONSHIP BETWEEN WEIGHT LOSS AND REGENERATION
# ============================================================

panel_1c_data <- tradeoff_no_igf %>%
  filter(Week == 6)

panel_1c <- ggplot(
  panel_1c_data,
  aes(
    x = Perc.WeightLoss,
    y = Perc.Regeneration
  )
) +
  geom_point(
    shape = 21,
    fill = "white",
    color = "black",
    size = 3,
    stroke = 0.8,
    alpha = 0.9
  ) +
  geom_smooth(
    method = "lm",
    formula = y ~ x,
    se = TRUE,
    color = "black",
    fill = "grey75",
    linewidth = 1.1
  ) +
  labs(
    title = "Weight loss constrains regeneration",
    subtitle = "Vehicle indviduals at Week 6",
    x = "Percent body-mass loss",
    y = "Percent regeneration"
  ) +
  theme_classic(base_size = 12) +
  theme(
    plot.title = element_text(
      face = "bold",
      size = 14
    ),
    plot.subtitle = element_text(
      face = "bold",
      size = 11
    ),
    axis.title = element_text(
      face = "bold"
    ),
    axis.text = element_text(
      face = "bold"
    )
  )


# ============================================================
# MODEL FOR PANEL D
# CHANGE IN WEIGHT-LOSS/REGENERATION RELATIONSHIP OVER TIME
# ============================================================

model_tradeoff_no_igf <- lmerTest::lmer(
  Perc.Regeneration ~
    Week * Perc.WeightLoss +
    (1 | Animal_ID),
  data = tradeoff_no_igf
)

tradeoff_slopes <- emtrends(
  model_tradeoff_no_igf,
  ~ Week,
  var = "Perc.WeightLoss",
  at = list(
    Week = c(2, 4, 6, 8)
  )
)

tradeoff_slopes_df <- as.data.frame(
  summary(
    tradeoff_slopes,
    infer = TRUE
  )
) %>%
  mutate(
    sig_zero = case_when(
      p.value < 0.001 ~ "***",
      p.value < 0.01  ~ "**",
      p.value < 0.05  ~ "*",
      p.value < 0.10  ~ "†",
      TRUE            ~ ""
    )
  )


# ============================================================
# PANEL D
# ENERGETIC CONSTRAINT STRENGTHENS OVER TIME
# ============================================================



# Estimate the slope of regeneration vs. weight loss
# at Weeks 2, 4, 6, and 8
tradeoff_slopes <- emtrends(
  model_tradeoff_no_igf,
  ~ Week,
  var = "Perc.WeightLoss",
  at = list(
    Week = c(2, 4, 6, 8)
  )
)

# Convert to dataframe and add significance symbols
tradeoff_slopes_df <- as.data.frame(
  summary(
    tradeoff_slopes,
    infer = TRUE
  )
) %>%
  mutate(
    sig_zero = case_when(
      p.value < 0.001 ~ "***",
      p.value < 0.01  ~ "**",
      p.value < 0.05  ~ "*",
      p.value < 0.10  ~ "†",
      TRUE            ~ ""
    ),

    # Position significance symbols below each confidence interval
    sig_y = lower.CL - 0.04
  )


# Make Panel D
panel_1d <- ggplot(
  tradeoff_slopes_df,
  aes(
    x = Week,
    y = Perc.WeightLoss.trend
  )
) +

  # Zero = no relationship between weight loss and regeneration
  geom_hline(
    yintercept = 0,
    linetype = "dashed",
    color = "black",
    linewidth = 0.7
  ) +

  # 95% confidence intervals
  geom_errorbar(
    aes(
      ymin = lower.CL,
      ymax = upper.CL
    ),
    width = 0.18,
    linewidth = 0.8,
    color = "grey30"
  ) +

  # Estimated slopes
  geom_point(
    size = 4.5,
    color = "grey30"
  ) +

  # Significance symbols
  geom_text(
    aes(
      y = sig_y,
      label = sig_zero
    ),
    size = 5,
    fontface = "bold",
    color = "grey30"
  ) +

  # X axis
  scale_x_continuous(
    breaks = c(2, 4, 6, 8)
  ) +

  # Y axis extended downward so Week 8 significance
  # symbols fit below the confidence interval
  scale_y_continuous(
    breaks = seq(-0.5, 0.1, by = 0.1),
    limits = c(-0.50, 0.15),
    labels = scales::number_format(
      accuracy = 0.1
    ),
    expand = expansion(
      mult = c(0.02, 0.05)
    )
  ) +

  labs(
    title = "Energetic constraint strengthened over time",
    subtitle = "Change in regeneration per 1% increase in weight loss",
    x = "Week",
    y = "Slope of Regeneration\nvs. Weight Loss"
  ) +

  theme_classic(
    base_size = 12
  ) +

  theme(
    plot.title = element_text(
      face = "bold",
      size = 14
    ),

    plot.subtitle = element_text(
      face = "bold",
      size = 11
    ),

    axis.title = element_text(
      face = "bold"
    ),

    axis.text = element_text(
      face = "bold"
    ),

    axis.title.y = element_text(
      margin = margin(r = 10)
    )
  )

panel_1d

# ============================================================
# COMBINE PANELS A-D
# ============================================================

figure_1 <-
  (panel_1a + panel_1b) /
  (panel_1c + panel_1d) +
  plot_layout(
    heights = c(1, 1)
  ) +
  plot_annotation(
    title = "Dietary Restriction Establishes an Energetic Constraint",
    tag_levels = "A"
  ) &
  theme(
    plot.tag = element_text(
      face = "bold",
      size = 20
    )
  )

figure_1

# ============================================================
# SAVE FIGURE
# ============================================================

ggsave(
  filename = "dietary_restriction_energetic_constraint_multiplot.png",
  plot = figure_1,
  width = 11,
  height = 8,
  dpi = 600
)

11.1.2 Popout box for Figure 1

# ============================================================
# PANEL C
# WEEK 6 MODEL-BASED POPOUT
# ============================================================

panel_1c_data <- tradeoff_no_igf %>%
  filter(Week == 6)

# Prediction grid across observed Week 6 weight-loss values
panel_1c_predictions <- data.frame(
  Week = 6,
  Perc.WeightLoss = seq(
    min(panel_1c_data$Perc.WeightLoss, na.rm = TRUE),
    max(panel_1c_data$Perc.WeightLoss, na.rm = TRUE),
    length.out = 100
  )
)

# Population-level predictions from SAME model used in Panel D
panel_1c_predictions$Predicted_Regeneration <- predict(
  model_tradeoff_no_igf,
  newdata = panel_1c_predictions,
  re.form = NA
)


panel_1c <- ggplot() +

  # Individual observed animals
  geom_point(
    data = panel_1c_data,
    aes(
      x = Perc.WeightLoss,
      y = Perc.Regeneration
    ),
    shape = 21,
    fill = "white",
    color = "black",
    size = 2.7,
    stroke = 0.7,
    alpha = 0.75
  ) +

  # Model-predicted Week 6 relationship
  geom_line(
    data = panel_1c_predictions,
    aes(
      x = Perc.WeightLoss,
      y = Predicted_Regeneration
    ),
    color = "black",
    linewidth = 1.3
  ) +

  labs(
    title = "Relationship at Week 6",
    subtitle = "Model prediction with individual observations",
    x = "Percent body-mass loss",
    y = "Percent regeneration"
  ) +

  theme_classic(base_size = 12) +

  theme(
    # Entire pop-out is grey90
    plot.background = element_rect(
      fill = "grey90",
      color = NA
    ),

    panel.background = element_rect(
      fill = "grey90",
      color = NA
    ),

    plot.title = element_text(
      face = "bold",
      size = 14
    ),

    plot.subtitle = element_text(
      face = "bold",
      size = 10
    ),

    axis.title = element_text(
      face = "bold"
    ),

    axis.text = element_text(
      face = "bold"
    ),

    plot.margin = margin(
      8, 8, 8, 8
    )
  )

panel_1c

ggsave(
  filename = "Fig1_panelC_week6_popout.png",
  plot = panel_1c,
  width = 3.8,
  height = 5,
  units = "in",
  dpi = 600,
  bg = "white"
)
# ============================================================
# PANEL D
# ENERGETIC CONSTRAINT STRENGTHENS OVER TIME
# ============================================================

# Estimate the slope of regeneration vs. weight loss
# at Weeks 2, 4, 6, and 8
tradeoff_slopes <- emtrends(
  model_tradeoff_no_igf,
  ~ Week,
  var = "Perc.WeightLoss",
  at = list(
    Week = c(2, 4, 6, 8)
  )
)

# Convert to dataframe and add significance symbols
tradeoff_slopes_df <- as.data.frame(
  summary(
    tradeoff_slopes,
    infer = TRUE
  )
) %>%
  mutate(
    sig_zero = case_when(
      p.value < 0.001 ~ "***",
      p.value < 0.01  ~ "**",
      p.value < 0.05  ~ "*",
      p.value < 0.10  ~ "†",
      TRUE            ~ ""
    ),

    # Position significance symbols below each confidence interval
    sig_y = lower.CL - 0.04
  )


# Make Panel D
panel_1d <- ggplot(
  tradeoff_slopes_df,
  aes(
    x = Week,
    y = Perc.WeightLoss.trend
  )
) +

  # ----------------------------------------------------------
  # Highlight Week 6
  # ----------------------------------------------------------
  annotate(
    "rect",
    xmin = 5.55,
    xmax = 6.45,
    ymin = -Inf,
    ymax = Inf,
    fill = "grey90",
    color = NA
  ) +

  # Zero = no relationship between weight loss and regeneration
  geom_hline(
    yintercept = 0,
    linetype = "dashed",
    color = "black",
    linewidth = 0.7
  ) +

  # 95% confidence intervals
  geom_errorbar(
    aes(
      ymin = lower.CL,
      ymax = upper.CL
    ),
    width = 0.18,
    linewidth = 0.8,
    color = "grey30"
  ) +

  # Estimated slopes
  geom_point(
    size = 4.5,
    color = "grey30"
  ) +

  # Significance symbols
  geom_text(
    aes(
      y = sig_y,
      label = sig_zero
    ),
    size = 5,
    fontface = "bold",
    color = "grey30"
  ) +

  # X axis
  scale_x_continuous(
    breaks = c(2, 4, 6, 8)
  ) +

  # Y axis
  scale_y_continuous(
    breaks = seq(-0.5, 0.1, by = 0.1),
    limits = c(-0.50, 0.15),
    labels = scales::number_format(
      accuracy = 0.1
    ),
    expand = expansion(
      mult = c(0.02, 0.05)
    )
  ) +

  labs(
    title = "Energetic constraint strengthened over time",
    subtitle = "Change in regeneration per 1% increase in weight loss",
    x = "Week",
    y = "Slope of Regeneration\nvs. Weight Loss"
  ) +

  theme_classic(
    base_size = 12
  ) +

  theme(
    plot.title = element_text(
      face = "bold",
      size = 14
    ),

    plot.subtitle = element_text(
      face = "bold",
      size = 11
    ),

    axis.title = element_text(
      face = "bold"
    ),

    axis.text = element_text(
      face = "bold"
    ),

    axis.title.y = element_text(
      margin = margin(r = 10)
    ),

    # Transparent plotting area
    panel.background = element_rect(
      fill = "transparent",
      color = NA
    ),

    # Transparent area surrounding plot
    plot.background = element_rect(
      fill = "transparent",
      color = NA
    )
  )

panel_1d_tagged <- panel_1d +
  labs(tag = "C") +
  theme(
    plot.tag = element_text(
      face = "bold",
      size = 18
    ),
    plot.tag.position = c(0, 1)
  )

panel_1d_tagged

ggsave(
  filename = "Fig1_panelD_constraint_slopes.png",
  plot = panel_1d_tagged,
  width = 6,
  height = 4,
  units = "in",
  dpi = 600,
  bg = "transparent"
)

11.2 Figure 2B: cumulative change in slope

cumulative_slope_change <- week_specific_slopes_df %>%
  filter(Week %in% c(2, 8)) %>%
  select(
    IGF_type,
    Week,
    Perc.WeightLoss.trend
  ) %>%
  pivot_wider(
    names_from = Week,
    values_from = Perc.WeightLoss.trend,
    names_prefix = "Week_"
  ) %>%
  mutate(
    change = Week_8 - Week_2
  )

plot_cumulative_effect <- ggplot(
  cumulative_slope_change,
  aes(
    x = IGF_type,
    y = change,
    fill = IGF_type
  )
) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  geom_col(width = 0.65) +
  scale_fill_manual(values = igf_colors) +
  scale_x_discrete(
    labels = c(
      "No_IGF" = "None",
      "IGF1" = "IGF1",
      "IGF2" = "IGF2"
    )
  ) +
  labs(
    title = "Cumulative Effects",
    x = "Supplementation",
    y = expression(Delta ~ slope ~ relative ~ to ~ control),
    fill = "IGF supplementation"
  ) +
  theme_classic() +
  theme(legend.position = "none")

plot_cumulative_effect

ggsave(
  filename = "IGF2_means.png",
  plot = plot_cumulative_effect,
  width = 2.5,
  height = 2,
  dpi = 600
)

11.3 Figure 2A: week-specific rescue slopes

# Prepare plotting data
rescue_plot_data <- week_specific_slopes_df %>%
  mutate(
    significance = case_when(
      IGF_type == "IGF2" & Week == 4 ~ "\u2020",
      IGF_type == "IGF2" & Week == 6 ~ "**",
      IGF_type == "IGF2" & Week == 8 ~ "***",
      TRUE ~ ""
    )
  )

plot_rescue <- ggplot(
  rescue_plot_data,
  aes(
    x = Week,
    y = Perc.WeightLoss.trend,
    color = IGF_type
  )
) +

  # Background shading
annotate(
  "rect",
  xmin = 5.55,
  xmax = 6.45,
  ymin = -Inf,
  ymax = Inf,
  fill = "grey90",
  alpha = 0.7
) +

  # Zero line
  geom_hline(
    yintercept = 0,
    linetype = "dashed",
    color = "black"
  ) +

  # Confidence intervals
  geom_errorbar(
    aes(ymin = lower.CL, ymax = upper.CL),
    width = 0.15,
    position = position_dodge(width = 0.25),
    linewidth = 0.8
  ) +

  # Point estimates
  geom_point(
    size = 4,
    position = position_dodge(width = 0.25)
  ) +

  # Significance symbols
  geom_text(
    aes(
      y = upper.CL + 0.05,
      label = significance
    ),
    position = position_dodge(width = 0.25),
    size = 8,
    fontface = "bold",
    show.legend = FALSE
  ) +

  scale_color_manual(
    values = c(
      "IGF1" = "#E87722",
      "IGF2" = "#0C2340",
      "No_IGF" = "grey40"
    ),
    labels = c(
      "IGF1" = "IGF1",
      "IGF2" = "IGF2",
      "No_IGF" = "None"
    )
  ) +

labs(
  title = "IGF2 alters the energetic constraint on regeneration",
  subtitle = "Relationship between weight loss and regeneration across time",
  x = "Week",
  y = "Slope of Regeneration vs. Weight Loss"
) +
  labs(
  color = "IGF Type"
)+

  theme_classic(base_size = 12) +
  theme(
    axis.title = element_text(face = "bold"),
    axis.text = element_text(face = "bold"),
    plot.title = element_text(face = "bold"),
    plot.subtitle = element_text(face = "bold", size = 10),
    legend.position = "bottom",
    legend.title = element_text(face = "bold")
  ) +
  theme(
  plot.background = element_rect(
    fill = "transparent",
    color = NA
  ),
  panel.background = element_rect(
    fill = "transparent",
    color = NA
  )
)

plot_rescue

ggsave(
  filename = "IGF2_energetic_rescue_figure.png",
  plot = plot_rescue,
  width = 7,
  height = 4,
  dpi = 600, 
  bg="transparent"
)
# ============================================================
# WEEK 6 POPOUT
# Raw data + regression lines by supplementation
# ============================================================

week6_popout_data <- data %>%
  filter(
    Week == 6,
    IGF_type %in% c("No_IGF", "IGF1", "IGF2"),
    !is.na(Perc.WeightLoss),
    !is.na(Perc.Regeneration)
  ) %>%
  mutate(
    Supplementation = case_when(
      IGF_type == "No_IGF" ~ "None",
      IGF_type == "IGF1"   ~ "IGF1",
      IGF_type == "IGF2"   ~ "IGF2"
    ),
    Supplementation = factor(
      Supplementation,
      levels = c("None", "IGF1", "IGF2")
    )
  )


week6_popout <- ggplot(
  week6_popout_data,
  aes(
    x = Perc.WeightLoss,
    y = Perc.Regeneration,
    color = Supplementation
  )
) +

  # Individual animals
  geom_point(
    size = 2.5,
    alpha = 0.65
  ) +

  # Raw regression lines
  geom_smooth(
    method = "lm",
    formula = y ~ x,
    se = FALSE,
    linewidth = 1.2
  ) +

  # Match colors in main figure
  scale_color_manual(
    values = c(
      "None" = "grey40",
      "IGF1" = "#E66101",
      "IGF2" = "#0B1F3A"
    )
  ) +

  labs(
    title = "Week 6 Pop-Out:",
    x = "Percent body-mass loss",
    y = "Percent regeneration"
  ) +

  theme_classic(base_size = 11) +

  theme(
    # Light gray background across entire pop-out
    plot.background = element_rect(
      fill = "grey90",
      color = "grey90",
      linewidth = 0.8
    ),

    # Keep plotting area gray as well
    panel.background = element_rect(
      fill = "grey90",
      color = NA
    ),

    plot.title = element_text(
      face = "bold",
      size = 13,
      hjust = 0.5,
      margin = margin(b = 2)
    ),

    plot.subtitle = element_text(
      face = "plain",
      size = 10.5,
      hjust = 0.5,
      margin = margin(b = 7)
    ),

    axis.title = element_text(
      face = "bold",
      size = 10
    ),

    axis.text = element_text(
      size = 9
    ),

    legend.position = "none",

    # Gives the gray box some breathing room
    plot.margin = margin(
      t = 10,
      r = 10,
      b = 10,
      l = 10
    )
  )

week6_popout

ggsave(
  filename = "week6_popout.png",
  plot = week6_popout,
  width = 3,
  height = 4,
  units = "in",
  dpi = 600,
  bg = "grey90"
)
model_igf_split <- lmer(
  Perc.Regeneration ~ Week * Perc.WeightLoss * IGF_type +
    (1 | Animal_ID),
  data = data
)

# ============================================================
# WEEK 6 MODEL-BASED POPOUT
# Individual observations + predictions from Panel B model
# ============================================================

week6_points <- data %>%
  filter(
    Week == 6,
    IGF_type %in% c("No_IGF", "IGF1", "IGF2"),
    !is.na(Perc.WeightLoss),
    !is.na(Perc.Regeneration)
  ) %>%
  mutate(
    Supplementation = case_when(
      IGF_type == "No_IGF" ~ "None",
      IGF_type == "IGF1"   ~ "IGF1",
      IGF_type == "IGF2"   ~ "IGF2"
    ),
    Supplementation = factor(
      Supplementation,
      levels = c("None", "IGF1", "IGF2")
    )
  )


# Create weight-loss values spanning the observed range
prediction_grid <- expand.grid(
  Week = 6,
  Perc.WeightLoss = seq(
    min(week6_points$Perc.WeightLoss, na.rm = TRUE),
    max(week6_points$Perc.WeightLoss, na.rm = TRUE),
    length.out = 100
  ),
  IGF_type = c("No_IGF", "IGF1", "IGF2")
)

# Make population-level predictions from the SAME model as Panel B
prediction_grid$Predicted_Regeneration <- predict(
  model_igf_split,
  newdata = prediction_grid,
  re.form = NA
)

prediction_grid <- prediction_grid %>%
  mutate(
    Supplementation = case_when(
      IGF_type == "No_IGF" ~ "None",
      IGF_type == "IGF1"   ~ "IGF1",
      IGF_type == "IGF2"   ~ "IGF2"
    ),
    Supplementation = factor(
      Supplementation,
      levels = c("None", "IGF1", "IGF2")
    )
  )

week6_popout <- ggplot() +

  # Raw individual observations
  geom_point(
    data = week6_points,
    aes(
      x = Perc.WeightLoss,
      y = Perc.Regeneration,
      color = Supplementation
    ),
    size = 2.3,
    alpha = 0.55
  ) +

  # MODEL-PREDICTED relationships
  geom_line(
    data = prediction_grid,
    aes(
      x = Perc.WeightLoss,
      y = Predicted_Regeneration,
      color = Supplementation
    ),
    linewidth = 1.3
  ) +

  scale_color_manual(
    values = c(
      "None" = "grey40",
      "IGF1" = "#E66101",
      "IGF2" = "#0B1F3A"
    )
  ) +

  labs(
    title = "Week 6 Pop-out",
    subtitle = "Model-Predicted Relationships",
    x = "Percent Body-Mass Loss",
    y = "Percent Regeneration"
  ) +

  theme_classic(base_size = 11) +

  theme(
    plot.background = element_rect(
      fill = "grey92",
      color = NA
    ),
    panel.background = element_rect(
      fill = "grey92",
      color = NA
    ),

    plot.title = element_text(
      face = "bold",
      size = 11,
      hjust = 0.5
    ),

    axis.title = element_text(
      face = "bold",
      size = 9
    ),

    axis.text = element_text(
      size = 8
    ),

    legend.position = "none",

    plot.margin = margin(8, 8, 8, 8)
  )

week6_popout

ggsave(
  filename = "week6_popout.png",
  plot = week6_popout,
  width = 3,
  height = 4,
  units = "in",
  dpi = 600,
  bg = "grey90"
)
# ============================================================
# WEIGHT LOSS THROUGHOUT REGENERATION
# Individual trajectories + treatment means
# ============================================================

# ============================================================
# PRE-AUTOTOMY PERIOD
# Include ALL negative weeks from diet dataset
# ============================================================

pre_autotomy <- diet %>%
  filter(
    Week < 0,
    !is.na(Perc_WL)
  ) %>%
  dplyr::select(
    Animal_ID,
    Week,
    Perc.WeightLoss = Perc_WL
  ) %>%
  
  # Add IGF assignment from regeneration dataset
  left_join(
    data %>%
      dplyr::select(Animal_ID, IGF_type) %>%
      distinct(),
    by = "Animal_ID"
  )


# ============================================================
# REGENERATION PERIOD
# ============================================================

regen_weight <- data %>%
  filter(
    Week >= 0,
    IGF_type %in% c("No_IGF", "IGF1", "IGF2"),
    !is.na(Perc.WeightLoss)
  ) %>%
  dplyr::select(
    Animal_ID,
    Week,
    Perc.WeightLoss,
    IGF_type
  )


# ============================================================
# COMBINE PRE-AUTOTOMY + REGENERATION
# ============================================================

weight_trajectory_data <- bind_rows(
  pre_autotomy,
  regen_weight
) %>%
  filter(
    IGF_type %in% c("No_IGF", "IGF1", "IGF2")
  ) %>%
  mutate(
    Supplementation = case_when(
      IGF_type == "No_IGF" ~ "None",
      IGF_type == "IGF1"   ~ "IGF1",
      IGF_type == "IGF2"   ~ "IGF2"
    ),
    Supplementation = factor(
      Supplementation,
      levels = c("None", "IGF1", "IGF2")
    )
  )


body_mass_plot <- ggplot(
  weight_trajectory_data,
  aes(
    x = Week,
    y = Perc.WeightLoss,
    color = Supplementation
  )
) +

  # Individual animal trajectories
  geom_line(
    aes(group = Animal_ID),
    alpha = 0.12,
    linewidth = 0.5
  ) +

  # Treatment means
  stat_summary(
    aes(group = Supplementation),
    fun = mean,
    geom = "line",
    linewidth = 1.2
  ) +

  # Mean points
  stat_summary(
    aes(group = Supplementation),
    fun = mean,
    geom = "point",
    size = 3
  ) +

  # 0 = original body mass
  geom_hline(
    yintercept = 0,
    linetype = "dashed",
    linewidth = 0.7,
    color = "grey40"
  ) +

  scale_color_manual(
    values = c(
      "None" = "grey40",
      "IGF1" = "#E66101",
      "IGF2" = "#0B1F3A"
    )
  ) +

scale_x_continuous(
  breaks = c(-4, -2, -1, 2, 4, 6, 8),
  labels = c("-4", "-2", "-1", "2", "4", "6", "8")
) +
  coord_cartesian(
  ylim = c(-30, 50)
)+
  geom_vline(
  xintercept = 0,
  linetype = "dotted",
  linewidth = 0.7,
  color = "grey40"
) +

  labs(
    title = "Weight loss remains stable",
    subtitle = "Individual trajectories and treatment means",
    x = "Week",
    y = "Percent Weight Loss",
    color = "Supplementation"
  ) +

 annotate(
  "label",
  x = 1.5,
  y = 47,
  label = "Autotomy",
  fontface = "bold",
  size = 5,
  color = "#0B1F3A",       # text color
  fill = "white",          # box fill
  label.size = 0.7,        # border thickness
  label.padding = unit(0.25, "lines"),
  hjust = 1
) +
  theme_classic(base_size = 12) +

  theme(
    plot.title = element_text(
      face = "bold",
      size = 14
    ),

    plot.subtitle = element_text(
      size = 10.5,
      margin = margin(b = 7)
    ),

    axis.title = element_text(
      face = "bold",
      size = 11
    ),

    axis.text = element_text(
      size = 10
    ),

    legend.position = "bottom",
    legend.title = element_text(face = "bold"),

    plot.margin = margin(
      t = 8,
      r = 8,
      b = 8,
      l = 8
    )
  )

body_mass_plot

ggsave(
  filename = "Fig2_stableWL.png",
  plot = body_mass_plot,
  width = 6,
  height = 4.5,
  units = "in",
  dpi = 600,
  bg = "white"
)

#12. Initial Investment Modeling

# Extract energetic condition immediately before tail loss
initial_condition <- diet %>%
  filter(Week == -1) %>%
  dplyr::select(
    Animal_ID,
    Initial_WeightLoss = Perc_WL
  ) %>%
  distinct()

summary(initial_condition$Initial_WeightLoss)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
## -17.028   1.857   6.059   6.101  11.211  27.965       7
nrow(initial_condition)
## [1] 76
regen_initial <- data %>%
  filter(Week %in% c(2, 4, 6, 8)) %>%
  left_join(initial_condition, by = "Animal_ID")

# Check the join
sum(is.na(regen_initial$Initial_WeightLoss))
## [1] 28
length(unique(regen_initial$Animal_ID))
## [1] 76
model_initial_condition <- lmer(
  Perc.Regeneration ~ Week * Initial_WeightLoss  +
    (1 | Animal_ID),
  data = regen_initial
)

summary(model_initial_condition)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Perc.Regeneration ~ Week * Initial_WeightLoss + (1 | Animal_ID)
##    Data: regen_initial
## 
## REML criterion at convergence: 1404.4
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.85976 -0.54156  0.06924  0.60861  2.75144 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 20.68    4.548   
##  Residual              28.17    5.308   
## Number of obs: 214, groups:  Animal_ID, 64
## 
## Fixed effects:
##                          Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)              -6.39141    1.32520 190.40917  -4.823 2.89e-06 ***
## Week                      4.68879    0.21895 161.85585  21.415  < 2e-16 ***
## Initial_WeightLoss        0.01689    0.15246 190.36370   0.111    0.912    
## Week:Initial_WeightLoss  -0.03805    0.02530 162.91730  -1.504    0.134    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Week   Int_WL
## Week        -0.750              
## Intl_WghtLs -0.627  0.465       
## Wk:Intl_WgL  0.463 -0.612 -0.748
anova(model_initial_condition)
## Type III Analysis of Variance Table with Satterthwaite's method
##                          Sum Sq Mean Sq NumDF  DenDF  F value Pr(>F)    
## Week                    12919.6 12919.6     1 161.86 458.5876 <2e-16 ***
## Initial_WeightLoss          0.3     0.3     1 190.36   0.0123 0.9119    
## Week:Initial_WeightLoss    63.7    63.7     1 162.92   2.2623 0.1345    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
initial_slopes <- emtrends(
  model_initial_condition,
  ~ Week,
  var = "Initial_WeightLoss",
  at = list(Week = c(2, 4, 6, 8))
)

summary(initial_slopes, infer = TRUE)
##  Week Initial_WeightLoss.trend    SE    df lower.CL upper.CL t.ratio p.value
##     2                  -0.0592 0.119 112.4   -0.296  0.17738  -0.496  0.6210
##     4                  -0.1353 0.102  62.9   -0.339  0.06846  -1.327  0.1893
##     6                  -0.2114 0.108  74.9   -0.427  0.00377  -1.957  0.0540
##     8                  -0.2875 0.134 141.7   -0.553 -0.02176  -2.139  0.0342
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95

The relationship between energetic condition and regeneration emerged over time and differed among treatment groups. Although animals largely maintained body mass during the regeneration period, vehicle-treated individuals with greater body-mass loss exhibited progressively lower cumulative regeneration, with this relationship becoming statistically detectable by week 6 and strengthening through week 8. In contrast, this relationship was altered by IGF2 treatment. Body-mass loss measured immediately before tail loss did not independently predict subsequent regeneration after accounting for original tail length, indicating that the later relationship between energetic condition and regeneration was not simply explained by differences in body condition at the onset of regeneration.

13. After IGF injection/autotomy, did IGF2 animals change in mass differently over time than IGF1 or Low Vehicle animals?

# ============================================================
# DOES IGF2 ALTER MASS FOLLOWING AUTOTOMY / INJECTION?
# Low-diet animals only
# ============================================================

library(dplyr)
library(lmerTest)
library(emmeans)

# Get each animal's mass immediately before autotomy/injection
pre_autotomy_mass <- diet %>%
  filter(Week == -1) %>%
  dplyr::select(
    Animal_ID,
    PreAutotomy_Mass = Mass
  ) %>%
  distinct()


# Post-autotomy observations for the three low-diet groups
post_igf_mass <- data %>%
  filter(
    Treatment %in% c("DR", "DR.IGF1", "DR.IGF2"),
    Week %in% c(2, 4, 6, 8),
    !is.na(Mass)
  ) %>%
  left_join(
    pre_autotomy_mass,
    by = "Animal_ID"
  ) %>%
  mutate(
    # Absolute change from immediately before autotomy
    Mass_Change = Mass - PreAutotomy_Mass,

    # Percent change from immediately before autotomy
    Percent_Mass_Change =
      ((Mass - PreAutotomy_Mass) / PreAutotomy_Mass) * 100,

    Supplementation = case_when(
      Treatment == "DR"      ~ "Vehicle",
      Treatment == "DR.IGF1" ~ "IGF1",
      Treatment == "DR.IGF2" ~ "IGF2"
    ),

    Supplementation = factor(
      Supplementation,
      levels = c("Vehicle", "IGF1", "IGF2")
    )
  )

post_igf_mass %>%
  group_by(Supplementation, Week) %>%
  summarise(
    n = n(),
    mean_change = mean(Percent_Mass_Change, na.rm = TRUE),
    sd_change = sd(Percent_Mass_Change, na.rm = TRUE),
    .groups = "drop"
  )
## # A tibble: 12 × 5
##    Supplementation  Week     n mean_change sd_change
##    <fct>           <int> <int>       <dbl>     <dbl>
##  1 Vehicle             2    13      -10.3       4.90
##  2 Vehicle             4    13      -14.9       6.63
##  3 Vehicle             6    10      -16.9       8.90
##  4 Vehicle             8     7      -10.3       8.29
##  5 IGF1                2    16       -9.00      7.40
##  6 IGF1                4    14      -10.9       5.76
##  7 IGF1                6    12       -9.17      5.85
##  8 IGF1                8    11       -6.04      8.63
##  9 IGF2                2    17      -11.9       7.04
## 10 IGF2                4    17      -15.9       6.82
## 11 IGF2                6    16      -16.0       9.22
## 12 IGF2                8    11       -7.86      6.78
model_post_igf_mass <- lmer(
  Percent_Mass_Change ~
    Week * Supplementation +
    (1 | Animal_ID),
  data = post_igf_mass
)

summary(model_post_igf_mass)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Percent_Mass_Change ~ Week * Supplementation + (1 | Animal_ID)
##    Data: post_igf_mass
## 
## REML criterion at convergence: 1025.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.4911 -0.4713 -0.0123  0.5277  3.8203 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Animal_ID (Intercept) 29.63    5.444   
##  Residual              28.39    5.328   
## Number of obs: 157, groups:  Animal_ID, 46
## 
## Fixed effects:
##                          Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)              -11.1939     2.4593 123.4563  -4.552 1.26e-05 ***
## Week                      -0.4673     0.4015 115.4443  -1.164    0.247    
## SupplementationIGF1        0.2233     3.2914 122.5491   0.068    0.946    
## SupplementationIGF2       -3.5735     3.2464 122.4341  -1.101    0.273    
## Week:SupplementationIGF1   0.8454     0.5271 115.8798   1.604    0.112    
## Week:SupplementationIGF2   0.7516     0.5175 113.9692   1.452    0.149    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Week   SpIGF1 SpIGF2 W:SIGF1
## Week        -0.712                             
## SpplmntIGF1 -0.747  0.532                      
## SpplmntIGF2 -0.758  0.540  0.566               
## Wk:SpplIGF1  0.543 -0.762 -0.705 -0.411        
## Wk:SpplIGF2  0.553 -0.776 -0.413 -0.712  0.591
anova(
  model_post_igf_mass,
  type = 3
)
## Type III Analysis of Variance Table with Satterthwaite's method
##                      Sum Sq Mean Sq NumDF  DenDF F value Pr(>F)
## Week                  2.808   2.808     1 114.73  0.0989 0.7537
## Supplementation      54.332  27.166     2 122.03  0.9570 0.3869
## Week:Supplementation 84.117  42.058     2 114.64  1.4816 0.2316
mass_slopes <- emtrends(
  model_post_igf_mass,
  ~ Supplementation,
  var = "Week"
)

mass_slopes
##  Supplementation Week.trend    SE  df lower.CL upper.CL
##  Vehicle             -0.467 0.402 115   -1.264    0.330
##  IGF1                 0.378 0.342 116   -0.300    1.056
##  IGF2                 0.284 0.327 112   -0.363    0.932
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95
pairs(
  mass_slopes,
  adjust = "tukey"
)
##  contrast       estimate    SE  df t.ratio p.value
##  Vehicle - IGF1  -0.8454 0.528 116  -1.601  0.2495
##  Vehicle - IGF2  -0.7516 0.518 114  -1.450  0.3189
##  IGF1 - IGF2      0.0937 0.473 114   0.198  0.9786
## 
## Degrees-of-freedom method: kenward-roger 
## P value adjustment: tukey method for comparing a family of 3 estimates
mass_slopes <- emtrends(
  model_post_igf_mass,
  ~ Supplementation,
  var = "Week"
)

mass_slopes
##  Supplementation Week.trend    SE  df lower.CL upper.CL
##  Vehicle             -0.467 0.402 115   -1.264    0.330
##  IGF1                 0.378 0.342 116   -0.300    1.056
##  IGF2                 0.284 0.327 112   -0.363    0.932
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95
pairs(mass_slopes, adjust = "tukey")
##  contrast       estimate    SE  df t.ratio p.value
##  Vehicle - IGF1  -0.8454 0.528 116  -1.601  0.2495
##  Vehicle - IGF2  -0.7516 0.518 114  -1.450  0.3189
##  IGF1 - IGF2      0.0937 0.473 114   0.198  0.9786
## 
## Degrees-of-freedom method: kenward-roger 
## P value adjustment: tukey method for comparing a family of 3 estimates
mass_slopes
##  Supplementation Week.trend    SE  df lower.CL upper.CL
##  Vehicle             -0.467 0.402 115   -1.264    0.330
##  IGF1                 0.378 0.342 116   -0.300    1.056
##  IGF2                 0.284 0.327 112   -0.363    0.932
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95
test(mass_slopes)
##  Supplementation Week.trend    SE  df t.ratio p.value
##  Vehicle             -0.467 0.402 115  -1.162  0.2478
##  IGF1                 0.378 0.342 116   1.105  0.2716
##  IGF2                 0.284 0.327 112   0.870  0.3862
## 
## Degrees-of-freedom method: kenward-roger

14. Session Information

Including session information documents the versions of R and all attached packages used to render the analysis.

sessionInfo()
## R version 4.5.2 (2025-10-31)
## Platform: aarch64-apple-darwin20
## Running under: macOS Sequoia 15.7.4
## 
## Matrix products: default
## BLAS:   /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.1
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## time zone: America/New_York
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] lmerTest_3.2-1      hrbrthemes_0.9.3    viridis_0.6.5      
##  [4] viridisLite_0.4.3   plotly_4.12.0       patchwork_1.3.2    
##  [7] ggplot2_4.0.2       Rmisc_1.5.1         plyr_1.8.9         
## [10] lattice_0.22-7      ggeffects_2.3.2     performance_0.16.0 
## [13] emmeans_2.0.2       multcomp_1.4-30     TH.data_1.1-5      
## [16] MASS_7.3-65         survival_3.8-3      mvtnorm_1.3-5      
## [19] bestNormalize_1.9.2 MuMIn_1.48.19       lmtest_0.9-40      
## [22] zoo_1.8-15          car_3.1-5           carData_3.0-6      
## [25] nlme_3.1-168        lme4_2.0-6          Matrix_1.7-4       
## [28] tidyr_1.3.2         dplyr_1.2.1        
## 
## loaded via a namespace (and not attached):
##   [1] RColorBrewer_1.1-3  rstudioapi_0.18.0   jsonlite_2.0.0     
##   [4] datawizard_1.3.1    magrittr_2.0.4      estimability_1.5.1 
##   [7] farver_2.1.2        nloptr_2.2.1        rmarkdown_2.30     
##  [10] ragg_1.5.2          vctrs_0.7.1         minqa_1.2.8        
##  [13] base64enc_0.1-6     butcher_0.4.0       htmltools_0.5.9    
##  [16] forcats_1.0.1       progress_1.2.3      haven_2.5.5        
##  [19] broom_1.0.12        Formula_1.2-5       sass_0.4.10        
##  [22] parallelly_1.47.0   bslib_0.10.0        htmlwidgets_1.6.4  
##  [25] pbkrtest_0.5.5      sandwich_3.1-1      lubridate_1.9.5    
##  [28] cachem_1.1.0        lifecycle_1.0.5     iterators_1.0.14   
##  [31] pkgconfig_2.0.3     sjlabelled_1.2.0    R6_2.6.1           
##  [34] fastmap_1.2.0       rbibutils_2.4.1     future_1.70.0      
##  [37] digest_0.6.39       numDeriv_2016.8-1.1 colorspace_2.1-2   
##  [40] textshaping_1.0.5   Hmisc_5.2-5         labeling_0.4.3     
##  [43] timechange_0.4.0    mgcv_1.9-3          httr_1.4.8         
##  [46] abind_1.4-8         compiler_4.5.2      rngtools_1.5.2     
##  [49] withr_3.0.2         doParallel_1.0.17   htmlTable_2.5.0    
##  [52] S7_0.2.1            backports_1.5.0     lava_1.9.0         
##  [55] tools_4.5.2         foreign_0.8-90      otel_0.2.0         
##  [58] future.apply_1.20.2 nnet_7.3-20         glue_1.8.1         
##  [61] grid_4.5.2          checkmate_2.3.4     cluster_2.1.8.2    
##  [64] generics_0.1.4      recipes_1.3.2       gtable_0.3.6       
##  [67] class_7.3-23        data.table_1.18.2.1 hms_1.1.4          
##  [70] utf8_1.2.6          stringr_1.6.0       foreach_1.5.2      
##  [73] pillar_1.11.1       splines_4.5.2       tidyselect_1.2.1   
##  [76] knitr_1.51          reformulas_0.4.4    gridExtra_2.3      
##  [79] stats4_4.5.2        xfun_0.61           hardhat_1.4.3      
##  [82] timeDate_4052.112   stringi_1.8.7       lazyeval_0.2.2     
##  [85] yaml_2.3.12         boot_1.3-32         evaluate_1.0.5     
##  [88] codetools_0.2-20    tibble_3.3.1        cli_3.6.5          
##  [91] rpart_4.1.24        xtable_1.8-8        systemfonts_1.3.2  
##  [94] Rdpack_2.6.6        jquerylib_0.1.4     Rcpp_1.1.2         
##  [97] globals_0.19.1      coda_0.19-4.1       parallel_4.5.2     
## [100] gower_1.0.2         prettyunits_1.2.0   doRNG_1.8.6.3      
## [103] listenv_0.10.1      ipred_0.9-15        scales_1.4.0       
## [106] prodlim_2026.03.11  insight_1.5.1       purrr_1.2.2        
## [109] crayon_1.5.3        rlang_1.3.0