I have this ggplot2:

    ggplot(DATABASE) + 
  geom_point(aes(x=GLUCOSE, y=BLOODPRESSURE), shape=1, size=1, alpha=0.8) + 
   theme(legend.position='none') +
    geom_smooth(aes(x=GLUCOSE, y=BLOODPRESSURE), method=lm, se=TRUE, color="black")

How can I add a p-value for three different regression models, where the first one is unadjusted, the second one adjusts for the variables "AGE" and "SEX", and the third one adjusts for the variables "AGE", "SEX", "CHOLESTEROL", and "BMI"?

In the actual case I have, the regression models have many more variables, and I will need to re-use them in other graphs and tables as well, so the best would be if there is a way to first define the regression models (i.e. which variables are adjusted for in each model). Is there a way I can do this?

Best regards

Edit: someone (@IRTFM?) closed the question with a link to an answer that didn't answer my original question. However, I would like to contribute with the actual answer, in case others would be interested:

### Linear regression for models 1-3
figure_2_model_1 = lm(GLUCOSE ~ BLOODPRESSURE, data=DATABASE)
summary(figure_2_model_1)
figure_2_model_1_p_values <- as.numeric(summary(figure_2_model_1)$coefficients[2,4])
figure_2_model_1_p_values <- round(figure_2_model_1_p_values, 3)
figure_2_model_1_p_values <- ifelse(figure_2_model_1_p_values < 0.001, "Model 1: P<0.001", paste0("Model 1: P=", figure_2_model_1_p_values))

figure_2_model_2 = lm(GLUCOSE ~ BLOODPRESSURE + AGE + SEX, data=DATABASE)
summary(figure_2_model_2)
figure_2_model_2_p_values <- as.numeric(summary(figure_2_model_2)$coefficients[2,4])
figure_2_model_2_p_values <- round(figure_2_model_2_p_values, 3)
figure_2_model_2_p_values <- ifelse(figure_2_model_2_p_values < 0.001, "Model 2: P<0.001", paste0("Model 2: P=", figure_2_model_2_p_values))

figure_2_model_3 = lm(GLUCOSE ~ BLOODPRESSURE + AGE + SEX + CHOLESTEROL + BMI, data=DATABASE)
summary(figure_2_model_3)
figure_2_model_3_p_values <- as.numeric(summary(figure_2_model_3)$coefficients[2,4])
figure_2_model_3_p_values <- round(figure_2_model_3_p_values, 3)
figure_2_model_3_p_values <- ifelse(figure_2_model_3_p_values < 0.001, "Model 3: P<0.001", paste0("Model 3: P=", figure_2_model_3_p_values))


figure_2_p_values <- c(figure_2_model_1_p_values, figure_2_model_2_p_values, figure_2_model_3_p_values)

figure_2_model_names <- c("Model 1: P", "Model 2: P", "Model 3: P")
var.labels = c(figure_2_p_values="P-values for figure 2")
figure_2_p_values_print <- as.data.frame(figure_2_p_values)
figure_2_p_values_print


### Additional code to ggplot:
  annotate("text", x = 120, y = c(85,80,75,70), label = c(figure_2_model_1_p_values, figure_2_model_2_p_values, figure_2_model_3_p_values, figure_2_model_4_p_values), size = 4)