Fit (Mantle)

fit_model <- function(d) {
  fit <- lm(OPS ~ I(Age - 30) + I((Age - 30)^2), data = d)
  b <- coef(fit)
  Age_max <- 30 - b[2] / b[3] / 2
  Max <- b[1] - b[2] ^ 2 / b[3] / 4
  list(fit = fit, Age_max = Age_max, Max = Max)
}
F2 <- fit_model(Mantle)
F2 |>
  pluck("fit") |>
  coef()

##     (Intercept)     I(Age - 30) I((Age - 30)^2) 
##     1.043134189    -0.022883024    -0.003868915

• What we really want though:

c(F2$Age_max, F2$Max)

## I(Age - 30) (Intercept) 
##    27.04271     1.07697